1/*
    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(ineq_r,
   43	[
   44	    ineq/4,
   45	    ineq_one/4,
   46	    ineq_one_n_n_0/1,
   47	    ineq_one_n_p_0/1,
   48	    ineq_one_s_n_0/1,
   49	    ineq_one_s_p_0/1
   50	]).   51:- use_module(bv_r,
   52	[
   53	    backsubst/3,
   54	    backsubst_delta/4,
   55	    basis_add/2,
   56	    dec_step/2,
   57	    deref/2,
   58	    determine_active_dec/1,
   59	    determine_active_inc/1,
   60	    export_binding/1,
   61	    get_or_add_class/2,
   62	    inc_step/2,
   63	    lb/3,
   64	    pivot_a/4,
   65	    rcbl_status/6,
   66	    reconsider/1,
   67	    same_class/2,
   68	    solve/1,
   69	    ub/3,
   70	    unconstrained/4,
   71	    var_intern/3,
   72	    var_with_def_intern/4
   73	]).   74:- use_module(store_r,
   75	[
   76	    add_linear_11/3,
   77	    add_linear_ff/5,
   78	    normalize_scalar/2
   79	]).   80
   81% ineq(H,I,Nf,Strictness)
   82%
   83% Solves the inequality Nf < 0 or Nf =< 0 where Nf is in normal form
   84% and H and I are the homogene and inhomogene parts of Nf.
   85
   86ineq([],I,_,Strictness) :- ineq_ground(Strictness,I).
   87ineq([v(K,[X^1])|Tail],I,Lin,Strictness) :-
   88	ineq_cases(Tail,I,Lin,Strictness,X,K).
   89
   90ineq_cases([],I,_,Strictness,X,K) :-	% K*X + I < 0 or K*X + I =< 0
   91	ineq_one(Strictness,X,K,I).
   92ineq_cases([_|_],_,Lin,Strictness,_,_) :-
   93	deref(Lin,Lind),	% Id+Hd =< 0
   94	Lind = [Inhom,_|Hom],
   95	ineq_more(Hom,Inhom,Lind,Strictness).
   96
   97% ineq_ground(Strictness,I)
   98%
   99% Checks whether a grounded inequality I < 0 or I =< 0 is satisfied.
  100
  101ineq_ground(strict,I) :- I < -1.0e-10.		% I < 0
  102ineq_ground(nonstrict,I) :- I < 1.0e-10.	% I =< 0
  103
  104% ineq_one(Strictness,X,K,I)
  105%
  106% Solves the inequality K*X + I < 0 or K*X + I =< 0
  107
  108ineq_one(strict,X,K,I) :-
  109	(   K > 1.0e-10 % K > 0.0
  110	->  (   I >= -1.0e-10, % I =:= 0.0
  111		I =< 1.0e-10
  112	    ->  ineq_one_s_p_0(X)	% K*X < 0, K > 0 => X < 0
  113	    ;   Inhom is I/K,
  114		ineq_one_s_p_i(X,Inhom)	% K*X + I < 0, K > 0 => X + I/K < 0
  115	    )
  116	;   (   I >= -1.0e-10, % I =:= 0.0
  117		I =< 1.0e-10
  118	    ->  ineq_one_s_n_0(X)	% K*X < 0, K < 0 => -X < 0
  119	    ;   Inhom is -I/K,
  120		ineq_one_s_n_i(X,Inhom)	% K*X + I < 0, K < 0 => -X - I/K < 0
  121	    )
  122	).
  123ineq_one(nonstrict,X,K,I) :-
  124	(   K > 1.0e-10 % K > 0.0
  125	->  (   I >= -1.0e-10,	% I =:= 0
  126		I =< 1.0e-10
  127	    ->  ineq_one_n_p_0(X)	% K*X =< 0, K > 0 => X =< 0
  128	    ;   Inhom is I/K,
  129		ineq_one_n_p_i(X,Inhom)	% K*X + I =< 0, K > 0 => X + I/K =< 0
  130	    )
  131	;   (   I >= -1.0e-10,	% I =:= 0
  132		I =< 1.0e-10
  133	    ->  ineq_one_n_n_0(X)	% K*X =< 0, K < 0 => -X =< 0
  134	    ;   Inhom is -I/K,
  135		ineq_one_n_n_i(X,Inhom)	% K*X + I =< 0, K < 0 => -X - I/K =< 0
  136	    )
  137	).
  138
  139% --------------------------- strict ----------------------------
  140
  141% ineq_one_s_p_0(X)
  142%
  143% Solves the inequality X < 0
  144
  145ineq_one_s_p_0(X) :-
  146	get_attr(X,clpqr_itf,Att),
  147	arg(4,Att,lin([Ix,_|OrdX])),
  148	!, % old variable, this is deref
  149	(   \+ arg(1,Att,clpr)
  150	->  throw(error(permission_error('mix CLP(Q) variables with',
  151		'CLP(R) variables:',X),context(_)))
  152	;   ineq_one_old_s_p_0(OrdX,X,Ix)
  153	).
  154ineq_one_s_p_0(X) :-	% new variable, nothing depends on it
  155	var_intern(t_u(0.0),X,1). % put a strict inactive upperbound on the variable
  156
  157% ineq_one_s_n_0(X)
  158%
  159% Solves the inequality X > 0
  160
  161ineq_one_s_n_0(X) :-
  162	get_attr(X,clpqr_itf,Att),
  163	arg(4,Att,lin([Ix,_|OrdX])),
  164	!,
  165	(   \+ arg(1,Att,clpr)
  166	->  throw(error(permission_error('mix CLP(Q) variables with',
  167		'CLP(R) variables:',X),context(_)))
  168	;   ineq_one_old_s_n_0(OrdX,X,Ix)
  169	).
  170ineq_one_s_n_0(X) :-
  171	var_intern(t_l(0.0),X,2). % puts a strict inactive lowerbound on the variable
  172
  173% ineq_one_s_p_i(X,I)
  174%
  175% Solves the inequality X < -I
  176
  177ineq_one_s_p_i(X,I) :-
  178	get_attr(X,clpqr_itf,Att),
  179	arg(4,Att,lin([Ix,_|OrdX])),
  180	!,
  181	(   \+ arg(1,Att,clpr)
  182	->  throw(error(permission_error('mix CLP(Q) variables with',
  183		'CLP(R) variables:',X),context(_)))
  184	;   ineq_one_old_s_p_i(OrdX,I,X,Ix)
  185	).
  186ineq_one_s_p_i(X,I) :-
  187	Bound is -I,
  188	var_intern(t_u(Bound),X,1). % puts a strict inactive upperbound on the variable
  189
  190% ineq_one_s_n_i(X,I)
  191%
  192% Solves the inequality X > I
  193
  194ineq_one_s_n_i(X,I) :-
  195	get_attr(X,clpqr_itf,Att),
  196	arg(4,Att,lin([Ix,_|OrdX])),
  197	!,
  198	(   \+ arg(1,Att,clpr)
  199	->  throw(error(permission_error('mix CLP(Q) variables with',
  200		'CLP(R) variables:',X),context(_)))
  201	;   ineq_one_old_s_n_i(OrdX,I,X,Ix)
  202	).
  203ineq_one_s_n_i(X,I) :- var_intern(t_l(I),X,2). % puts a strict inactive lowerbound on the variable
  204
  205% ineq_one_old_s_p_0(Hom,X,Inhom)
  206%
  207% Solves the inequality X < 0 where X has linear equation Hom + Inhom
  208
  209ineq_one_old_s_p_0([],_,Ix) :- Ix < -1.0e-10. % X = I: Ix < 0
  210ineq_one_old_s_p_0([l(Y*Ky,_)|Tail],X,Ix) :-
  211	(   Tail = [] % X = K*Y + I
  212	->  Bound is -Ix/Ky,
  213	    update_indep(strict,Y,Ky,Bound)	% X < 0, X = K*Y + I => Y < -I/K or Y > -I/K (depending on K)
  214	;   Tail = [_|_]
  215	->  get_attr(X,clpqr_itf,Att),
  216	    arg(2,Att,type(Type)),
  217	    arg(3,Att,strictness(Old)),
  218	    arg(4,Att,lin(Lin)),
  219	    udus(Type,X,Lin,0.0,Old)	% update strict upperbound
  220	).
  221
  222% ineq_one_old_s_p_0(Hom,X,Inhom)
  223%
  224% Solves the inequality X > 0 where X has linear equation Hom + Inhom
  225
  226ineq_one_old_s_n_0([],_,Ix) :- Ix > 1.0e-10. % X = I: Ix > 0
  227ineq_one_old_s_n_0([l(Y*Ky,_)|Tail], X, Ix) :-
  228	(   Tail = []	% X = K*Y + I
  229	->  Coeff is -Ky,
  230	    Bound is Ix/Coeff,
  231	    update_indep(strict,Y,Coeff,Bound)
  232	;   Tail = [_|_]
  233	->  get_attr(X,clpqr_itf,Att),
  234	    arg(2,Att,type(Type)),
  235	    arg(3,Att,strictness(Old)),
  236	    arg(4,Att,lin(Lin)),
  237	    udls(Type,X,Lin,0.0,Old)	% update strict lowerbound
  238	).
  239
  240% ineq_one_old_s_p_i(Hom,C,X,Inhom)
  241%
  242% Solves the inequality X + C < 0 where X has linear equation Hom + Inhom
  243
  244ineq_one_old_s_p_i([],I,_,Ix) :- Ix + I < -1.0e-10. % X = I
  245ineq_one_old_s_p_i([l(Y*Ky,_)|Tail],I,X,Ix) :-
  246	(   Tail = []	% X = K*Y + I
  247	->  Bound is -(Ix + I)/Ky,
  248	    update_indep(strict,Y,Ky,Bound)
  249	;   Tail = [_|_]
  250	->  Bound is -I,
  251	    get_attr(X,clpqr_itf,Att),
  252	    arg(2,Att,type(Type)),
  253	    arg(3,Att,strictness(Old)),
  254	    arg(4,Att,lin(Lin)),
  255	    udus(Type,X,Lin,Bound,Old)	% update strict upperbound
  256	).
  257
  258% ineq_one_old_s_n_i(Hom,C,X,Inhom)
  259%
  260% Solves the inequality X  - C > 0 where X has linear equation Hom + Inhom
  261
  262ineq_one_old_s_n_i([],I,_,Ix) :- -Ix + I < -1.0e-10. % X = I
  263ineq_one_old_s_n_i([l(Y*Ky,_)|Tail],I,X,Ix) :-
  264	(   Tail = []	% X = K*Y + I
  265	->  Coeff is -Ky,
  266	    Bound is (Ix - I)/Coeff,
  267	    update_indep(strict,Y,Coeff,Bound)
  268	;   Tail = [_|_]
  269	->  get_attr(X,clpqr_itf,Att),
  270	    arg(2,Att,type(Type)),
  271	    arg(3,Att,strictness(Old)),
  272	    arg(4,Att,lin(Lin)),
  273	    udls(Type,X,Lin,I,Old)	% update strict lowerbound
  274	).
  275
  276% -------------------------- nonstrict --------------------------
  277
  278% ineq_one_n_p_0(X)
  279%
  280% Solves the inequality X =< 0
  281
  282ineq_one_n_p_0(X) :-
  283	get_attr(X,clpqr_itf,Att),
  284	arg(4,Att,lin([Ix,_|OrdX])),
  285	!,	% old variable, this is deref
  286	(   \+ arg(1,Att,clpr)
  287	->  throw(error(permission_error('mix CLP(Q) variables with',
  288		'CLP(R) variables:',X),context(_)))
  289	;   ineq_one_old_n_p_0(OrdX,X,Ix)
  290	).
  291ineq_one_n_p_0(X) :-	% new variable, nothing depends on it
  292	var_intern(t_u(0.0),X,0).	% nonstrict upperbound
  293
  294% ineq_one_n_n_0(X)
  295%
  296% Solves the inequality X >= 0
  297
  298ineq_one_n_n_0(X) :-
  299	get_attr(X,clpqr_itf,Att),
  300	arg(4,Att,lin([Ix,_|OrdX])),
  301	!,
  302	(   \+ arg(1,Att,clpr)
  303	->  throw(error(permission_error('mix CLP(Q) variables with',
  304		'CLP(R) variables:',X),context(_)))
  305	;   ineq_one_old_n_n_0(OrdX,X,Ix)
  306	).
  307ineq_one_n_n_0(X) :-
  308	var_intern(t_l(0.0),X,0).	% nonstrict lowerbound
  309
  310% ineq_one_n_p_i(X,I)
  311%
  312% Solves the inequality X =< -I
  313
  314ineq_one_n_p_i(X,I) :-
  315	get_attr(X,clpqr_itf,Att),
  316	arg(4,Att,lin([Ix,_|OrdX])),
  317	!,
  318	(   \+ arg(1,Att,clpr)
  319	->  throw(error(permission_error('mix CLP(Q) variables with',
  320		'CLP(R) variables:',X),context(_)))
  321	;   ineq_one_old_n_p_i(OrdX,I,X,Ix)
  322	).
  323ineq_one_n_p_i(X,I) :-
  324	Bound is -I,
  325	var_intern(t_u(Bound),X,0).	% nonstrict upperbound
  326
  327% ineq_one_n_n_i(X,I)
  328%
  329% Solves the inequality X >= I
  330
  331ineq_one_n_n_i(X,I) :-
  332	get_attr(X,clpqr_itf,Att),
  333	arg(4,Att,lin([Ix,_|OrdX])),
  334	!,
  335	(   \+ arg(1,Att,clpr)
  336	->  throw(error(permission_error('mix CLP(Q) variables with',
  337		'CLP(R) variables:',X),context(_)))
  338	;   ineq_one_old_n_n_i(OrdX,I,X,Ix)
  339	).
  340ineq_one_n_n_i(X,I) :-
  341	var_intern(t_l(I),X,0).	% nonstrict lowerbound
  342
  343% ineq_one_old_n_p_0(Hom,X,Inhom)
  344%
  345% Solves the inequality X =< 0 where X has linear equation Hom + Inhom
  346
  347ineq_one_old_n_p_0([],_,Ix) :- Ix < 1.0e-10. % X =I
  348ineq_one_old_n_p_0([l(Y*Ky,_)|Tail],X,Ix) :-
  349	(   Tail = []	%  X = K*Y + I
  350	->  Bound is -Ix/Ky,
  351	    update_indep(nonstrict,Y,Ky,Bound)
  352	;   Tail = [_|_]
  353	->  get_attr(X,clpqr_itf,Att),
  354	    arg(2,Att,type(Type)),
  355	    arg(3,Att,strictness(Old)),
  356	    arg(4,Att,lin(Lin)),
  357	    udu(Type,X,Lin,0.0,Old)	% update nonstrict upperbound
  358	).
  359
  360% ineq_one_old_n_n_0(Hom,X,Inhom)
  361%
  362% Solves the inequality X >= 0 where X has linear equation Hom + Inhom
  363
  364ineq_one_old_n_n_0([],_,Ix) :- Ix > -1.0e-10.	% X = I
  365ineq_one_old_n_n_0([l(Y*Ky,_)|Tail], X, Ix) :-
  366	(   Tail = []	% X = K*Y + I
  367	->  Coeff is -Ky,
  368	    Bound is Ix/Coeff,
  369	    update_indep(nonstrict,Y,Coeff,Bound)
  370	;   Tail = [_|_]
  371	->  get_attr(X,clpqr_itf,Att),
  372	    arg(2,Att,type(Type)),
  373	    arg(3,Att,strictness(Old)),
  374	    arg(4,Att,lin(Lin)),
  375	    udl(Type,X,Lin,0.0,Old)	% update nonstrict lowerbound
  376	).
  377
  378% ineq_one_old_n_p_i(Hom,C,X,Inhom)
  379%
  380% Solves the inequality X  + C =< 0 where X has linear equation Hom + Inhom
  381
  382ineq_one_old_n_p_i([],I,_,Ix) :- Ix + I < 1.0e-10.	% X = I
  383ineq_one_old_n_p_i([l(Y*Ky,_)|Tail],I,X,Ix) :-
  384	(   Tail = []	% X = K*Y + I
  385	->  Bound is -(Ix + I)/Ky,
  386	    update_indep(nonstrict,Y,Ky,Bound)
  387	;   Tail = [_|_]
  388	->  Bound is -I,
  389	    get_attr(X,clpqr_itf,Att),
  390	    arg(2,Att,type(Type)),
  391	    arg(3,Att,strictness(Old)),
  392	    arg(4,Att,lin(Lin)),
  393	    udu(Type,X,Lin,Bound,Old)	% update nonstrict upperbound
  394	).
  395
  396% ineq_one_old_n_n_i(Hom,C,X,Inhom)
  397%
  398% Solves the inequality X  - C >= 0 where X has linear equation Hom + Inhom
  399
  400ineq_one_old_n_n_i([],I,_,Ix) :- -Ix + I < 1.0e-10. % X = I
  401ineq_one_old_n_n_i([l(Y*Ky,_)|Tail],I,X,Ix) :-
  402	(   Tail = []
  403	->  Coeff is -Ky,
  404	    Bound is (Ix - I)/Coeff,
  405	    update_indep(nonstrict,Y,Coeff,Bound)
  406	;   Tail = [_|_]
  407	->  get_attr(X,clpqr_itf,Att),
  408	    arg(2,Att,type(Type)),
  409	    arg(3,Att,strictness(Old)),
  410	    arg(4,Att,lin(Lin)),
  411	    udl(Type,X,Lin,I,Old)
  412	).
  413
  414% ---------------------------------------------------------------
  415
  416% ineq_more(Hom,Inhom,Lin,Strictness)
  417%
  418% Solves the inequality Lin < 0 or Lin =< 0 with Lin = Hom + Inhom
  419
  420ineq_more([],I,_,Strictness) :- ineq_ground(Strictness,I).	% I < 0 or I =< 0
  421ineq_more([l(X*K,_)|Tail],Id,Lind,Strictness) :-
  422	(   Tail = []
  423	->  % X*K < Id or X*K =< Id
  424	    % one var: update bound instead of slack introduction
  425	    get_or_add_class(X,_),	% makes sure X belongs to a class
  426	    Bound is -Id/K,
  427	    update_indep(Strictness,X,K,Bound)	% new bound
  428	;   Tail = [_|_]
  429	->  ineq_more(Strictness,Lind)
  430	).
  431
  432% ineq_more(Strictness,Lin)
  433%
  434% Solves the inequality Lin < 0 or Lin =< 0
  435
  436ineq_more(strict,Lind) :-
  437	(   unconstrained(Lind,U,K,Rest)
  438	->  % never fails, no implied value
  439	    % Lind < 0 => Rest < -K*U where U has no bounds
  440	    var_intern(t_l(0.0),S,2),	% create slack variable S
  441	    get_attr(S,clpqr_itf,AttS),
  442	    arg(5,AttS,order(OrdS)),
  443	    Ki is -1.0/K,
  444	    add_linear_ff(Rest,Ki,[0.0,0.0,l(S*1.0,OrdS)],Ki,LinU),	% U = (-1/K)*Rest + (-1/K)*S
  445	    LinU = [_,_|Hu],
  446	    get_or_add_class(U,Class),
  447	    same_class(Hu,Class),	% put all variables of new lin. eq. of U in the same class
  448	    get_attr(U,clpqr_itf,AttU),
  449	    arg(5,AttU,order(OrdU)),
  450	    arg(6,AttU,class(ClassU)),
  451	    backsubst(ClassU,OrdU,LinU)	% substitute U by new lin. eq. everywhere in the class
  452	;   var_with_def_intern(t_u(0.0),S,Lind,1),	% Lind < 0 => Lind = S with S < 0
  453	    basis_add(S,_),			% adds S to the basis
  454	    determine_active_dec(Lind),		% activate bounds
  455	    reconsider(S)			% reconsider basis
  456	).
  457ineq_more(nonstrict,Lind) :-
  458	(   unconstrained(Lind,U,K,Rest)
  459	->  % never fails, no implied value
  460	    % Lind =< 0 => Rest =< -K*U where U has no bounds
  461	    var_intern(t_l(0.0),S,0),	% create slack variable S
  462	    Ki is -1.0/K,
  463	    get_attr(S,clpqr_itf,AttS),
  464	    arg(5,AttS,order(OrdS)),
  465	    add_linear_ff(Rest,Ki,[0.0,0.0,l(S*1.0,OrdS)],Ki,LinU),	% U = (-1K)*Rest + (-1/K)*S
  466	    LinU = [_,_|Hu],
  467	    get_or_add_class(U,Class),
  468	    same_class(Hu,Class),	% put all variables of new lin. eq of U in the same class
  469	    get_attr(U,clpqr_itf,AttU),
  470	    arg(5,AttU,order(OrdU)),
  471	    arg(6,AttU,class(ClassU)),
  472	    backsubst(ClassU,OrdU,LinU)	% substitute U by new lin. eq. everywhere in the class
  473	;   % all variables are constrained
  474	    var_with_def_intern(t_u(0.0),S,Lind,0),	% Lind =< 0 => Lind = S with S =< 0
  475	    basis_add(S,_),				% adds S to the basis
  476	    determine_active_dec(Lind),
  477	    reconsider(S)
  478	).
  479
  480
  481% update_indep(Strictness,X,K,Bound)
  482%
  483% Updates the bound of independent variable X where X < Bound or X =< Bound
  484% or X > Bound or X >= Bound, depending on Strictness and K.
  485
  486update_indep(strict,X,K,Bound) :-
  487	get_attr(X,clpqr_itf,Att),
  488	arg(2,Att,type(Type)),
  489	arg(3,Att,strictness(Old)),
  490	arg(4,Att,lin(Lin)),
  491	(   K < -1.0e-10
  492	->  uils(Type,X,Lin,Bound,Old)	% update independent lowerbound strict
  493	;   uius(Type,X,Lin,Bound,Old)	% update independent upperbound strict
  494	).
  495update_indep(nonstrict,X,K,Bound) :-
  496	get_attr(X,clpqr_itf,Att),
  497	arg(2,Att,type(Type)),
  498	arg(3,Att,strictness(Old)),
  499	arg(4,Att,lin(Lin)),
  500	(   K < -1.0e-10
  501	->  uil(Type,X,Lin,Bound,Old)	% update independent lowerbound nonstrict
  502	;   uiu(Type,X,Lin,Bound,Old)	% update independent upperbound nonstrict
  503	).
  504
  505
  506% ---------------------------------------------------------------------------------------
  507
  508%
  509% Update a bound on a var xi
  510%
  511%   a) independent variable
  512%
  513%	a1) update inactive bound: done
  514%
  515%	a2) update active bound:
  516%	    Determine [lu]b including most constraining row R
  517%	      If we are within: done
  518%	    else pivot(R,xi) and introduce bound via (b)
  519%
  520%	a3) introduce a bound on an unconstrained var:
  521%	    All vars that depend on xi are unconstrained (invariant) ->
  522%	      the bound cannot invalidate any Lhs
  523%
  524%   b) dependent variable
  525%
  526%	repair upper or lower (maybe just swap with an unconstrained var from Rhs)
  527%
  528
  529%
  530% Sign = 1,0,-1 means inside,at,outside
  531%
  532
  533% Read following predicates as update (dependent/independent) (lowerbound/upperbound) (strict)
  534
  535% udl(Type,X,Lin,Bound,Strict)
  536%
  537% Updates lower bound of dependent variable X with linear equation
  538% Lin that had type Type and strictness Strict, to the new non-strict
  539% bound Bound.
  540
  541udl(t_none,X,Lin,Bound,_Sold) :-
  542	get_attr(X,clpqr_itf,AttX),
  543	arg(5,AttX,order(Ord)),
  544	setarg(2,AttX,type(t_l(Bound))),
  545	setarg(3,AttX,strictness(0)),
  546	(   unconstrained(Lin,Uc,Kuc,Rest)
  547	->  % X = Lin => -1/K*Rest + 1/K*X = U where U has no bounds
  548	    Ki is -1.0/Kuc,
  549	    add_linear_ff(Rest,Ki,[0.0,0.0,l(X* -1.0,Ord)],Ki,LinU),
  550	    get_attr(Uc,clpqr_itf,AttU),
  551	    arg(5,AttU,order(OrdU)),
  552	    arg(6,AttU,class(Class)),
  553	    backsubst(Class,OrdU,LinU)
  554	;   % no unconstrained variables in Lin: make X part of basis and reconsider
  555	    basis_add(X,_),
  556	    determine_active_inc(Lin),
  557	    reconsider(X)
  558	).
  559udl(t_l(L),X,Lin,Bound,Sold) :-
  560	TestBL is Bound - L,
  561	(   TestBL < -1.0e-10
  562	->  true	% new bound is smaller than old one: keep old
  563	;   TestBL > 1.0e-10
  564	->  % new bound is larger than old one: use new and reconsider basis
  565	    Strict is Sold /\ 1,
  566	    get_attr(X,clpqr_itf,Att),
  567	    setarg(2,Att,type(t_l(Bound))),
  568	    setarg(3,Att,strictness(Strict)),
  569	    reconsider_lower(X,Lin,Bound)	% makes sure that Lin still satisfies lowerbound Bound
  570	;   true	% new bound is equal to old one, new one is nonstrict: keep old
  571	).
  572
  573udl(t_u(U),X,Lin,Bound,_Sold) :-
  574	TestUB is U - Bound,
  575	(   TestUB < -1.0e-10
  576	->  fail	% new bound is larger than upperbound: fail
  577	;   TestUB > 1.0e-10
  578	->  % new bound is smaller than upperbound: add new and reconsider basis
  579	    get_attr(X,clpqr_itf,Att),
  580	    setarg(2,Att,type(t_lu(Bound,U))),
  581	    reconsider_lower(X,Lin,Bound)	% makes sure that Lin still satisfies lowerbound Bound
  582	;   solve_bound(Lin,Bound)	% new bound is equal to upperbound: solve
  583	).
  584udl(t_lu(L,U),X,Lin,Bound,Sold) :-
  585	TestBL is Bound - L,
  586	(   TestBL < -1.0e-10
  587	->  true	% smaller than lowerbound: keep
  588	;   TestBL > 1.0e-10
  589	->  % larger than lowerbound: check upperbound
  590	    TestUB is U - Bound,
  591	    (   TestUB < -1.0e-10
  592	    ->  fail	% larger than upperbound: fail
  593	    ;   TestUB > 1.0e-10
  594	    ->  % smaller than upperbound: use new and reconsider basis
  595		Strict is Sold /\ 1,
  596		get_attr(X,clpqr_itf,Att),
  597		setarg(2,Att,type(t_lu(Bound,U))),
  598		setarg(3,Att,strictness(Strict)),
  599		reconsider_lower(X,Lin,Bound)
  600	    ;   % equal to upperbound: if strictness matches => solve
  601		Sold /\ 1 =:= 0,
  602		solve_bound(Lin,Bound)
  603	    )
  604	;   true	% equal to lowerbound and nonstrict: keep
  605	).
  606
  607% udls(Type,X,Lin,Bound,Strict)
  608%
  609% Updates lower bound of dependent variable X with linear equation
  610% Lin that had type Type and strictness Strict, to the new strict
  611% bound Bound.
  612
  613udls(t_none,X,Lin,Bound,_Sold) :-
  614	get_attr(X,clpqr_itf,AttX),
  615	arg(5,AttX,order(Ord)),
  616	setarg(2,AttX,type(t_l(Bound))),
  617	setarg(3,AttX,strictness(2)),
  618	(   unconstrained(Lin,Uc,Kuc,Rest)
  619	->  % X = Lin => U = -1/K*Rest + 1/K*X with U an unconstrained variable
  620	    Ki is -1.0/Kuc,
  621	    add_linear_ff(Rest,Ki,[0.0,0.0,l(X* -1.0,Ord)],Ki,LinU),
  622	    get_attr(Uc,clpqr_itf,AttU),
  623	    arg(5,AttU,order(OrdU)),
  624	    arg(6,AttU,class(Class)),
  625	    backsubst(Class,OrdU,LinU)
  626	;   % no unconstrained variables: add X to basis and reconsider basis
  627	    basis_add(X,_),
  628	    determine_active_inc(Lin),
  629	    reconsider(X)
  630	).
  631udls(t_l(L),X,Lin,Bound,Sold) :-
  632	TestBL is Bound - L,
  633	(   TestBL < -1.0e-10
  634	->  true	% smaller than lowerbound: keep
  635	;   TestBL > 1.0e-10
  636	->  % larger than lowerbound: use new and reconsider basis
  637	    Strict is Sold \/ 2,
  638	    get_attr(X,clpqr_itf,Att),
  639	    setarg(2,Att,type(t_l(Bound))),
  640	    setarg(3,Att,strictness(Strict)),
  641	    reconsider_lower(X,Lin,Bound)
  642	;   % equal to lowerbound: check strictness
  643	    Strict is Sold \/ 2,
  644	    get_attr(X,clpqr_itf,Att),
  645	    setarg(3,Att,strictness(Strict))
  646	).
  647udls(t_u(U),X,Lin,Bound,Sold) :-
  648	U - Bound > 1.0e-10,	% smaller than upperbound: set new bound
  649	Strict is Sold \/ 2,
  650	get_attr(X,clpqr_itf,Att),
  651	setarg(2,Att,type(t_lu(Bound,U))),
  652	setarg(3,Att,strictness(Strict)),
  653	reconsider_lower(X,Lin,Bound).
  654udls(t_lu(L,U),X,Lin,Bound,Sold) :-
  655	TestBL is Bound - L,
  656	(   TestBL < -1.0e-10
  657	->  true	% smaller than lowerbound: keep
  658	;   TestBL > 1.0e-10
  659	->  % larger than lowerbound: check upperbound and possibly use new and reconsider basis
  660	    U - Bound > 1.0e-10,
  661	    Strict is Sold \/ 2,
  662	    get_attr(X,clpqr_itf,Att),
  663	    setarg(2,Att,type(t_lu(Bound,U))),
  664	    setarg(3,Att,strictness(Strict)),
  665	    reconsider_lower(X,Lin,Bound)
  666	;   % equal to lowerbound: put new strictness
  667	    Strict is Sold \/ 2,
  668	    get_attr(X,clpqr_itf,Att),
  669	    setarg(3,Att,strictness(Strict))
  670	).
  671
  672% udu(Type,X,Lin,Bound,Strict)
  673%
  674% Updates upper bound of dependent variable X with linear equation
  675% Lin that had type Type and strictness Strict, to the new non-strict
  676% bound Bound.
  677
  678udu(t_none,X,Lin,Bound,_Sold) :-
  679	get_attr(X,clpqr_itf,AttX),
  680	arg(5,AttX,order(Ord)),
  681	setarg(2,AttX,type(t_u(Bound))),
  682	setarg(3,AttX,strictness(0)),
  683	(   unconstrained(Lin,Uc,Kuc,Rest)
  684	->  % X = Lin => U = -1/K*Rest + 1/K*X with U an unconstrained variable
  685	    Ki is -1.0/Kuc,
  686	    add_linear_ff(Rest,Ki,[0.0,0.0,l(X* -1.0,Ord)],Ki,LinU),
  687	    get_attr(Uc,clpqr_itf,AttU),
  688	    arg(5,AttU,order(OrdU)),
  689	    arg(6,AttU,class(Class)),
  690	    backsubst(Class,OrdU,LinU)
  691	;   % no unconstrained variables: add X to basis and reconsider basis
  692	    basis_add(X,_),
  693	    determine_active_dec(Lin),	% try to lower R
  694	    reconsider(X)
  695	).
  696udu(t_u(U),X,Lin,Bound,Sold) :-
  697	TestUB is U - Bound,
  698	(   TestUB < -1.0e-10
  699	->  true	% larger than upperbound: keep
  700	;   TestUB > 1.0e-10
  701	->  % smaller than upperbound: update and reconsider basis
  702	    Strict is Sold /\ 2,
  703	    get_attr(X,clpqr_itf,Att),
  704	    setarg(2,Att,type(t_u(Bound))),
  705	    setarg(3,Att,strictness(Strict)),
  706	    reconsider_upper(X,Lin,Bound)
  707	;   true	% equal to upperbound and nonstrict: keep
  708	).
  709udu(t_l(L),X,Lin,Bound,_Sold) :-
  710	TestBL is Bound - L,
  711	(   TestBL < -1.0e-10
  712	->  fail	% smaller than lowerbound: fail
  713	;   TestBL > 1.0e-10
  714	->  % larger than lowerbound: use new and reconsider basis
  715	    get_attr(X,clpqr_itf,Att),
  716	    setarg(2,Att,type(t_lu(L,Bound))),
  717	    reconsider_upper(X,Lin,Bound)
  718	;   solve_bound(Lin,Bound)	% equal to lowerbound: solve
  719	).
  720udu(t_lu(L,U),X,Lin,Bound,Sold) :-
  721	TestUB is U - Bound,
  722	(   TestUB < -1.0e-10
  723	->  true	% larger than upperbound: keep
  724	;   TestUB > 1.0e-10
  725	->  % smaller than upperbound: check lowerbound
  726	    TestBL is Bound - L,
  727	    (   TestBL < -1.0e-10
  728	    ->  fail	% smaller than lowerbound: fail
  729	    ;   TestBL > 1.0e-10
  730	    ->  % larger than lowerbound: update and reconsider basis
  731		Strict is Sold /\ 2,
  732		get_attr(X,clpqr_itf,Att),
  733		setarg(2,Att,type(t_lu(L,Bound))),
  734		setarg(3,Att,strictness(Strict)),
  735		reconsider_upper(X,Lin,Bound)
  736	    ;   % equal to lowerbound: check strictness and possibly solve
  737		Sold /\ 2 =:= 0,
  738		solve_bound(Lin,Bound)
  739	    )
  740	;   true	% equal to upperbound and nonstrict: keep
  741	).
  742
  743% udus(Type,X,Lin,Bound,Strict)
  744%
  745% Updates upper bound of dependent variable X with linear equation
  746% Lin that had type Type and strictness Strict, to the new strict
  747% bound Bound.
  748
  749udus(t_none,X,Lin,Bound,_Sold) :-
  750	get_attr(X,clpqr_itf,AttX),
  751	arg(5,AttX,order(Ord)),
  752	setarg(2,AttX,type(t_u(Bound))),
  753	setarg(3,AttX,strictness(1)),
  754	(   unconstrained(Lin,Uc,Kuc,Rest)
  755	->   % X = Lin => U = -1/K*Rest + 1/K*X with U an unconstrained variable
  756	    Ki is -1.0/Kuc,
  757	    add_linear_ff(Rest,Ki,[0.0,0.0,l(X* -1.0,Ord)],Ki,LinU),
  758	    get_attr(Uc,clpqr_itf,AttU),
  759	    arg(5,AttU,order(OrdU)),
  760	    arg(6,AttU,class(Class)),
  761	    backsubst(Class,OrdU,LinU)
  762	;   % no unconstrained variables: add X to basis and reconsider basis
  763	    basis_add(X,_),
  764	    determine_active_dec(Lin),
  765	    reconsider(X)
  766	).
  767udus(t_u(U),X,Lin,Bound,Sold) :-
  768	TestUB is U - Bound,
  769	(   TestUB < -1.0e-10
  770	->  true	% larger than upperbound: keep
  771	;   TestUB > 1.0e-10
  772	->  % smaller than upperbound: update bound and reconsider basis
  773	    Strict is Sold \/ 1,
  774	    get_attr(X,clpqr_itf,Att),
  775	    setarg(2,Att,type(t_u(Bound))),
  776	    setarg(3,Att,strictness(Strict)),
  777	    reconsider_upper(X,Lin,Bound)
  778	;   % equal to upperbound: set new strictness
  779	    Strict is Sold \/ 1,
  780	    get_attr(X,clpqr_itf,Att),
  781	    setarg(3,Att,strictness(Strict))
  782	).
  783udus(t_l(L),X,Lin,Bound,Sold) :-
  784	Bound - L > 1.0e-10,	% larger than lowerbound: update and reconsider basis
  785	Strict is Sold \/ 1,
  786	get_attr(X,clpqr_itf,Att),
  787	setarg(2,Att,type(t_lu(L,Bound))),
  788	setarg(3,Att,strictness(Strict)),
  789	reconsider_upper(X,Lin,Bound).
  790udus(t_lu(L,U),X,Lin,Bound,Sold) :-
  791	TestUB is U - Bound,
  792	(   TestUB < -1.0e-10
  793	->  true	% larger than upperbound: keep
  794	;   TestUB > 1.0e-10
  795	->  % smaller than upperbound: check lowerbound, possibly update and reconsider basis
  796	    Bound - L > 1.0e-10,
  797	    Strict is Sold \/ 1,
  798	    get_attr(X,clpqr_itf,Att),
  799	    setarg(2,Att,type(t_lu(L,Bound))),
  800	    setarg(3,Att,strictness(Strict)),
  801	    reconsider_upper(X,Lin,Bound)
  802	;   % equal to upperbound: update strictness
  803	    Strict is Sold \/ 1,
  804	    get_attr(X,clpqr_itf,Att),
  805	    setarg(3,Att,strictness(Strict))
  806	).
  807
  808% uiu(Type,X,Lin,Bound,Strict)
  809%
  810% Updates upper bound of independent variable X with linear equation
  811% Lin that had type Type and strictness Strict, to the new non-strict
  812% bound Bound.
  813
  814uiu(t_none,X,_Lin,Bound,_) :-	% X had no bounds
  815	get_attr(X,clpqr_itf,Att),
  816	setarg(2,Att,type(t_u(Bound))),
  817	setarg(3,Att,strictness(0)).
  818uiu(t_u(U),X,_Lin,Bound,Sold) :-
  819	TestUB is U - Bound,
  820	(   TestUB < -1.0e-10
  821	->  true	% larger than upperbound: keep
  822	;   TestUB > 1.0e-10
  823	->  % smaller than upperbound: update.
  824	    Strict is Sold /\ 2,	% update strictness: strictness of lowerbound is kept,
  825					% strictness of upperbound is set to non-strict
  826	    get_attr(X,clpqr_itf,Att),
  827	    setarg(2,Att,type(t_u(Bound))),
  828	    setarg(3,Att,strictness(Strict))
  829	;   true	% equal to upperbound and nonstrict: keep
  830	).
  831uiu(t_l(L),X,Lin,Bound,_Sold) :-
  832	TestBL is Bound - L,
  833	(   TestBL < -1.0e-10
  834	->  fail	% Lowerbound was smaller than new upperbound: fail
  835	;   TestBL > 1.0e-10
  836	->   % Upperbound is larger than lowerbound: store new bound
  837	    get_attr(X,clpqr_itf,Att),
  838	    setarg(2,Att,type(t_lu(L,Bound)))
  839	;   solve_bound(Lin,Bound) % Lowerbound was equal to new upperbound: solve
  840	).
  841uiu(t_L(L),X,Lin,Bound,_Sold) :-
  842	TestBL is Bound - L,
  843	(   TestBL < -1.0e-10
  844	->  fail	% Same as for t_l
  845	;   TestBL > 1.0e-10
  846	->  % Same as for t_l (new bound becomes t_Lu)
  847	    get_attr(X,clpqr_itf,Att),
  848	    setarg(2,Att,type(t_Lu(L,Bound)))
  849	;   solve_bound(Lin,Bound)	% Same as for t_l
  850	).
  851uiu(t_lu(L,U),X,Lin,Bound,Sold) :-
  852	TestUB is U - Bound,
  853	(   TestUB < -1.0e-10
  854	->  true	% Upperbound was smaller than new bound: keep
  855	;   TestUB > 1.0e-10
  856	->  TestBL is Bound - L,	% Upperbound was larger than new bound: check lowerbound
  857	    (   TestBL < -1.0e-10
  858	    ->	fail	% Lowerbound was larger than new bound: fail
  859	    ;   TestBL > 1.0e-10
  860	    ->  % Lowerbound was smaller than new bound: store new bound
  861		Strict is Sold /\ 2,
  862		get_attr(X,clpqr_itf,Att),
  863		setarg(2,Att,type(t_lu(L,Bound))),
  864		setarg(3,Att,strictness(Strict))
  865	    ;	% Lowerbound was equal to new bound: solve
  866		Sold /\ 2 =:= 0,	% Only solve when strictness matches
  867		solve_bound(Lin,Bound)
  868	    )
  869	;   true	% Upperbound was equal to new bound and new bound non-strict: keep
  870	).
  871uiu(t_Lu(L,U),X,Lin,Bound,Sold) :-	% See t_lu case
  872	TestUB is U - Bound,
  873	(   TestUB < -1.0e-10
  874	->  true
  875	;   TestUB > 1.0e-10
  876	->  TestBL is Bound - L,
  877	    (   TestBL < -1.0e-10
  878	    ->  fail
  879	    ;   TestBL > 1.0e-10
  880	    ->  Strict is Sold /\ 2,
  881		get_attr(X,clpqr_itf,Att),
  882		setarg(2,Att,type(t_Lu(L,Bound))),
  883		setarg(3,Att,strictness(Strict))
  884	    ;   Sold /\ 2 =:= 0,
  885		solve_bound(Lin,Bound)
  886	    )
  887	;   true
  888	).
  889uiu(t_U(U),X,_Lin,Bound,Sold) :-
  890	TestUB is U - Bound,
  891	(   TestUB < -1.0e-10
  892	->  true	% larger than upperbound: keep
  893	;   TestUB > 1.0e-10
  894	->  % smaller than active upperbound: check how much active upperbound can be lowered.
  895	    % if enough, just lower bound, otherwise update the bound, make X dependent and reconsider basis
  896	    Strict is Sold /\ 2,
  897	    (   get_attr(X,clpqr_itf,Att),
  898		arg(5,Att,order(OrdX)),
  899		arg(6,Att,class(ClassX)),
  900		lb(ClassX,OrdX,Vlb-Vb-Lb),
  901		Bound - (Lb + U) < 1.0e-10
  902	    ->  get_attr(X,clpqr_itf,Att2), % changed?
  903		setarg(2,Att2,type(t_U(Bound))),
  904		setarg(3,Att2,strictness(Strict)),
  905		pivot_a(Vlb,X,Vb,t_u(Bound)),
  906		reconsider(X)
  907	    ;   get_attr(X,clpqr_itf,Att),
  908		arg(5,Att,order(OrdX)),
  909		arg(6,Att,class(ClassX)),
  910		setarg(2,Att,type(t_U(Bound))),
  911		setarg(3,Att,strictness(Strict)),
  912		Delta is Bound - U,
  913		backsubst_delta(ClassX,OrdX,X,Delta)
  914	    )
  915	;   true	% equal to upperbound and non-strict: keep
  916	).
  917uiu(t_lU(L,U),X,Lin,Bound,Sold) :-
  918	TestUB is U - Bound,
  919	(   TestUB < -1.0e-10
  920	->  true	% larger than upperbound: keep
  921	;   TestUB > 1.0e-10
  922	->  TestBL is Bound-L,
  923	    (   TestBL < -1.0e-10
  924	    ->  fail	% smaller than lowerbound: fail
  925	    ;   TestBL > 1.0e-10
  926	    ->  % larger than lowerbound: see t_U case for rest
  927		Strict is Sold /\ 2,
  928		(   get_attr(X,clpqr_itf,Att),
  929		    arg(5,Att,order(OrdX)),
  930		    arg(6,Att,class(ClassX)),
  931		    lb(ClassX,OrdX,Vlb-Vb-Lb),
  932		    Bound - (Lb + U) < 1.0e-10
  933		->  get_attr(X,clpqr_itf,Att2), % changed?
  934		    setarg(2,Att2,type(t_lU(L,Bound))),
  935		    setarg(3,Att2,strictness(Strict)),
  936		    pivot_a(Vlb,X,Vb,t_lu(L,Bound)),
  937		    reconsider(X)
  938		;   get_attr(X,clpqr_itf,Att),
  939		    arg(5,Att,order(OrdX)),
  940		    arg(6,Att,class(ClassX)),
  941		    setarg(2,Att,type(t_lU(L,Bound))),
  942		    setarg(3,Att,strictness(Strict)),
  943		    Delta is Bound - U,
  944		    backsubst_delta(ClassX,OrdX,X,Delta)
  945		)
  946	    ;	% equal to lowerbound: check strictness and solve
  947		Sold /\ 2 =:= 0,
  948		solve_bound(Lin,Bound)
  949	    )
  950	;   true	% equal to upperbound and non-strict: keep
  951			% smaller than upperbound: check lowerbound
  952	).
  953
  954% uius(Type,X,Lin,Bound,Strict)
  955%
  956% Updates upper bound of independent variable X with linear equation
  957% Lin that had type Type and strictness Strict, to the new strict
  958% bound Bound. (see also uiu/5)
  959
  960uius(t_none,X,_Lin,Bound,_Sold) :-
  961	get_attr(X,clpqr_itf,Att),
  962	setarg(2,Att,type(t_u(Bound))),
  963	setarg(3,Att,strictness(1)).
  964uius(t_u(U),X,_Lin,Bound,Sold) :-
  965	TestUB is U - Bound,
  966	(   TestUB < -1.0e-10
  967	->  true
  968	;   TestUB > 1.0e-10
  969	->  Strict is Sold \/ 1,
  970	    get_attr(X,clpqr_itf,Att),
  971	    setarg(2,Att,type(t_u(Bound))),
  972	    setarg(3,Att,strictness(Strict))
  973	;   Strict is Sold \/ 1,
  974	    get_attr(X,clpqr_itf,Att),
  975	    setarg(3,Att,strictness(Strict))
  976	).
  977uius(t_l(L),X,_Lin,Bound,Sold) :-
  978	Bound - L > 1.0e-10,
  979	Strict is Sold \/ 1,
  980	get_attr(X,clpqr_itf,Att),
  981	setarg(2,Att,type(t_lu(L,Bound))),
  982	setarg(3,Att,strictness(Strict)).
  983uius(t_L(L),X,_Lin,Bound,Sold) :-
  984	Bound - L > 1.0e-10,
  985	Strict is Sold \/ 1,
  986	get_attr(X,clpqr_itf,Att),
  987	setarg(2,Att,type(t_Lu(L,Bound))),
  988	setarg(3,Att,strictness(Strict)).
  989uius(t_lu(L,U),X,_Lin,Bound,Sold) :-
  990	TestUB is U - Bound,
  991	(   TestUB < -1.0e-10
  992	->  true
  993	;   TestUB > 1.0e-10
  994	->  Bound - L > 1.0e-10,
  995	    Strict is Sold \/ 1,
  996	    get_attr(X,clpqr_itf,Att),
  997	    setarg(2,Att,type(t_lu(L,Bound))),
  998	    setarg(3,Att,strictness(Strict))
  999	;   Strict is Sold \/ 1,
 1000	    get_attr(X,clpqr_itf,Att),
 1001	    setarg(3,Att,strictness(Strict))
 1002	).
 1003uius(t_Lu(L,U),X,_Lin,Bound,Sold) :-
 1004	TestUB is U - Bound,
 1005	(   TestUB < -1.0e-10
 1006	->  true
 1007	;   TestUB > 1.0e-10
 1008	->  Bound - L > 1.0e-10,
 1009	    Strict is Sold \/ 1,
 1010	    get_attr(X,clpqr_itf,Att),
 1011	    setarg(2,Att,type(t_Lu(L,Bound))),
 1012	    setarg(3,Att,strictness(Strict))
 1013	;   Strict is Sold \/ 1,
 1014	    get_attr(X,clpqr_itf,Att),
 1015	    setarg(3,Att,strictness(Strict))
 1016	).
 1017uius(t_U(U),X,_Lin,Bound,Sold) :-
 1018	TestUB is U - Bound,
 1019	(   TestUB < -1.0e-10
 1020	->  true
 1021	;   TestUB > 1.0e-10
 1022	->  Strict is Sold \/ 1,
 1023	    (   get_attr(X,clpqr_itf,Att),
 1024		arg(5,Att,order(OrdX)),
 1025		arg(6,Att,class(ClassX)),
 1026		lb(ClassX,OrdX,Vlb-Vb-Lb),
 1027		Bound - (Lb + U) < 1.0e-10
 1028	    ->  get_attr(X,clpqr_itf,Att2), % changed?
 1029		setarg(2,Att2,type(t_U(Bound))),
 1030		setarg(3,Att2,strictness(Strict)),
 1031		pivot_a(Vlb,X,Vb,t_u(Bound)),
 1032		reconsider(X)
 1033	    ;   get_attr(X,clpqr_itf,Att),
 1034		arg(5,Att,order(OrdX)),
 1035		arg(6,Att,class(ClassX)),
 1036		setarg(2,Att,type(t_U(Bound))),
 1037		setarg(3,Att,strictness(Strict)),
 1038		Delta is Bound - U,
 1039		backsubst_delta(ClassX,OrdX,X,Delta)
 1040	    )
 1041	;   Strict is Sold \/ 1,
 1042	    get_attr(X,clpqr_itf,Att),
 1043	    setarg(3,Att,strictness(Strict))
 1044	).
 1045uius(t_lU(L,U),X,_Lin,Bound,Sold) :-
 1046	TestUB is U - Bound,
 1047	(   TestUB < -1.0e-10
 1048	->  true
 1049	;   TestUB > 1.0e-10
 1050	->  Bound - L > 1.0e-10,
 1051	    Strict is Sold \/ 1,
 1052	    (   get_attr(X,clpqr_itf,Att),
 1053		arg(5,Att,order(OrdX)),
 1054		arg(6,Att,class(ClassX)),
 1055		lb(ClassX,OrdX,Vlb-Vb-Lb),
 1056		Bound - (Lb + U) < 1.0e-10
 1057	    ->  get_attr(X,clpqr_itf,Att2), % changed?
 1058		setarg(2,Att2,type(t_lU(L,Bound))),
 1059		setarg(3,Att2,strictness(Strict)),
 1060		pivot_a(Vlb,X,Vb,t_lu(L,Bound)),
 1061		reconsider(X)
 1062	    ;	get_attr(X,clpqr_itf,Att),
 1063		arg(5,Att,order(OrdX)),
 1064		arg(6,Att,class(ClassX)),
 1065		setarg(2,Att,type(t_lU(L,Bound))),
 1066		setarg(3,Att,strictness(Strict)),
 1067		Delta is Bound - U,
 1068		backsubst_delta(ClassX,OrdX,X,Delta)
 1069	    )
 1070	;   Strict is Sold \/ 1,
 1071	    get_attr(X,clpqr_itf,Att),
 1072	    setarg(3,Att,strictness(Strict))
 1073	).
 1074
 1075% uil(Type,X,Lin,Bound,Strict)
 1076%
 1077% Updates lower bound of independent variable X with linear equation
 1078% Lin that had type Type and strictness Strict, to the new non-strict
 1079% bound Bound. (see also uiu/5)
 1080
 1081
 1082uil(t_none,X,_Lin,Bound,_Sold) :-
 1083	get_attr(X,clpqr_itf,Att),
 1084	setarg(2,Att,type(t_l(Bound))),
 1085	setarg(3,Att,strictness(0)).
 1086uil(t_l(L),X,_Lin,Bound,Sold) :-
 1087	TestBL is Bound - L,
 1088	(   TestBL < -1.0e-10
 1089	->  true
 1090	;   TestBL > 1.0e-10
 1091	->  Strict is Sold /\ 1,
 1092	    get_attr(X,clpqr_itf,Att),
 1093	    setarg(2,Att,type(t_l(Bound))),
 1094	    setarg(3,Att,strictness(Strict))
 1095	;   true
 1096	).
 1097uil(t_u(U),X,Lin,Bound,_Sold) :-
 1098	TestUB is U - Bound,
 1099	(   TestUB < -1.0e-10
 1100	->  fail
 1101	;   TestUB > 1.0e-10
 1102	->  get_attr(X,clpqr_itf,Att),
 1103	    setarg(2,Att,type(t_lu(Bound,U)))
 1104	;   solve_bound(Lin,Bound)
 1105	).
 1106uil(t_U(U),X,Lin,Bound,_Sold) :-
 1107	TestUB is U - Bound,
 1108	(   TestUB < -1.0e-10
 1109	->  fail
 1110	;   TestUB > 1.0e-10
 1111	->  get_attr(X,clpqr_itf,Att),
 1112	    setarg(2,Att,type(t_lU(Bound,U)))
 1113	;   solve_bound(Lin,Bound)
 1114	).
 1115uil(t_lu(L,U),X,Lin,Bound,Sold) :-
 1116	TestBL is Bound - L,
 1117	(   TestBL < -1.0e-10
 1118	->  true
 1119	;   TestBL > 1.0e-10
 1120	->  TestUB is U - Bound,
 1121	    (   TestUB < -1.0e-10
 1122	    ->  fail
 1123	    ;   TestUB > 1.0e-10
 1124	    ->  Strict is Sold /\ 1,
 1125		get_attr(X,clpqr_itf,Att),
 1126		setarg(2,Att,type(t_lu(Bound,U))),
 1127		setarg(3,Att,strictness(Strict))
 1128	    ;   Sold /\ 1 =:= 0,
 1129		solve_bound(Lin,Bound)
 1130	    )
 1131	;   true
 1132	).
 1133uil(t_lU(L,U),X,Lin,Bound,Sold) :-
 1134	TestBL is Bound - L,
 1135	(   TestBL < -1.0e-10
 1136	->  true
 1137	;   TestBL > 1.0e-10
 1138	->  TestUB is U - Bound,
 1139	    (   TestUB < -1.0e-10
 1140	    ->  fail
 1141	    ;   TestUB > 1.0e-10
 1142	    ->  Strict is Sold /\ 1,
 1143		get_attr(X,clpqr_itf,Att),
 1144		setarg(2,Att,type(t_lU(Bound,U))),
 1145		setarg(3,Att,strictness(Strict))
 1146	    ;   Sold /\ 1 =:= 0,
 1147		solve_bound(Lin,Bound)
 1148	    )
 1149	;   true
 1150	).
 1151uil(t_L(L),X,_Lin,Bound,Sold) :-
 1152	TestBL is Bound - L,
 1153	(   TestBL < -1.0e-10
 1154	->  true
 1155	;   TestBL > 1.0e-10
 1156	->  Strict is Sold /\ 1,
 1157	    (   get_attr(X,clpqr_itf,Att),
 1158		arg(5,Att,order(OrdX)),
 1159		arg(6,Att,class(ClassX)),
 1160		ub(ClassX,OrdX,Vub-Vb-Ub),
 1161		Bound - (Ub + L) > -1.0e-10
 1162	    ->  get_attr(X,clpqr_itf,Att2), % changed?
 1163		setarg(2,Att2,type(t_L(Bound))),
 1164		setarg(3,Att2,strictness(Strict)),
 1165		pivot_a(Vub,X,Vb,t_l(Bound)),
 1166		reconsider(X)
 1167	    ;   get_attr(X,clpqr_itf,Att),
 1168		arg(5,Att,order(OrdX)),
 1169		arg(6,Att,class(ClassX)),
 1170		setarg(2,Att,type(t_L(Bound))),
 1171		setarg(3,Att,strictness(Strict)),
 1172		Delta is Bound - L,
 1173		backsubst_delta(ClassX,OrdX,X,Delta)
 1174	    )
 1175	;   true
 1176	).
 1177uil(t_Lu(L,U),X,Lin,Bound,Sold) :-
 1178	TestBL is Bound - L,
 1179	(   TestBL < -1.0e-10
 1180	->  true
 1181	;   TestBL > 1.0e-10
 1182	->  TestUB is U - Bound,
 1183	    (   TestUB < -1.0e-10
 1184	    ->  fail
 1185	    ;   TestUB > 1.0e-10
 1186	    ->  Strict is Sold /\ 1,
 1187		(   get_attr(X,clpqr_itf,Att),
 1188		    arg(5,Att,order(OrdX)),
 1189		    arg(6,Att,class(ClassX)),
 1190		    ub(ClassX,OrdX,Vub-Vb-Ub),
 1191		    Bound - (Ub + L) > -1.0e-10
 1192		->  get_attr(X,clpqr_itf,Att2), % changed?
 1193		    setarg(2,Att2,t_Lu(Bound,U)),
 1194		    setarg(3,Att2,strictness(Strict)),
 1195		    pivot_a(Vub,X,Vb,t_lu(Bound,U)),
 1196		    reconsider(X)
 1197		;   get_attr(X,clpqr_itf,Att),
 1198		    arg(5,Att,order(OrdX)),
 1199		    arg(6,Att,class(ClassX)),
 1200		    setarg(2,Att,type(t_Lu(Bound,U))),
 1201		    setarg(3,Att,strictness(Strict)),
 1202		    Delta is Bound - L,
 1203		    backsubst_delta(ClassX,OrdX,X,Delta)
 1204		)
 1205	    ;	Sold /\ 1 =:= 0,
 1206		solve_bound(Lin,Bound)
 1207	    )
 1208	;   true
 1209	).
 1210
 1211% uils(Type,X,Lin,Bound,Strict)
 1212%
 1213% Updates lower bound of independent variable X with linear equation
 1214% Lin that had type Type and strictness Strict, to the new strict
 1215% bound Bound. (see also uiu/5)
 1216
 1217uils(t_none,X,_Lin,Bound,_Sold) :-
 1218	get_attr(X,clpqr_itf,Att),
 1219	setarg(2,Att,type(t_l(Bound))),
 1220	setarg(3,Att,strictness(2)).
 1221uils(t_l(L),X,_Lin,Bound,Sold) :-
 1222	TestBL is Bound - L,
 1223	(   TestBL < -1.0e-10
 1224	->  true
 1225	;   TestBL > 1.0e-10
 1226	->  Strict is Sold \/ 2,
 1227	    get_attr(X,clpqr_itf,Att),
 1228	    setarg(2,Att,type(t_l(Bound))),
 1229	    setarg(3,Att,strictness(Strict))
 1230	;   Strict is Sold \/ 2,
 1231	    get_attr(X,clpqr_itf,Att),
 1232	    setarg(3,Att,strictness(Strict))
 1233	).
 1234uils(t_u(U),X,_Lin,Bound,Sold) :-
 1235	U - Bound > 1.0e-10,
 1236	Strict is Sold \/ 2,
 1237	get_attr(X,clpqr_itf,Att),
 1238	setarg(2,Att,type(t_lu(Bound,U))),
 1239	setarg(3,Att,strictness(Strict)).
 1240uils(t_U(U),X,_Lin,Bound,Sold) :-
 1241	U - Bound > 1.0e-10,
 1242	Strict is Sold \/ 2,
 1243	get_attr(X,clpqr_itf,Att),
 1244	setarg(2,Att,type(t_lU(Bound,U))),
 1245	setarg(3,Att,strictness(Strict)).
 1246uils(t_lu(L,U),X,_Lin,Bound,Sold) :-
 1247	TestBL is Bound - L,
 1248	(   TestBL < -1.0e-10
 1249	->  true
 1250	;   TestBL > 1.0e-10
 1251	->  U - Bound > 1.0e-10,
 1252	    Strict is Sold \/ 2,
 1253	    get_attr(X,clpqr_itf,Att),
 1254	    setarg(2,Att,type(t_lu(Bound,U))),
 1255	    setarg(3,Att,strictness(Strict))
 1256	;   Strict is Sold \/ 2,
 1257	    get_attr(X,clpqr_itf,Att),
 1258	    setarg(3,Att,strictness(Strict))
 1259	).
 1260uils(t_lU(L,U),X,_Lin,Bound,Sold) :-
 1261	TestBL is Bound - L,
 1262	(   TestBL < -1.0e-10
 1263	->  true
 1264	;   TestBL > 1.0e-10
 1265	->  U - Bound > 1.0e-10,
 1266	    Strict is Sold \/ 2,
 1267	    get_attr(X,clpqr_itf,Att),
 1268	    setarg(2,Att,type(t_lU(Bound,U))),
 1269	    setarg(3,Att,strictness(Strict))
 1270	;   Strict is Sold \/ 2,
 1271	    get_attr(X,clpqr_itf,Att),
 1272	    setarg(3,Att,strictness(Strict))
 1273	).
 1274uils(t_L(L),X,_Lin,Bound,Sold) :-
 1275	TestBL is Bound - L,
 1276	(   TestBL < -1.0e-10
 1277	->  true
 1278	;   TestBL > 1.0e-10
 1279	->  Strict is Sold \/ 2,
 1280	    (   get_attr(X,clpqr_itf,Att),
 1281		arg(5,Att,order(OrdX)),
 1282		arg(6,Att,class(ClassX)),
 1283		ub(ClassX,OrdX,Vub-Vb-Ub),
 1284		Bound - (Ub + L) > -1.0e-10
 1285	    ->  get_attr(X,clpqr_itf,Att2), % changed?
 1286		setarg(2,Att2,type(t_L(Bound))),
 1287		setarg(3,Att2,strictness(Strict)),
 1288		pivot_a(Vub,X,Vb,t_l(Bound)),
 1289		reconsider(X)
 1290	    ;   get_attr(X,clpqr_itf,Att),
 1291		arg(5,Att,order(OrdX)),
 1292		arg(6,Att,class(ClassX)),
 1293		setarg(2,Att,type(t_L(Bound))),
 1294		setarg(3,Att,strictness(Strict)),
 1295		Delta is Bound - L,
 1296		backsubst_delta(ClassX,OrdX,X,Delta)
 1297	    )
 1298	;   Strict is Sold \/ 2,
 1299	    get_attr(X,clpqr_itf,Att),
 1300	    setarg(3,Att,strictness(Strict))
 1301	).
 1302uils(t_Lu(L,U),X,_Lin,Bound,Sold) :-
 1303	TestBL is Bound - L,
 1304	(   TestBL < -1.0e-10
 1305	->  true
 1306	;   TestBL > 1.0e-10
 1307	->  U - Bound > 1.0e-10,
 1308	    Strict is Sold \/ 2,
 1309	    (   get_attr(X,clpqr_itf,Att),
 1310		arg(5,Att,order(OrdX)),
 1311		arg(6,Att,class(ClassX)),
 1312		ub(ClassX,OrdX,Vub-Vb-Ub),
 1313		Bound - (Ub + L) > -1.0e-10
 1314	    ->  get_attr(X,clpqr_itf,Att2), % changed?
 1315		setarg(2,Att2,type(t_Lu(Bound,U))),
 1316		setarg(3,Att2,strictness(Strict)),
 1317		pivot_a(Vub,X,Vb,t_lu(Bound,U)),
 1318		reconsider(X)
 1319	    ;   get_attr(X,clpqr_itf,Att),
 1320		arg(5,Att,order(OrdX)),
 1321		arg(6,Att,class(ClassX)),
 1322		setarg(2,Att,type(t_Lu(Bound,U))),
 1323		setarg(3,Att,strictness(Strict)),
 1324		Delta is Bound - L,
 1325		backsubst_delta(ClassX,OrdX,X,Delta)
 1326	    )
 1327	;   Strict is Sold \/ 2,
 1328	    get_attr(X,clpqr_itf,Att),
 1329	    setarg(3,Att,strictness(Strict))
 1330	).
 1331
 1332% reconsider_upper(X,Lin,U)
 1333%
 1334% Checks if the upperbound of X which is U, satisfies the bounds
 1335% of the variables in Lin: let R be the sum of all the bounds on
 1336% the variables in Lin, and I be the inhomogene part of Lin, then
 1337% upperbound U should be larger or equal to R + I (R may contain
 1338% lowerbounds).
 1339% See also rcb/3 in bv.pl
 1340
 1341reconsider_upper(X,[I,R|H],U) :-
 1342	R + I - U > -1.0e-10,	% violation
 1343	!,
 1344	dec_step(H,Status),	% we want to decrement R
 1345	rcbl_status(Status,X,[],Binds,[],u(U)),
 1346	export_binding(Binds).
 1347reconsider_upper( _, _, _).
 1348
 1349% reconsider_lower(X,Lin,L)
 1350%
 1351% Checks if the lowerbound of X which is L, satisfies the bounds
 1352% of the variables in Lin: let R be the sum of all the bounds on
 1353% the variables in Lin, and I be the inhomogene part of Lin, then
 1354% lowerbound L should be smaller or equal to R + I (R may contain
 1355% upperbounds).
 1356% See also rcb/3 in bv.pl
 1357
 1358reconsider_lower(X,[I,R|H],L) :-
 1359	R + I - L < 1.0e-10,	% violation
 1360	!,
 1361	inc_step(H,Status),	% we want to increment R
 1362	rcbl_status(Status,X,[],Binds,[],l(L)),
 1363	export_binding(Binds).
 1364reconsider_lower(_,_,_).
 1365
 1366%
 1367% lin is dereferenced
 1368%
 1369
 1370% solve_bound(Lin,Bound)
 1371%
 1372% Solves the linear equation Lin - Bound = 0
 1373% Lin is the linear equation of X, a variable whose bounds have narrowed to value Bound
 1374
 1375solve_bound(Lin,Bound) :-
 1376	Bound >= -1.0e-10,
 1377	Bound =< 1.0e-10,
 1378	!,
 1379	solve(Lin).
 1380solve_bound(Lin,Bound) :-
 1381	Nb is -Bound,
 1382	normalize_scalar(Nb,Nbs),
 1383	add_linear_11(Nbs,Lin,Eq),
 1384	solve(Eq)