1/* $Id$ 2 3 Part of CLP(R) (Constraint Logic Programming over Reals) 4 5 Author: Leslie De Koninck 6 E-mail: Leslie.DeKoninck@cs.kuleuven.be 7 WWW: http://www.swi-prolog.org 8 http://www.ai.univie.ac.at/cgi-bin/tr-online?number+95-09 9 Copyright (C): 2004, K.U. Leuven and 10 1992-1995, Austrian Research Institute for 11 Artificial Intelligence (OFAI), 12 Vienna, Austria 13 14 This software is part of Leslie De Koninck's master thesis, supervised 15 by Bart Demoen and daily advisor Tom Schrijvers. It is based on CLP(Q,R) 16 by Christian Holzbaur for SICStus Prolog and distributed under the 17 license details below with permission from all mentioned authors. 18 19 This program is free software; you can redistribute it and/or 20 modify it under the terms of the GNU General Public License 21 as published by the Free Software Foundation; either version 2 22 of the License, or (at your option) any later version. 23 24 This program is distributed in the hope that it will be useful, 25 but WITHOUT ANY WARRANTY; without even the implied warranty of 26 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 27 GNU General Public License for more details. 28 29 You should have received a copy of the GNU Lesser General Public 30 License along with this library; if not, write to the Free Software 31 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 32 33 As a special exception, if you link this library with other files, 34 compiled with a Free Software compiler, to produce an executable, this 35 library does not by itself cause the resulting executable to be covered 36 by the GNU General Public License. This exception does not however 37 invalidate any other reasons why the executable file might be covered by 38 the GNU General Public License. 39*/ 40 41 42:- module(clpcd_fourmotz, 43 [ 44 fm_elim/4 45 ]). 46 47:- use_module(library(clpcd/class)). 48:- use_module(library(clpcd/project)). 49:- use_module(library(clpcd/indep)). 50:- use_module(library(clpcd/redund)). 51:- use_module(library(clpcd/store)). 52:- use_module(library(clpcd/solve)). 53:- use_module(library(clpcd/bv)). 54:- use_module(library(clpcd/domain_ops)). 55 56fm_elim(CLP,Vs,Target,Pivots) :- 57 prefilter(Vs,Vsf), 58 fm_elim_int(Vsf, CLP, Target, Pivots). 59 60% prefilter(Vars,Res) 61% 62% filters out target variables and variables that do not occur in bounded linear equations. 63% Stores that the variables in Res are to be kept independent. 64 65prefilter([],[]). 66prefilter([V|Vs],Res) :- 67 ( get_attr(V,clpcd_itf,Att), 68 arg(9,Att,n), 69 occurs(V) 70 -> % V is a nontarget variable that occurs in a bounded linear equation 71 Res = [V|Tail], 72 setarg(10,Att,keep_indep), 73 prefilter(Vs,Tail) 74 ; prefilter(Vs,Res) 75 ). 76 77% 78% the target variables are marked with an attribute, and we get a list 79% of them as an argument too 80% 81fm_elim_int([], _, _, Pivots) :- % done 82 unkeep(Pivots). 83fm_elim_int(Vs, CLP, Target, Pivots) :- 84 Vs = [_|_], 85 ( best(CLP, Vs, Best, Rest) 86 -> occurences(Best,Occ), 87 elim_min(CLP, Best, Occ, Target, Pivots, NewPivots) 88 ; % give up 89 NewPivots = Pivots, 90 Rest = [] 91 ), 92 fm_elim_int(Rest, CLP, Target, NewPivots). 93 94% best(Vs,Best,Rest) 95% 96% Finds the variable with the best result (lowest Delta) in fm_cp_filter 97% and returns the other variables in Rest. 98 99best(CLP, Vs, Best, Rest) :- 100 findall(Delta-N,fm_cp_filter(CLP, Vs, Delta, N),Deltas), 101 keysort(Deltas,[_-N|_]), 102 select_nth(Vs,N,Best,Rest). 103 104% fm_cp_filter(Vs,Delta,N) 105% 106% For an indepenent variable V in Vs, which is the N'th element in Vs, 107% find how many inequalities are generated when this variable is eliminated. 108% Note that target variables and variables that only occur in unbounded equations 109% should have been removed from Vs via prefilter/2 110 111fm_cp_filter(CLP, Vs, Delta, N) :- 112 length(Vs,Len), % Len = number of variables in Vs 113 mem(Vs,X,Vst), % Selects a variable X in Vs, Vst is the list of elements after X in Vs 114 get_attr(X,clpcd_itf,Att), 115 arg(4,Att,lin(Lin)), 116 arg(5,Att,order(OrdX)), 117 arg(9,Att,n), % no target variable 118 indep(CLP, Lin, OrdX), % X is an independent variable 119 occurences(X,Occ), 120 Occ = [_|_], 121 cp_card(Occ,0,Lnew), 122 length(Occ,Locc), 123 Delta is Lnew-Locc, 124 length(Vst,Vstl), 125 N is Len-Vstl. % X is the Nth element in Vs 126 127% mem(Xs,X,XsT) 128% 129% If X is a member of Xs, XsT is the list of elements after X in Xs. 130 131mem([X|Xs],X,Xs). 132mem([_|Ys],X,Xs) :- mem(Ys,X,Xs). 133 134% select_nth(List,N,Nth,Others) 135% 136% Selects the N th element of List, stores it in Nth and returns the rest of the list in Others. 137 138select_nth(List,N,Nth,Others) :- 139 select_nth(List,1,N,Nth,Others). 140 141select_nth([X|Xs],N,N,X,Xs) :- !. 142select_nth([Y|Ys],M,N,X,[Y|Xs]) :- 143 M1 is M+1, 144 select_nth(Ys,M1,N,X,Xs). 145 146% 147% fm_detach + reverse_pivot introduce indep t_none, which 148% invalidates the invariants 149% 150elim_min(CLP, V, Occ, Target, Pivots, NewPivots) :- 151 crossproduct(Occ, CLP, New, []), 152 activate_crossproduct(New, CLP), 153 reverse_pivot(Pivots, CLP), 154 fm_detach(Occ, CLP), 155 allvars(V,All), 156 redundancy_vars(All, CLP), % only for New \== [] 157 make_target_indep(Target,NewPivots), 158 drop_dep(All). 159 160% 161% restore NF by reverse pivoting 162% 163reverse_pivot([], _). 164reverse_pivot([I:D|Ps], CLP) :- 165 get_attr(D,clpcd_itf,AttD), 166 arg(2,AttD,type(Dt)), 167 setarg(11,AttD,n), % no longer 168 get_attr(I,clpcd_itf,AttI), 169 arg(2,AttI,type(It)), 170 arg(5,AttI,order(OrdI)), 171 arg(6,AttI,class(ClI)), 172 pivot(CLP, D, ClI, OrdI, Dt, It), 173 reverse_pivot(Ps, CLP). 174 175% unkeep(Pivots) 176% 177% 178 179unkeep([]). 180unkeep([_:D|Ps]) :- 181 get_attr(D,clpcd_itf,Att), 182 setarg(11,Att,n), 183 drop_dep_one(D), 184 unkeep(Ps). 185 186 187% 188% All we drop are bounds 189% 190fm_detach([], _). 191fm_detach([V:_|Vs], CLP) :- 192 detach_bounds(CLP, V), 193 fm_detach(Vs, CLP). 194 195% activate_crossproduct(Lst) 196% 197% For each inequality Lin =< 0 (or Lin < 0) in Lst, a new variable is created: 198% Var = Lin and Var =< 0 (or Var < 0). Var is added to the basis. 199 200activate_crossproduct([], _). 201activate_crossproduct([lez(Strict,Lin)|News], CLP) :- 202 var_with_def_intern(t_u(0), CLP, Var, Lin, Strict), 203 % Var belongs to same class as elements in Lin 204 basis_add(Var,_), 205 activate_crossproduct(News, CLP). 206 207% ------------------------------------------------------------------------------ 208 209% crossproduct(Lst,Res,ResTail) 210% 211% See crossproduct/4 212% This predicate each time puts the next element of Lst as First in crossproduct/4 213% and lets the rest be Next. 214 215crossproduct([], _) --> []. 216crossproduct([A|As], CLP) --> 217 crossproduct(As, CLP, CLP, A), 218 crossproduct(As, CLP). 219 220% crossproduct(Next,First,Res,ResTail) 221% 222% Eliminates a variable in linear equations First + Next and stores the generated 223% inequalities in Res. 224% Let's say A:K1 = First and B:K2 = first equation in Next. 225% A = ... + K1*V + ... 226% B = ... + K2*V + ... 227% Let K = -K2/K1 228% then K*A + B = ... + 0*V + ... 229% from the bounds of A and B, via cross_lower/7 and cross_upper/7, new inequalities 230% are generated. Then the same is done for B:K2 = next element in Next. 231 232crossproduct([], CLP, CLP, _) --> []. 233crossproduct([B:Kb|Bs], CLP, CLP, A:Ka) --> 234 { 235 get_attr(A,clpcd_itf,AttA), 236 arg(2,AttA,type(Ta)), 237 arg(3,AttA,strictness(Sa)), 238 arg(4,AttA,lin(LinA)), 239 get_attr(B,clpcd_itf,AttB), 240 arg(2,AttB,type(Tb)), 241 arg(3,AttB,strictness(Sb)), 242 arg(4,AttB,lin(LinB)), 243 div_d(CLP, -Kb, Ka, K), 244 add_linear_f1(CLP, LinA, K, LinB, Lin) % Lin doesn't contain the target variable anymore 245 }, 246 ( { K > 0 } % K > 0: signs were opposite 247 -> { Strict is Sa \/ Sb }, 248 cross_lower(CLP, Ta, Tb, K, Lin, Strict), 249 cross_upper(CLP, Ta, Tb, K, Lin, Strict) 250 ; % La =< A =< Ua -> -Ua =< -A =< -La 251 { 252 flip(Ta,Taf), 253 flip_strict(Sa,Saf), 254 Strict is Saf \/ Sb 255 }, 256 cross_lower(CLP, Taf, Tb, K, Lin, Strict), 257 cross_upper(CLP, Taf, Tb, K, Lin, Strict) 258 ), 259 crossproduct(Bs, CLP, CLP, A:Ka). 260 261% cross_lower(Ta,Tb,K,Lin,Strict,Res,ResTail) 262% 263% Generates a constraint following from the bounds of A and B. 264% When A = LinA and B = LinB then Lin = K*LinA + LinB. Ta is the type 265% of A and Tb is the type of B. Strict is the union of the strictness 266% of A and B. If K is negative, then Ta should have been flipped (flip/2). 267% The idea is that if La =< A =< Ua and Lb =< B =< Ub (=< can also be <) 268% then if K is positive, K*La + Lb =< K*A + B =< K*Ua + Ub. 269% if K is negative, K*Ua + Lb =< K*A + B =< K*La + Ub. 270% This predicate handles the first inequality and adds it to Res in the form 271% lez(Sl,Lhs) meaning K*La + Lb - (K*A + B) =< 0 or K*Ua + Lb - (K*A + B) =< 0 272% with Sl being the strictness and Lhs the lefthandside of the equation. 273% See also cross_upper/7 274 275cross_lower(CLP, Ta, Tb, K, Lin, Strict) --> 276 { 277 lower(Ta,La), 278 lower(Tb,Lb), 279 !, 280 eval_d(CLP, K*La+Lb, L), 281 normalize_scalar(L,Ln), 282 add_linear_f1(CLP, Lin, -1, Ln, Lhs), 283 Sl is Strict >> 1 % normalize to upper bound 284 }, 285 [ lez(Sl,Lhs) ]. 286cross_lower(_, _, _, _, _, _) --> []. 287 288% cross_upper(Ta,Tb,K,Lin,Strict,Res,ResTail) 289% 290% See cross_lower/7 291% This predicate handles the second inequality: 292% -(K*Ua + Ub) + K*A + B =< 0 or -(K*La + Ub) + K*A + B =< 0 293 294cross_upper(CLP, Ta, Tb, K, Lin, Strict) --> 295 { 296 upper(Ta,Ua), 297 upper(Tb,Ub), 298 !, 299 eval_d(CLP, -(K*Ua+Ub), U), 300 normalize_scalar(U,Un), 301 add_linear_11(CLP, Un, Lin, Lhs), 302 Su is Strict /\ 1 % normalize to upper bound 303 }, 304 [ lez(Su,Lhs) ]. 305cross_upper(_, _, _, _, _, _) --> []. 306 307% lower(Type,Lowerbound) 308% 309% Returns the lowerbound of type Type if it has one. 310% E.g. if type = t_l(L) then Lowerbound is L, 311% if type = t_lU(L,U) then Lowerbound is L, 312% if type = t_u(U) then fails 313 314lower(t_l(L),L). 315lower(t_lu(L,_),L). 316lower(t_L(L),L). 317lower(t_Lu(L,_),L). 318lower(t_lU(L,_),L). 319 320% upper(Type,Upperbound) 321% 322% Returns the upperbound of type Type if it has one. 323% See lower/2 324 325upper(t_u(U),U). 326upper(t_lu(_,U),U). 327upper(t_U(U),U). 328upper(t_Lu(_,U),U). 329upper(t_lU(_,U),U). 330 331% flip(Type,FlippedType) 332% 333% Flips the lower and upperbound, so the old lowerbound becomes the new upperbound and 334% vice versa. 335 336flip(t_l(X),t_u(X)). 337flip(t_u(X),t_l(X)). 338flip(t_lu(X,Y),t_lu(Y,X)). 339flip(t_L(X),t_u(X)). 340flip(t_U(X),t_l(X)). 341flip(t_lU(X,Y),t_lu(Y,X)). 342flip(t_Lu(X,Y),t_lu(Y,X)). 343 344% flip_strict(Strict,FlippedStrict) 345% 346% Does what flip/2 does, but for the strictness. 347 348flip_strict(0,0). 349flip_strict(1,2). 350flip_strict(2,1). 351flip_strict(3,3). 352 353% cp_card(Lst,CountIn,CountOut) 354% 355% Counts the number of bounds that may generate an inequality in 356% crossproduct/3 357 358cp_card([],Ci,Ci). 359cp_card([A|As],Ci,Co) :- 360 cp_card(As,A,Ci,Cii), 361 cp_card(As,Cii,Co). 362 363% cp_card(Next,First,CountIn,CountOut) 364% 365% Counts the number of bounds that may generate an inequality in 366% crossproduct/4. 367 368cp_card([],_,Ci,Ci). 369cp_card([B:Kb|Bs],A:Ka,Ci,Co) :- 370 get_attr(A,clpcd_itf,AttA), 371 arg(2,AttA,type(Ta)), 372 get_attr(B,clpcd_itf,AttB), 373 arg(2,AttB,type(Tb)), 374 ( sign(Ka) =\= sign(Kb) 375 -> cp_card_lower(Ta,Tb,Ci,Cii), 376 cp_card_upper(Ta,Tb,Cii,Ciii) 377 ; flip(Ta,Taf), 378 cp_card_lower(Taf,Tb,Ci,Cii), 379 cp_card_upper(Taf,Tb,Cii,Ciii) 380 ), 381 cp_card(Bs,A:Ka,Ciii,Co). 382 383% cp_card_lower(TypeA,TypeB,SIn,SOut) 384% 385% SOut = SIn + 1 if both TypeA and TypeB have a lowerbound. 386 387cp_card_lower(Ta,Tb,Si,So) :- 388 lower(Ta,_), 389 lower(Tb,_), 390 !, 391 So is Si+1. 392cp_card_lower(_,_,Si,Si). 393 394% cp_card_upper(TypeA,TypeB,SIn,SOut) 395% 396% SOut = SIn + 1 if both TypeA and TypeB have an upperbound. 397 398cp_card_upper(Ta,Tb,Si,So) :- 399 upper(Ta,_), 400 upper(Tb,_), 401 !, 402 So is Si+1. 403cp_card_upper(_,_,Si,Si). 404 405% ------------------------------------------------------------------------------ 406 407% occurences(V,Occ) 408% 409% Returns in Occ the occurrences of variable V in the linear equations of dependent variables 410% with bound =\= t_none in the form of D:K where D is a dependent variable and K is the scalar 411% of V in the linear equation of D. 412 413occurences(V,Occ) :- 414 get_attr(V,clpcd_itf,Att), 415 arg(5,Att,order(OrdV)), 416 arg(6,Att,class(C)), 417 class_allvars(C,All), 418 occurences(All,OrdV,Occ). 419 420% occurences(De,OrdV,Occ) 421% 422% Returns in Occ the occurrences of variable V with order OrdV in the linear equations of 423% dependent variables De with bound =\= t_none in the form of D:K where D is a dependent 424% variable and K is the scalar of V in the linear equation of D. 425 426occurences(De,_,[]) :- 427 var(De), 428 !. 429occurences([D|De],OrdV,Occ) :- 430 ( get_attr(D,clpcd_itf,Att), 431 arg(2,Att,type(Type)), 432 arg(4,Att,lin(Lin)), 433 occ_type_filter(Type), 434 nf_coeff_of(Lin,OrdV,K) 435 -> Occ = [D:K|Occt], 436 occurences(De,OrdV,Occt) 437 ; occurences(De,OrdV,Occ) 438 ). 439 440% occ_type_filter(Type) 441% 442% Succeeds when Type is any other type than t_none. Is used in occurences/3 and occurs/2 443 444occ_type_filter(t_l(_)). 445occ_type_filter(t_u(_)). 446occ_type_filter(t_lu(_,_)). 447occ_type_filter(t_L(_)). 448occ_type_filter(t_U(_)). 449occ_type_filter(t_lU(_,_)). 450occ_type_filter(t_Lu(_,_)). 451 452% occurs(V) 453% 454% Checks whether variable V occurs in a linear equation of a dependent variable with a bound 455% =\= t_none. 456 457occurs(V) :- 458 get_attr(V,clpcd_itf,Att), 459 arg(5,Att,order(OrdV)), 460 arg(6,Att,class(C)), 461 class_allvars(C,All), 462 occurs(All,OrdV). 463 464% occurs(De,OrdV) 465% 466% Checks whether variable V with order OrdV occurs in a linear equation of any dependent variable 467% in De with a bound =\= t_none. 468 469occurs(De,_) :- 470 var(De), 471 !, 472 fail. 473occurs([D|De],OrdV) :- 474 ( get_attr(D,clpcd_itf,Att), 475 arg(2,Att,type(Type)), 476 arg(4,Att,lin(Lin)), 477 occ_type_filter(Type), 478 nf_coeff_of(Lin,OrdV,_) 479 -> true 480 ; occurs(De,OrdV) 481 )