Expressive Applications of Constraint Logic Programming

Abstract

This introduction to the Constraint Logic Programming languageCLP([real]) uses applications to provide insight tothe language's strengths. An overview of CLP([real])is followed by a discussion of three applications that illustratethe language's unifying treatment both of numeric and symboliccomputing and of engineering analysis and synthesis problems.Another discussion dissects the interpreter's constraint solverand clarifies how a problem's search space can be restricteddeclaratively. The final example is an extended description ofthe construction of a network of interpreters, which can be usedto distributively solve a set of linear equations. This extensionrequires no modification of the CLP([real]) interpreterand points out the benefits of revisiting established algorithmsvis-a-vis CLP([real]).

Artificial Intelligence Review 11: 427452, 1997. 427

c[n0d] 1997 Kluwer Academic Publishers. Printed in the Netherlands.Expressive Applications of Constraint Logic ProgrammingWILLIAM B. DAYDepartment of Computer Science and Engineering, Auburn University, AL 36849-5347Abstract. This introduction to the Constraint Logic Programming language CLP(R) uses

applications to provide insight to the languages strengths. An overview of CLP(R) is followed

by a discussion of three applications that illustrate the languages unifying treatment both of

numeric and symbolic computing and of engineering analysis and synthesis problems. Another

discussion dissects the interpreters constraint solver and clarifies how a problems search space

can be restricted declaratively. The final example is an extended description of the construction

of a network of interpreters, which can be used to distributively solve a set of linear equations.

This extension requires no modification of the CLP(R) interpreter and points out the benefits

of revisiting established algorithms vis-a-vis CLP(R).Key words: constraint logic programming, symbolic and numeric computing, declarative

programming, constrained search, linear equations, distributed problem-solving, constraint

solver.1. IntroductionConstraint Logic Programming (CLP) refer to the class of logic-based languages in which the expressive power has been extended by generalizing

unification. These languages have been used successfully in various applications that involve both symbolic and numeric computations. In symbolic

calculations a CLP languages ability to declare constraints provides a clear

and concise means for restricting a problems search space. In numeric calculations a CLP language can be an effective tool for distributively solving

an enormous set of linear equations. This paper is an introduction to CLP

languages and how they provide new insights into solving familiar problems.Section 2 of this paper is an overview of CLP, including the motivation

and design goals of the specific language CLP(R). Short examples are also

given in this section to illustrate new and useful features of the language.

The CLP(R) architecture is added to complete a description of the systems

basic structure. Section 3 presents three application areas of CLP(R): (1)

financial analysis, (2) problems from electrical engineering, and (3) analysis

and partial synthesis of truss structures. These two sections are included not

for their originality, but to offer a simple introduction to CLP languages and

to show the practical diversity of characteristic applications.428 WILLIAM B. DAYIn [JM94] CLP applications are divided into two general classes: (1)

combinatorial search problems, and (2) transparent representation problems.

Sections 4 and 5 of this paper deal with two new applications of CLP(R),

one from each class. An example of a combinatorial search problem appears

in Section 4. There the Magic Square of Primes is introduced and a comparison is made between the efficiency of a Prolog programs solution and

a CLP(R) programs solution. In particular, a detailed explanation of the

constraint solver of CLP(R) is described using the explicit constraints of the

program for the Magic Square of Primes. The constraints of this example are

symbolic. Numerically, CLP(R) reduces the number of comparisons needed

in naive search from 268 to fewer than 105.In Section 5 the role of CLP(R) in a distributed environment is discussed.

It is assumed that a copy of the interpreter is available to each node in

the network. In this section a classical numeric constraint problem, finite

difference approximations, is solved distributively on a network of CLP(R)

interpreters. This form of distribution is often ignored by traditional finite

difference algorithms although the CLP representation is transparent. Section

6 contains the conclusions of this paper.2. Basics of Constraint Logic ProgrammingPredicates are the building blocks of all logic programming languages. For

example, p(X, Y, 5) is a predicate with three arguments, the first two of which

are variables (as indicated by capital letters X and Y). From predicates are

constructed rules. For example,p[n28]X ; Y [n29] : r[n28]X [n29]; s[n28]Y [n29]:This rule means that predicate p(X, Y), is true IF (:) r(X) AND (,) s(Y) are

true. In the special case that no predicates follow the IF, it is omitted and the

rule is a fact; e.g., q(16). All facts and rules are terminated by a period.A logic program consists of facts and rules. Questions are posed by a user,

and the languages interpreter attempts a proof by resolution; i.e., the question

is denied and a contradiction is found. When the resolution is successful, the

variables of the question are typically bound to specific values. At each step of

resolution, the interpreter examines the rules and facts of the program in order

from the first to the last and resolves each rule from left to right, beginning

after the IF.Some common built-in predicates are included by the interpreter. For example, X is Y + 2 is an assignment to X of the value of Y added to 2. However,

X=Y is a boolean comparison of the values of variables X and Y.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 429The merger of constraint satisfaction with logic programming to form

CLP languages caused much excitement in computer science because of

its conglomeration of such diverse areas as numerical analysis, artificial

intelligence, operations research, logic, mathematics, and formal languages

(Cohen 1990). The combination of logic programming with constraint satisfaction is a natural extension of traditional logic programming languages

since the fundamental unification process can be viewed as a special case of

constraint satisfaction; i.e. unification of f(X, a) with f(b, Y) is equivalent to

solving f(X, a)=f(b, Y).CLP is a class of logic programming languages that are modeled on rulebased constraint programming. The specific instant CLP(R) uses the real

numbers (R) as its structure of computation. In particular, the constraints are

relations over real arithmetic terms with uninterpreted functors (Jaffar and

Michaylor 1987). Since all CLP languages are soundly based on a single

framework of formal semantics there are no restrictions in CLP(R) to the

Herbrand Universe, to unification, or to allowed equations.There are several features of CLP(R) that are noteworthy. The most important features are briefly described here and illustrated more fully later. The

output of CLP(R) programs is symbolic; thus, it is not only possible to answer

queries with responses like X=6 or Y=bill, but one may also have responses

like X*X = Y+4.CLP(R) expands the use of declarative programming. For example, in using

an equation of the form X = Y + 4, it is no longer necessary to write different

rules, depending on whether X or Y is instantiated (has a known value) at the

time that the equation is invoked. This equation is merely a constraint, and

CLP (R) is able to deal with all four possible instantiation cases: neither X

nor Y is ground (i.e., bound to a constant), X is ground, Y is ground, or both

X and Y are ground. From the programmers view, only a simple statement

of the equation is needed. CLP(R) then makes the appropriate interpretation

during execution.CLP(R) is able to run Prolog programs with no changes. Indeed, one of

the early design goals of CLP(R) was not only to execute Prolog programs

but to have them run as efficiently as normal. This was done in CLP(R) by

separating the interpreter into an inference engine and a constraint solver. In

running a Prolog program, there should be no need for the constraint solver.

Consequently, the inference engine is essentially a Prolog interpreter.A second fundamental design goal of CLP(R) dictated that only linear

constraints would be considered by the constraint solver. This limitation was

ameliorated by defining delay and wake-up conditions. Basically, simplification of a non-linear constraint is delayed until a sufficient number of the

variables are instantiated to make the constraint linear.430 WILLIAM B. DAYA third design goal required CLP(R) to deal with constraints incrementally.

This means that in the execution of a query it is not necessary to verify the

solvability of an entire set of constraints when a single new constraint is added.

Rather it is necessary only to verify that the new constraint is compatible with

the old set of constraints.Other design goals of CLP(R) included choosing appropriate algorithms for

testing the satisfiability of systems of constraints, demanding that the system

be able to simplify a system off constraints, and finding a canonical form

to represent internally the constraints. These latter objectives deal with the

very practical issues of constructing an efficient implementation. Additional

implementation details can be found in Jaffar et al. (1992).2.1. Constraints versus AssignmentsHere an example is given to illustrate the increased power of expression

available to CLP(R). Considr defining the fibonacci series in Prolog and in

CLP(R), each with two facts and one rule.in Prolog in CLP(R)fib(0,1). fib(0,1).fib(0,1). fib(1,1).fib(N,R): fib(N,R1 + R2):N1 is N 1, N>=2.fib(N1,R1), fib(N 1,R1)N2 is N 2, fib(N 2,R2).fib(N2,R1),R is R1 +R2.The predicate fib uses the first argument for the index N and the second

argument for the value of the Nth Fibonacci number.Both of these programs can handle the question ?-fib(10,F) by answering

F=89, but only the CLP(R) program can answer the question ?-fib(N,89)

with N=10. This example illustrates several features of CLP(R). First, the

second program is more declarative: to find the Nth Fibonacci number if

N is larger than one, add the values of the (N-1)st and (N-2)nd Fibonacci

numbers. Secondly, the CLP(R) definition allows a reversal of roles of the

input variable and the output variable. Although this is common with

many other Prolog predicates (compare, the standard definitions of member

and append), it is not the original intention here. Finally, the CLP(R)version

allows arithmetic operations to occur in the arguments of a predicate. This is

a radical notion only for Prolog since most procedural languages generally

provide such facilities.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 4312.2. The CLP(R) InterpreterArithmetic terms in CLP(R) are constructed from real constants, variables,

and these interpreted functors: +, , *, /, sin, cos, tan, pow (the power

function), min, and max. The CLP(R) interpreter is composed of five parts:(1) the inference engine, (2) the interface between the inference engine and

the constraint solver, (3) the equation solver, (4) the inequality solver, and (5)

the output module. A complete description of the parts is given in Jaffar et

al. (1992) and is summarized here.The constraint solver itself is composed of the equation solver and the

inequality solver. The inference engine is basically a Prolog interpreter, which

controls the execution of the derivation steps and updates variable bindings.

The usual stack mechanisms and symbol tables are maintained here. The

interface modules primary purpose is to smooth the transition between the

inference engine and the pieces of the constraint solver. The interface evaluates complex arithmetic expressions and transforms constraints into standard

forms.Constraints that are received by the constraint solver occur in one of the

three forms: linear equations, linear inequalities, or non-linear equations.

Linear equations are sent to the equation solver while linear inequalities are

sent to the inequality solver. Non-linear equations are delayed until so many of

its variables are ground that the equations become linear. The output module

finally converts the internal representations of constraints to a canonical

output form.As progress is made with the linear constraints, the involved variables

become ground. It is then possible to check if these variables are involved in

any of the delayed non-linear constraints. After making appropriate adjustments to the non-linear constraints, a queue is formed of all non-linear constraints that have become linear. These linear constraints are then passed to

either the equation solver or the inequality solver before proceeding with the

next derivation step of the proof. Chaining effects in which one linearized

constraint causes another non-linear constraint to become linear are possible.A simple bookkeeping device is used to determine when backtracking in

the constraint solver is necessary. In general, backtracking in the constraint

solver always corresponds to a backtracking point in the inference engine but

not vice versa.Linear equations are stored in parametric form. Thus, a constraint like X

= Y, which involves program variables X and Y may be represented as these

two parametric equations: X = t and Y = t. In general, a program variable will

be converted to this parametric form: X = b + c1t2 + ::: + cktk. A detailed

example of these parametric equations is given in Section 4.432 WILLIAM B. DAY3. Applications of CLP(R)CLP(R) has been used for an ever increasing number of applications. Three

examples are cited in this section to display this diversity: (1) options trading

on the stock market, (2) electrical engineering circuit analysis, and (3) truss

analysis and synthesis of structures in civil engineering.3.1. CLP(R) in Financial AnalysisOptions are used to express confidence levels in the intrinsic value of an asset.

For example, one may make a call (put) option to buy (sell) a stock before

some specific date. Should the stock drop very low before the specified date,

the owner of the call option would typically purchase the stock at its low

point. In the worst case, the stock increases monotonically, and the owner

must buy the stock at a higher price when the specified date arrives.Two features of options trading have been cited (Lassez 1987) that make it

an appropriate application of combining symbolic and numeric computation.

They are that traders typically combine (1) their personal, expert heuristics

with (2) numeric valuation functions in decision-making.The indifference of rules in CLP(R) to dictating which variables are input

and which variables are output is an attractive characteristic to option trading

since it is not necessary to write additional rules for reversing the direction

of computation. Furthermore, in CLP(R) it is possible to unite symbolic output and constraints on goals to obtain what if analyses. For example, the

following query describes a straddle, which occurs when a call option and a

put option are both purchased, (Lassez 1987)? PosVal = V1 + V2, Call=4, Put=3, Int=0.05, Ex=50, PosVal

>=3,payoff(call, buy, S, Call, , Int, Ex, , V1),payoff(put, buy, S, , Put, Int, Ex, , V2).Here S is the stocks current price, Int is the current risk-free interest rate, and

Ex is the exercise price. The underscore in the argument of predicate payoff

denotes a dont care value for the variable. This answer is displayed by the

interpreter:PosVal = 42.65 S, 0< S <= 39.65,PosVal = S 57.35, S>=60.35.The solution indicates that a profit exists only when the stocks price lies

within a window.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 4333.2. CLP(R) in Electrical EngineeringA second area of applications using CLP(R) has been with some electrical

engineering problems (Heintze 1986). Three categories of these problems are

lumped circuits, digital filtering networks, and electronmagnetic fields.Circuit problems include both analysis and synthesis of analogue circuits,

in which the constraints are either local (e.g., Ohms law for a resistor) or

global (e.g., Kirchoffs law in nodal analysis). Specific examples of these

problems are analysis of steady-state RLC circuits and diode components in

RLC circuits. Bipolar junction transistors for amplifier and logic circuits have

also been analyzed and designed.Consider the DC circuit with two resistors R1 and R2 in parallel and currents

I1 and I2 flowing through the resistors, respectively. Then the simple ruleresistor(V,I,R) : V = I*R.allows numerous queries to be posed. The answers will depend on which

variables are instantiated. For example,? V=5, R1=50, R2=25, resistor(V,I1,R1), resistor(V,I2,R2),

I=I1+I2.For this query the interpreter will respond with these values for the currents:

I1=0.1, I2=0.2, I=0.3. Current I is the current that flows through the power

source of 5 volts. Changing the first line of this query to R1=50, R2=25,

I1=10 will yield V=500, I2=20, I=30. Finally, with the first line of the query

set to R1=50, I1=10, the answer involves constraints and is this: V=5000,

R2*I2=500, I=10+I2.In the area of digital signal flow, the simulations of a linear shift-invariant

digital system with the independent variable time have been successfully

modeled in terms of linear signal-flow graphs. In these graphs signal values

at the nodes at successive time steps are simulated by summing constraints

for each node and by collecting branch equations. In CLP(R) one needs only

to declaratively state the equations while the interpreter decides the particular

method for solving the equations.The third category of electrical engineering problems involves numerical,

approximate solutions for the pertinent partial differential equations of electromagnetic fields. A finite difference approximation is the common approach

for these problems. The examples of this category illustrate the elegance of

CLP(R) in using finite difference methods by its concise representations and

its use of the constraint collection mechanism. Section 5 discusses in detail a

distributed solution for finite difference methods.434 WILLIAM B. DAY3.3. CLP(R) in Truss StructuresAnalysis and synthesis problems are common to all areas of engineering. Recently Lakmazaheri and Rasdorf (1989) has shown that constraint

logic programming is suitable for the analysis and partial synthesis of truss

structures. Typically these two problem areas of analysis and synthesis are

considered antithetical in engineering; however, a CLP(R) representation is

able to use only one representation for the two. A particular query dictates

the appropriate class.This work on truss structures has shown that the technique of local propagation is inadequate for satisfying simultaneous constraints when the associated

constraint graph contains cycles. But in the constraint-based approach the

constraints can be viewed as relations rather than functions; hence, there is

no explicit commitment to a specific variable being either input or output.The constraints associated with structural components are classified by

these categories: (1) truss components, which satisfy a two or three dimensional form off the vector equation KD=F where K is an n by n (n=2 or 3)

stiffness matrix, D is the n by 1 displacement vector at the two end nodes and

F is the n by 1 force vector at the two end nodes; (2) support components of

the type pin support, rollerX support, or rollerY support; and (3) load components, each of which is associated with a force vector and a displacement

vector.The position and displacement vectors are functions, and a nodes force

history is a list. Nodes in two dimensions are terms of this form:n(p(Px, Py), d(Dx,Dy), H)where Px and Py are the x and y coordinates of the position vector p, Dx and

Dy are the x and y coordinates of the displacement vector d,and H is the node

force history list.One can then describe general rules for the local nodal behaviors. This is

analogous to Ohms law in electrical circuits. For example, this ruleA valid node exists if its coordinates are valid.is represented by this CLP(R) rulenode(n(p(Px,Py), , )): ordinate(Px),ordinate(Py).Since the behavior of a structural model is derived from the behavior of its

structural components, rules can be added that insure the satisfiability of the

constraints associated with the structural components.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 435Four examples are presented in Lakmacaheri and Rasdorf (1989). The

first example illustrates the capability of the constraint logic program to

perform numerical analysis of a statically indeterminate truss structure. The

second example shows the capability of the program to perform symbolic

computation for structural analysis. The third example gives a partial synthesis

of a simple truss structure using the same program. Finally, the last example

demonstrates the applicability of the program for solving more complex

problems.The main contribution of this work is the development of a general and

declarative constraint-based formulation for the analysis and partial synthesis

of truss structures using the framework of the constraint logic programming

methodology. A single and uniform framework for modeling engineering

design problems facilitates the integration of the different design activities of

synthesis, analysis, and evaluation.4. The Magic Square of PrimesThe Magic Square of Primes is perfect for illustrating how CLP(R) stores

and uses constraints which are linear equations. Here is a condensed version

of the Prime Magic Square program. It consists of twenty-six facts (the

first twenty-six prime numbers), eight constraints (which require the rows,

columns, and diagonals of a three-by-three matrix to have a constant sum),

and eight requirements that the eight unknowns be prime numbers. Only the

(3,2)-element is known a priori; it is one. the [::: ] notation is used for a

list.p(2). p(3). p(5). p(7). p(11). p(13).p(17). p(19). p(23). p(29). p(31). p(37).p(41). p(43). p(47). p(53). p(59). p(61).p(67). p(71). p(73). p(79). p(83). p(89).p(97). p(101).ms([A1,A2,A3], [[B1,B2,B3],[C1,1,C3]) :Sum=A1+A2+A3, Sum=B1+B2+B3,Sum=C1+1+C3, Sum=A1+B1+C1,Sum=A2+B2+1, Sum=A3+B3+C3,Sum=A1+B2+C3, Sum=A3+B2+C1,p(A1), p(A2), p(A3), p(B1),l p(B2), p(B3), p(C1), p(C3).As the constraints are encountered in the rule, they are converted to

CLP(R)s internal format, using parametric variables. In the first constraint,436 WILLIAM B. DAYSum=A1+A2+A3, all four variables are initially uninstantiated. Thus, in a

left-to-right order the constraint solver assigns Sum to parametric variable ti,

A1 to the new variable t2, A2 to the new variable t3, and A4 to the new variable t4. The equation is then rewritten in the parametric variables in a normal

form which puts the lowest numbered parametric variable on the left-hand

side. At the end of this standardization step, the substitutions and the equation

would be as follows:Sum = t1 t1 = t2 + t3 + t4A1 = t2A2 = t3A3 = t4Since parametric variables are defined in terms of mathematical equations,

there is no need to use capital letters as is the case for CLP(R) variables.Now the constraint solver encounters the second constraint: Sum=B1+

B2+B3. Since the variable Sum has already been numbered, the parametric

variable numbering continues with B1, B2, and B3. They are assigned t5, t6

and t7, respectively. The parametric form of this constraint is initially t1 = t5

+ t6 + t7. But t1 is substituted, and the resulting equations are solved for t2,

the lowest numbered subscript. This leads to the following structure:Substitutions: Sum=t1,A1=t2,A2=t3,A3=t4,B1=t5, B2=t6, B3=t7Equations: t2 = t3 t4 + t5 + t6 + t7t1 = t5 + t6 + t7Independent Variables: ti, i = 3, ::: , 7.After converting the third constraint (t1 = t8 + 1+ t9), the following structure

is obtained:Substituations: Sum=t1,A1=t2,A2=t3,A3=t4,B1=t5, B2=t6, B3=t7, C1=t8,C2=t9Equations: t5 = t6 t8 + 1 + it9t2 = t3 t4 + t8 + 1+ t9t1 = t8 + 1+ t9Independent Variables: t3, t4, t6, t7, t8, t9.Continuing in this manner, one finds after converting the fourth constraint,

four equations for the four dependent, parametric variables t3, t5, t2, and t1.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 437Addition of the fifth constraint adds t4 to the list of dependent, parametric

variables. At the sixth step of this processes there is no change; i.e., the sixth

constraint is redundant. This is a surprising and delightful development since

it makes the problem even more under-determined than was assumed: nine

original variables (including Sum) and only seven linearly independent equations. After two more steps, one arrives at this final structure:Substitutions: same as in the third stepEquations: t7 = 4/3 t8 2/3 t9 + 1/3t6 = 1/3 t8 + 1/3 t9 + 1/3t4 = 1/3t8 + 2/3 t9 + 2/3t3 = 2/3 t8 + 2/3 t9 1/3t5 = 2/3 t8 + 4/3 t9 + 1/3t2 = 2/3 t8 1/3 t9 + 2/3t1 = t8 + t9 + 1Independent Variables: t8, t9.It is at this point that checks begin to determine if each of the eight original

variables (excluding Sum) is a prime number. Because of the stored, parametric form of the constraints, it is not difficult to see how CLP(R) eliminates

most of the billions of choices.Consider the initial case where the interpreter guesses A1 = 2 and A2 =2. Internally, this requires t2 and t3 = 2. Hence, the constraint solver simultaneously solves these two equations:2 = 2/3 t8 1/3 t9 + 2/32 = 2/3 t8 + 2/3 t9 1/3finds t9 =1, t8 = 5/2, and ultimately A3 = 1/2/. Up to this point, the interpreter

has verified that A1 is a prime and that A2 is a prime. When it now tries o

confirm that A3 is a prime, it fails. In backtracking, it will be unnecessary to

check any other combinations which have A1 and A2 set to 2. This is a huge

savings since a naive search would be required to check all 26 possibilities

for the other six original variables; this is roughly 309 million cases.The savings can be even more impressive if one considers only the case A1

= 2. Here naive search requires checking 267 cases, about 8 trillion. In reality,

CLP(R) checks fewer than fifty of these before abandoning the quest for a

solution with A1 = 2. One can check these numbers with the trace facility

of CLP(R), but it is more enlightening to see exactly what the interpreter is

doing in these test cases.438 WILLIAM B. DAYTable 1.A2 A3 B1 B2 B3 C1 C32 1/23 15 2 3 3 3 47 3 5 411 5 913 617 819 923 11 2129 1431 1537 1841 2043 2147 23 4553 2659 29 5761 3067 3371 3573 3679 3983 41 8189 4497 48101 50In general, if A1 is set to two and A2 is set to the constant K, then one

obtains the following solutions from the six equations which relate the eight

parametric variables:t9=C2=K 1 t8=C1=(K+3)/2t4=A3=(K 1)/2 t5=B1=K 2t6=B2=(K+1/2 t7=B3=3Ultimately, all twenty-six prime facts must be compared with an assigned

value of one of the original variables in order to show that the value is not

a prime. For the cases of A1 = 2, the required assignments are shown in

Table 1.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 439The number of comparisons that must be made for integer N is P in the

case that N is not a prime. Here P is 26. In 19 of the 26 rows this is the case.

In 5 of the 26 rows there are 2P comparisons since during backtracking it is

necessary to scan the remaining list of prime facts; e.g., in the fifth row

A2=11, A3+5, and B1=9.During Forward Execution, it is verified that A3=5 is a prime after 3

comparisons, and it is verified that A3=9 is not a prime after 26 comparisons.

But on back-tracking it is necessary to go through the other 23 comparisons

for A3 to be assured that there is no other fact or rule for p(5). Thus, 52=2P

comparisons are needed. One sees in this case that 37P comparisons are need

to abandon A1=2 as a possible solution.In the worst case for any choice of A1, the constraint Solver would have to

make 6P comparisons on each of the 26 rows of table like that shown above.

Thus 156P times the 26 choices for A1 exhausts the possible solutions. This

number is approximately a hundred thousand, and is well below 268, the

number of comparisons required in naive search, which a Prolog program

must use. This is caused by the inability of Prolog to handle the code given

above; Prolog must check all combinations.Since this problem actually has two solutions the comparisons given here

are the total number required to find all solutions. Empirically, CLP(R) finds

both solutions in less than a minute while Prolog had not found the first

solution after thirty minutes.5. Distributed CLP(R) in Banded MatricesThe example in this section examines solving a familiar numeric problem with

a distributed form of CLP(R). This solution is predicated on the expressive

ability of the CLP(R) to answer queries symbolically. Specifically, CLP(R)

solves under-constrained matrix equations by returning a reduced set of constraints. Although this method of solving matrix equations can be used with

other programming languages, it requires construction of a mechanism for

solving and recording the solutions of under-constrained problems. In CLP(R)

this construction is gratis.Matrix operations are often the heart of numerical problem solving, and

matrix inversion is one of the most fundamental of these operations. Beyond

inverting a general matrix, banded matrices, which have the nonzero elements

lying within a narrow band of the diagonal, have been studied for both

theoretical and practical benefits. The numerical example presented below is

that of a banded matrix. The bandwidth is defined as the maximum over all

rows of the difference in indices between the first and last nonzero elements

of the row. It is common to assume that the order of the matrix (the larger of440 WILLIAM B. DAYnumber of rows and the number of columns) is orders of magnitude larger

than the bandwidth. We do not make his assumption, although such banded

matrices are more amenable to decomposition algorithms.Consider the problem of solving thirty equations with thirty unknowns.

Assume the nonzero elements for this system form a tri-diagonal matrix.

Then subdividing the system into thirds yields these three subsystems: (1) the

first ten equations, using the first eleven variables; (2) the middle ten equations

using twelve variables, number ten through number twenty-one; (3) the last

ten equations using the last eleven variables. Variables ten and eleven appear

in the first two subsystems, and variables twenty and twenty-one appear in

the last two subsystems. All other variables appear in only one of the three

divisions. Since each of the three parts is underconstrained, it is possible to

solve each partially by expressing the unknowns in terms of the four shared

variables.Thus variables one through ten can all be expressed in terms of variable eleven in the first division. Variables eleven through twenty can all be

expressed in terms of variables ten and twenty-one in the second division.

Variables twenty-one through thirty can all be expressed in terms of variable

twenty in the third division. Furthermore, these three reductions can be done

concurrently.After completing the three pieces, it is possible to use four particular

equations to construct the final solution. One of these equations (variable ten

in term of variable eleven) comes from the first set. Two equations (variables

ten and twenty-one in terms of variables eleven and twenty) come from the

second set. One equation (variable twenty-one in terms of variable twenty)

comes from the third set. Finally, solving these four equations in four variables

provides the basis for substitutions needed for the solution. In general, N

equations solved on k processors will concurrently solve k subdivisions of

N/k equations in N/k + (1 or 2) variables and then produce 2k 2 equations in

2k 2 variables to be solved.This section begins with an introduction to the Parallel Virtual Machine

system, which has been used in the implementation to avoid writing de novo

code for datagram sockets. Implementation details are described including a

discussion of the parabolic, second-order partial differential equation (PDE)

that is used as the primary application and a discussion of the architecture

of the system, the drones, and the queen. Experimental results conclude this

variation of CLP(R).5.1. The Parallel Virtual Machine (PVM) SystemThe PVM system allows the designer of a distributed system to avoid dealing

with communication details of remote processors (Sunderam 1990). TheEXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 441designer is thereby allowed to concentrate on constructing a model that

appears to be running on only one processor. This is accomplished through a

library of fewer than two dozen functions that coordinate daemons running

on remote machines and enable programs to send and received messages

between remote processors. The programmers code can ignore opening

sockets between processors, designing a communication protocol, and similar

issues.For example, consider a user who is logged into one machine (elm) and

wishes to use three other machines (peach, aspen and spruce). It will be necessary to have each of the programs which will run concurrently to be resident

on the correct machine, but the entire system will be initially linked by using

the following command:elm: pvmd HOSTFILE &Here the executable file pvmd begins running a pvm daemon on elm and, one

at a time, opens a socket and initiates another pvm daemon to run on each of

the remote machines listed in HOSTFILE. The character & in the command

line puts the pvmd process in the background. At this point the main program

residing on elm can be started: the client will initiate server programs on

the remote machines. This is done with the initiate function in the client

program.Basically each transmission of data between machines require three function calls by each of the participants. The message sender uses (1) initsend(1)

to prepare for packing a new message, (2) putstring(message) to load the

message into the sockets buffer, and (3) snd(name, inum, type) to transmit

the message. In function snd, the name and inum variables are the name and

instance number of the process to which to send the message and type is the

user-defined type of the message.At the receiving end, one needs (1) rcv(type) to remove a message from

the processs receiving queue, (2) rcvinfo(length, type, proc, inum) to obtain

information about the last message, and (3) getstring(buffer) to unload the

message from the sockets buffer. Here length is the length of the message in

bytes, type has been described, and proc and inum hold the name and instance

number of the sending process.The PVM system is much more general than has been described with this

example. It can be used on a heterogeneous network, in which the individual

machine architectures may be one of twenty different choices; e.g. SUN4,

BBN Butterfly TC2000, Intel IPSC/2, NeXT, etc. Other features of PVM can

be found in (Sunderam 1990).442 WILLIAM B. DAY5.2. A Parabolic Second-Order Partial Differential EquationIn Fox et al., (1988) are examples where a finite difference matrix has order

N = O(m2) and band-width w = O(m). These solution techniques apply only

if the given problem is completely determined; i.e., one in which N equations

and N unknowns have a unique solution. The common decompositions and

algorithms discussed in Chu and George (1987) and Fox et al., (1988) cannot

deal with underdetermined linear systems. But the expressiveness of CLP(R)

allows such a solution since answers may be symbolic. Consequently, an

answer, expressed as a constraint (equation) is perfectly acceptable in CLP(R).

Furthermore, the distributed CLP(R) solution maintains a balanced work load

and does not required either the size of the matrix (N) or the number of

available processors (k) to be special; e.g., a perfect square.The first example tested on the distributed CLP(R) system is one from

numerical analysis: a one dimensional, second-order, parabolic PDE:[n14]@2u=@x2 = c@u=@t:The solution of this equation describes the unsteady-state heat flow in a bar

with @ thermal conductivity, c heat capacity and density. Typical auxiliary

conditions are these:Initial Condition: u(x,0) = 2.0Boundary Conditions: u(0,t) = 0.0u(L,t) = 5.0.The CrankNicolson method of solving this equation uses finite difference

approximations for the two partial derivatives, and arrives at a set of N linear

algebraic equations in N unknowns. In the case that [n14] [n01]t = c ([n01]x)2, one

obtains these equations (Gerald 1970), p. 225):4*x1 x2 = 2.0xx 1 + 4* xi+1 = 4.0, for i = 2..N 1,xN 1 + 4 * xN = 12.0In matrix form this system is represented by AX = B, where the constant

matrix A is tri-diagonal, with 4s along the diagonal and 1s on both the first

sub-diagonal and the first super-diagonal.Consider subdividing the system into two parts. The result will put these

two adjacent equations into separate parts:(0 or 1) * xk 1 + 4 * x k + (0 or 1) * xk+1 = ck(0 or 1) * x k + 4* xk+1 + (0 or 1) * xk+2 = ck+1.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 443Since the first subdivision contains k equations and (k +1) variables, it is

under-determined. With a rank of k, this subset of k equations can be solved

by expressing all variables xj (for j = 1..k) in terms of xk+1. Similarly, in

the second subdivision the number of variables is one more than the number

of equations, and the variables xm (for m=(k+1)..N) can all be expressed in

terms on xk.The complete set of equations is then solved by solving simultaneously

two equations: (1) the equation from the first subdivision which expressions

xk in terms of xk+1 and (2) the equation from the second subdivision which

expresses xk+1 in terms of xk. Thus, we have replaced solving N equations

in N unknowns with solving three sets of equations. One set of equations

has k equations and k+1 unknowns; the second set of equations has (N k)

equations and (N k)+1 unknowns; the third set of equations has two equations

and two unknowns. Furthermore, the first two sets of equations can be solved

simultaneously. The set of two equations can only be solved after solving

both of the larger sets.This idea can be repeated with k subdivisions rather than two. It is assumed

without loss of generality that these subdivisions are equally sized. The final

result will require solving (k+1) sets of equations:(1) two sets (the first and last subdivisions) will each have N/k equations and

(N/k + 1) unknowns;(2) (k 2) sets (all the middle subdivisions) will each have N/k) equations

and (N/k + 2) unknowns;(3) one set of (2k 2) equations and (2k 2) unknowns.All k sets of type (1) or type (2) can be solved concurrently. The one set of

type (3) must wait for the completion of the other sets. Note that type (1)

sets each contribute one equation and one unknown to the type (3) set, while

type (2) sets each contribute two equations and two unknowns to the type (3)

set.With k processors available, the time required to solve the set of N equations

is the sum of (1) the longest time to solve a set of N/k (tri-diagonal) equations;(2) the time to solve a set of 2k 2 equations; and (3) the communication time.

It is possible that one processor acts as the group manager by coordinating

the cooperative effort and by acting as one of the k type (1) or type (2)

system solvers. The communication time consists of k 1 unique orders (one

to each of the other equation solvers) and k 1 unique solution responses

(each of which consists of one or two equations). The group manager must

finally direct the appropriately solved unknowns of the type (3) system to

their corresponding processors so that all final substitutions can be made, but

we exclude this last elaboration here.444 WILLIAM B. DAYAn upper bound on the size of k can be found as follows. Since the type (1)

and type (2) systems are tri-diagonal systems and since the type (3) system is

not, we do not want the type (3) system to be larger than the typical type (1)

or type (2) system. Thus, we require that2k 2 < N/k.For N=100, we consider only subdivisions with two to seven pieces. For N =

1000, we need 2k < 22. In general, k < O(pN). This theory can be easily

generalized to a banded matrix with a fixed bandwidth.5.3. The Systems ArchitectureThe architecture for the entire system is shown in Figure 1. In that figure,

it is assumed that the problem under consideration is being solved on three

machines, one of which (Queen) is used as the interface to the user. The other

two machines, Drone1 and Drone2, service the Queen.The user initiates the system with the command pvmd HOSTFILE &, as

discussed earlier. This causes each of the pvm daemons to begin and creates

the network communications link (the dotted arc) from the Queen to each

drone. During execution, this can be the most time-consuming part of this

implementation. The user sees the links being made each time the name of a

drone appears on the Queens screen.The HOSTFILE contains a network name for each machine that is to

be used in the network, beginning with the name of the Queens machine.

For example,elm.csepeach.csespruce.cse#hickory, csewould be all that PVM needs in a HOSTFILE to use machine elm as a network

Queen (listed first) and machines peach and spruce as drones. A third machine,

hickory, is currently not to be included in the network, as indicated by the #

flag.After the pvm daemons are running, the user starts the execution of the

client problem on the Queens machine. A complete listing of the program

client.c and all other code discussed in this paper can be found in Day (1993).

The clients C program must register itself with the PVM. Then each server

program on the drones must be started and registered with the PVM. The

initiation of the servers is done from client.c. This registration procedure[n14]EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 445Figure 1. Network architecture.allows servers to respond to more than one client and allows more than one

process to be active on the same machine.As shown in Figure 1, messages are sent from the Queen to Drone1 and

Drone2 through the pvm daemons. Each message directs a drone to solve a

specific piece of a problem. For example, if a set of 50 equations is split in

half, the Queen may order Drone1 to solve Equations 1 through 25 and order

Drone2 to solve Equations 26 through 50. Each drone accesses a C shell,

driver.sh, to organize the necessary work. The results of the collection of446 WILLIAM B. DAYprograms that comprise the driver-package are written to a file, answers.The

contents of this file are then returned to the Queen via PVM. When the Queen

has received all the drones responses, the assembled solutions are converted

into a single-predicate CLP(R) program. The program is run on the Queens

CLP(R) interpreter.5.4. A Drones ArchitectureFigure 2 is representation of the collection of programs that are accessed

by the server.c program and are labeled the driver-package. From server.c,

a typical command system(driver.sh 21 40 60) is issued. This command

requests execution of the C shell driver.sh with the three input arguments

Nlow = 21, Nhi = 40 and Nall = 60. These three arguments instruct this

process to solve a piece (equations 21 through 40) of the system with 60

equations.First, driver.sq creates another C shell to invoke and record the responses

of the CLP(R) interpreter. This is necessary because it is not possible to pass

variables (the values of Nlow, Nhi, and Nall) dynamically to the interpreter.

However, only seven lines of this shell differ from one case to another.

Therefore, driver.sq calls a C program, instantiate.c, to read a template shell,

template.sq, modify the appropriate seven lines, and write the instantiated

copy to the new shell, instance.x.y.z. The suffix (.x.y.z) refers to the values of

Nlow, Nhi, and Nall, respectively.Second, driver.sh invokes the new shall, instance.x.y.z. This second shell

invokes the CLP(R) interpreter and captures in file session.x.y.z all output

that normally appears on a terminal screen. The CLP(R) program that is

actually invoked is shown in the Appendix. Most of the lines in the new file

session.x.y.z are useless. There is only one (very long line) that holds the

solution. In this line the entire solution is written as a list, one variable per list

component. Some re-formating of the file session.x.y.z is done using common

Unix filters and sed instructions. The result in store in file format.x.y.z.Third, driver.sh invokes the C program, coder.c, to encode the answer in a

form that is appropriate for transmission to the Queen process. This consists

of extracting only two of the answer components and combining them into an

equation, if this server is actually working on the first subdivision (Nlow=1)

or the last subdivision (Nhi=Nall) of the total problem. In all other cases, four

of the answer components are needed to construct two equations to return to

the Queen.The coder.c uses format.x.y.z as its input file and stores its final answers

in answer.x. Intermediate results are held in the file temp.x. After these three

actions driver.sh deletes all temporary files and shells, except for answer.x.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 447Figure 2. The driver package.The code of server.c passes the contents f answer.x through the pvm daemon

to the Queen. The server.c program then terminates.5.5. The Queens ProgramThe Queens software consists solely of the client.c program. This program

determines which remote machine listed in the HOSTFILE will actually be

used in problem-solving by requesting from the user the total number of448 WILLIAM B. DAYequations (Nall) and the number of remote machines (Nprocess). The first

Nprocess machines listed in the HOSTFILE are used. On each of these

machines a copy of server.c is ordered to begin execution with the PVM

command initiateM(server,host[i]). A message is then sent to each server.c, giving the exact portion of the problem to solve. Although these request

messages are sent in order (of the HOSTFILE listing), the results that are

returned to the Queen are accepted on a first-arrival basis. All the replies are

written to one file whose first line is of this form:solve(T20,T21,T40,T41):and whose typical reply message is this form:T20 = 1.4332 + 4.2273 * T21.Thus, a CLP(R) program (consisting of a single predicate) is constructed

in the file composite. The last invocation of the CLP(R) interpreter is done

using composite as input. Finally, one returns the pvm daemon running on the

local machine (elm) to the foreground and stops the local and all remote pvm

daemons with a Control-C.5.6. Results and GeneralizationsThis section discusses three examples to illustrate the systems power. First

consider the case of solving twenty-four equations by the CrankNicolson

algorithm. On a single machine, this is an easy problem for CLP(R)(the code

appears in the Appendix) with the following solution vector:[0.98203, 1.71281, 1.92305, 1.97938, 1.99448, 1.99852,1.9996, 1.99989 1.99997, 1.99999, 2, 2, 2, 2, 2.0001,2.00004, 2.00016, 2.0006, 2.00222, 2.00829, 2.03093,2.11543, 2.43078, 3.6077].In contrast, by solving the same problem on three machines, broken into

three divisions of eight equations each, the following set of constraints is

returned by the three subdivisions:T8 = 59360 + 40545 * T16 10864 * T17,T9 = 15904 + 10864 * T16 2911 * T17,T8 = 1.464 + 0.267949 * T9,T17 = 1.46425 + 0.267949 * T16.EXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 449By prefixing these constraints with solve([T8, T9, T16, T17]): one has a

CLP(R) program. After running the interpreter on this program, this solution

vector is found:[1.99989, 1.99997, 2.00004, 2.00016].These four components are exactly the same as the corresponding four components of the 24-vector given above. Notice also that the order of the equations

implies that the responses to the Queens requests (Drone1, Drone2, Drone3)

was (Drone2, Drone1, Drone3).The second example is one for 200 equations solved using 20 machines,

each with ten equations. It is not possible to compare this solution with the onemachine solution since the interpreter itself cannot handle so many equations.

This is caused by a size limit in the constraint solver of the interpreter.

However, the exact solution of this problem for more than a dozen equations

is known: all components of the N-vector should be 2, with the exception of

the first dozen and the last dozen. This can be confirmed empirically for sizes

between 20 and 50. The composite program that is assembled by the Queen

contains thirty-eight equations and eighteen unknown [T10, T11, T20, T21,

::: T180, T181, T190, T191]. Each components value is 2 except that T10 =1.99999 and T191 = 2.00001.The last example has not actually been executed, but one can imagine that

it would require only minor code modifications to be successful. Suppose we

have 20 machines available for solving 10,000 equations with a tri-diagonal

matrix (like the CrankNicolson algorithm). The new twist is that one is not

allowed to request more than ten equations at a time to be solved on any one

processor. The method of subdividing can still solve this problem as follows:1. Order the first 200 equations to be solve on the twenty machines, ten

equations per machine. This will be repeated 49 times for a total of 50

waves.2. The first and last waves will each return thirty-eight equations and the

other forty-eight waves will return forty equations. The total is 1998

equations.3. Solve the 1998 equations in ten waves of 200 equations per wave. (The

last wave has only 198 equations).4. This second round of waves (step 3) will return 398 equations.5. A third round of waves will have 200 equations on the first wave and 198

equations on the second wave.6. Finally the third round of waves (step 5) will return 79 equations, which

can be solved using only eight of the machines.450 WILLIAM B. DAYEach machine has been required to solve at most sixty-three sets of ten

equations. Although this illustrates the variable granularity that is possible

with this system, the primary focus of this hypothetical example is that the

system can account for the limited computing capacity of machines.6. ConclusionsCLP languages have motivated new ways to think about old problems. Annexing a constraint solver to a logic programming inference engine has provided

a declarative way to search efficiently. The constraint solver of CLP(R)isable

to selectively limit large search spaces through its internal, canonical representation of constraints. Restricted searching in the Magic Square of Primes is

an simple illustration of this hiding of (detailed) information. Other examples

of the constrain-and-generate paradigm can be found in van Hentenryck

(1989) and in Jaffar and Maher (1994).The increased transparency in problem representations that is offered by

CLP languages has also been promoted in this paper. This has lead applications researchers both to unify disparate regimes (as in the three examples of

Section 3) and to reconsider fundamental algorithms of mathematics (as in

Section 5). For expert systems, CLP languages provide a lucid combination of

symbolic and numeric calculations and erect a platform on which to perform

hypothetical scenarios. In engineering problems, CLP languages unite the

areas of analysis and synthesis for a programmer. This promotes abstraction

of a problem and moves further toward declarative programming. Finally,

through the ability of CLP languages to use constraints in both queries and

answers, it is possible to re-examine classic solution techniques (as shown in

Section 5) by the set of linear equations, which are partitioned for distributed

processing.

/* This CLP(R) program solves a particular first-order parabolic PDEs using the */

/* CrankNicolson method */go(S):

write( Enter the lower index value, the upper index value, and),

write( the total number of equations.). nl,write( Nlow = ), read(Nlow), nl,write( Nhi = ), read(Nhi), nl,write( Ntotal =), read(Nall), nl,choose direction(S,Nlow,Nhi,Nall).AppendixEXPRESSIVE APPLICATIONS OF CONSTRAINT LOGIC PROGRAMMING 451choose direction(S,Nlow,Nhi,Nall) :Nhi < Nall, Index = Nlow,

solve forward(S,Nlow,Index,Nhi,Nall).

choose direction(S,Nlow,Nhi,Nall):Nlow=1,Nhi=Nall,Index=Nlow,

solve forward(S,Nlow,Index,Nhi,Nall).

choose direction(S,Nlow,Nhi,Nall):Nlow>1, Nhi=Nall,Index=Nhi,

solve backward(S,Nlow,Index,Nhi,Nall).solve forward([X,Y,ZkT],Nlow,I,Nhi,Nall):I=1, 4.0*X Y = 2.0,

solve forward([X,Y,ZjT], Nlow,2,NHi,Nall).

solve forward([X,Y,ZjT],Nlow,I,Nhi,Nall) :I > 1, I < Nall, I>= Nlow, I < Nhi,X + 4.0*Y Z = 4.0,

solve forward([Y,ZjT],Nlow,I+1,Nhi,Nall).

solve forward([X,Y,Zj[]],Nlow,I,Nhi,Nall):I > 1, I < Nall, I >= Nlow, I = Nhi,X + 4.0*Y Z = 4.0.

solve forward([X,Yj[]],Nlow,I,Nhi,Nall):I=Nall, X + 4.0*Y = 12.0.solve backward([X,Y,ZjT],Nlow,I,Nhi,Nall):I=Nhi,4.0*X Y = 12.0,

solve backward([X,Y,ZjT],Nlow,I 1,Nhi,Nall).

solve backward([X,Y,ZjT],Nlow,I,Nhi,Nall):I > Nlow, I < Nhi, X + 4.0*Y Z = 4.0,

solve backward([Y,ZjT],Nlow,I 1,Nhi,Nall).

solve backward([X,Y,Zj[]],Nlow,I,Nhi,Nall):I=Nlow, X + 4.0*Y Z = 4.0

Chu, E. & George, A. (1987) Gaussian Elimination with Partial and Load Balancing on aMultiprocessor. Parallel Computing 5: 65.Cohen, J. (1990). Constraint Logic Programming Languages. Communications of the ACM,5268.Day, W. (1993). Planning Model Implementation in CLP. Final Technical Report, ContractNo. C/UB-1754-B, Rome Laboratory, Griffiss AFB, NY.

Fox, G. C., Johnson, M. A., Lyzenga, G. A., Otto, S. W., Salmon, J. K. & Walker, D. W. (1988).Solving Problems on Concurrent Processors. Vol. I, Prentice Hall, Englewood Cliffs, NJ.

Gerald, Curtis F. (1970). Applied Numerical Analysis. Addison-Wesley.Heintze, N., Michaylov, S. & Stuckey, P (1986). CLP(R) and Some Electrical EngineeringProblem. Proceedings of The 4th Int. Conf. on Logic Programming, J-L.Lassez(ed.).MIT

Press.Jaffar, J. & Michaylov, S. (1986). Methodology and Implementation of a CLP System. Proceedings of The 4th Int. Conf. on Logic Programming, J-L. Lassez (ed.), MIT Press.References452 WILLIAM B. DAYJaffar, J., Michaylov, S., Stuckey, P. & Yap, R. H. C. (1992). The CLP(R) Language andSystem. ACM Trans. Programming Languages 14: 339395.

Jaffar, J. & Maher, M. J. (1994). Constraint Logic Programming: A Survey. J. Logic Programming 19, 20: 503581.Lakmazaheri, S. & Rasdorf, W. J. (1989). Constraint Logic Programming for the Analysis

and Partial Syntheses of Truss Structures. Artifical Intelligence for Engineering Design,

Analysis, & Manufacturing 8: 157, 173.Lassez, C., MacAloon, K. & Yap, R. (1987) Constraint Logic Programming and OptionTrading. IEEE Expert 2: 4250.Sunderam, V. (1990), PVM: A Framework for Parallel Distributed Computing. Concurrency:Practice & Experience 2: 315339.

van Hentenryck, P. (1989). Constraint Satisfaction in Logic Programming. MIT Press.