clpBNR.pl -- clpBNR: Constraint Logic Programming over Continuous Domain of Reals

CLP(BNR) (library(clpbnr), henceforth just clpBNR) is a CLP over the domain of real numbers extended with ±∞. Since integers are a proper subset of reals, and booleans (0 or 1) a subset of integers, these "sub-domains" are also supported.

Since the set of real numbers is continuous it's not possible to represent an aribitray real number, e.g., π in the finite resources of a computer. So clpBNR uses intervals to represent the domain of a numeric variable. A real variable X has a domain of (L,U) if L ≤ X ≤ U where L and U are numeric values which can be finitely represented, e.g., floats, integers or rationals.

The use of intervals (and interval arithmetic) provides guarantees of completeness and correctness - unlike floating point arithmetic - by sacrificing some precision since calulations using floating point domain bounds will be outward rounded.

Finiteness is guaranteed since intervals can only get narrower over the course of a computation. Certainty is only guaranteed if there are no solutions (i.e., query fails) - final interval values may contain 0, 1, or many solutions. When this occurs, the application can further constrain the solution, e.g., by testing specific (point) values in the domain, or by making use of some external knowledge of the problem being solved.

More extensive documentation and many examples are provided in A Guide to CLP(BNR) (HTML version included with this pack in directory docs/).

Documentation for exported predicates follows. The "custom" types include:

interval : a variable with a clpBNR attribute
numeric : an interval or a number
*_list : a list of *
*_List : a * or a list of *

clpStatistics is det

Resets clpBNR statistics - always succeeds.

clpBNR collects a number of "operational measurements" on a per-thread basis and combines them with some system statistics for subsequent querying. clpBNR measurements include:

`narrowingOps`	number of interval primitives called
`narrowingFails`	number of interval primitive failures
`node_count`	number of nodes in `clpBNR` constraint network
`max_iterations`	maximum number of iterations before throttling occurs (`max/limit`

System statistics included in clpStatistics:

`userTime`	from `statistics:cputime`
`gcTime`	from `statistics:garbage_collection.Time`
`globalStack`	from `statistics:globalused/statistics:global`
`trailStack`	from `statistics:trailused/statistics:trail`
`localStack`	from `statistics:localused/statistics:local`
`inferences`	from `statistics:inferences`

clpStatistic(?S) is nondet

Succeeds if S unifies with a clpStatistic value; otherwise fails. On backtracking all values that unify with S will be generated. Examples:

?- clpStatistics, X::real, {X**4-4*X**3+4*X**2-4*X+3==0}, clpStatistic(narrowingOps(Ops)).
Ops = 2245,
X::real(-1.509169756145379, 4.18727500493995).

?- clpStatistics, X::real, {X**4-4*X**3+4*X**2-4*X+3==0}, clpStatistic(S).
S = userTime(0.02277600000000035),
X::real(-1.509169756145379, 4.18727500493995) ;
S = gcTime(0.0),
X::real(-1.509169756145379, 4.18727500493995) ;
S = globalStack(43696/524256),
X::real(-1.509169756145379, 4.18727500493995) ;
S = trailStack(664/133096),
X::real(-1.509169756145379, 4.18727500493995) ;
S = localStack(1864/118648),
X::real(-1.509169756145379, 4.18727500493995) ;
S = inferences(86215),
X::real(-1.509169756145379, 4.18727500493995) ;
S = narrowingOps(2245),
X::real(-1.509169756145379, 4.18727500493995) ;
S = narrowingFails(0),
X::real(-1.509169756145379, 4.18727500493995) ;
S = node_count(9),
X::real(-1.509169756145379, 4.18727500493995) ;
S = max_iterations(2245/3000),
X::real(-1.509169756145379, 4.18727500493995).

clpStatistics(?Ss) is semidet

Succeeds if Ss unifies with a list of clpStatistic's values; otherwise fails. Example:

?- clpStatistics, X::real, {X**4-4*X**3+4*X**2-4*X+3==0}, clpStatistics(Ss).
Ss = [userTime(0.023398999999999504), gcTime(0.001), globalStack(19216/131040), trailStack(1296/133096), localStack(2184/118648), inferences(82961), narrowingOps(2245), narrowingFails(0), node_count(9), max_iterations(2245/3000)],
X::real(-1.509169756145379, 4.18727500493995).

list(?X:list) is semidet

Succeeds if X is a list; otherwise fails. Note: not equivalent to is_list/1 but executes in O(1) time. This filter is provided for historical compatability.

interval(?X:interval) is semidet

Succeeds if X is an interval, i.e., a variable with a clpBNR attribute; otherwise fails.

interval_degree(?X:numeric, ?N:integer) is semidet

Succeeds if X is numeric and N = number of clpBNR constraints on X; otherwise fails. If X is a number, N = 0. Examples:

?- {X==Y+1}, interval_degree(X,N).
N = 1,
X::real(-1.0Inf, 1.0Inf),
Y::real(-1.0Inf, 1.0Inf).

?- interval_degree(42,N).
N = 0.

nb_setbounds(?X:interval, +Bs:number_list) is semidet

Succeeds if X is an interval and can be narrowed to the bounds Bs = [L,U]; otherwise fails. On backtracking, this value is not undone.

Caution: this predicate is non-logical and intended for specialized use case, e.g., some branch-and-bound algorithms (narrow to current solution, then backtrack to next solution).

range(?X, ?Bs:number_list) is semidet

Succeeds if X is numeric and Bs unifies with a list containing the lower and upper bound of X; otherwise fails. If X is a logic variable range(X,[2,3]) is equivalent to X::real(2,3). If X is a number the lower and upper bounds are the same. Examples:

?- X::integer(1,10), range(X,Bs).
Bs = [1, 10],
X::integer(1, 10).

?- range(42,Bs).
Bs = [42, 42].

?- range(X,[2,3]).
X::real(2, 3).

domain(?X:interval, ?Dom) is semidet

Succeeds if X is an interval and Dom unifies with the domain of X; otherwise fails. Dom is a compound term with a functor specifying the type (real or integer) and the arguments specifying the bounds. If X has a domain of integer(0,1), Dom will be boolean. Examples:

?- range(X,[2,3]), domain(X,Dom).
Dom = real(2, 3),
X::real(2, 3).

?- X::integer(0,1),domain(X,Dom).
Dom = boolean,
X::boolean.

?- domain(X,Dom).
false.

Note: unlike range/2, domain/2 will not change X.

delta(?X:numeric, ?W:number) is semidet

Succeeds if X is numeric and W unifies with the width of X (upper bound-lowerbound); otherwise fails. Examples:

?- X:: real(1r2,5r3),delta(X,D).
D = 7r6,
X::real(0.5, 1.6666666666666667).

?- delta(42,W).
W = 0.

delta is also available as an arithmetic function:

?- X::real(1r2,pi), W is delta(X).
W = 2.6415926535897936,
X::real(0.5, 3.1415926535897936).

midpoint(?X:numeric, ?M:number) is semidet

Succeeds if X is numeric and M unifies with the midpoint of X; otherwise fails. Examples:

?- X:: real(1r2,5r3), midpoint(X,M).
M = 13r12,
X::real(0.5, 1.6666666666666667).

?- midpoint(42,M).
M = 42.

midpoint is also available as an arithmetic function:

?- X::real(1r2,pi), M is midpoint(X).
M = 1.8207963267948968,
X::real(0.5, 3.1415926535897936).

median(?X:numeric, ?M:float) is semidet

Succeeds if X is numeric and M unifies with the median of X; otherwise fails. The median is 0 if the domain of X contains 0; otherwise it is the floating point value which divides the interval into two sub-domains each containing approximately equal numbers of floating point values. Examples:

?- X:: real(1r2,5r3), median(X,M).
M = 0.9128709291752769,
X::real(0.5, 1.6666666666666667).

?- median(42,M).
M = 42.0.

median is also available as an arithmetic function:

?- X::real(1r2,pi), M is median(X).
M = 1.2533141373155003,
X::real(0.5, 3.1415926535897936).

lower_bound(?X:numeric) is semidet

Succeeds if X is numeric and unifies with the lower bound of its domain. Examples:

?- X::integer(1,10),lower_bound(X).
X = 1.

?- X = 42, lower_bound(X).
X = 42.

Note that lower_bound will unify X with a number on success, but it may fail if this value is inconsistent with current constraints.

upper_bound(?X:numeric) is semidet

Succeeds if X is numeric and unifies with the upper bound of its domain. Examples:

?- X::integer(1,10),upper_bound(X).
X = 10.

?- X = 42, upper_bound(X).
X = 42.

Note that upper_bound will unify X with a number on success, but it may fail if this value is inconsistent with current constraints.

::(-X:numeric_List, ?Dom) is semidet

Succeeds if variable X has domain Dom; otherwise fails. If Dom, or either bound of Dom, is a variable (or missing), it will be unified with the default value depending on its type. Default domains are real(-1.0e+16, 1.0e+16) and integer(-72057594037927936, 72057594037927935). Examples:

?- X::real(-pi/2,pi/2).
X::real(-1.5707963267948968, 1.5707963267948968).

?- X::real, Y::integer.
X::real(-1.0e+16, 1.0e+16),
Y::integer(-72057594037927936, 72057594037927935).

?- Y::integer(1,_), Y::Dom.
Dom = integer(1, 72057594037927935),
Y::integer(1, 72057594037927935).

?- B::boolean.
B::boolean.

?- 42::Dom.
false.

Note that bounds can be defined using arithmetic expressions.

Alternatively, the first argument may be a list of variables:

?- [B1,B2,B3]::boolean.
B1::boolean,
B2::boolean,
B3::boolean.

?- length(Vs,3), Vs::real(-1,1).
Vs = [_A, _B, _C],
_A::real(-1, 1),
_B::real(-1, 1),
_C::real(-1, 1).

{+Constraints} is semidet

Succeeds if Constraints is a sequence of one or more boolean expressions (typically equalities and inequalities) defining a set of valid and consistent constraints; otherwise fails. The ',' binary operator is used to syntactically separate individual constraints and has semantics and (similar to the use of ',' in clause bodies). Arithmetic expressions are expressed using the same set of operators used in functional Prolog arithmetic (listed below) with addition of boolean operators to support that sub-domain of reals.

Table of supported interval relations:

`+ - * /`	arithmetic
`**`	includes real exponent, odd/even integer
`abs`	absolute value
`sqrt`	positive square root
`min max`	binary min/max
`== is <> =\= =< >= < >`	comparison (`is` and `=\=` synonyms for `==` and `<>`)
`<=`	included (one way narrowing)
`and or nand nor xor -> ,`	boolean (`,` synonym for `and`)
`- ~`	unary negate and not (boolean)
`exp log`	exp/ln
`sin asin cos acos tan atan`	trig functions
`integer`	must be an integer value

clpBNR defines the following additional operators for use in constraint expressions:

`op(200, fy, ~)`	boolean 'not'
`op(500, yfx, and)`	boolean 'and'
`op(500, yfx, or)`	boolean 'or'
`op(500, yfx, nand)`	boolean 'nand'
`op(500, yfx, nor)`	boolean 'nor'

Note that the comparison operators <>, =\=, '<' and '>' are unsound (due to incompleteness) over the real domain but sound over the integer domain. Strict inequality (<> and =\=) is disallowed for type real (will be converted to type integer) but < and > are supported for reals since they may be useful for things like branch and bound searches (with caution). The boolean functions are restricted to type 'boolean', constraining their argument types accordingly. Some examples:

?- {X == Y+1, Y >= 1}.
X::real(2, 1.0Inf),
Y::real(1, 1.0Inf).

?- {X == cos(X)}.
X:: 0.73908513321516... .

?- X::real, {X**4-4*X**3+4*X**2-4*X+3==0}.
X::real(-1.509169756145379, 4.18727500493995).

?- {A or B, C and D}.
C = D, D = 1,
A::boolean,
B::boolean.

Note that any variable in a constraint expression with no domain will be assigned the most general value consistent with the operator types, e.g., real(-1.0Inf,1.0Inf), boolean, etc.

print_interval(?T) is semidet

Succeeds printing term T with interval expanded to domains and vars labelled V*; fails if output fails. Provided for historical compatibility, use system output facilities instead. Example:

?- X::real,print_interval(f(X)),X=42.
f(V0::real(-1.0e+16,1.0e+16))
X = 42.

deprecated: - use format/2

print_interval(+Stream, ?T) is semidet

Same as print_interval with output to a stream. It uses format/3 so extended stream options, e.g., atom(A), are supported.

deprecated: - use format/3

enumerate(?Term:term_List) is nondet

Succeeds non-deterministically by enumerating values of any Term. Enumerating is defined as follows:

enumerating a list enumerates each item in the list (in that order).
enumerating an interval of type integer sets the interval to each value in the domain starting from the lower bound.

enumeration of any other term succeeds with no other effect. Examples:

?- X::integer(1,2), enumerate(X).
X = 1 ;
X = 2.

?- Is=[X,Y], Is::integer(1,2), enumerate(Is).
Is = [1, 1],
X = Y, Y = 1 ;
Is = [1, 2],
X = 1,
Y = 2 ;
Is = [2, 1],
X = 2,
Y = 1 ;
Is = [2, 2],
X = Y, Y = 2.

?- X::real, enumerate(X).
X::real(-1.0e+16, 1.0e+16).

?- enumerate(sam).
true.

?- B::boolean, enumerate([42,-1.0,B,Z]).
B = 0 ;
B = 1.

small(+X:numeric_List) is semidet

Succeeds if the width of the domain of Numeric is less than the value defined by the environment flag clpBNR_default_precision which is a positive integer specifying number of digits; otherwise fails. For example, a clpBNR_default_precision value of 6 (the default) defines a domain width limit of 1e-7. Numbers have a domain width of 0, so they are always "small".

If Numeric is a list of numerics, all elements of the list must be "small". Examples:

?- X::real, small(X).
false.

?- X::real(-1e-10,1e-10), small(X).
X::real(-1.0000000000000002e-10, 1.0000000000000002e-10).

?- X::real(-1e-10,1e-10), small([X,42]).
X::real(-1.0000000000000002e-10, 1.0000000000000002e-10).

Note that this is really only useful for real intervals; integer intervals are not small until they become point values.

small(+Numeric, +Precision) is semidet

Succeeds if the width of the domain of Numeric is less than the value defined by P (digits of precision); otherwise fails. As with small/1, Numeric can be a single numeric or list of numerics.

See also: - small/1

global_maximum(+Exp, ?Z:numeric) is semidet

The maximum allowable width for the generated maximum is determined by the current default precision (environment flag clpBNR_default_precision). Example:

?- X::real(0,3r4*pi), global_maximum(X*sin(4*X),Z).
X:: 1.994...,
Z:: 1.97918... .

Note that intervals in the expression may not narrow significantly if more than one maximum can found using the the initial domains. In such cases, additional "searching", e.g., using solve/1, may be necessary.

See also: - global_maximize/2

global_maximum(+Exp, ?Z:numeric, +Precision:integer) is semidet

Same as global_maximum/2 with additional argument defining precision (overrides environment flag clpBNR_default_precision). Example:

?- X::real(0,3r4*pi), global_maximum(X*sin(4*X),Z,4).
X:: 1.99...,
Z:: 1.979... .

See also: - global_maximum/2

global_minimum(+Exp, ?Z:numeric) is semidet

Succeeds if Z unifies with the global minimum of (evaluated) expression Exp; otherwise fails. This is analogous to global_maximum/2 for finding minima. See global_maximum/2 for more details. Example:

?- X::real(0,1r2*pi), global_minimum(X*sin(4*X),Z).
X:: 1.228...,
Z:: -1.203617... .

See also: - global_maximum/2, global_minimize/2

global_minimum(+Exp, ?Z:numeric, +Precision:integer) is semidet

Same as global_minimum/2 with additional argument defining precision (overrides environment flag clpBNR_default_precision). Example:

?- X::real(0,1r2*pi),global_minimum(X*sin(4*X),Z,4).
X:: 1.23...,
Z:: -1.203... .

See also: - global_minimum/2

global_maximize(+Exp, ?Z:numeric) is semidet

Succeeds if Z unifies with the global maximum of (evaluated) expression Exp; otherwise fails. Exp must be an actual expression which can be evaluated as part of the search for the global optima, not an interval equal to the evaluated expression. Any solutions will be constrained to be at the global maximum (Z), and the intervals in Exp will be narrowed to the maximizers found for the global maximum value. If the global maximum satsfies multiple sets of maximizers, they will be lost. Finding a single set of maximizers is often sufficient for many practical problems.

The maximum allowable width for the generated maximum is determined by the current default precision (environment flag clpBNR_default_precision). Example:

?- X::real(0,3r4*pi), global_maximize(X*sin(4*X),Z).
X:: 1.994666...,
Z:: 1.97918... .

See also: - global_maximum/2

global_maximize(+Exp, ?Z:numeric, +Precision:integer) is semidet

Same as global_maximize/2 with additional argument defining precision (overrides environment flag clpBNR_default_precision). Example:

?- X::real(0,3r4*pi), global_maximize(X*sin(4*X),Z,4).
X:: 1.9947...,
Z:: 1.979... .

See also: - global_maximum/3

global_minimize(+Exp, ?Z:numeric) is semidet

Succeeds if Z unifies with the global minimum of (evaluated) expression Exp; otherwise fails. This is analogous to global_maximize/2 for finding minima and a single set of minimizers. See global_maximize/2 for more details. Example:

?- X::real(0,1r2*pi), global_minimize(X*sin(4*X),Z).
X:: 1.228295...,
Z:: -1.203617... .

See also: - global_maximize/2

global_minimize(+Exp, ?Z:numeric, +Precision:integer) is semidet

Same as global_minimize/2 with additional argument defining precision (overrides environment flag clpBNR_default_precision). Example:

?- X::real(0,1r2*pi),global_minimize(X*sin(4*X),Z,4).
X:: 1.2283...,
Z:: -1.203... .

See also: - global_minimize/2

solve(X:numeric_List) is nondet

Succeeds if a solution can be found for all values in X where the resultant domain of any value is narrower than the limit specified by the default precision (number of digits as defined by the environment flag clpBNR_default_precision); otherwise fails. This is done by splitting any intervals in round robin order of their widths until all domains are smaller than the required limit. Splitting can only be done at points not in the solution space (unlike splitsolve/1); this avoids the splitting a single solution range into multiple solutions (although this can still occur for other reasons). Other solutions can be generated on backtracking. Examples:

?- X::real, {17*X**256+35*X**17-99*X==0}, solve(X).
X:: 0.0000000000000000... ;
X:: 1.005027892894011... .

?- [X,Y]::real, {X+Y==1,X-Y==1}, solve([X,Y]).
X:: 1.0000000000000...,
Y::real(-4.96269692007445e-14, 4.96269692007445e-14).

The two main use cases for solve/1 are a) to separate multiple solutions in a within domain (or set of domains), and b) to overcome the well known dependancy issue when using interval arithmetic. (In clpfd terminology, solve/1 is a labelling predicate.)

See also: - splitsolve/1

solve(X:numeric_List, Precision:integer) is nondet

Same as solve/1 with precision defined by Precision.

See also: - solve/1

splitsolve(X:numeric_List) is nondet

Succeeds if a solution can be found for all values in X where the resultant domain of any value if narrower than the limit specified by the default precision (number of digits as defined by the environment flag clpBNR_default_precision); otherwise fails. This is done by splitting any intervals in order of their widths until all domains are smaller than the required limit. Other solutions can be generated on backtracking.

Normally solve/1 is a better choice but this predicate can be used when solve/1 cannot find a suitable non-solution value to use to split an interval (or intervals). This predicate is also less computationally expensive, but may result in many solutions being produced for a single wider interval. (This is why solve/1 splits on non-solutions.)

See also: - solve/1

splitsolve(X:numeric_List, Precision:integer) is nondet

Same as splitsolve/1 with precision defined by Precision.

See also: - splitsolve/1

absolve(X:numeric_List) is semidet

Succeeds if a if X is a numeric (or list of numeric); otherwise fails. absolve is intended solely to trim up the boundaries of what is essentially a single (non-point) solution to a problem. The strategy used is to work in from the edges of the interval ("nibbling away") at subdomains which are inconsistent until you cannot go farther, then reduce step size and resume nibbling. In this case, the the environment flag clpBNR_default_precision is used to specify the number of step size reductions to apply; the initial step size is half the interval width.

absolve can be used to further narrow intervals after solve to "sharpen" the result; example:

?- X::real,{X**4-4*X**3+4*X**2-4*X+3==0},solve(X).
X:: 1.000000... ;
X:: 1.0000000... ;
X:: 3.000000... ;
X:: 3.000000... ;
false.

?- X::real,{X**4-4*X**3+4*X**2-4*X+3==0},solve(X),absolve(X).
X:: 1.00000000... ;
X:: 3.00000000... ;
false.

absolve(X:numeric_List, Precision:integer) is nondet

Same as absolve/1 with precision defined by Precision.

See also: - absolve/1

partial_derivative(+Exp, -X, ?Drv) is semidet

Suucceds if the (symbolic) partial derivative of Exp with respect to variable X is Drv; otherwise fails. The syntax of Exp is determined by normal clpBNR (and Prolog) arithmetic syntax. Examples:

?- partial_derivative(X**2,X,Drv).
Drv = 2*X.

?- partial_derivative(X/Y,X,Drv).
Drv = 1/Y.

?- partial_derivative(X/Y,Y,Drv).
Drv = -1*X/Y**2.

?- partial_derivative(max(X,Y),Y,Drv).
false.

This predicate can be used in generating additional constraints, e.g., local optima with a gradient of 0, or in constructing meta-contractors like the Taylor series contractor described in the User Guide.

watch(+X:interval_List, +Action:atom) is semidet

Succeeds if X is an interval and Action is an atom; otherwise fails. If successful, and Action is not none, a watchpoint is placed on X. Watchpoints are only "actioned" when the debug topic clpBNR is enabled. If Action = log a debug message is printed when the interval doman narrows. If Action = trace the debugger is invoked. If Action = none the watchpoint is removed.

trace_clpBNR(?B:boolean) is semidet

Succeeds if B can be unified with the current value of the clpBNR trace flag or if the trace flag can be set to B (true or false); otherwise fails. If the trace flag is true and the clpBNR debug topic is enabled, a trace of the fixed point iteration is displayed.