1/*
    2
    3    Part of CLP(Q,R) (Constraint Logic Programming over Rationals and 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): 2006, K.U. Leuven and
   10		   1992-1995, Austrian Research Institute for
   11		              Artificial Intelligence (OFAI),
   12			      Vienna, Austria
   13
   14    This software is based on CLP(Q,R) by Christian Holzbaur for SICStus
   15    Prolog and distributed under the license details below with permission from
   16    all mentioned authors.
   17
   18    This program is free software; you can redistribute it and/or
   19    modify it under the terms of the GNU General Public License
   20    as published by the Free Software Foundation; either version 2
   21    of the License, or (at your option) any later version.
   22
   23    This program is distributed in the hope that it will be useful,
   24    but WITHOUT ANY WARRANTY; without even the implied warranty of
   25    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   26    GNU General Public License for more details.
   27
   28    You should have received a copy of the GNU Lesser General Public
   29    License along with this library; if not, write to the Free Software
   30    Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
   31
   32    As a special exception, if you link this library with other files,
   33    compiled with a Free Software compiler, to produce an executable, this
   34    library does not by itself cause the resulting executable to be covered
   35    by the GNU General Public License. This exception does not however
   36    invalidate any other reasons why the executable file might be covered by
   37    the GNU General Public License.
   38*/
   39
   40
   41% attribute = t(CLP,type(_),strictness(_),lin(_),order(_),class(_),forward(_),
   42%		nonzero,target,keep_indep,keep)
   43
   44:- module(clpcd_itf,
   45	[
   46	    dump_linear/3,
   47	    dump_nonzero/3
   48	]).   49
   50:- use_module(library(apply)).   51:- use_module(library(neck)).   52:- use_module(library(clpcd/class)).   53:- use_module(library(clpcd/domain_ops)).   54:- use_module(library(clpcd/indep)).   55:- use_module(library(clpcd/bv)).   56:- use_module(library(clpcd/nf)).   57:- use_module(library(clpcd/store)).   58:- use_module(library(clpcd/solve)).   59:- init_expansors.   60
   61:- multifile attribute_goals//1.   62
   63% ----------------------------------- dump ------------------------------------
   64
   65% dump_strict(FilteredStrictness,Nonstrict,Strict,Res)
   66%
   67% Unifies Res with either Nonstrict or Strict depending on FilteredStrictness.
   68% FilteredStrictness is the component of strictness related to the bound: 0
   69% means nonstrict, 1 means strict upperbound, 2 means strict lowerbound,
   70% 3 is filtered out to either 1 or 2.
   71
   72dump_strict(0,Result,_,Result).
   73dump_strict(1,_,Result,Result).
   74dump_strict(2,_,Result,Result).
   75
   76% dump_nz(CLP,V,H,I,Dump,DumpTail)
   77%
   78% Returns in Dump a representation of the nonzero constraint of variable V
   79% which has linear
   80% equation H + I.
   81
   82dump_nz(CLP,_,H,I) -->
   83	{
   84	    H = [l(_*K,_)|_],
   85	    div_d(CLP, 1, K, Kr),
   86	    eval_d(CLP, -Kr*I, I1),
   87	    mult_hom(H,CLP,Kr,H1),
   88	    nf2sum(H1,CLP,0,Sum)
   89	},
   90	[Sum =\= I1].
   91
   92% dump_var(Type,Var,I,H,Dump,DumpTail)
   93%
   94% Returns in Dump a representation of the linear constraint on variable
   95% Var which has linear equation H + I and has type Type.
   96
   97dump_v(t_none,CLP,V,I,H) -->
   98	!,
   99	(   {
  100		H = [l(W*K,_)],
  101		V == W,
  102                compare_d(CLP, =, 0, I),
  103                compare_d(CLP, =, K, 1)
  104	    }
  105	->  % indep var
  106	    []
  107	;   {nf2sum(H,CLP,I,Sum)},
  108	    [V = Sum]
  109	).
  110dump_v(t_L(L),CLP,V,I,H) -->
  111	!,
  112	dump_v(t_l(L),CLP,V,I,H).
  113% case lowerbound: V >= L or V > L
  114% say V >= L, and V = K*V1 + ... + I, then K*V1 + ... + I >= L
  115% and K*V1 + ... >= L-I and V1 + .../K = (L-I)/K
  116dump_v(t_l(L),CLP,V,I,H) -->
  117	!,
  118	{
  119	    H = [l(_*K,_)|_], % avoid 1 >= 0
  120	    get_attr(V,clpcd_itf,Att),
  121	    arg(3,Att,strictness(Strict)),
  122	    Sm is Strict /\ 2,
  123            div_d(CLP, 1, K, Kr),
  124	    eval_d(CLP, Kr*(L - I), Li),
  125	    mult_hom(H,CLP,Kr,H1),
  126	    nf2sum(H1,CLP,0,Sum),
  127	    (   compare_d(CLP, <, 0, K)
  128	    ->	dump_strict(Sm,Sum >= Li,Sum > Li,Result)
  129	    ;   dump_strict(Sm,Sum =< Li,Sum < Li,Result)
  130	    )
  131	},
  132	[Result].
  133dump_v(t_U(U),C,V,I,H) -->
  134	!,
  135	dump_v(t_u(U),C,V,I,H).
  136dump_v(t_u(U),CLP,V,I,H) -->
  137	!,
  138	{
  139	    H = [l(_*K,_)|_], % avoid 0 =< 1
  140	    get_attr(V,clpcd_itf,Att),
  141	    arg(3,Att,strictness(Strict)),
  142	    Sm is Strict /\ 1,
  143	    div_d(CLP, 1, K, Kr),
  144	    eval_d(CLP, Kr*(U-I), Ui),
  145	    mult_hom(H,CLP,Kr,H1),
  146	    nf2sum(H1,CLP,0,Sum),
  147	    (   compare_d(CLP, <, 0, K)
  148	    ->	dump_strict(Sm,Sum =< Ui,Sum < Ui,Result)
  149	    ;   dump_strict(Sm,Sum >= Ui,Sum > Ui,Result)
  150	    )
  151	},
  152	[Result].
  153dump_v(t_Lu(L,U),C,V,I,H) -->
  154	!,
  155	dump_v(t_l(L),C,V,I,H),
  156	dump_v(t_u(U),C,V,I,H).
  157dump_v(t_lU(L,U),C,V,I,H) -->
  158	!,
  159	dump_v(t_l(L),C,V,I,H),
  160	dump_v(t_u(U),C,V,I,H).
  161dump_v(t_lu(L,U),C,V,I,H) -->
  162	!,
  163	dump_v(t_l(L),C,V,I,H),
  164	dump_v(t_U(U),C,V,I,H).
  165dump_v(T,_,V,I,H) --> % should not happen
  166	[V:T:I+H].
  167
  168
  169% nf2sum(Lin,CLP,Sofar,Term)
  170%
  171% Transforms a linear expression into a sum
  172% (e.g. the expression [5,_,[l(X*2,OrdX),l(Y*-1,OrdY)]] gets transformed into 5 + 2*X - Y)
  173
  174nf2sum([],_,I,I).
  175nf2sum([X|Xs],CLP,I,Sum) :-
  176	(   compare_d(CLP, =, 0, I)
  177	->  X = l(Var*K,_),
  178 	    (   % K =:= 1.0
  179                compare_d(CLP, =, K, 1)
  180	    ->  hom2sum(Xs,CLP,Var,Sum)
  181	    ;   % K =:= -1.0
  182                compare_d(CLP, =, K, -1)
  183	    ->  hom2sum(Xs,CLP,-Var,Sum)
  184	    ;	hom2sum(Xs,CLP,K*Var,Sum)
  185	    )
  186	;   hom2sum([X|Xs],CLP,I,Sum)
  187 	).
  188
  189% hom2sum(Hom,_,Sofar,Term)
  190%
  191% Transforms a linear expression into a sum
  192% this predicate handles all but the first term
  193% (the first term does not need a concatenation symbol + or -)
  194% see also nf2sum/3
  195
  196hom2sum([],_,Term,Term).
  197hom2sum([l(Var*K,_)|Cs],CLP,Sofar,Term) :-
  198	(   % K =:= 1.0
  199            compare_d(CLP, =, K, 1)
  200	->  Next = Sofar + Var
  201	;   % K =:= -1.0
  202            compare_d(CLP, =, K, -1)
  203	->  Next = Sofar - Var
  204	;   % K < 0.0
  205	    compare_d(CLP, <, K, 0)
  206	->  eval_d(CLP, -K, Ka),
  207	    Next = Sofar - Ka*Var
  208	;   Next = Sofar + K*Var
  209	),
  210	hom2sum(Cs,CLP,Next,Term).
  211
  212dump_linear(V) -->
  213	{
  214	    get_attr(V,clpcd_itf,Att),
  215	    arg(1,Att,CLP),
  216	    arg(2,Att,type(Type)),
  217	    arg(4,Att,lin(Lin)),
  218	    !,
  219	    Lin = [I,_|H]
  220	},
  221	(   {
  222		Type=t_none
  223	    ;	arg(9,Att,n)
  224	    }
  225	->  []
  226	;   dump_v(t_none,CLP,V,I,H)
  227	),
  228	(   {
  229		Type=t_none,
  230		arg(9,Att,n) % attribute should not have changed by dump_v...
  231	    }
  232	->  % nonzero produces such
  233	    []
  234	;   dump_v(Type,CLP,V,I,H)
  235	).
  236dump_linear(_) --> [].
  237
  238dump_nonzero(V) -->
  239	{
  240	    get_attr(V,clpcd_itf,Att),
  241	    arg(1,Att,CLP),
  242	    arg(4,Att,lin(Lin)),
  243	    arg(8,Att,nonzero),
  244	    !,
  245	    Lin = [I,_|H]
  246	},
  247	dump_nz(CLP,V,H,I).
  248dump_nonzero(_) --> [].
  249
  250attr_unify_hook(t(CLP,n,n,n,n,n,n,n,_,_,_),Y) :-
  251	!,
  252	(   get_attr(Y,clpcd_itf,AttY),
  253	    \+ arg(1,AttY,CLP)
  254	->  throw(error(permission_error('mix CLP(CD) variables with',
  255		'CLP(CD) variables of other subdomain:',Y),context(_)))
  256	;   true
  257	).
  258attr_unify_hook(t(CLP,Ty,St,Li,Or,Cl,_,No,_,_,_),Y) :-
  259	(   get_attr(Y,clpcd_itf,AttY),
  260	    \+ arg(1,AttY,CLP)
  261	->  throw(error(permission_error('mix CLP(CD) variables with',
  262		'CLP(CD) variables of other subdomain:',Y),context(_)))
  263	;   true
  264	),
  265	do_checks(CLP,Y,Ty,St,Li,Or,Cl,No,Later),
  266	maplist(call,Later).
  267
  268do_checks(CLP,Y,Ty,St,Li,Or,Cl,No,Later) :-
  269    numbers_only(CLP,Y),
  270    verify_nonzero(No,CLP,Y),
  271    verify_type(Ty,CLP,St,Y,Later,[]),
  272    verify_lin(Or,CLP,Cl,Li,Y),
  273    maplist(call,Later).
  274
  275% verify_nonzero(Nonzero,CLP,Y)
  276%
  277% if Nonzero = nonzero, then verify that Y is not zero
  278% (if possible, otherwise set Y to be nonzero)
  279
  280verify_nonzero(nonzero,CLP,Y) :-
  281    (   var(Y)
  282    ->  (   get_attr(Y,clpcd_itf,Att)
  283	->  setarg(8,Att,nonzero)
  284	;   put_attr(Y,clpcd_itf,t(CLP,n,n,n,n,n,n,nonzero,n,n,n))
  285	)
  286    ;   compare_d(CLP, \=, 0, Y)
  287    ).
  288verify_nonzero(n,_,_). % X is not nonzero
  289
  290% verify_type(type(Type),strictness(Strict),Y,[OL|OLT],OLT)
  291%
  292% if possible verifies whether Y satisfies the type and strictness of X
  293% if not possible to verify, then returns the constraints that follow from
  294% the type and strictness
  295
  296verify_type(type(Type),CLP,strictness(Strict),Y) -->
  297    verify_type2(CLP,Y,Type,Strict).
  298verify_type(n,_,n,_) --> [].
  299
  300verify_type2(C,Y,TypeX,StrictX) -->
  301    {var(Y)},
  302    !,
  303    verify_type_var(TypeX,C,Y,StrictX).
  304verify_type2(C,Y,TypeX,StrictX) -->
  305    {verify_type_nonvar(TypeX,C,Y,StrictX)}.
  306
  307% verify_type_nonvar(Type,CLP,Nonvar,Strictness)
  308%
  309% verifies whether the type and strictness are satisfied with the Nonvar
  310
  311verify_type_nonvar(t_none,_,_,_).
  312verify_type_nonvar(t_l(L),C,Value,S) :- ilb(S,C,L,Value).
  313verify_type_nonvar(t_u(U),C,Value,S) :- iub(S,C,U,Value).
  314verify_type_nonvar(t_lu(L,U),C,Value,S) :-
  315    ilb(S,C,L,Value),
  316    iub(S,C,U,Value).
  317verify_type_nonvar(t_L(L),C,Value,S) :- ilb(S,C,L,Value).
  318verify_type_nonvar(t_U(U),C,Value,S) :- iub(S,C,U,Value).
  319verify_type_nonvar(t_Lu(L,U),C,Value,S) :-
  320    ilb(S,C,L,Value),
  321    iub(S,C,U,Value).
  322verify_type_nonvar(t_lU(L,U),C,Value,S) :-
  323    ilb(S,C,L,Value),
  324    iub(S,C,U,Value).
  325
  326% ilb(Strict,CLP,Lower,Value) & iub(Strict,CLP,Upper,Value)
  327%
  328% check whether Value is satisfiable with the given lower/upper bound and
  329% strictness.
  330% strictness is encoded as follows:
  331% 2 = strict lower bound
  332% 1 = strict upper bound
  333% 3 = strict lower and upper bound
  334% 0 = no strict bounds
  335
  336ilb(S,C,L,V) :-
  337    between(0, 3, S),
  338    ( S /\ 2 =:= 0
  339    ->Op = (=<) % non-strict
  340    ; Op = (<)  % strict
  341    ),
  342    neck,
  343    compare_d(C, Op, L, V).
  344
  345iub(S,C,U,V) :-
  346    between(0, 3, S),
  347    ( S /\ 1 =:= 0
  348    ->Op = (=<) % non-strict
  349    ; Op = (<)  % strict
  350    ),
  351    neck,
  352    compare_d(C, Op, V, U).
  353
  354%
  355% Running some goals after X=Y simplifies the coding. It should be possible
  356% to run the goals here and taking care not to put_atts/2 on X ...
  357%
  358
  359% verify_type_var(Type,Var,Strictness,[OutList|OutListTail],OutListTail)
  360%
  361% returns the inequalities following from a type and strictness satisfaction
  362% test with Var
  363
  364verify_type_var(t_none,_,_,_) --> [].
  365verify_type_var(t_l(L),C,Y,S) --> llb(S,C,L,Y).
  366verify_type_var(t_u(U),C,Y,S) --> lub(S,C,U,Y).
  367verify_type_var(t_lu(L,U),C,Y,S) -->
  368    llb(S,C,L,Y),
  369    lub(S,C,U,Y).
  370verify_type_var(t_L(L),C,Y,S) --> llb(S,C,L,Y).
  371verify_type_var(t_U(U),C,Y,S) --> lub(S,C,U,Y).
  372verify_type_var(t_Lu(L,U),C,Y,S) -->
  373    llb(S,C,L,Y),
  374    lub(S,C,U,Y).
  375verify_type_var(t_lU(L,U),C,Y,S) -->
  376    llb(S,C,L,Y),
  377    lub(S,C,U,Y).
  378
  379% llb(Strict,Lower,Value,[OL|OLT],OLT) and lub(Strict,Upper,Value,[OL|OLT],OLT)
  380%
  381% returns the inequalities following from the lower and upper bounds and the
  382% strictness see also lb and ub
  383llb(S,C,L,V) -->
  384    {S /\ 2 =:= 0 },
  385    !,
  386    [C:{L =< V}].
  387llb(_,C,L,V) -->
  388    [C:{L < V}].
  389
  390lub(S,C,U,V) -->
  391    {S /\ 1 =:= 0 },
  392     !,
  393     [C:{V =< U}].
  394lub(_,C,U,V) -->
  395    [C:{V < U}].
  396
  397%
  398% We used to drop X from the class/basis to avoid trouble with subsequent
  399% put_atts/2 on X. Now we could let these dead but harmless updates happen.
  400% In R however, exported bindings might conflict, e.g. 0 \== 0.0
  401%
  402% If X is indep and we do _not_ solve for it, we are in deep shit
  403% because the ordering is violated.
  404%
  405verify_lin(order(OrdX),CLP,class(Class),lin(LinX),Y) :-
  406    !,
  407    (   indep(CLP,LinX,OrdX)
  408    ->  detach_bounds_vlv(CLP,OrdX,LinX,Class,Y,NewLinX),
  409	% if there were bounds, they are requeued already
  410	class_drop(Class,Y),
  411	nf(-Y, CLP, NfY),
  412	deref(CLP,NfY,LinY),
  413	add_linear_11(CLP, NewLinX, LinY, Lind),
  414	(   nf_coeff_of(Lind,OrdX,_)
  415	->  % X is element of Lind
  416	    solve_ord_x(CLP, Lind, OrdX, Class)
  417	;   solve(CLP, Lind)	% X is gone, can safely solve Lind
  418	)
  419    ;   class_drop(Class,Y),
  420	nf(-Y, CLP, NfY),
  421	deref(CLP,NfY,LinY),
  422	add_linear_11(CLP, LinX, LinY, Lind),
  423	solve(CLP, Lind)
  424    ).
  425verify_lin(_,_,_,_,_)