<* PRAGMA LL *> MODULEA; IMPORT Egraph, Equiv, RedundantSolve AS NonLinearSolve, RTVal; IMPORT RefSeq, Word; EXCEPTION Unsolvable; TYPE Type = { Any, Pair, Num, Text, Null }; (* There is a flat partial order on types, with "Type.Any" as bottom. *) VAR (* CONST *) PlusOp, TimesOp, SinOp, CosOp, AtanOp, ExpOp, PairOp := Init(NEW(Var)); zero := RTVal.FromReal(0.0); emptyString := RTVal.FromText(""); TYPE EC = Private; (* The type "EC" represents an equivalence class with an extra field to record if the class contains a pair expression. *) REVEAL Private = Egraph.T BRANDED "JunoSolve.Private" OBJECT cons: EC; OVERRIDES union := Merge; END; Var = Public BRANDED "JunoSolve.Var" OBJECT type: Type; index: INTEGER; marked: BOOLEAN; uses: UseList; next, availLink: Var; END; TYPE UseList = REF RECORD c: Constraint; mask := 0; next, availLink: UseList END; JunoSolve
cons node
is a node whose car
is PairOp
. It represents a pair
expression (x, y)
.
The call v.union(w)
where v
is of type EC
requires that v.root #
w.root
.
The following are invariants on valid EC
's and Var
's:
P1: If an equivalence class whose root is rt
contains a cons node, then
rt.cons
points to a cons node in rt
's class. Conversely, if rt
's
class does not contain a cons node, then rt.cons = NIL
.
P2: The type
field of a Var
records its type; it is only valid for a
variable that is the root of its equivalence class.
P3: For any variable rt
that is the root of its equivalence class,
rt.known => rt.type # Type.Any
. (The converse is not true, since a
REAL or TEXT constraint can set the type of an unknown.)
P4: The index
field of a Var
v
reflects the index in the global array
numeric_vars
of the numeric unknown represented by v
. It is -1
if
v
does not correspond to any variable in the array.
P5: The uses
and next
fields of a Var
are relevant only in
NumericSolve()
, and they are valid only for Var
's that are the
roots of their equivalence classes.
TYPE Args = ARRAY [0..2] OF Var; ConType = { Equal, Cons, Plus, Times, Atan, Sin, Cos, Exp, Real, Text }; REVEAL Constraint = BRANDED "JunoSolve.Constraint" OBJECT type: ConType; arg: Args; hintCnt: CARDINAL; mask: CARDINAL; nls_con: NonLinearSolve.Constraint; unknown: Var; nls_con_mask: CARDINAL; next: Constraint; availLink: Constraint END; (* These fields are only valid during the procedure "NumericSolve". *) VAR mu := NEW(MUTEX); <* LL >= { mu } *> unifyList := NEW(RefSeq.T).init(sizeHint := 100); numeric_vars := NEW(Vars, 60); numeric_vals := NEW(REF ARRAY OF RTVal.Real, 60); numeric_con := NEW(REF ARRAY OF NonLinearSolve.Constraint, 30);The
unify
list (u1, v1, u2, v2, ...)
contains all pairs (ui, vi)
of
cons nodes that are equivalent, but have not yet been merged.
All Var
's in numeric_vars
are roots of their equivalence classes.
PROCEDURENew (known := FALSE; val: RTVal.T := NIL): Var = BEGIN RETURN NewPair(NIL, NIL, known, val) END New; VAR varAvail, varInUse: Var := NIL; conAvail, conInUse: Constraint := NIL; useAvail, useInUse: UseList := NIL; PROCEDURENewPair (car, cdr: Var; known := FALSE; val: RTVal.T := NIL): Var= VAR res: Var; BEGIN IF varAvail # NIL THEN res := varAvail; varAvail := varAvail.availLink ELSE res := NEW(Var) END; res.availLink := varInUse; varInUse := res; res.known := known; res.val := val; res.car := car; res.cdr := cdr; RETURN Init(res) END NewPair; PROCEDURENewCon (type: ConType): Constraint = VAR res: Constraint; BEGIN IF conAvail # NIL THEN res := conAvail; conAvail := conAvail.availLink ELSE res := NEW(Constraint) END; res.availLink := conInUse; conInUse := res; res.type := type; RETURN res END NewCon; PROCEDURENewEqual (x, y: Var): Constraint = VAR res := NewCon(ConType.Equal); BEGIN res.arg[0] := x; res.arg[1] := y; RETURN res END NewEqual; PROCEDURENewCons (x, y, z: Var): Constraint = VAR res := NewCon(ConType.Cons); BEGIN res.arg[0] := x; res.arg[1] := y; res.arg[2] := z; RETURN res END NewCons; PROCEDURENewPlus (x, y, z: Var): Constraint = VAR res := NewCon(ConType.Plus); BEGIN res.arg[0] := x; res.arg[1] := y; res.arg[2] := z; RETURN res END NewPlus; PROCEDURENewTimes (x, y, z: Var): Constraint = VAR res := NewCon(ConType.Times); BEGIN res.arg[0] := x; res.arg[1] := y; res.arg[2] := z; RETURN res END NewTimes; PROCEDURENewAtan (x, y, z: Var): Constraint = VAR res := NewCon(ConType.Atan); BEGIN res.arg[0] := x; res.arg[1] := y; res.arg[2] := z; RETURN res END NewAtan; PROCEDURENewSin (x, y: Var): Constraint = VAR res := NewCon(ConType.Sin); BEGIN res.arg[0] := x; res.arg[1] := y; RETURN res END NewSin; PROCEDURENewCos (x, y: Var): Constraint = VAR res := NewCon(ConType.Cos); BEGIN res.arg[0] := x; res.arg[1] := y; RETURN res END NewCos; PROCEDURENewExp (x, y: Var): Constraint = VAR res := NewCon(ConType.Exp); BEGIN res.arg[0] := x; res.arg[1] := y; RETURN res END NewExp; PROCEDURENewReal (x: Var): Constraint = VAR res := NewCon(ConType.Real); BEGIN res.arg[0] := x; RETURN res END NewReal; PROCEDURENewText (x: Var): Constraint = VAR res := NewCon(ConType.Text); BEGIN res.arg[0] := x; RETURN res END NewText; PROCEDURENewUse (con: Constraint; next: UseList): UseList = VAR res: UseList; BEGIN IF useAvail # NIL THEN res := useAvail; useAvail := useAvail.availLink; res.c := con; res.mask := 0; res.next := next ELSE res := NEW(UseList, c := con, next := next) END; res.availLink := useInUse; useInUse := res; RETURN res END NewUse; PROCEDUREDispose () = BEGIN VAR l := varInUse; BEGIN IF l # NIL THEN WHILE l.availLink # NIL DO l := l.availLink END; l.availLink := varAvail; varAvail := varInUse; varInUse := NIL END END; VAR l := conInUse; BEGIN IF l # NIL THEN WHILE l.availLink # NIL DO l := l.availLink END; l.availLink := conAvail; conAvail := conInUse; conInUse := NIL END END; VAR l := useInUse; BEGIN IF l # NIL THEN WHILE l.availLink # NIL DO l := l.availLink END; l.availLink := useAvail; useAvail := useInUse; useInUse := NIL END END END Dispose; CONST NoIndex = -1; PROCEDUREInit (v: Var): Var =
Modifyv
into a variablev'
such thatv'.root = v'
and such that P1 - P5 hold ofv'
. Returnv'
.
BEGIN EVAL Egraph.T.init(v); v.index := NoIndex; v.marked := FALSE; v.uses := NIL; v.next := NIL; IF v.known THEN TYPECASE v.val OF <* NOWARN *> | RTVal.Null => v.type := Type.Null | RTVal.Number => v.type := Type.Num | RTVal.Text => v.type := Type.Text | RTVal.Pair => v.type := Type.Pair END ELSE v.type := Type.Any END; IF v.car = PairOp THEN v.cons := v ELSE v.cons := NIL END; RETURN v END Init; PROCEDUREP (READONLY c: ARRAY OF Constraint): BOOLEAN = BEGIN LOCK mu DO IF NUMBER(numeric_vars^) < 3 * NUMBER(c) THEN numeric_vars := NEW(Vars, 3 * NUMBER(c)); END; TRY EVAL unifyList.init(); ConstructECs(c); UnifyClose(); NumericSolve(c); ConstructSoln(c); EXCEPT Equiv.Forbidden, Unsolvable => RETURN FALSE END; RETURN TRUE END END P; PROCEDUREEtpLogConstraint (<*UNUSED*>type: INTEGER) =
Types are: 0 = Equal; 1 = CONS; 2 = +; 3 = *; 4 = REAL; 5 = TEXT; 6 = ATAN; 7 = SIN; 8 = COS; 9 = EXP.
BEGIN END EtpLogConstraint; PROCEDUREConstructECs (READONLY c: ARRAY OF Constraint) RAISES {Equiv.Forbidden} = VAR con: Constraint; BEGIN FOR i := FIRST(c) TO LAST(c) DO con := c[i]; CASE con.type OF <* NOWARN *> | ConType.Equal => EtpLogConstraint(0); WITH r1 = con.arg[0].root, r2 = con.arg[1].root DO IF r1 # r2 THEN EVAL r1.union(r2) END END | ConType.Cons => EtpLogConstraint(1); VAR pair := List3(PairOp, con.arg[1], con.arg[2]); BEGIN pair.type := Type.Pair; EVAL (con.arg[0].root).union(pair) END | ConType.Plus => EtpLogConstraint(2); VAR plus := List3(PlusOp, con.arg[1], con.arg[2]); BEGIN SetFuncType(plus, Type.Num); EVAL (con.arg[0].root).union(plus) END | ConType.Times => EtpLogConstraint(3); VAR times := List3(TimesOp, con.arg[1], con.arg[2]); BEGIN SetFuncType(times, Type.Num); EVAL (con.arg[0].root).union(times) END | ConType.Atan => EtpLogConstraint(6); VAR atan := List3(AtanOp, con.arg[1], con.arg[2]); BEGIN SetFuncType(atan, Type.Num); EVAL (con.arg[0].root).union(atan) END | ConType.Sin => EtpLogConstraint(7); VAR sin := List2(SinOp, con.arg[1]); BEGIN SetFuncType(sin, Type.Num); EVAL (con.arg[0].root).union(sin) END | ConType.Cos => EtpLogConstraint(8); VAR cos := List2(CosOp, con.arg[1]); BEGIN SetFuncType(cos, Type.Num); EVAL (con.arg[0].root).union(cos) END | ConType.Exp => EtpLogConstraint(9); VAR exp := List2(ExpOp, con.arg[1]); BEGIN SetFuncType(exp, Type.Num); EVAL (con.arg[0].root).union(exp) END | ConType.Real => EtpLogConstraint(4); SetType(con.arg[0], Type.Num) | ConType.Text => EtpLogConstraint(5); SetType(con.arg[0], Type.Text) END END; END ConstructECs; PROCEDURESetFuncType (v: Var; t: Type) RAISES {Equiv.Forbidden} =
Sets the type ofv
tot
. Verifies that the type of each argument inv.cdr
is at mostt
in the flat partial order withType.Any
as bottom, and promotes the type of each argument with typeType.Any
tot
.Requires
v
to be a root. RaisesEquiv.Forbidden
if some argument inv.cdr
has a type that is not at mostt
.
BEGIN <* ASSERT v.root = v *> v.type := t; VAR curr: Var := v.cdr; BEGIN WHILE curr # NIL DO SetType(curr.car, t); curr := curr.cdr END END END SetFuncType; PROCEDURESetType (v: Var; t: Type) RAISES {Equiv.Forbidden} =
Sets the type of the class containingv
tot
, or raisesEquiv.Forbidden
if the current type ofv
's class is different fromt
and greater than bottom.
VAR rt: Var := v.root; BEGIN IF rt.type = Type.Any THEN rt.type := t ELSIF rt.type # t THEN RAISE Equiv.Forbidden END; END SetType; PROCEDUREList3 (x, y, z: Egraph.T): Var = BEGIN RETURN NewPair(car := x, cdr := NewPair(car := y, cdr := NewPair(car := z, cdr := NIL))) END List3; PROCEDUREList2 (x, y: Egraph.T): Var = BEGIN RETURN NewPair(car := x, cdr := NewPair(car := y, cdr := NIL)) END List2; PROCEDUREMerge (x: EC; y: Equiv.T): Equiv.T RAISES {Equiv.Forbidden} =
By the precondition onEquiv.T.union()
,x
andy
must be the roots of their equivalence classes. They will also be of typeJunoSolve.Var
. Letrx
andry
be the results of narrowingx
andy
, respectively, to the typeJunoSolve.Var
. Letroot
be the root of their resulting combined equivalence class. ThenMerge
either raises the exception or establishes the following post-conditions:Q1:
root.type
= the meet ofrx.type
andry.type
in the flat partial order. If the meet is undefined, raiseEquiv.Forbidden
.Q2:
root.known = rx.known OR ry.known
.root.val = rx.val
ifrx.known
, androot.val = ry.val
ifry.known
. If both are known, but their values are different, raiseEquiv.Forbidden
. IfNOT root.known
, thenroot.val
isrx.val
orry.val
for whichever of those two values is non-NIL, or either if both are non-NIL.Q3:
root.cons
is selected fromrx.cons
andry.cons
so that it isNIL
only if both of them are. If both are non-NIL, then the valuesrx.cons
andry.cons
are added to the global unify list.Q4: If
root.cons # NIL
androot.val
is a pair, then the car and cdr ofroot.val
have been propagated to the arguments of the cons. Ifroot.known
androot.cons # NIL
, raiseEquiv.Forbidden
if the known value is not a pair.
<* LL >= { mu } *> VAR rx: Var := x; ry: Var := y; propagateKnowns := rx.known # ry.known AND (rx.cons = NIL) = rx.known AND (ry.cons = NIL) = ry.known; propagateHints := (rx.val # NIL) # (ry.val # NIL) AND (rx.cons = NIL) = (rx.val # NIL) AND (ry.cons = NIL) = (ry.val # NIL); (* Work is necessary to establish Q4 iff "propagateKnowns OR propagateHints". *) BEGIN <* ASSERT rx # ry *> VAR root: Var := Egraph.T.union(rx, ry); nonroot: Var; BEGIN IF root = rx THEN nonroot := ry ELSE nonroot := rx END; (* Establish Q1 *) IF root.type = Type.Any THEN root.type := nonroot.type ELSIF nonroot.type # Type.Any AND root.type # nonroot.type THEN RAISE Equiv.Forbidden END; (* Establish Q2 *) IF nonroot.known THEN IF NOT root.known THEN root.known := TRUE; root.val := nonroot.val ELSIF NOT RTVal.Equal(rx.val, ry.val) THEN RAISE Equiv.Forbidden END; ELSIF NOT root.known THEN IF nonroot.val # NIL THEN root.val := nonroot.val END END; (* Establish Q3 *) IF nonroot.cons # NIL THEN IF root.cons = NIL THEN root.cons := nonroot.cons; ELSE unifyList.addhi(rx.cons); unifyList.addhi(ry.cons) END END; (* Establish Q4 *) IF propagateHints OR propagateKnowns THEN TYPECASE root.val OF | RTVal.Pair (p) => Fix(root.cons.cdr.car, root.known, p.car); Fix(root.cons.cdr.cdr.car, root.known, p.cdr) ELSE IF root.known THEN RAISE Equiv.Forbidden END END; END; RETURN root END; END Merge; PROCEDUREFix (v: Var; known: BOOLEAN; val: RTVal.T) RAISES {Equiv.Forbidden} = VAR t := New(known, val); BEGIN EVAL (v.root).union(t) END Fix; PROCEDUREUnifyClose () RAISES {Equiv.Forbidden} = BEGIN WHILE unifyList.size() # 0 DO VAR y: EC := unifyList.remhi(); x: EC := unifyList.remhi(); BEGIN <* ASSERT x.car = PairOp AND y.car = PairOp *> x := x.cdr; y := y.cdr; WHILE x # NIL DO VAR xa := x.car.root; ya := y.car.root; BEGIN IF xa # ya THEN EVAL xa.union(ya) END END; x := x.cdr; y := y.cdr END; <* ASSERT y = NIL *> END END END UnifyClose; PROCEDURENumericSolve (READONLY c: ARRAY OF Constraint) RAISES { Unsolvable } =
By a {\em variable}, we mean aVar
that is the root of its equivalence class and appears in a numeric constraint inc
. We say a variablev
is {\em labeled} ifv.index # NoIndex
andnumeric_vars[v.index] = v
. We say a variablev
{\em occurs} in a constraintcon
if it is unlabeled andcon
contains aVar
whose root isv
.Q1: For each constraint
con
,con.hintCnt
is the number of distinct unlabeled variables occurring incon
, andcon.mask
is2_b2b1b0
, wherebi
== ``con.arg[i]
is unlabeled''.Q2: There are two linked lists of constraints whose heads are
ghostReady
andtrueReady
. Constraints on both lists havehintCnt <= 1
, and have validnls_con
andunknown
fields set by Functional(). Functional constraints withhintCnt
equal to 1 are on theghostReady
queue, while non-functional constraints withhintCnt
equal to 1 are on thetrueReady
queue. All constraints withhintCnt
equal to 0 appear on one or the other of the two queues.Q3: For any variable
v
,v.uses
is a list that contains an entryul
for each constraint in whichv
occurs. The valueul.c
is the constraint, andul.mask
is2_b2b1b0
wherebi
==(ul.c.arg[i] = v)
.Q4: The hint for an unlabeled variable is either
NIL
or a numeric value.Q5: The variable
var_cnt
is the number of non-constant variables. The variablenumeric_cnt
is the number of numeric constraints.Q6: Every constant has been labeled by a number between
next_const + 1
andLAST(numeric_vars^)
.Q7: The variable
hvl
is the list of hinted unlabeled variables (possibly containing some variables which have become labeled), anduvl
is the list of unhinted, unlabeled variables (also possible containing some variables which have become labeled).Q8: numeric_cnt <= NUMBER(numeric_con).
Q9: All knowns (constants) are labeled, and all unknowns are unlabeled.
R1: The variable
lo
is the number of labeled ``true'' variables, each of which has been labeled with a number between0
andlo - 1
.R2: The expression
var_cnt - hi
is the number of labeled ghost variables, each of which has been labeled by a number betweenhi
andvar_cnt - 1
.R3: The variable
c_lo
is the number of ghost constraints; these constraints are stored in the firstc_lo
entries ofnumeric_con
. The expressionnumeric_cnt - c_hi
is the number of true constraints; these constraints are stored in the last entries ofnumeric_con
.
DURING VARIABLE PROCESSING: numeric_vars[] __ ________ |..| | | |..| = unused | True | |__| | Vars | |________| lo -> |........| |........| numeric_con[] hi -> |________| _____________ | | | | | Ghost | | Ghost | | Vars | | Constraints | |________| |_____________| var_cnt -> |........| |.............| <- c_lo |........| |.............| next_const -> |________| |_____________| <- c_hi | | | | | Consts | | True | |________| | Constraints | |_____________| |.............| <- numeric_cnt |.............| |_____________| AFTER VARIABLE PROCESSING: numeric_vars[] numeric_vals[] ________ ________ | | | | | True | | True | | Vars | | Vals | | | | | numeric_con[] hi -> |________| |________| _____________ lo -> | | | | | | | | | | | Ghost | | | | | | Constraints | | Ghost | | Ghost | | | | Vars | | Vals | | | |________| |________| |_____________| <- c_hi var_cnt -> |........| |........| | | <- c_lo |........| |........| | | next_const -> |________| |________| | True | | | | | | Constraints | | Consts | | Consts | |_____________| |________| |________| |.............| <- numeric_cnt |.............| |_____________|
VAR lo, hi, c_lo, c_hi: CARDINAL; var_cnt: CARDINAL := 0; (* total # of variables *) numeric_cnt: CARDINAL := 0; (* total # of numeric constraints *) next_const := LAST(numeric_vars^); (* index of next known value *) ghostReady, trueReady: Constraint := NIL; hvl, uvl: Var := NIL; PROCEDURE NumericArgCnt(c: Constraint): CARDINAL = (* Return the number of variables in constraint "c" if it is a numeric constraint; 0 otherwise. *) BEGIN CASE c.type OF ConType.Plus, ConType.Times, ConType.Atan => RETURN 3 | ConType.Sin, ConType.Cos, ConType.Exp => RETURN 2 ELSE RETURN 0 END END NumericArgCnt; PROCEDURE AddConToQueue(con: Constraint) = BEGIN IF Functional(con) THEN con.next := ghostReady; ghostReady := con ELSE con.next := trueReady; trueReady := con END END AddConToQueue; PROCEDURE UpdateUses(l: UseList) = BEGIN WHILE l # NIL DO DEC(l.c.mask, l.mask); DEC(l.c.hintCnt); IF l.c.hintCnt = 1 THEN AddConToQueue(l.c) END; l := l.next END END UpdateUses; (* NumericSolve *) BEGIN FOR i := FIRST(c) TO LAST(c) DO VAR con := c[i]; cnt := NumericArgCnt(con); BEGIN IF cnt > 0 THEN INC(numeric_cnt); con.hintCnt := 0; con.mask := 0; FOR j := 0 TO cnt - 1 DO VAR v: Var := con.arg[j].root; BEGIN IF NOT v.known THEN (* "v" is an unknown *) IF v.uses = NIL THEN (* first time seeing "v" *) INC(var_cnt); (* disregard any non-numeric hint *) IF v.val # NIL AND NOT ISTYPE(v.val, RTVal.Number) THEN v.val := NIL END; IF v.val = NIL THEN v.next := uvl; uvl := v (* unhinted *) ELSE v.next := hvl; hvl := v (* hinted *) END END; (* add "con" to "v.uses" if necessary *) IF v.uses = NIL OR v.uses.c # con THEN v.uses := NewUse(con, next := v.uses); INC(con.hintCnt) END; (* update the variable and constraint masks *) VAR bit := Word.LeftShift(1, j); BEGIN v.uses.mask := Word.Or(v.uses.mask, bit); con.mask := Word.Or(con.mask, bit) END ELSIF v.index = NoIndex THEN numeric_vars[next_const] := v; v.index := next_const; DEC(next_const) END END END; (* FOR *) (* Add to "ghostReady" or "trueReady" queue if necessary *) IF con.hintCnt <= 1 THEN AddConToQueue(con) END END END END; (* Q1 - Q7, Q9 *) IF numeric_cnt = 0 THEN RETURN END; IF numeric_cnt > NUMBER(numeric_con^) THEN numeric_con := NEW(REF ARRAY OF NonLinearSolve.Constraint, MAX(numeric_cnt, 2 * NUMBER(numeric_con^))) END; (* Q1 - Q8, Q9 *) lo := 0; hi := var_cnt; c_lo := 0; c_hi := numeric_cnt; LOOP (* Q1 - Q8, R1 - R3 *) IF ghostReady # NIL THEN VAR con := ghostReady; BEGIN <* ASSERT con.hintCnt <= 1 *> ghostReady := ghostReady.next; IF con.hintCnt = 0 THEN (* add "con" to the "trueReady" queue *) con.next := trueReady; trueReady := con ELSE (* label "con.unknown" as a ghost variable *) DEC(hi); numeric_vars[hi] := con.unknown; con.unknown.index := hi; UpdateUses(con.unknown.uses); (* finish constructing "con.nls_con" *) con.nls_con.arg[0] := hi; (* make "con.nls_con" a ghost constraint *) numeric_con[c_lo] := con.nls_con; INC(c_lo) END END ELSIF hvl # NIL THEN IF hvl.index = NoIndex THEN (* make "hvl" a true variable *) numeric_vars[lo] := hvl; hvl.index := lo; INC(lo); UpdateUses(hvl.uses) END; hvl := hvl.next ELSIF trueReady # NIL THEN VAR con := trueReady; BEGIN <* ASSERT con.hintCnt <= 1 *> trueReady := trueReady.next; IF con.hintCnt = 1 THEN (* make "con.unknown" a true variable *) numeric_vars[lo] := con.unknown; con.unknown.index := lo; INC(lo); UpdateUses(con.unknown.uses) END; (* finish constructing "con.nls_con" from "con.nls_con_mask" *) VAR i := 0; BEGIN WHILE con.nls_con_mask > 0 DO IF Word.And(con.nls_con_mask, 2_1) = 2_1 THEN con.nls_con.arg[i] := con.unknown.index END; INC(i); con.nls_con_mask := Word.RightShift(con.nls_con_mask, 1) END END; (* make "con.nls_con" a true constraint *) DEC(c_hi); numeric_con[c_hi] := con.nls_con END ELSIF uvl # NIL THEN IF uvl.index = NoIndex THEN (* make "uvl" a true variable *) numeric_vars[lo] := uvl; uvl.index := lo; INC(lo); UpdateUses(uvl.uses) END; uvl := uvl.next ELSE EXIT END END; (* Solve *) InitVals(numeric_vals, numeric_vars, lo, next_const); WITH constraints = SUBARRAY(numeric_con^, 0, numeric_cnt) DO IF NOT NonLinearSolve.P(lo, var_cnt, numeric_vals^, constraints) THEN RAISE Unsolvable END END; NonLinearSolve.Dispose() END NumericSolve; PROCEDUREInitVals (VAR numeric_vals: REF ARRAY OF RTVal.Real; numeric_vars: Vars; lo, next_const: INTEGER) = (* Fill in "numeric_vals" from "numeric_vars" *) BEGIN IF NUMBER(numeric_vals^) < NUMBER(numeric_vars^) THEN numeric_vals := NEW(REF ARRAY OF RTVal.Real, NUMBER(numeric_vars^)) END; (* fill in ``true'' variable values *) FOR i := 0 TO lo - 1 DO TYPECASE numeric_vars[i].val OF | RTVal.Null => numeric_vals[i] := 0.0 | RTVal.Number (r) => numeric_vals[i] := r.val ELSE numeric_vals[i] := 0.0 END END; (* fill in constant (known) values *) FOR i := next_const + 1 TO LAST(numeric_vars^) DO TYPECASE numeric_vars[i].val OF <* NOWARN *> | RTVal.Number (r) => numeric_vals[i] := r.val END END END InitVals; CONST NoArg = FIRST(INTEGER); PROCEDUREFunctional (con: Constraint): BOOLEAN =
Requirescon.hintCnt <= 1
, i.e., that there is at most one unknown in the constraint.If
con.mask
# 0 and it indicates thatcon
is functional in the single unknown variable, then setcon.unknown
to that unknown variable, setcon.nls_con
to a constraint that solves it, and return TRUE.Otherwise, set
con.nls_con
to aNonLinearSolve.Constraint
clone ofcon
, setcon.unknown
to NIL, and return FALSE.In either case, the arguments in
con.nls_con
corresponding to the unknown are unset.
CONST LoBit = ARRAY [1..7] OF [0..2]{0, 1, 0, 2, 0, 1, 0}; BEGIN <* ASSERT con.hintCnt <= 1 *> IF con.mask = 0 THEN con.unknown := NIL ELSE con.unknown := con.arg[LoBit[con.mask]].root END; CASE con.type OF <* NOWARN *> ConType.Plus => CASE con.mask OF <* NOWARN *> 2_001 => con.nls_con := SetArgs(NonLinearSolve.NewPlus(), con, 1, 2); con.nls_con_mask := 2_001 | 2_010 => con.nls_con := SetArgs(NonLinearSolve.NewMinus(), con, 0, 2); con.nls_con_mask := 2_001 | 2_100 => con.nls_con := SetArgs(NonLinearSolve.NewMinus(), con, 0, 1); con.nls_con_mask := 2_001 | 2_110 => con.nls_con := SetArgs(NonLinearSolve.NewHalve(), con, 0); con.nls_con_mask := 2_001 | 2_000, 2_101, 2_011, 2_111 => con.nls_con := SetArgs(NonLinearSolve.NewPlus(), con, 1, 2, 0); con.nls_con_mask := con.mask; RETURN FALSE END | ConType.Times => IF con.mask = 2_001 THEN con.nls_con := SetArgs(NonLinearSolve.NewTimes(), con, 1, 2); con.nls_con_mask := 2_001 ELSE con.nls_con := SetArgs(NonLinearSolve.NewTimes(), con, 1, 2, 0); con.nls_con_mask := con.mask; RETURN FALSE END | ConType.Atan => CASE con.mask OF <* NOWARN *> | 2_001 => con.nls_con := SetArgs(NonLinearSolve.NewAtan(), con, 1, 2); con.nls_con_mask := 2_001 | 2_010 => con.nls_con := SetArgs(NonLinearSolve.NewMultTan(), con, 2, 0); con.nls_con_mask := 2_001 | 2_000, 2_100, 2_011, 2_101, 2_110, 2_111 => con.nls_con := SetArgs(NonLinearSolve.NewAtan(), con, 1, 2, 0); con.nls_con_mask := con.mask; RETURN FALSE END | ConType.Sin => IF con.mask = 2_01 THEN con.nls_con := SetArgs(NonLinearSolve.NewSin(), con, 1); con.nls_con_mask := 2_01 ELSE con.nls_con := SetArgs(NonLinearSolve.NewSin(), con, 1, NoArg, 0); con.nls_con_mask := con.mask; RETURN FALSE END | ConType.Cos => IF con.mask = 2_01 THEN con.nls_con := SetArgs(NonLinearSolve.NewCos(), con, 1); con.nls_con_mask := 2_01 ELSE con.nls_con := SetArgs(NonLinearSolve.NewCos(), con, 1, NoArg, 0); con.nls_con_mask := con.mask; RETURN FALSE END | ConType.Exp => IF con.mask = 2_01 THEN con.nls_con := SetArgs(NonLinearSolve.NewExp(), con, 1); con.nls_con_mask := 2_01 ELSE con.nls_con := SetArgs(NonLinearSolve.NewExp(), con, 1, NoArg, 0); con.nls_con_mask := con.mask; RETURN FALSE END END; RETURN TRUE END Functional; PROCEDURESetArgs ( nlc: NonLinearSolve.Constraint; con: Constraint; a1: CARDINAL; a2, a0: INTEGER := NoArg) : NonLinearSolve.Constraint = BEGIN nlc.arg[1] := NARROW(con.arg[a1].root, Var).index; IF a2 # NoArg THEN nlc.arg[2] := NARROW(con.arg[a2].root, Var).index END; IF a0 # NoArg THEN nlc.arg[0] := NARROW(con.arg[a0].root, Var).index END; RETURN nlc END SetArgs; PROCEDUREConstructSoln (READONLY c: ARRAY OF Constraint) RAISES {Unsolvable} =
Set theval
field of each unknown inc
to the value of its root.
BEGIN FOR i := FIRST(c) TO LAST(c) DO VAR argCnt: CARDINAL; BEGIN CASE c[i].type OF <* NOWARN *> ConType.Cons, ConType.Plus, ConType.Times, ConType.Atan => argCnt := 3 | ConType.Sin, ConType.Cos, ConType.Exp, ConType.Equal => argCnt := 2 | ConType.Real, ConType.Text => argCnt := 1 END; FOR j := 0 TO argCnt - 1 DO VAR arg := c[i].arg[j]; BEGIN IF NOT arg.known THEN arg.val := Value(arg.root) END END END END END END ConstructSoln; PROCEDUREValue (v: Var): RTVal.T RAISES {Unsolvable} =
Return the value of v
, which must be the root of its equivalence class.
BEGIN <* ASSERT v = v.root *> IF NOT v.known THEN IF v.cons # NIL THEN IF v.marked THEN RAISE Unsolvable ELSE v.marked := TRUE; v.val := RTVal.FromPair( car := Value(v.cons.cdr.car.root), cdr := Value(v.cons.cdr.cdr.car.root)) END ELSIF v.index # NoIndex THEN v.val := RTVal.FromReal(numeric_vals[v.index]) ELSE (* "v" is an unknown var that is not equivalent to a cons nor involved in any numeric constraints. It may not have a type, or it may have a non-trivial type due to a REAL or TEXT constraint. Assign "v" a valid value for its type if it does not already have one. *) CASE v.type OF <* NOWARN *> | Type.Any => IF v.val = NIL THEN v.val := RTVal.nil END | Type.Num => IF v.val = NIL OR NOT ISTYPE(v.val, RTVal.Number) THEN v.val := zero END | Type.Text => IF v.val = NIL OR NOT ISTYPE(v.val, RTVal.Text) THEN v.val := emptyString END END END; v.known := TRUE; END; RETURN v.val END Value; BEGIN END JunoSolve.