23
31print_interval(Term) :-
32 printed_term_(Term,Out),
33 format('~w',[Out]). 34
35print_interval(Stream,Term) :-
36 printed_term_(Term,Out),
37 format(Stream,'~w',[Out]). 38
39printed_term_(Term,Out) :-
40 copy_term_nat(Term,Out), 41 term_variables(Term,TVars),
42 term_variables(Out,OVars),
43 name_vars_(TVars,OVars,0). 44
45name_vars_([],[],_).
46name_vars_([TVar|TVars],[OVar|OVars],N) :-
47 (domain(TVar,Dom)
48 -> OVar = (V::Dom)
49 ; OVar = V
50 ),
51 number_chars(N,C),atom_chars(V,['V'|C]),
52 N1 is N+1,
53 name_vars_(TVars,OVars,N1).
54
58attr_portray_hook(interval(Type,Val,_Node,_Flags),_Int) :-
59 interval_domain_(Type,Val,Dom),
60 format('~w',Dom). 61
65user:portray('$clpBNR...'(Out)) :- 66 string(Out),
67 format('~W...',[Out,[quoted(false)]]). 68
74user:expand_answer(Bindings,Bindings) :- 75 current_prolog_flag(clpBNR_verbose,Verbose), 76 add_names_(Bindings,Verbose).
77
79add_names_([],_).
80add_names_([Var|Bindings],Verbose) :- var(Var), !,
81 add_names_(['' = Var|Bindings],Verbose).
82add_names_([Name = Var|Bindings],Verbose) :-
83 (get_interval_flags_(Var,Flags)
84 -> set_interval_flags_(Var,[name(Name,Verbose)|Flags]), 85 (Verbose == false -> reset_interval_nodelist_(Var) ; true) 86 ; (compound(Var) 87 -> term_variables(Var,Vars),
88 add_names_(Vars,Verbose) 89 ; true 90 )
91 ),
92 add_names_(Bindings,Verbose).
93
97:- if(current_prolog_flag(clpBNR_swish, true)). 98
100:- use_module(library(http/html_write)). 101:- multifile(term_html:portray//2). 102
103term_html:portray('$clpBNR...'(Out),_) --> 104 { string(Out) },
105 html(span(class('pl-float'), [Out, '...'])).
106
108:- multifile(swish_trace:post_context/1). 109
110swish_trace:post_context(Dict) :-
111 expand_answer(Dict.bindings,_).
112
113:- endif. % current_prolog_flag(clpBNR_swish, true)
114
115%
116% Constructs goals to build interval(X)
117%
118attribute_goals(X) -->
119 {(interval_object(X, Type, Val, Nodes)
120 -> (get_interval_flags_(X,Flags), memberchk(name(_,false),Flags), var(Nodes)
121 -> intValue_(Val,Type,Dom), % non-verbose and empty nodelist => "pretty" value
122 Goals = [X::Dom]
123 ; interval_domain_(Type, Val, Dom), % full copy
124 constraints_(Nodes,X,Cs), % reverse map to list of original constraints
125 to_comma_exp_(Cs,X::Dom,Goals)
126 )
127 ; Goals = []
128 )
129 },
130 list_dcg_(Goals). 131
132list_dcg_([]) --> [].
133list_dcg_([H|T]) --> [H], list_dcg_(T).
134
135constraints_([Node],_,[]) :- var(Node), !. 136constraints_([node(Op,P,_,Args)|Nodes],X,[C|Cs]) :-
137 var(P), 138 term_variables(Args,[V|_]), 139 V==X,
140 map_constraint_op_(Op,Args,C),
141 constraints_(Nodes,X,Cs).
142constraints_([_|Nodes],X,Cs) :-
143 constraints_(Nodes,X,Cs).
144
145to_comma_exp_([],Decl,[Decl]).
146to_comma_exp_([N|NCs],Decl,[(Decl,{Cs})]) :-
147 to_comma_cexp_(NCs,N,Cs).
148
149to_comma_cexp_([],In,In).
150to_comma_cexp_([C|Cs],In,Out) :-
151 to_comma_cexp_(Cs,(In,C),Out).
152
163intValue_((0,1),integer,boolean). 164intValue_((L,H),real,'$clpBNR...'(Out)) :- 165 zero_float_(L),
166 zero_float_(H), !,
167 format(chars(Zero),"~16f",0.0),
168 string_chars(Out,[' '|Zero]). 169intValue_((L,H),real,'$clpBNR...'(Out)) :- 170 float_chars_(L,LC),
171 float_chars_(H,HC),
172 sign_chars_(LC,HC,ULC,UHC,Codes/Match),
173 leading_zero(ULC,UHC,ZLC),
174 matching_(ZLC,UHC,Match,0,Dec,MLen),
175 MLen>Dec+1, !, 176 string_codes(Out,Codes).
177intValue_((L,H),Type,Dom) :- 178 rational_fraction_(L,FL),
179 rational_fraction_(H,FH),
180 interval_domain_(Type, (FL,FH), Dom).
181
182zero_float_(B) :- float(B), float_class(B,C), (C=zero ; C=subnormal).
183
184float_chars_(B,Cs) :-
185 rational_fraction_(B,F),
186 float(F), float_class(F, normal),
187 format(chars(Cs),'~16f',F), 188 length(Cs,Len), Len=<32. 189
190rational_fraction_(B,FB) :-
191 rational(B,_,D), D \= 1, !, 192 FB is float(B).
193rational_fraction_(F,0.0) :-
194 float(F),
195 float_class(F,subnormal), !.
196rational_fraction_(B,B).
197
198sign_chars_(['-'|LC],['-'|HC],HC,LC,[' ','-'|T]/T) :- !. 199sign_chars_(LC,HC,LC,HC,[' '|T]/T). 200
201leading_zero(['9'|ULC],['1','0'|_],['0','9'|ULC]) :- !.
202leading_zero(ULC,_,ULC).
203
206matching_([C,'.'|LCs],[C,'.'|HCs],[C,'.'|Cs],N,Dec,Nout) :- !, 207 digit_(C,_),
208 Dec is N+1,
209 N1 is N+2,
210 matching_(LCs,HCs,Cs,N1,Dec,Nout).
211matching_([LC,'.',LC1|LCs],[HC,'.',HC1|HCs],[HC,'.'|Cs],N,Dec,Nout) :- !, 212 digit_match_(LC,LC1,HC,HC1),
213 Dec is N+1,
214 N1 is N+2,
215 matching_([LC1|LCs],[HC1|HCs],Cs,N1,Dec,Nout).
216matching_([C|LCs],[C|HCs],[C|Cs],N,Dec,Nout) :- !, 217 digit_(C,_),
218 N1 is N+1,
219 matching_(LCs,HCs,Cs,N1,Dec,Nout).
220matching_([LC,LC1|LCs],[HC,HC1|HCs],[HC|Cs],N,Dec,Nout) :- 221 digit_match_(LC,LC1,HC,HC1), !,
222 N1 is N+1,
223 matching_([LC1|LCs],[HC1|HCs],Cs,N1,Dec,Nout).
224matching_(_LCs,_HCs,[],N,Dec,N) :- 225 nonvar(Dec).
226
227digit_match_(LC,LC1,HC,HC1) :- 228 digit_(LC,L), digit_(LC1,L1), digit_(HC,H), digit_(HC1,H1),
229 H is (L+1) mod 10,
230 L1 >= 5, H1 < 5.
231
233digit_(DC,D) :- atom_number(DC,D), integer(D), 0=<D,D=<9.
234
238enumerate(X) :-
239 interval_object(X, integer, (L,H), _), !,
240 between(L,H,X). 241enumerate(X) :- list(X), !, 242 enumerate_list_(X).
243enumerate(_X). 244
245enumerate_list_([]).
246enumerate_list_([X|Xs]) :-
247 (integer(X) -> true ; enumerate(X)), 248 enumerate_list_(Xs).
249
254small(V) :-
255 current_prolog_flag(clpBNR_default_precision,P),
256 small(V,P).
257
258small(V,P) :-
259 Err is 10.0**(-P),
260 small_(V,Err).
261
262small_(V,Err) :-
263 getValue(V,(L,H)), !,
264 chk_small(L,H,Err).
265small_(L,Err) :-
266 list(L),
267 small_list_(L,Err).
268
269small_list_([],_).
270small_list_([X|Xs],Err) :-
271 small_(X,Err),
272 small_list_(Xs,Err).
273
274chk_small(L,H,Err) :- 275 W is H-L,
276 ( W < Err 277 -> true
278 ; ErrH is Err*sqrt((L**2+H**2)/2), 279 W < ErrH, 280 ErrH \= 1.0Inf 281 ).
282
294global_maximum(Exp,Z) :-
295 global_minimum(-Exp,NZ),
296 constrain_(Z== -NZ).
297global_maximum(Exp,Z,P) :-
298 global_minimum(-Exp,NZ,P),
299 constrain_(Z== -NZ).
300
301global_minimum(Exp,Z) :-
302 current_prolog_flag(clpBNR_default_precision,P),
303 global_minimum(Exp,Z,P).
304global_minimum(Exp,Z,_P) :- ground(Exp), !,
305 Z is Exp.
306global_minimum(Exp,Z,P) :-
307 global_optimum_(Exp,Z,P,false).
308
310global_maximize(Exp,Z) :-
311 global_minimize(-Exp,NZ),
312 constrain_(Z== -NZ).
313global_maximize(Exp,Z,P) :-
314 global_minimize(-Exp,NZ,P),
315 constrain_(Z== -NZ).
316
317global_minimize(Exp,Z) :-
318 current_prolog_flag(clpBNR_default_precision,P),
319 global_minimize(Exp,Z,P).
320global_minimize(Exp,Z,_P) :- ground(Exp), !,
321 Z is Exp.
322global_minimize(Exp,Z,P) :-
323 global_optimum_(Exp,Z,P,true).
324
325global_optimum_(Exp,Z,P,BindVars) :-
326 term_variables(Exp,Xs), 327 {Z==Exp}, 328 box_([Z|Xs],[(Zl,Zh)|XVs]), 329 iterate_MS(Z,Xs,P,Zl-(Zh,XVs),_ZTree,BindVars). 330
331iterate_MS(Z,Xs,P,Zl-(Zh,XVs),ZTree,BindVars) :-
332 continue_MS(Zl,Zh,P,Xs,XVs,False), 333 widest_MS(Xs,XVs,Xf,XfMid), !, 334 eval_MS(False,Z,Xs,XVs,Xf=< ::(XfMid,XfMid),V1), 335 tree_insert(ZTree,V1,ZTree1), 336 eval_MS(False,Z,Xs,XVs,Xf>= ::(XfMid,XfMid),V2), 337 tree_insert(ZTree1,V2,ZTree2), 338 select_min(ZTree2,NxtY,ZTreeY), 339 iterate_MS(Z,Xs,P,NxtY,ZTreeY,BindVars). % Step 13.
340iterate_MS(Z,Xs,_P,Zl-(Zh,XVs),_ZTree,BindVars) :- % termination criteria (Step 12.) met
341 Z::real(Zl,Zh),
342 (BindVars -> optimize_vars_(Xs,XVs) ; true). 343
344optimize_vars_([],[]).
345optimize_vars_([X|Xs],[(Xl,Xh)|XVs]) :-
346 X::real(Xl,Xh),
347 optimize_vars_(Xs,XVs).
348
349continue_MS(Zl,Zh,P,_,_,_) :- 350 Err is 10.0**(-P),
351 \+chk_small(Zl,Zh,Err), 352 !.
353continue_MS(_,_,P,Xs,XVs,_) :- 354 build_box_MS(Xs,XVs,T/T),
355 SErr is 10.0**(-min(6,P+2)), 356 simplesolveall_(Xs,SErr), 357 !, fail. 358continue_MS(_Zl,_Zh,_,_,_XVs,false). 359
361eval_MS(False,_,_,_,_,[]) :- False == false, !. 362eval_MS(_,Z,Xs,XVs,XConstraint,_FV) :- 363 nb_setval('clpBNR:eval_MS',[]), 364 buildConstraint_(XConstraint,T/T,Agenda),
365 build_box_MS(Xs,XVs,Agenda),
366 box_([Z|Xs],[(Zl,Zh)|NXVs]), 367 nb_setval('clpBNR:eval_MS',Zl-(Zh,NXVs)), 368 fail. 369eval_MS(_,_,_,_,_,FV) :-
370 nb_getval('clpBNR:eval_MS',FV). 371
372sandbox:safe_global_variable('clpBNR:eval_MS').
373
374build_box_MS([],[],Agenda) :-
375 stable_(Agenda).
376build_box_MS([X|Xs],[XNew|XVs],Agenda) :-
377 getValue(X,XCurrent),
378 updateValue_(XCurrent,XNew,X,1,Agenda,NextAgenda), !, 379 build_box_MS(Xs,XVs,NextAgenda).
380
381tree_insert(Tree,[],Tree) :-!.
382tree_insert(MT,Data,Ntree) :- var(MT), !,
383 Ntree = n(_L,_R,Data).
384tree_insert(n(L,R,K-Data), NK-NData, Ntree) :- -1 is cmpr(NK,K), !, 385 tree_insert(L, NK-NData, NewL),
386 Ntree = n(NewL,R,K-Data).
387tree_insert(n(L,R,K-Data), NK-NData, Ntree) :- 1 is cmpr(NK,K), !, 388 tree_insert(R, NK-NData, NewR),
389 Ntree = n(L,NewR,K-Data).
390tree_insert(n(L,R,K-(U,Data)), K-(NU,NData), Ntree) :- -1 is cmpr(NU,U), !, 391 tree_insert(L, K-(NU,NData), NewL),
392 Ntree = n(NewL,R,K-(U,Data)).
393tree_insert(n(L,R,Data), NData, n(L,NewR,Data)) :- 394 Data = NData
395 -> NewR = R 396 ; tree_insert(R, NData, NewR). 397
398select_min(n(L,R,Min),Min,R) :- var(L), !,
399 nonvar(Min). 400select_min(n(L,R,KData),Min,n(NewL,R,KData)) :-
401 select_min(L,Min,NewL).
402
404box_([],[]).
405box_([X|Xs],[R|Rs]) :-
406 getValue(X,R),
407 box_(Xs,Rs).
408
410widest_MS([X|Xs],[XV|XVs],Xf,XfMid) :-
411 widest1_MS(Xs,XVs,XV,X,Xf,XfMid).
412
413widest1_MS([],[],(L,H),Xf,Xf,XfMid) :-
414 midpoint_(L,H,XfMid), !, 415 -2 is cmpr(L,XfMid) + cmpr(XfMid,H). 416widest1_MS([X|Xs],[(L,H)|XVs],(L0,H0),_X0,Xf,XfMid) :-
417 1 is cmpr((H-L),(H0-L0)), !, 418 widest1_MS(Xs,XVs,(L,H),X,Xf,XfMid).
419widest1_MS([_X|Xs],[_XV|XVs],W0,X0,Xf,XfMid) :-
420 widest1_MS(Xs,XVs,W0,X0,Xf,XfMid).
421
426
431splitsolve(X) :-
432 current_prolog_flag(clpBNR_default_precision,P),
433 splitsolve(X,P).
434
435splitsolve(X,_P) :-
436 number(X), !. 437splitsolve(X,P) :-
438 interval(X), !, 439 Err is 10.0**(-P),
440 simplesolveall_([X],Err).
441splitsolve(X,P) :- 442 flatten_list(X,[],XF), 443 Err is 10.0**(-P),
444 simplesolveall_(XF,Err).
445
447flatten_list(V,Tail,[V|Tail]) :- var(V), !.
448flatten_list([],Tail,Tail) :- !.
449flatten_list([H|T], Tail, List) :- !,
450 flatten_list(H, HTail, List),
451 flatten_list(T, Tail, HTail).
452flatten_list(N,Tail,[N|Tail]).
453
454simplesolveall_(Xs,Err) :-
455 select_wide_(Xs,_,X),
456 interval_object(X, Type, V, _),
457 split_choices_(Type,X,V,Err,Choices),
458 !,
459 xpsolve_choice(Choices), 460 simplesolveall_(Xs,Err).
461simplesolveall_(_Xs,_Err). 462
463select_wide_([X],_,X) :- !. 464select_wide_([X1,X2|Xs],D1,X) :- 465 (var(D1) -> delta(X1,D1) ; true),
466 delta(X2,D2),
467 (D1 >= D2
468 -> select_wide_([X1|Xs],D1,X)
469 ; select_wide_([X2|Xs],D2,X)
470 ).
471
472split_choices_(real,X,V,Err,split(X,::(SPt,SPt))) :- !,
473 splitinterval_real_(V,SPt,Err). 474split_choices_(integer,X,V,_Err,Cons) :-
475 splitinterval_(integer,X,V,_,Cons). 476
497
498absolve( X ):-
499 current_prolog_flag(clpBNR_default_precision,P),
500 absolve(X,P),!.
501
502absolve(X, _ ):- number(X), !.
503absolve(X, Limit ):- interval_object(X, Type, Val, _), !, 504 delta_(Type,Val,Delta),
505 506 (not(not(lower_bound(X))) -> true ; absolve_l(X,Type,Delta,1,Limit)),
507 (not(not(upper_bound(X))) -> true ; absolve_r(X,Type,Delta,1,Limit)).
508
509absolve([],_). 510absolve([X|Xs],Lim):- absolve(X,Lim),!, absolve(Xs,Lim).
511
512delta_(integer,(L,U),D) :- D is U div 2 - L div 2.
513delta_(real, (L,U),D) :- D is U/2 - L/2.
514
515absolve_l(X, Type, DL, NL, Limit):- NL<Limit, 516 getValue(X,(LB1,UB1)),
517 trim_point_(NL,NL1,Type,Limit,DL,DL1), 518 Split is LB1 + DL1,
519 LB1 < Split, Split < UB1, 520 not(constrain_(X=< ::(Split,Split))),!,
521 constrain_(X>= ::(Split,Split)), 522 absolve_l(X,Type, DL1, NL1, Limit).
523absolve_l(_,_,_,_,_). 524
525absolve_r(X, Type, DU, NU, Limit):- NU<Limit, 526 getValue(X,(LB1,UB1)),
527 trim_point_(NU,NU1,Type,Limit,DU,DU1), 528 Split is UB1 - DU1,
529 LB1 < Split, Split < UB1, 530 not(constrain_(X>= ::(Split,Split))),!,
531 constrain_(X=< ::(Split,Split)), 532 absolve_r(X,Type, DU1, NU1,Limit).
533absolve_r(_,_,_,_,_). 534
535trim_point_( N,N, _Type, _Limit, Delta, Delta).
536trim_point_( N,M, integer, Limit, Delta, Result):- N<Limit,N1 is N + 1,
537 D is Delta div 2,
538 trim_point_(N1,M, integer, Limit,D, Result).
539trim_point_( N,M, real, Limit, Delta, Result):- N<Limit,N1 is N + 1,
540 D is Delta/2,
541 trim_point_(N1,M,real, Limit,D, Result).
542
546solve(X) :-
547 current_prolog_flag(clpBNR_default_precision,P),
548 solve(X,P).
549
550solve(X,P) :-
551 interval(X), !, 552 solve([X],P).
553solve(X,_P) :-
554 number(X), !. 555solve(X,P) :- 556 Err is 10.0**(-(P+1)), 557 xpsolveall_(X,Err).
558
559xpsolveall_([],_Err) :- !.
560xpsolveall_(Xs,Err) :-
561 xpsolve_each_(Xs,Us,Err), 562 xpsolveall_(Us,Err). 563
564xpsolve_each_([],[],_Err).
565xpsolve_each_([X|Xs],[X|Us],Err) :-
566 interval_object(X,Type,V,_), 567 splitinterval_(Type,X,V,Err,Choices), 568 !,
569 xpsolve_choice(Choices), 570 xpsolve_each_(Xs,Us,Err). 571xpsolve_each_([X|Xs],Us,Err) :- 572 X==[], !, 573 xpsolve_each_(Xs,Us,Err). 574xpsolve_each_([X|Xs],[U|Us],Err) :-
575 list(X), !, 576 xpsolve_each_(X,U,Err), 577 xpsolve_each_(Xs,Us,Err). 578xpsolve_each_([_X|Xs],Us,Err) :-
579 xpsolve_each_(Xs,Us,Err). 580
581xpsolve_choice(split(X,SPt)) :- constrain_(X =< SPt). 582xpsolve_choice(split(X,SPt)) :- constrain_(SPt =< X).
583xpsolve_choice(split_integer(X,SPt)) :- constrain_(X =< SPt).
584xpsolve_choice(split_integer(X,SPt)) :- constrain_(SPt < X).
585xpsolve_choice(enumerate(X)) :- enumerate(X).
586
590splitinterval_(real,X,V,Err,split(X,::(SPt,SPt))) :- 591 splitinterval_real_(V,Pt,Err), 592 split_real_(X,V,Pt,Err,SPt).
593splitinterval_(integer,X,V,_,Cons) :-
594 V = (L,H),
595 ( H-L < 15
596 -> Cons = enumerate(X) 597 ; splitinterval_integer_(V,Pt), 598 Cons = split_integer(X,Pt) 599 ).
602
604split_real_(X,_,Pt,_,Pt) :- 605 X = Pt -> fail ; !. 606split_real_(X,(L,_H),Pt,Err,NPt) :- 607 split_real_lo(X,L,Pt,NPt,Err), !.
608split_real_(X,(_L,H),Pt,Err,NPt) :- 609 split_real_hi(X,Pt,H,NPt,Err).
610
611split_real_lo(X,L,Pt,NPt,Err) :- 612 splitinterval_real_((L,Pt),SPt,Err),
613 (X\=SPt -> NPt=SPt ; split_real_lo(X,L,SPt,NPt,Err)).
614
615split_real_hi(X,Pt,H,NPt,Err) :- 616 splitinterval_real_((Pt,H),SPt,Err),
617 (X\=SPt -> NPt=SPt ; split_real_hi(X,SPt,H,NPt,Err)).
618
622splitinterval_integer_((L,H),0) :-
623 L < 0, H > 0, !. 624splitinterval_integer_((-1.0Inf,H),Pt) :- !, 625 finite_interval(integer, (Pt,_)), 626 H > Pt.
627splitinterval_integer_((L,1.0Inf), Pt) :- !, 628 finite_interval(integer, (_,Pt)), 629 L < Pt.
630splitinterval_integer_((L,H),Pt) :- 631 Pt is (L div 2) + (H div 2). 632
633splitinterval_real_((L,H),0,E) :- 634 L < 0, H > 0, !, 635 (H-L) > E. 636splitinterval_real_((-1.0Inf,H),Pt,_) :- 637 !, 638 Pt is float(H*10-1), 639 Pt > -1.0Inf. 640splitinterval_real_((L,1.0Inf),Pt,_) :- 641 !, 642 Pt is float(L*10+1), 643 Pt < 1.0Inf. 644splitinterval_real_((L,H),Pt,E) :- 645 \+ chk_small(L,H,E), 646 splitMean_(L,H,Pt), !,
647 L < Pt,Pt < H. 648
651
653splitMean_(L,H,M) :- L >=0, !, 654 (L=:=0 -> M is min(H/2, sqrt( H)) ; M is sqrt(L)*sqrt(H)).
655splitMean_(L,H,M) :- 656 (H=:=0 -> M is max(L/2,-sqrt(-L)) ; M is -sqrt(-L)*sqrt(-H)).
657
658
671partial_derivative(F,X,DFX) :-
672 pd(F,X,DD),
673 simplify(DD,DFX).
674 675 676
677pd(F,X,DD) :-
678 ( var(F)
679 -> ( F==X -> DD=1 ; DD=0)
680 ; ( atomic(F)
681 -> DD=0
682 ; pd_f(F,X,DD)
683 )
684 ).
685
686:- style_check(-singleton). 687
688pd_f(-U,X,DX) :- 689 pd(U,X,DU),
690 pd_minus(DU,DX).
691
692pd_f(U+V,X,DX) :- 693 pd(U,X,DU), pd(V,X,DV),
694 pd_add(DU,DV,DX).
695
696pd_f(U-V,X,DX) :- 697 pd(U,X,DU), pd(V,X,DV),
698 pd_sub(DU,DV,DX).
699
700pd_f(U*V,X,DX) :- 701 pd(U,X,DU), pd(V,X,DV),
702 pd_mul(V,DU,VDU),
703 pd_mul(U,DV,UDV),
704 pd_add(VDU,UDV,DX).
705
706pd_f(U/V,X,DX) :- 707 pd(U,X,DU), pd(V,X,DV),
708 pd_mul(V,DU,VDU),
709 pd_mul(U,DV,UDV),
710 pd_sub(VDU,UDV,DXN),
711 pd_pow(V,2,DXD),
712 pd_div(DXN,DXD,DX).
713
714pd_f(U**N,X,DX) :- 715 pd(U,X,DU),
716 pd_sub(N,1,N1), 717 pd_pow(U,N1,UN1),
718 pd_mul(N,UN1,DX1),
719 pd_mul(DX1,DU,DX).
720
721pd_f(log(U),X,DX) :- 722 pd(U,X,DU),
723 pd_div(DU,U,DX).
724
725pd_f(exp(U),X,DX) :- 726 pd(U,X,DU),
727 pd_exp(U,ExpU),
728 pd_mul(ExpU,DU,DX).
729
730pd_f(sqrt(U),X,DX) :- 731 pd(U,X,DU),
732 pd_sqrt(U,SqrtU),
733 pd_mul(2,SqrtU,DXD),
734 pd_div(DU,DXD,DX).
735
736pd_f(sin(U),X,DX) :- 737 pd(U,X,DU),
738 pd_cos(U,CosU),
739 pd_mul(CosU,DU,DX).
740
741pd_f(cos(U),X,DX) :- 742 pd(U,X,DU),
743 pd_sin(U,SinU),
744 pd_mul(SinU,DU,DSU),
745 pd_minus(DSU,DX).
746
747pd_f(tan(U),X,DX) :- 748 pd(U,X,DU),
749 pd_cos(U,CosU),
750 pd_pow(CosU,2,DXD),
751 pd_div(DU,DXD,DX).
752
753pd_f(asin(U),X,DX) :- 754 pd(U,X,DU),
755 pd_pow(U,2,USQ),
756 pd_sub(1,USQ,USQ1),
757 pd_sqrt(USQ1,DXD),
758 pd_div(DU,DXD,DX).
759
760pd_f(acos(U),X,DX) :- 761 pd(asin(U),X,DasinX),
762 pd_minus(DasinX,DX).
763
764pd_f(atan(U),X,DX) :- 765 pd(U,X,DU),
766 pd_pow(U,2,USQ),
767 pd_add(1,USQ,USQ1),
768 pd_div(DU,USQ1,DX).
769
770
773pd_minus(DU,DX) :- ground(DU) -> DX is -DU ; DX = -DU.
774
775pd_add(DU,DV,DV) :- DU==0, !.
776pd_add(DU,DV,DU) :- DV==0, !.
777pd_add(DU,DV,DX) :- ground((DU,DV)) -> DX is DU+DV ; DX = DU+DV.
778
779pd_sub(DU,DV,DX) :- DU==0, !, pd_minus(DV,DX).
780pd_sub(DU,DV,DU) :- DV==0, !.
781pd_sub(DU,DV,DX) :- ground((DU,DV)) -> DX is DU-DV ; DX = DU-DV.
782
783pd_mul(DU,DV,0) :- DU==0, !.
784pd_mul(DU,DV,0) :- DV==0, !.
785pd_mul(DU,DV,DV) :- DU==1, !.
786pd_mul(DU,DV,DU) :- DV==1, !.
787pd_mul(DU,DV,DX) :- DU== -1, !, pd_minus(DV,DX).
788pd_mul(DU,DV,DX) :- DV== -1, !, pd_minus(DU,DX).
789pd_mul(DU,DV,DX) :- ground((DU,DV)) -> DX is DU*DV ; DX = DU*DV.
790
791pd_div(DU,DV,0) :- DU==0, !.
792pd_div(DU,DV,0) :- DV==0, !, fail.
793pd_div(DU,DV,DU) :- DV==1, !.
794pd_div(DU,DV,DU) :- DV== -1, !, pd_minus(DU,DX).
795pd_div(DU,DV,DX) :- ground((DU,DV)) -> DX is DU/DV ; DX = DU/DV.
796
797pd_pow(DU,DV,0) :- DU==0, !.
798pd_pow(DU,DV,1) :- DV==0, !.
799pd_pow(DU,DV,1) :- DU==1, !.
800pd_pow(DU,DV,DU) :- DV==1, !.
801pd_pow(DU,DV,DX) :- ground((DU,DV)) -> DX is DU**DV ; DX = DU**DV.
802
803pd_exp(DU,DX) :- ground(DU) -> DX is exp(DU) ; DX = exp(DU).
804
805pd_sqrt(DU,DX) :- ground(DU) -> DX is sqrt(DU) ; DX = sqrt(DU).
806
807pd_sin(DU,DX) :- ground(DU) -> DX is sin(DU) ; DX = sin(DU).
808
809pd_cos(DU,DX) :- ground(DU) -> DX is cos(DU) ; DX = cos(DU).
810
811:- style_check(+singleton).