| 1 | | % (c) 2020-2026 Lehrstuhl fuer Softwaretechnik und Programmiersprachen, |
| 2 | | % Heinrich Heine Universitaet Duesseldorf |
| 3 | | % This software is licenced under EPL 1.0 (http://www.eclipse.org/org/documents/epl-v10.html) |
| 4 | | |
| 5 | | :- module(kernel_reals, [construct_real/2, construct_negative_real/2, |
| 6 | | construct_real_number/2, |
| 7 | | is_real/1, is_real/2, is_real_wf/2, |
| 8 | | is_float/1, is_float_wf/2, |
| 9 | | is_not_real/1, is_not_float/1, |
| 10 | | is_ground_real/1, |
| 11 | | is_largest_positive_float/1, is_smallest_positive_float/1, |
| 12 | | is_next_larger_float/2, is_next_smaller_float/2, |
| 13 | | convert_int_to_real/2, |
| 14 | | real_floor/2, real_ceiling/2, real_truncate/2, |
| 15 | | real_addition_wf/4, real_subtraction_wf/4, |
| 16 | | real_multiplication_wf/4, real_division_wf/5, |
| 17 | | real_unary_minus_wf/3, |
| 18 | | real_absolute_value_wf/3, |
| 19 | | real_sign_wf/3, |
| 20 | | real_square_root_wf/4, |
| 21 | | real_round_wf/3, |
| 22 | | real_unop_wf/4, real_unop_wf/5, |
| 23 | | real_binop_wf/6, |
| 24 | | real_power_of_wf/5, |
| 25 | | real_less_than_wf/3, real_less_than_equal_wf/3, |
| 26 | | real_comp_wf/5, |
| 27 | | real_maximum_of_set/4, real_minimum_of_set/4, |
| 28 | | enumerate_real_wf/3, |
| 29 | | use_clpfd_real_solver/0, do_not_double_check_solution/0, |
| 30 | | post_real_equal_expr_wf/4 |
| 31 | | ]). |
| 32 | | |
| 33 | | % Reals/Floats in ProB are represented by terms of the form term(floating(Number)) |
| 34 | | % in future a proper wrapper such as real(Number) might be created |
| 35 | | |
| 36 | | :- meta_predicate real_comp_wf(2,-,-,-,-). |
| 37 | | :- meta_predicate real_comp_float_blocking(2,-,-,-). |
| 38 | | |
| 39 | | :- use_module(module_information,[module_info/2]). |
| 40 | | :- module_info(group,kernel). |
| 41 | | :- module_info(description,'This module provides (external) functions to manipulate B reals and floats.'). |
| 42 | | |
| 43 | | :- use_module(error_manager). |
| 44 | | :- use_module(self_check). |
| 45 | | :- use_module(library(lists)). |
| 46 | | :- use_module(debug). |
| 47 | | |
| 48 | | :- use_module(kernel_objects,[exhaustive_kernel_check_wf/2,exhaustive_kernel_check_wf/3]). |
| 49 | | :- use_module(extension('counter/counter'), |
| 50 | | [next_smaller_abs_float/2, next_smaller_float/2,next_larger_float/2, next_float_twoards/3, |
| 51 | | smallest_float/1, largest_float/1, smallest_abs_float/1]). |
| 52 | | |
| 53 | | |
| 54 | | :- if((current_prolog_flag(dialect, sicstus), |
| 55 | | current_prolog_flag(version_data, sicstus(4,VN,_,_,_)), VN>=10)). % 4.10 or higher |
| 56 | | can_use_clpfd_real_solver. |
| 57 | | use_clpfd_real_solver :- get_preference(solver_for_reals,S), (S=precise_float_solver ; S=real_solver). |
| 58 | | :- use_module(library(clpfd),[fd_min/2, fd_max/2, fd_batch/1, maximum/2, minimum/2, |
| 59 | | '#<=>'/2, '$>='/2, '$=<'/2, '$='/2, '#='/2, labeling/2]). |
| 60 | | |
| 61 | | post_real_equal_expr_wf(X,Y,Span,WF) :- |
| 62 | | catch('$='(X,Y), |
| 63 | | error(evaluation_error(ERR),_), |
| 64 | | process_evaluation_error(ERR,'='(X,Y),Span,WF)). |
| 65 | | :- else. |
| 66 | | can_use_clpfd_real_solver :- fail. |
| 67 | | use_clpfd_real_solver :- fail. |
| 68 | | post_real_equal_expr_wf(X,Y,Span,WF) :- add_internal_error('Not available:',post_real_equal_expr_wf(X,Y,Span,WF)),fail. |
| 69 | | :- endif. |
| 70 | | |
| 71 | | |
| 72 | | % used to construct a Value from real(Atom) AST node |
| 73 | | construct_real(Atom,term(floating(Float))) :- |
| 74 | | atom_codes(Atom,C), |
| 75 | | number_codes(Nr,C), |
| 76 | | Float is float(Nr). % make sure we store a float; not required for AST literals |
| 77 | | |
| 78 | | construct_real_number(Atom,Nr) :- |
| 79 | | atom_codes(Atom,C), |
| 80 | | number_codes(Nr,C). |
| 81 | | |
| 82 | | construct_negative_real(Atom,term(floating(Float))) :- |
| 83 | | atom_codes(Atom,C), |
| 84 | | number_codes(Nr,C), |
| 85 | | Float is -float(Nr). |
| 86 | | |
| 87 | | is_real_wf(X,_WF) :- is_real(X,_). |
| 88 | | is_real(X) :- is_real(X,_). |
| 89 | | |
| 90 | | is_real(term(floating(Nr)),Nr). |
| 91 | | % TO DO: other formats |
| 92 | | |
| 93 | | is_float_wf(X,_WF) :- is_float(X). |
| 94 | | |
| 95 | | is_float(term(floating(_))). |
| 96 | | |
| 97 | | is_not_real(_) :- fail. |
| 98 | | is_not_float(_) :- fail. |
| 99 | | |
| 100 | | is_ground_real(term(floating(Nr))) :- number(Nr). |
| 101 | | |
| 102 | | is_largest_positive_float(term(floating(Nr))) :- largest_float(Nr). |
| 103 | | is_smallest_positive_float(term(floating(Nr))) :- smallest_abs_float(Nr). |
| 104 | | |
| 105 | | is_next_larger_float(term(floating(Nr)),term(floating(Next))) :- |
| 106 | | block_next_larger_float(Nr,Next). |
| 107 | | :- block block_next_larger_float(-,-). |
| 108 | | block_next_larger_float(Nr,Next) :- nonvar(Nr),!,next_larger_float(Nr,Next). |
| 109 | | block_next_larger_float(Nr,Next) :- next_smaller_float(Next,Nr). |
| 110 | | |
| 111 | | is_next_smaller_float(term(floating(Nr)),term(floating(Next))) :- |
| 112 | | block_next_larger_float(Next,Nr). |
| 113 | | |
| 114 | | |
| 115 | | % convert an integer to a real, corresponds to the real(.) function in Atelier-B; convert_real as AST |
| 116 | | convert_int_to_real(int(X),term(floating(R))) :- |
| 117 | | convert_int_to_real_aux(X,R). |
| 118 | | |
| 119 | | convert_int_to_real_aux(Int,Real) :- number(Int),!, |
| 120 | | Real is float(Int). |
| 121 | | convert_int_to_real_aux(Int,Real) :- |
| 122 | | use_clpfd_real_solver,!, |
| 123 | | '$='(Real,float(Int)). |
| 124 | | convert_int_to_real_aux(Int,Real) :- block_convert_int_to_real_aux(Int,Real). |
| 125 | | |
| 126 | | :- block block_convert_int_to_real_aux(-,-). |
| 127 | | block_convert_int_to_real_aux(Int,Real) :- number(Int),!, |
| 128 | | Real is float(Int). |
| 129 | | block_convert_int_to_real_aux(Int,Real) :- Int1 is floor(Real), |
| 130 | | Real is float(Int1), |
| 131 | | Int1=Int. |
| 132 | | |
| 133 | | |
| 134 | | % convert_int_floor as AST, Atelier-B floor(.) |
| 135 | | % The value is the greatest integer less or equal to X |
| 136 | | real_floor(term(floating(R)),int(X)) :- real_floor_aux(R,X). |
| 137 | | |
| 138 | | real_floor_aux(Real,Int) :- |
| 139 | | use_clpfd_real_solver,!, % + check if we use CLP(FD) for integers? |
| 140 | | '#='(Int,floor(Real)). |
| 141 | | real_floor_aux(Real,Int) :- block_real_floor_aux(Real,Int). |
| 142 | | |
| 143 | | :- block block_real_floor_aux(-,?). |
| 144 | | block_real_floor_aux(Real,Int) :- Int is floor(Real). |
| 145 | | |
| 146 | | |
| 147 | | % convert_int_ceiling as AST, Atelier-B ceiling(.) |
| 148 | | % The value is the least integer greater or equal to X. |
| 149 | | real_ceiling(term(floating(R)),int(X)) :- |
| 150 | | real_ceiling_aux(R,X). |
| 151 | | real_ceiling_aux(Real,Int) :- |
| 152 | | use_clpfd_real_solver,!, % + check if we use CLP(FD) for integers? |
| 153 | | '#='(Int,ceiling(Real)). |
| 154 | | real_ceiling_aux(Real,Int) :- block_real_ceiling_aux(Real,Int). |
| 155 | | |
| 156 | | :- block block_real_ceiling_aux(-,?). |
| 157 | | block_real_ceiling_aux(Real,Int) :- Int is ceiling(Real). |
| 158 | | |
| 159 | | |
| 160 | | % The value is the closest integer between X and 0 |
| 161 | | real_truncate(term(floating(R)),int(X)) :- |
| 162 | | real_truncate_aux(R,X). |
| 163 | | real_truncate_aux(Real,Int) :- |
| 164 | | use_clpfd_real_solver,!, % + check if we use CLP(FD) for integers? |
| 165 | | '#='(Int,truncate(Real)). |
| 166 | | real_truncate_aux(Real,Int) :- block_real_truncate_aux(Real,Int). |
| 167 | | |
| 168 | | :- block block_real_truncate_aux(-,?). |
| 169 | | block_real_truncate_aux(Real,Int) :- Int is truncate(Real). |
| 170 | | |
| 171 | | |
| 172 | | |
| 173 | | :- assert_must_succeed(exhaustive_kernel_check_wf([commutative],kernel_reals:real_addition_wf(term(floating(1.0)),term(floating(2.0)),term(floating(3.0)),WF),WF)). |
| 174 | | :- assert_must_succeed(exhaustive_kernel_check_wf([commutative],kernel_reals:real_addition_wf(term(floating(0.0)),term(floating(2.0)),term(floating(2.0)),WF),WF)). |
| 175 | | :- assert_must_succeed(exhaustive_kernel_check_wf([commutative],kernel_reals:real_addition_wf(term(floating(-1.0)),term(floating(1.0)),term(floating(0.0)),WF),WF)). |
| 176 | | |
| 177 | | real_addition_wf(term(floating(X)),term(floating(Y)),term(floating(R)),WF) :- |
| 178 | | real_add_wf_aux(X,Y,R,WF). |
| 179 | | |
| 180 | | %real_add_wf_aux(X,Y,R,WF) :- X==Y,!, '$='(R,2.0*X), % TODO: check if this works reliably, probably not |
| 181 | | % (do_not_double_check_solution -> true ; real_add_wf_aux2(X,Y,R,WF)). |
| 182 | | real_add_wf_aux(X,Y,R,WF) :- |
| 183 | | use_clpfd_real_solver,!, |
| 184 | | '$='(R,X+Y), |
| 185 | | (do_not_double_check_solution -> true ; real_add_wf_aux2(X,Y,R,WF)). |
| 186 | | real_add_wf_aux(X,Y,R,WF) :- block_real_add_wf_aux(X,Y,R,WF). |
| 187 | | |
| 188 | | :- block block_real_add_wf_aux(-,-,?,?), block_real_add_wf_aux(?,-,-,?), block_real_add_wf_aux(-,?,-,?). |
| 189 | | block_real_add_wf_aux(X,Y,R,WF) :- |
| 190 | | (nonvar(X) |
| 191 | | -> (nonvar(Y) |
| 192 | | -> safe_is(R,X+Y,WF) |
| 193 | | ; allow_constraint_propagation(Kind), |
| 194 | | safe_is(YY,R-X,WF), |
| 195 | | confirm_solution_for_expr(YY,X+SOL,SOL,R,Kind,Solutions,WF) |
| 196 | | -> set_solution_for_expr(Y,Solutions,WF) |
| 197 | | ; real_add_wf_aux2(X,Y,R,WF) % delay until Y is known |
| 198 | | ) |
| 199 | | ; allow_constraint_propagation(Kind), |
| 200 | | safe_is(XX,R-Y,WF), |
| 201 | | confirm_solution_for_expr(XX,SOL+Y,SOL,R,Kind,Solutions,WF) |
| 202 | | -> set_solution_for_expr(X,Solutions,WF) |
| 203 | | ; real_add_wf_aux2(X,Y,R,WF) % delay until X is known |
| 204 | | ). |
| 205 | | :- block real_add_wf_aux2(-,?,?,?), real_add_wf_aux2(?,-,?,?). |
| 206 | | real_add_wf_aux2(X,Y,R,WF) :- safe_is(R,X+Y,WF). |
| 207 | | |
| 208 | | |
| 209 | | % confirm_solution_for_expr(XX,ExprX,X,R : check if XX=X is a solution for R = ExprX, where ExprX contains variable X |
| 210 | | % depending on solver setting : compute a list of candidate solutions around XX |
| 211 | | confirm_solution_for_expr(XX,_,_,_,no_checking,Solutions,_) :- !, Solutions=XX. |
| 212 | | confirm_solution_for_expr(XX,ExprX,X,R,simple_checking,Solutions,WF) :- !, |
| 213 | | check_sol_direct(XX,ExprX,X,R,WF), Solutions=XX. % single solution |
| 214 | | confirm_solution_for_expr(XX,ExprX,X,R,check_unique_solution,Solutions,WF) :- !, |
| 215 | | check_sol_direct(XX,ExprX,X,R,WF), Solutions=XX, % single solution |
| 216 | | % Now check if there are other solutions: |
| 217 | | next_larger_float(XX,XN), |
| 218 | | (X=XN,safe_is(R,ExprX,WF) |
| 219 | | -> debug_format(19,'Reject ambiguous larger solution {~w,~w} for ~w = ~w~n',[XX,XN,ExprX,R]),fail |
| 220 | | ; true |
| 221 | | ), |
| 222 | | next_smaller_float(XX,XP), |
| 223 | | (X=XP,safe_is(R,ExprX,WF) |
| 224 | | -> debug_format(19,'Reject ambiguous smaller solution {~w,~w} for ~w = ~w~n',[XX,XP,ExprX,R]),fail |
| 225 | | ; true). |
| 226 | | confirm_solution_for_expr(XX,ExprX,X,R,check_small_interval_solution,Interval,WF) :- |
| 227 | | get_solutions_around(XX,ExprX,X,R,Interval,WF). |
| 228 | | |
| 229 | | :- use_module(kernel_waitflags,[get_wait_flag/4]). |
| 230 | | % set the solution based on result of confirm_solution_for_expr |
| 231 | | set_solution_for_expr(X,XX,_WF) :- number(XX),!,X=XX. % deterministic constraint propagation found precise solution |
| 232 | | set_solution_for_expr(X,[XX],_WF) :- !,X=XX. % ditto |
| 233 | | set_solution_for_expr(_,[],_WF) :- !, fail. % no solutions |
| 234 | | set_solution_for_expr(X,Interval,WF) :- |
| 235 | | length(Interval,Len), % we found an interval list of possible solutions |
| 236 | | get_wait_flag(Len,kernel_reals,WF,LWF), |
| 237 | | set_sol_aux(X,Interval,LWF). |
| 238 | | :- block set_sol_aux(-,?,-). |
| 239 | | set_sol_aux(X,Interval,_) :- |
| 240 | ? | member(X,Interval). |
| 241 | | |
| 242 | | |
| 243 | | % get a sorted list of solutions around the candidate XX |
| 244 | | get_solutions_around(XX,ExprX,X,R,Interval,WF) :- |
| 245 | | (check_sol_direct(XX,ExprX,X,R,WF) % XX is a solution |
| 246 | | -> get_preference(solver_strength,SS), |
| 247 | | Limit is 10+SS, % number of solutions we are willing to find |
| 248 | | next_larger_float(XX,XN), |
| 249 | | find_larger_bounds(XN,Limit,ExprX,X,R,WF,LargerSols,[]), |
| 250 | | next_smaller_float(XX,XP), |
| 251 | | find_lower_bounds(XP,Limit,ExprX,X,R,WF,Interval,[XX|LargerSols]) |
| 252 | | ; (X='X',debug_format(19,'Reject non-solution ~w for ~w = ~w~n',[XX,ExprX,R]),fail ; true), |
| 253 | | %next_smaller_float(XX,XP), (check_sol_direct(XP,ExprX,X,R,WF) -> write(xp(XP)),nl ; write(nxp(XP)),nl), |
| 254 | | %next_larger_float(XX,XN), (check_sol_direct(XN,ExprX,X,R,WF) -> write(xn(XN)),nl ; write(nxp(XN)),nl), |
| 255 | | Interval = [] % no floating point solutions; TODO: do we need to check smaller/larger values?? |
| 256 | | ). |
| 257 | | |
| 258 | | % e.g. for x+1.0 = 3.0 we have 2.0 as solution and the previous float 1.9999999999999998, nothing else |
| 259 | | % e.b. for x+2.0 = 3.0 we have four solutions {0.9999999999999998,0.9999999999999999,1.0,1.0000000000000002} |
| 260 | | % card({x|x + 2.999 = 3.0}) = 2049 |
| 261 | | % card({x|x + 2.9999 = 3.0}) = 32769 |
| 262 | | |
| 263 | | % increment floats until we hit the first non-solution or we reach the limit |
| 264 | | % TODO: try and compute this bound, rather then single-stepping |
| 265 | | % TODO: enumerate the first solution and then generate an enum warning if we reach limit |
| 266 | | find_larger_bounds(XP,Lim,ExprX,X,R,WF,[XP|Sols1],Sols2) :- |
| 267 | | check_sol_direct(XP,ExprX,X,R,WF),!, Lim>0, L1 is Lim-1, |
| 268 | | next_larger_float(XP,XP2), |
| 269 | | find_larger_bounds(XP2,L1,ExprX,X,R,WF,Sols1,Sols2). |
| 270 | | find_larger_bounds(_Bound,_Lim,_ExprX,_X,_R,_WF,Sols,Sols). |
| 271 | | |
| 272 | | % decrement floats until we hit the first non-solution or we reach the limit |
| 273 | | find_lower_bounds(XP,Lim,ExprX,X,R,WF,[XP|Sols1],Sols2) :- |
| 274 | | check_sol_direct(XP,ExprX,X,R,WF),!, Lim>0, L1 is Lim-1, |
| 275 | | next_smaller_float(XP,XP2), |
| 276 | | find_lower_bounds(XP2,L1,ExprX,X,R,WF,Sols1,Sols2). |
| 277 | | find_lower_bounds(_Bound,_Lim,_ExprX,_X,_R,_WF,Sols,Sols). |
| 278 | | |
| 279 | | check_sol_direct(XX,ExprX,X,R,WF) :- !, |
| 280 | | % (x + y = z <=> x = z - y) does not hold for all floats, e.g., when x very small and y very large (x+y)-y = 0 |
| 281 | | (\+((X=XX,safe_is(R,ExprX,WF))) -> fail ; true). |
| 282 | | |
| 283 | | |
| 284 | | :- use_module(probsrc(preferences), [get_preference/2]). |
| 285 | | |
| 286 | | allow_constraint_propagation(Kind) :- |
| 287 | | get_preference(solver_for_reals,Solver), |
| 288 | | get_solver_kind(Solver,Kind). |
| 289 | | |
| 290 | | get_solver_kind(real_solver,no_checking). % like CLP(R): do not check propagated solution |
| 291 | | get_solver_kind(float_solver,simple_checking). % check propagated solution, but do not check uniqueness |
| 292 | | %get_solver_kind(precise_float_solver,check_unique_solution). |
| 293 | | get_solver_kind(precise_float_solver,check_small_interval_solution). % perform limited constraint propagation, fail directly if no solution exists |
| 294 | | |
| 295 | | do_not_double_check_solution :- get_preference(solver_for_reals,real_solver). |
| 296 | | |
| 297 | | % we could have various other options: |
| 298 | | % aggressive: do it like CLP(R), which can generate non-solutions |
| 299 | | % | ?- {X+10.0e10 = 1.0e-9}, write(sol(X)),nl, {X+10.0e10 = 1.0e-9}. |
| 300 | | % prints sol(-1.0E+11) but then fails |
| 301 | | % conservative: do the precision check and also check that there are no other solutions (not done yet) |
| 302 | | % e.g., X+10000000000.0 = 10000000000.0 & (X=0.000000000001 or X=2.0) is FALSE with the float_solver, but TRUE with none |
| 303 | | |
| 304 | | |
| 305 | | :- assert_must_succeed(exhaustive_kernel_check_wf(kernel_reals:real_subtraction_wf(term(floating(3.0)),term(floating(2.0)),term(floating(1.0)),WF),WF)). |
| 306 | | real_subtraction_wf(term(floating(X)),term(floating(Y)),term(floating(R)),WF) :- |
| 307 | | real_sub_wf_aux(X,Y,R,WF). |
| 308 | | |
| 309 | | real_sub_wf_aux(X,Y,R,WF) :- |
| 310 | | use_clpfd_real_solver,!, |
| 311 | | '$='(R,X-Y), |
| 312 | | (do_not_double_check_solution -> true ; real_sub_wf_aux2(X,Y,R,WF)). |
| 313 | | real_sub_wf_aux(X,Y,R,WF) :- block_real_sub_wf_aux(X,Y,R,WF). |
| 314 | | |
| 315 | | :- block block_real_sub_wf_aux(-,-,?,?), block_real_sub_wf_aux(?,-,-,?), block_real_sub_wf_aux(-,?,-,?). |
| 316 | | block_real_sub_wf_aux(X,Y,R,WF) :- |
| 317 | | (nonvar(X) |
| 318 | | -> (nonvar(Y) |
| 319 | ? | -> safe_is(R,X-Y,WF) |
| 320 | | ; allow_constraint_propagation(Kind), |
| 321 | | safe_is(YY,X-R,WF), |
| 322 | | confirm_solution_for_expr(YY,X-SOL,SOL,R,Kind,Solutions,WF) |
| 323 | | -> set_solution_for_expr(Y,Solutions,WF) % deterministic constraint propagation found precise solution |
| 324 | | ; real_sub_wf_aux2(X,Y,R,WF) % delay until Y is known |
| 325 | | ) |
| 326 | | ; allow_constraint_propagation(Kind), |
| 327 | | safe_is(XX,R+Y,WF), |
| 328 | | confirm_solution_for_expr(XX,SOL-Y,SOL,R,Kind,Solutions,WF) |
| 329 | | -> set_solution_for_expr(X,Solutions,WF) % deterministic constraint propagation found precise solution |
| 330 | | ; real_sub_wf_aux2(X,Y,R,WF) % delay until X is known |
| 331 | | ). |
| 332 | | |
| 333 | | :- block real_sub_wf_aux2(-,?,?,?), real_sub_wf_aux2(?,-,?,?). |
| 334 | | real_sub_wf_aux2(X,Y,R,WF) :- safe_is(R,X-Y,WF). |
| 335 | | |
| 336 | | |
| 337 | | :- assert_must_succeed(exhaustive_kernel_check_wf([commutative],kernel_reals:real_multiplication_wf(term(floating(3.0)),term(floating(2.0)),term(floating(6.0)),WF),WF)). |
| 338 | | real_multiplication_wf(term(floating(X)),term(floating(Y)),term(floating(R)),WF) :- |
| 339 | | real_mul_wf_aux0(X,Y,R,WF). |
| 340 | | |
| 341 | | real_mul_wf_aux0(X,Y,R,WF) :- |
| 342 | | use_clpfd_real_solver,!, |
| 343 | | (X==Y -> '$='(R,X*X) ; % causes segmentation faults up to SICS 4.10.2-beta2 |
| 344 | | '$='(R,X*Y),propagate_square_inverse(X,Y,R,WF)), %tools_printing:print_term_summary(mul(R,X,Y)),nl, |
| 345 | | (do_not_double_check_solution -> true ; real_mul_wf_aux2(X,Y,R,WF)). |
| 346 | | real_mul_wf_aux0(X,Y,R,WF) :- block_real_mul_wf_aux0(X,Y,R,WF). |
| 347 | | |
| 348 | | % | ?- 2.0 $= X*0.0. -> clpfd:(X in fsup..fsup) ? i.e., CLPFD does not detect multiplication by 0.0 |
| 349 | | % CLPFD also does not take square root |
| 350 | | :- block propagate_square_inverse(-,-,-,?). |
| 351 | | propagate_square_inverse(X,Y,R,_WF) :- (X==0.0 ; Y==0.0), !, R=0.0. |
| 352 | | propagate_square_inverse(X,Y,R,WF) :- var(X),X==Y, !, % X*X = R -> compute sqrt |
| 353 | | %'$='(R,X^2.0). % makes test 2430 fail for c=card({y|y*y=x}) & x = 2.0E-309 & c=14 -> c=1 |
| 354 | | %'$='(R,X*X). % x * x = 1.0 -> segmentation fault |
| 355 | | kernel_square_inverse(X,R,WF). |
| 356 | | propagate_square_inverse(X,Y,R,_WF) :- nonvar(R), R \= 0.0, !, |
| 357 | | dif(X,0.0), dif(Y,0.0). |
| 358 | | propagate_square_inverse(_,_,_,_). |
| 359 | | |
| 360 | | :- block block_real_mul_wf_aux0(-,-,-,?). |
| 361 | | block_real_mul_wf_aux0(X,Y,R,_WF) :- (X==0.0 ; Y==0.0), !, R=0.0. |
| 362 | | block_real_mul_wf_aux0(X,Y,R,WF) :- var(X),X==Y, !, % X*X = R -> compute sqrt |
| 363 | | kernel_square_inverse(X,R,WF). |
| 364 | | block_real_mul_wf_aux0(X,Y,R,WF) :- real_mul_wf_aux(X,Y,R,WF). |
| 365 | | |
| 366 | | |
| 367 | | :- block real_mul_wf_aux(-,-,?,?), real_mul_wf_aux(?,-,-,?), real_mul_wf_aux(-,?,-,?). |
| 368 | | real_mul_wf_aux(X,Y,R,_WF) :- (X==0.0 ; Y==0.0), !, R=0.0. |
| 369 | | real_mul_wf_aux(X,Y,R,WF) :- |
| 370 | | (nonvar(X) |
| 371 | | -> (nonvar(Y) |
| 372 | ? | -> safe_is(R,X*Y,WF) |
| 373 | | ; allow_constraint_propagation(Kind), |
| 374 | | safe_is(YY,R/X,WF), |
| 375 | | confirm_solution_for_expr(YY,X*SOL,SOL,R,Kind,Solutions,WF) |
| 376 | | -> set_solution_for_expr(Y,Solutions,WF) |
| 377 | | ; real_mul_wf_aux2(X,Y,R,WF) % delay until Y is known |
| 378 | | ) |
| 379 | | ; allow_constraint_propagation(Kind), |
| 380 | | safe_is(XX,R/Y,WF), |
| 381 | | confirm_solution_for_expr(XX,SOL*Y,SOL,R,Kind,Solutions,WF) |
| 382 | | -> set_solution_for_expr(X,Solutions,WF) |
| 383 | | ; real_mul_wf_aux2(X,Y,R,WF) % delay until X is known |
| 384 | | ). |
| 385 | | |
| 386 | | :- block real_mul_wf_aux2(-,?,?,?), real_mul_wf_aux2(?,-,?,?). |
| 387 | | real_mul_wf_aux2(X,Y,R,WF) :- safe_is(R,X*Y,WF). |
| 388 | | |
| 389 | | % compute sqrt, R is known |
| 390 | | kernel_square_inverse(X,R,WF) :- |
| 391 | | (R < 0.0 -> fail |
| 392 | | ; allow_constraint_propagation(Kind), |
| 393 | | safe_is(XX,sqrt(R),WF), |
| 394 | | confirm_solution_for_expr(XX,SOL*SOL,SOL,R,Kind,SolutionsPos,WF), |
| 395 | | (XX=0.0 -> Solutions = SolutionsPos |
| 396 | | ; NX is -XX, |
| 397 | | confirm_solution_for_expr(NX,SOL*SOL,SOL,R,Kind,SolutionsNeg,WF), |
| 398 | | append(SolutionsPos,SolutionsNeg,Solutions) |
| 399 | | ) |
| 400 | | -> set_solution_for_expr(X,Solutions,WF) |
| 401 | | ; real_mul_wf_aux(X,X,R,WF)). |
| 402 | | |
| 403 | | |
| 404 | | % used for 'RSQRT' external function |
| 405 | | :- block 'real_square_root_wf'(-,-,?,?). |
| 406 | | real_square_root_wf(term(floating(X)),term(floating(R)),Span,WF) :- |
| 407 | | real_sqrt_wf_aux(X,R,Span,WF). |
| 408 | | |
| 409 | | real_sqrt_wf_aux(X,R,Span,WF) :- |
| 410 | | use_clpfd_real_solver,!, |
| 411 | | '$='(R,sqrt(X)), |
| 412 | | (do_not_double_check_solution -> true ; real_unop_wf_aux('sqrt',X,R,Span,WF)). |
| 413 | | real_sqrt_wf_aux(X,R,Span,WF) :- block_real_sqrt_wf_aux(X,R,Span,WF). |
| 414 | | |
| 415 | | :- block block_real_sqrt_wf_aux(-,-,?,?). |
| 416 | | block_real_sqrt_wf_aux(X,R,_Span,WF) :- var(X), |
| 417 | | allow_constraint_propagation(Kind), |
| 418 | | get_preference(find_abort_values,false), |
| 419 | | (R >= 0.0 |
| 420 | | -> safe_is(XX,R*R,WF), |
| 421 | | confirm_solution_for_expr(XX,sqrt(SOL),SOL,R,Kind,Solutions,WF) |
| 422 | | ; Solutions = []), % RSQRT always returns a positive number |
| 423 | | !, |
| 424 | | set_solution_for_expr(X,Solutions,WF). |
| 425 | | block_real_sqrt_wf_aux(X,R,Span,WF) :- |
| 426 | | real_unop_wf_aux('sqrt',X,R,Span,WF). |
| 427 | | |
| 428 | | |
| 429 | | :- assert_must_succeed(exhaustive_kernel_check_wf(kernel_reals:real_division_wf(term(floating(6.0)),term(floating(2.0)),term(floating(3.0)),unknown,WF),WF)). |
| 430 | | real_division_wf(term(floating(X)),term(floating(Y)),term(floating(R)),Span,WF) :- |
| 431 | | real_div_wf_aux(X,Y,R,Span,WF). |
| 432 | | |
| 433 | | real_div_wf_aux(X,Y,R,Span,WF) :- |
| 434 | | use_clpfd_real_solver,!, |
| 435 | | '$='(R,X / Y), |
| 436 | | (do_not_double_check_solution -> true ; check_div(X,Y,R,Span,WF)). |
| 437 | | real_div_wf_aux(X,Y,R,Span,WF) :- block_real_div_wf_aux(X,Y,R,Span,WF). |
| 438 | | |
| 439 | | :- block block_real_div_wf_aux(-,?,-,?,?), block_real_div_wf_aux(?,-,?,?,?). |
| 440 | ? | block_real_div_wf_aux(X,Y,R,Span,WF) :- nonvar(X),!,safe_is(R,X/Y,Span,WF). |
| 441 | | block_real_div_wf_aux(X,Y,R,Span,WF) :- Y \= 0.0, allow_constraint_propagation(Kind), |
| 442 | | safe_is(XX,R*Y,Span,WF), |
| 443 | | confirm_solution_for_expr(XX,SOL/Y,SOL,R,Kind,Solutions,WF), |
| 444 | | !, |
| 445 | | set_solution_for_expr(X,Solutions,WF). |
| 446 | | block_real_div_wf_aux(X,Y,R,Span,WF) :- check_div(X,Y,R,Span,WF). |
| 447 | | |
| 448 | | :- block check_div(-,?,?,?,?), check_div(?,-,?,?,?). |
| 449 | | check_div(X,Y,R,Span,WF) :- safe_is(R,X/Y,Span,WF). |
| 450 | | |
| 451 | | :- block 'real_unary_minus_wf'(-,-,?). |
| 452 | | real_unary_minus_wf(term(floating(X)),term(floating(R)),WF) :- |
| 453 | | real_um_wf_aux(X,R,WF). |
| 454 | | |
| 455 | | real_um_wf_aux(X,R,WF) :- |
| 456 | | use_clpfd_real_solver,!, |
| 457 | | '$='(R,-X), |
| 458 | | (do_not_double_check_solution -> true ; block_real_um_wf_aux2(X,R,WF)). |
| 459 | | real_um_wf_aux(X,R,WF) :- block_real_um_wf_aux2(X,R,WF). |
| 460 | | |
| 461 | | |
| 462 | | :- block block_real_um_wf_aux2(-,-,?). |
| 463 | | block_real_um_wf_aux2(X,R,_) :- nonvar(X),!,R is -X. |
| 464 | | block_real_um_wf_aux2(X,R,_) :- X is -R. |
| 465 | | |
| 466 | | % used for 'RABS' external function |
| 467 | | :- block 'real_absolute_value_wf'(-,-,?). |
| 468 | | real_absolute_value_wf(term(floating(X)),term(floating(R)),WF) :- |
| 469 | | real_abs_wf_aux(X,R,WF). |
| 470 | | real_abs_wf_aux(X,R,WF) :- |
| 471 | | use_clpfd_real_solver,!, |
| 472 | | '$='(R,abs(X)), |
| 473 | | (do_not_double_check_solution -> true ; when(nonvar(X), safe_is(R,abs(X),WF))). |
| 474 | | real_abs_wf_aux(X,R,WF) :- block_real_abs_wf_aux(X,R,WF). |
| 475 | | |
| 476 | | :- block block_real_abs_wf_aux(-,-,?). |
| 477 | | block_real_abs_wf_aux(X,R,_) :- nonvar(X),!, R is abs(X). |
| 478 | | block_real_abs_wf_aux(X,R,_) :- R = 0.0, !, X=0.0. |
| 479 | | block_real_abs_wf_aux(_,R,_) :- R < 0.0, !, fail. |
| 480 | | block_real_abs_wf_aux(X,R,WF) :- MR is -R, |
| 481 | | Solutions = [R,MR], |
| 482 | | set_solution_for_expr(X,Solutions,WF). |
| 483 | | |
| 484 | | |
| 485 | | % used for 'RSIGN' external function |
| 486 | | :- block 'real_sign_wf'(-,-,?). |
| 487 | | real_sign_wf(term(floating(X)),term(floating(R)),WF) :- |
| 488 | | real_sign_wf_aux(X,R,WF). |
| 489 | | |
| 490 | | real_sign_wf_aux(X,R,WF) :- |
| 491 | | use_clpfd_real_solver,!, |
| 492 | | '#='(IR,sign(X)), % Note: it is important to use sign(.) directly here with #=, |
| 493 | | % otherwise we get an expected 'numeric argument' error |
| 494 | | convert_int_to_real_aux(IR,R), |
| 495 | | (do_not_double_check_solution -> true ; block_real_unop_wf(sign,X,R,unknown,WF)). |
| 496 | | real_sign_wf_aux(X,R,WF) :- block_real_unop_wf(sign,X,R,unknown,WF). |
| 497 | | |
| 498 | | |
| 499 | | % used for 'ROUND' external function |
| 500 | | real_round_wf(term(floating(X)),int(R),_WF) :- real_round_aux(X,R). |
| 501 | | |
| 502 | | real_round_aux(Real,Int) :- |
| 503 | | use_clpfd_real_solver,!, |
| 504 | | '#='(Int,round(Real)). |
| 505 | | real_round_aux(Real,Int) :- block_real_round_aux(Real,Int). |
| 506 | | |
| 507 | | :- block block_real_round_aux(-,?). |
| 508 | | block_real_round_aux(Real,Int) :- Int is round(Real). |
| 509 | | |
| 510 | | |
| 511 | | :- use_module(kernel_waitflags,[add_wd_error/3, add_wd_error_span/4]). |
| 512 | | % a version of is/2 which catches overflows |
| 513 | | safe_is(R,Expr,WF) :- |
| 514 | ? | safe_is(R,Expr,unknown,WF). |
| 515 | | safe_is(R,Expr,Span,WF) :- |
| 516 | ? | catch(R is Expr, |
| 517 | | error(evaluation_error(ERR),_), |
| 518 | | process_evaluation_error(ERR,Expr,Span,WF)). |
| 519 | | |
| 520 | | process_evaluation_error(float_overflow,Expr,Span,WF) :- !, |
| 521 | | add_wd_error_span('Float Overflow while computing:',Expr,Span,WF). |
| 522 | | process_evaluation_error(zero_divisor,Expr,Span,WF) :- !, |
| 523 | | add_wd_error_span('Division by zero while computing:',Expr,Span,WF). |
| 524 | | process_evaluation_error(undefined,Expr,Span,WF) :- !, |
| 525 | | add_wd_error_span('Arithmetic operator undefined while computing:',Expr,Span,WF). |
| 526 | | process_evaluation_error(_,Expr,Span,WF) :- |
| 527 | | add_wd_error_span('Unknown evaluation error while computing:',Expr,Span,WF). |
| 528 | | |
| 529 | | |
| 530 | | % call a Prolog unary artihmetic operator |
| 531 | | real_unop_wf(OP,X,R,WF) :- |
| 532 | | real_unop_wf(OP,X,R,unknown,WF). |
| 533 | | |
| 534 | | real_unop_wf(OP,RX,RR,Span,WF) :- |
| 535 | | get_real(RX,X,OP,Span,WF), get_real(RR,R,OP,Span,WF), |
| 536 | | real_unop_wf_aux(OP,X,R,Span,WF). |
| 537 | | |
| 538 | | unsupported_unary_operator(float_fractional_part). % RFRACTION |
| 539 | | |
| 540 | | real_unop_wf_aux(OP,X,R,Span,WF) :- |
| 541 | | use_clpfd_real_solver, |
| 542 | | \+ unsupported_unary_operator(OP),!, |
| 543 | | Expr =.. [OP,X], |
| 544 | | post_real_equal_expr_wf(R,Expr,Span,WF), % TODO: catch WD issues (RASIN(>1.0), ...) |
| 545 | | (do_not_double_check_solution -> true |
| 546 | | ; block_real_unop_wf(OP,X,R,Span,WF)). % relevant e.g. for {y|x = RSQRT(y)}=r & x=4.0 |
| 547 | | real_unop_wf_aux(OP,X,R,Span,WF) :- block_real_unop_wf(OP,X,R,Span,WF). |
| 548 | | |
| 549 | | :- block block_real_unop_wf(-,?,?,?,?), block_real_unop_wf(?,-,?,?,?). |
| 550 | | block_real_unop_wf(OP,X,R,Span,WF) :- |
| 551 | | Expr =.. [OP,X], |
| 552 | | safe_is(R,Expr,Span,WF). |
| 553 | | |
| 554 | | |
| 555 | | :- use_module(probsrc(tools_strings),[ajoin/2]). |
| 556 | | :- use_module(kernel_waitflags,[add_error_wf/5]). |
| 557 | | get_real(Var,Real,_,_,_) :- var(Var),!, Var=term(floating(Real)). |
| 558 | | get_real(term(F),Real,_,_,_) :- !, F=floating(Real). |
| 559 | | get_real(Other,_,OP,Span,WF) :- |
| 560 | | ajoin(['Argument for ',OP,' is not a real number:'],Msg), |
| 561 | | add_error_wf(kernel_reals,Msg,Other,Span,WF). |
| 562 | | |
| 563 | | % when called for ** operator (power_of_real) the exponent Y is a natural number represented as a float; |
| 564 | | % for RPOW Y can be any float |
| 565 | | real_power_of_wf(RX,RY,RR,Span,WF) :- |
| 566 | | get_real(RX,X,'exp',Span,WF), get_real(RY,Y,'exp',Span,WF), get_real(RR,R,'exp',Span,WF), |
| 567 | | real_power_of_aux(X,Y,R,Span,WF). |
| 568 | | |
| 569 | | :- if(fail). %can_use_clpfd_real_solver). % % 4.10 or higher |
| 570 | | % code still deactivated as i:10..20 & 10.0**i + 1.0 = 10.0**i from test 2430 crashes SICStus |
| 571 | | real_power_of_aux(X,Y,R,Span,WF) :- % write(rpow(X,Y,R)),nl,flush_output, |
| 572 | | real_binop_wf_aux('^',X,Y,R,Span,WF). |
| 573 | | :- else. |
| 574 | | :- block real_power_of_aux(-,-,?,?,?), real_power_of_aux(?,-,-,?,?), real_power_of_aux(-,?,-,?,?). |
| 575 | | % we currently require to know at least the exponent, TODO: in future we could also compute the log |
| 576 | | real_power_of_aux(X,Y,R,Span,WF) :- |
| 577 | | var(X),!, |
| 578 | | real_power_inverse(X,Y,R,Span,WF). |
| 579 | | real_power_of_aux(X,Y,Res,Span,WF) :- |
| 580 | | real_power_direct(X,Y,Res,Span,WF). |
| 581 | | |
| 582 | | % X is variable, Y and R are known |
| 583 | | real_power_inverse(X,Y,R,_Span,WF) :- Y=2.0,!, |
| 584 | | kernel_square_inverse(X,R,WF). |
| 585 | | real_power_inverse(X,Y,R,_Span,_WF) :- Y=1.0,!, X=R. |
| 586 | | real_power_inverse(X,Y,R,Span,WF) :- |
| 587 | | real_power_direct(X,Y,R,Span,WF). |
| 588 | | % TODO: add more exponents and compute root |
| 589 | | :- endif. |
| 590 | | |
| 591 | | :- block real_power_direct(-,?,?,?,?), real_power_direct(?,-,?,?,?). |
| 592 | | real_power_direct(X,Y,R,Span,WF) :- |
| 593 | | safe_is(R,exp(X,Y),Span,WF). |
| 594 | | |
| 595 | | % call a Prolog binary artihmetic operator: min, max, atan2, log |
| 596 | | real_binop_wf(OP,RX,RY,RR,Span,WF) :- |
| 597 | | get_real(RX,X,OP,Span,WF), get_real(RY,Y,OP,Span,WF), get_real(RR,R,OP,Span,WF), |
| 598 | | real_binop_wf_aux(OP,X,Y,R,Span,WF). |
| 599 | | |
| 600 | | % atan2, log are not supported by CLP(FD) |
| 601 | | supported_binary_operator(min). |
| 602 | | supported_binary_operator(max). |
| 603 | | supported_binary_operator('^'). |
| 604 | | |
| 605 | | real_binop_wf_aux(OP,X,Y,R,Span,WF) :- supported_binary_operator(OP), |
| 606 | | use_clpfd_real_solver,!, |
| 607 | | Expr =.. [OP,X,Y], |
| 608 | | post_real_equal_expr_wf(R,Expr,Span,WF), % TODO: double check necessary?? |
| 609 | | (do_not_double_check_solution -> true |
| 610 | | ; block_real_binop_wf(OP,X,Y,R,Span,WF)). % necessary for min/max ? |
| 611 | | real_binop_wf_aux(OP,X,Y,R,Span,WF) :- block_real_binop_wf(OP,X,Y,R,Span,WF). |
| 612 | | |
| 613 | | :- block block_real_binop_wf(-,?,?,?,?,?), block_real_binop_wf(?,-,?,?,?,?), block_real_binop_wf(?,?,-,?,?,?). |
| 614 | | block_real_binop_wf(OP,X,Y,R,Span,WF) :- |
| 615 | | Expr =.. [OP,X,Y], |
| 616 | | safe_is(R,Expr,Span,WF). |
| 617 | | |
| 618 | | real_less_than_wf(term(floating(X)),term(floating(Y)),_WF) :- |
| 619 | | clpfd_less_real(X,Y). |
| 620 | | real_less_than_equal_wf(term(floating(X)),term(floating(Y)),_WF) :- |
| 621 | | clpfd_leq_real(X,Y). |
| 622 | | |
| 623 | | % call a Prolog binary artihmetic operator, reified version. |
| 624 | | real_comp_wf(OP,term(floating(X)),term(floating(Y)),R,_WF) :- |
| 625 | | real_comp_float(OP,X,Y,R). |
| 626 | | |
| 627 | | clpfd_leq_real(X,Y) :- use_clpfd_real_solver,!,'$=<'(X,Y). |
| 628 | | clpfd_leq_real(X,Y) :- real_comp_float_blocking('=<',X,Y,pred_true). |
| 629 | | |
| 630 | | clpfd_less_real(X,Y) :- use_clpfd_real_solver,!,'$=<'(X,Y), dif(X,Y). |
| 631 | | clpfd_less_real(X,Y) :- real_comp_float_blocking('<',X,Y,pred_true). |
| 632 | | |
| 633 | | real_comp_float('=<',X,Y,R) :- use_clpfd_real_solver, !, |
| 634 | | '#<=>'('$=<'(X,Y),R01), |
| 635 | | b_interpreter_check:prop_pred_01(R,R01). |
| 636 | | % TODO: reify <, =:=, =\\= |
| 637 | | real_comp_float(OP,X,Y,R) :- real_comp_float2(OP,X,Y,R). |
| 638 | | |
| 639 | | :- block real_comp_float2(?,-,?,-), real_comp_float2(?,?,-,-). |
| 640 | | real_comp_float2('<',X,Y,R) :- R == pred_true,!, clpfd_less_real(X,Y). |
| 641 | | real_comp_float2('<',X,Y,R) :- R == pred_false,!, clpfd_leq_real(Y,X). |
| 642 | | %real_comp_float2('=<',X,Y,R) :- R == pred_true,!, clpfd_leq_real(X,Y). % not necessary, as reified |
| 643 | | %real_comp_float2('=<',X,Y,R) :- R == pred_false,!, clpfd_less_real(Y,X). |
| 644 | | real_comp_float2(OP,X,Y,R) :- real_comp_float_blocking(OP,X,Y,R). |
| 645 | | |
| 646 | | :- block real_comp_float_blocking(-,?,?,?), real_comp_float_blocking(?,-,?,?), real_comp_float_blocking(?,?,-,?). |
| 647 | | real_comp_float_blocking(OP,X,Y,R) :- |
| 648 | | (call(OP,X,Y) -> R=pred_true ; R=pred_false). |
| 649 | | |
| 650 | | % ----------------- |
| 651 | | |
| 652 | | :- use_module(probsrc(kernel_tools),[ground_value_check/2]). |
| 653 | | :- use_module(probsrc(custom_explicit_sets),[expand_and_convert_to_avl_set/4, expand_custom_set_to_list_wf/5, |
| 654 | | max_of_explicit_set_wf/3, min_of_explicit_set_wf/3]). |
| 655 | | |
| 656 | | :- block real_maximum_of_set(-,?,?,?). |
| 657 | | real_maximum_of_set([],_Res,Span,WF) :- !, |
| 658 | | add_wd_error_span('max applied to empty set of reals:',[],Span,WF). |
| 659 | | real_maximum_of_set(avl_set(A),Res,_Span,WF) :- !, max_of_explicit_set_wf(avl_set(A),Res,WF). |
| 660 | | real_maximum_of_set([RX|T],RRes,Span,WF) :- use_clpfd_real_solver, !, |
| 661 | | get_real(RX,X,maximum,Span,WF), |
| 662 | | get_real(RRes,Res,maximum,Span,WF), |
| 663 | | clpfd_leq_real(X,Res), |
| 664 | | expand_custom_set_to_list_wf(T,ET,_,real_maximum_of_set,WF), |
| 665 | | real_max_list(ET,[X],Res,Span,WF). |
| 666 | | real_maximum_of_set(Set,Res,Span,WF) :- |
| 667 | | ground_value_check(Set,Gr), |
| 668 | | rmax_ground_set(Gr,Set,Res,Span,WF). |
| 669 | | |
| 670 | | :- block real_max_list(-,?,?,?). |
| 671 | | real_max_list([],List,Res,_Span,_WF) :- clpfd:maximum(Res,List). % requires new CLP(FD) solver in 4.10 |
| 672 | | real_max_list([RX|T],Acc,Res,Span,WF) :- |
| 673 | | get_real(RX,X,maximum,Span,WF), |
| 674 | | clpfd_leq_real(X,Res), % not useful if list skeleton known |
| 675 | | real_max_list(T,[X|Acc],Res,Span,WF). |
| 676 | | |
| 677 | | :- block rmax_ground_set(-,?,?,?,?). |
| 678 | | rmax_ground_set(_,Set,Res,Span,WF) :- |
| 679 | | expand_and_convert_to_avl_set(Set,ESet,'RMAXIMUM','RMAXIMUM'), |
| 680 | | (ESet=empty |
| 681 | | -> add_wd_error_span('max applied to empty set of reals:',Set,Span,WF) |
| 682 | | ; max_of_explicit_set_wf(avl_set(ESet),Res,WF) |
| 683 | | ). |
| 684 | | |
| 685 | | :- block real_minimum_of_set(-,?,?,?). |
| 686 | | real_minimum_of_set([],_Res,Span,WF) :- !, |
| 687 | | add_wd_error_span('min applied to empty set of reals:',[],Span,WF). |
| 688 | | real_minimum_of_set(avl_set(A),Res,_Span,WF) :- !, min_of_explicit_set_wf(avl_set(A),Res,WF). |
| 689 | | real_minimum_of_set([RX|T],RRes,Span,WF) :- use_clpfd_real_solver, !, |
| 690 | | get_real(RX,X,minumum,Span,WF), |
| 691 | | get_real(RRes,Res,minumum,Span,WF), |
| 692 | | clpfd_leq_real(Res,X), |
| 693 | | expand_custom_set_to_list_wf(T,ET,_,real_minimum_of_set,WF), |
| 694 | | real_min_list(ET,[X],Res,Span,WF). |
| 695 | | real_minimum_of_set(Set,Res,Span,WF) :- |
| 696 | | ground_value_check(Set,Gr), |
| 697 | | rmin_ground_set(Gr,Set,Res,Span,WF). |
| 698 | | |
| 699 | | :- block real_min_list(-,?,?,?). |
| 700 | | real_min_list([],List,Res,_Span,_WF) :- clpfd:minimum(Res,List). % requires new CLP(FD) solver in 4.10 |
| 701 | | real_min_list([RX|T],Acc,Res,Span,WF) :- |
| 702 | | get_real(RX,X,minimum,Span,WF), |
| 703 | | clpfd_leq_real(Res,X), % not useful if list skeleton known |
| 704 | | real_min_list(T,[X|Acc],Res,Span,WF). |
| 705 | | |
| 706 | | :- block rmin_ground_set(-,?,?,?,?). |
| 707 | | rmin_ground_set(_,Set,Res,Span,WF) :- |
| 708 | | expand_and_convert_to_avl_set(Set,ESet,'RMINIMUM','RMINIMUM'), |
| 709 | | (ESet=empty |
| 710 | | -> add_wd_error_span('min applied to empty set of reals:',Set,Span,WF) |
| 711 | | ; min_of_explicit_set_wf(avl_set(ESet),Res,WF) |
| 712 | | ). |
| 713 | | |
| 714 | | |
| 715 | | % not used yet: useful for kernel_strings ? |
| 716 | | %:- block real_to_string(-,?). |
| 717 | | %real_to_string(term(floating(I)),S) :- real_to_string2(I,S). |
| 718 | | % |
| 719 | | %:- block real_to_string2(-,?). |
| 720 | | %real_to_string2(Num,Res) :- |
| 721 | | % number_codes(Num,C), |
| 722 | | % atom_codes(S,C), Res=string(S). |
| 723 | | |
| 724 | | % ------------------------- |
| 725 | | |
| 726 | | :- use_module(kernel_objects,[gen_enum_warning_wf/6]). |
| 727 | | :- use_module(library(random),[random/3]). |
| 728 | | |
| 729 | | |
| 730 | | enumerate_real_wf(term(floating(F)),EnumWarning,WF) :- |
| 731 | ? | enum_real(F,EnumWarning,WF). |
| 732 | | |
| 733 | | enum_real(F,_,_) :- number(F),!. |
| 734 | | enum_real(F,EnumWarning,WF) :- use_clpfd_real_solver,!, |
| 735 | | clpfd_enum_real(F,EnumWarning,WF). |
| 736 | | enum_real(F,EnumWarning,WF) :- |
| 737 | | gen_enum_warning_wf('REAL',inf,'"0.0","1.0",...',EnumWarning,unknown,WF), |
| 738 | | ( F = 0.0 ; F = 1.0 ; F = -1.0 |
| 739 | | ; random(0.0,1.0,F) |
| 740 | | ; random(-1.0,0.0,F) |
| 741 | | ; preferences:preference(maxint,MaxInt), random(1.0,MaxInt,F) |
| 742 | | ; preferences:preference(minint,MinInt), random(MinInt,-1.0,F) |
| 743 | | ). |
| 744 | | |
| 745 | | :- if(can_use_clpfd_real_solver). % 4.10 or higher |
| 746 | | clpfd_enum_real(F,EnumWarning,WF) :- |
| 747 | | fd_min(F,Min), fd_max(F,Max), % finf .. fsup |
| 748 | | (number(Min),number(Max) -> |
| 749 | | (preferences:preference(real_solver_precision,Prec), Prec > 0.0 |
| 750 | | -> gen_enum_warning_wf('REAL',Min..Max,precision(Prec),EnumWarning,unknown,WF), |
| 751 | | labeling([random, precision(Prec)],[F]) |
| 752 | | ; labeling([random],[F]) |
| 753 | | ) |
| 754 | | ; number(Min), preferences:preference(maxint,MaxInt), |
| 755 | | MaxR is max(float(MaxInt),Min+1.0) % should we use next float instead of + 1.0? or precision(.) preference? |
| 756 | | -> gen_enum_warning_wf('REAL',Min..fsup,Min..MaxR,EnumWarning,unknown,WF), |
| 757 | | '$=<'(F,MaxR), labeling([random, precision(0.5)],[F]) % no need to do exhaustive enumeration; use precision |
| 758 | | ; number(Max), preferences:preference(minint,MinInt), |
| 759 | | MinR is min(float(MinInt),Max-1.0) % maybe use previous float or use precision(.) preference |
| 760 | | -> gen_enum_warning_wf('REAL',finf..Max,MinR..Max,EnumWarning,unknown,WF), |
| 761 | | '$>='(F,MinR), labeling([random, precision(0.5)],[F]) % no need to do exhaustive enumeration; use precision |
| 762 | | ; preferences:preference(minint,MinInt), MinR is float(MinInt), |
| 763 | | preferences:preference(maxint,MaxInt), MaxR is float(MaxInt) |
| 764 | | -> gen_enum_warning_wf('REAL',finf..fsup,MinR..MaxR,EnumWarning,unknown,WF), |
| 765 | | fd_batch(['$>='(F,MinR), '$=<'(F,MaxR)]), |
| 766 | | labeling([random, precision(0.5)],[F]) % no need to do exhaustive enumeration; enum warning anyway |
| 767 | | ). |
| 768 | | :- else. |
| 769 | | clpfd_enum_real(F,_,_) :- add_internal_error('Requires SICStus 4.10:',clpfd_enum_real(F,_,_)). |
| 770 | | :- endif. |
| 771 | | |
| 772 | | /* |
| 773 | | |
| 774 | | A few reference tests with SICS 4.10 real support in CLP(FD) |
| 775 | | |
| 776 | | | ?- X $>= 1.0, fd_dom(X,D), fd_min(X,Min), fd_max(X,Max). |
| 777 | | D = 1.0 .. fsup, |
| 778 | | Min = 1.0, |
| 779 | | Max = fsup, |
| 780 | | clpfd:(X in1.0 .. fsup) ? |
| 781 | | |
| 782 | | | ?- X in 1.0 .. 1.000000000000001, labeling([random],[X]). |
| 783 | | X = 1.0000000000000004 ? ; |
| 784 | | X = 1.0000000000000009 ? ; |
| 785 | | X = 1.000000000000001 ? ; |
| 786 | | X = 1.0000000000000007 ? ; |
| 787 | | X = 1.0 ? ; |
| 788 | | X = 1.0000000000000002 ? ; |
| 789 | | no |
| 790 | | |
| 791 | | | ?- X in 1.0 .. 2.0, X + 0.5 $= 2.0, labeling([precision(0.0)],[X]), Delta is (X + 0.5) - 2.0, X + 0.5 =:= 2.0. |
| 792 | | X = 1.5, |
| 793 | | Delta = 0.0 ? ; |
| 794 | | X = 1.5000000000000002, |
| 795 | | Delta = 0.0 ? ; |
| 796 | | no |
| 797 | | |
| 798 | | % enumeration always starts with 1.0 here and not very random: |
| 799 | | | ?- X $>= -1.0, X $=< 3.0, labeling([random],[X]). |
| 800 | | X = 1.0 ? ; |
| 801 | | X = -5.551115123125783E-17 ? ; |
| 802 | | X = 0.49999999999999994 ? ; |
| 803 | | X = 0.24999999999999992 ? ; |
| 804 | | X = 0.3749999999999999 ? ; |
| 805 | | X = 0.4374999999999999 ? ; |
| 806 | | X = 0.4062499999999999 ? ; |
| 807 | | X = 0.3906249999999999 ? ; |
| 808 | | X = 0.3984374999999999 ? ; |
| 809 | | X = 0.4023437499999999 ? ; |
| 810 | | X = 0.4042968749999999 ? ; |
| 811 | | X = 0.4033203124999999 ? ; |
| 812 | | X = 0.4038085937499999 ? ; |
| 813 | | X = 0.4035644531249999 ? ; |
| 814 | | X = 0.4034423828124999 ? ; |
| 815 | | X = 0.4035034179687499 ? ; |
| 816 | | X = 0.4035339355468749 ? ; |
| 817 | | X = 0.4035186767578124 ? ; |
| 818 | | X = 0.40351104736328114 ? ; |
| 819 | | |
| 820 | | | ?- X $= 1.0 + Y, X = 2.0. |
| 821 | | X = 2.0, |
| 822 | | clpfd:(Y in0.9999999999999996 .. 1.0000000000000004) ? |
| 823 | | yes |
| 824 | | |
| 825 | | | ?- X $= 1.0 + Y, X = 2.0, labeling([random],[Y]), write(Y),nl, (X is 1.0 + Y -> true ; write(ko),nl,fail). |
| 826 | | 1.0 |
| 827 | | X = 2.0, |
| 828 | | Y = 1.0 ? ; |
| 829 | | 0.9999999999999998 |
| 830 | | ko |
| 831 | | 0.9999999999999997 |
| 832 | | ko |
| 833 | | 0.9999999999999996 |
| 834 | | ko |
| 835 | | 0.9999999999999999 |
| 836 | | X = 2.0, |
| 837 | | Y = 0.9999999999999999 ? ; |
| 838 | | 1.0000000000000004 |
| 839 | | ko |
| 840 | | 1.0000000000000002 |
| 841 | | X = 2.0, |
| 842 | | Y = 1.0000000000000002 ? ; |
| 843 | | no |
| 844 | | | ?- (X $>= 1.0) #<=> R, X $>= 1.1. |
| 845 | | R = 1, |
| 846 | | clpfd:(X in1.1 .. fsup) ? |
| 847 | | yes |
| 848 | | |
| 849 | | | ?- X in 11..12, Y $>= float(X). |
| 850 | | clpfd:(Y in 11.0..fsup), |
| 851 | | X in 11..12 ? |
| 852 | | yes |
| 853 | | | ?- Y in 11.1 .. 12.0, X #>= floor(Y). |
| 854 | | clpfd:(Y in 11.1..12.0), |
| 855 | | X in 11..sup ? |
| 856 | | yes |
| 857 | | | ?- Y in 11.1 .. 12.0, X #>= ceiling(Y). |
| 858 | | clpfd:(Y in 11.1..12.0), |
| 859 | | X in 12..sup ? |
| 860 | | | ?- X in 1.0 .. 2.5, Y $= sin(X). |
| 861 | | clpfd:(X in 1.0..2.5), |
| 862 | | clpfd:(Y in 0.5984721441039559..1.0) ? |
| 863 | | yes |
| 864 | | | ?- X in 1.0 .. 2.5, Y $= cos(X). |
| 865 | | clpfd:(X in 1.0..2.5), |
| 866 | | clpfd:(Y in-0.8011436155469343..0.54030230586814) ? |
| 867 | | yes |
| 868 | | | ?- X in 1.0 .. 2.5, Y $= sinh(X). |
| 869 | | clpfd:(X in 1.0..2.5), |
| 870 | | clpfd:(Y in 1.1752011936438012..6.050204481039793) ? |
| 871 | | yes |
| 872 | | | ?- X in 1.0 .. 2.5, Y $= tanh(X). |
| 873 | | clpfd:(X in 1.0..2.5), |
| 874 | | clpfd:(Y in 0.7615941559557647..0.9866142981514303) ? |
| 875 | | yes |
| 876 | | |
| 877 | | | ?- I in 1.0 .. 2.0, J in 1.5 .. 2.1, maximum(Max,[I,J]). |
| 878 | | clpfd:(I in 1.0..2.0), |
| 879 | | clpfd:(J in 1.5..2.1), |
| 880 | | clpfd:(Max in 1.5..2.1) ? |
| 881 | | yes |
| 882 | | */ |