3 // Cassowary Incremental Constraint Solver
4 // Original Smalltalk Implementation by Alan Borning
5 // This C++ Implementation by Greg J. Badros, <gjb@cs.washington.edu>
6 // http://www.cs.washington.edu/homes/gjb
7 // (C) 1998, 1999 Greg J. Badros and Alan Borning
8 // See ../LICENSE for legal details regarding this software
14 #include <cassowary/debug.h>
15 #include <cassowary/ClSimplexSolver.h>
16 #include <cassowary/ClErrors.h>
17 #include <cassowary/ClVariable.h>
18 #include <cassowary/ClPoint.h>
19 #include <cassowary/ClSlackVariable.h>
20 #include <cassowary/ClObjectiveVariable.h>
21 #include <cassowary/ClDummyVariable.h>
22 #include <cassowary/cl_auto_ptr.h>
31 #define CONFIG_H_INCLUDED
34 const char *szCassowaryVersion = VERSION;
36 // Need to delete all expressions
37 // and all slack and dummy variables
38 // See NewExpression -- all allocation is done in there
39 ClSimplexSolver::~ClSimplexSolver()
41 #ifdef CL_SOLVER_STATS
42 cerr << "_slackCounter == " << _slackCounter
43 << "\n_artificialCounter == " << _artificialCounter
44 << "\n_dummyCounter == " << _dummyCounter << endl;
45 cerr << "stayMinusErrorVars " << _stayMinusErrorVars.size() << ", "
46 << "stayPlusErrorVars " << _stayPlusErrorVars.size() << ", "
47 << "errorVars " << _errorVars.size() << ", "
48 << "markerVars " << _markerVars.size() << endl;
50 // Cannot print *this here, since local ClVariable-s may have been
54 // Add the constraint cn to the tableau
56 ClSimplexSolver::AddConstraint(ClConstraint *const pcn)
59 Tracer TRACER(__FUNCTION__);
60 cerr << "(" << *pcn << ")" << endl;
63 if (!pcn->FIsOkayForSimplexSolver()) {
64 throw ExCLTooDifficultSpecial("SimplexSolver cannot handle this constraint object");
67 if (pcn->IsStrictInequality()) {
68 // cannot handle strict inequalities
69 throw ExCLStrictInequalityNotAllowed();
72 if (pcn->ReadOnlyVars().size() > 0) {
73 // cannot handle read-only vars
74 throw ExCLReadOnlyNotAllowed();
77 if (pcn->IsEditConstraint())
79 ClEditConstraint *pcnEdit = dynamic_cast<ClEditConstraint *>(pcn);
80 const ClVariable &v = pcnEdit->variable();
81 if (!v.IsExternal() ||
82 (!FIsBasicVar(v) && !ColumnsHasKey(v)))
84 // we could try to make this case work,
85 // but it'd be unnecessarily inefficient --
86 // and probably easier for the client application
88 throw ExCLEditMisuse("(ExCLEditMisuse) Edit constraint on variable not in tableau.");
90 ClEditInfo *pcei = PEditInfoFromClv(v);
93 // we need to only add a partial _editInfoList entry for this
94 // edit constraint since the variable is already being edited.
95 // otherwise a more complete entry is added later in this function
96 _editInfoList.push_back(new ClEditInfo(v, NULL, clvNil, clvNil, 0));
101 ClVariable clvEplus, clvEminus;
102 Number prevEConstant;
103 ClLinearExpression *pexpr = NewExpression(pcn, /* output to: */
106 bool fAddedOkDirectly = false;
110 // If possible Add expr directly to the appropriate tableau by
111 // choosing a subject for expr (a variable to become basic) from
112 // among the current variables in expr. If this doesn't work use an
113 // artificial variable. After adding expr re-Optimize.
114 fAddedOkDirectly = TryAddingDirectly(*pexpr);
116 catch (ExCLRequiredFailure &error)
119 cerr << "could not Add directly -- caught ExCLRequiredFailure error" << endl;
121 RemoveConstraintInternal(pcn);
125 if (!fAddedOkDirectly)
126 { // could not Add directly
127 ExCLRequiredFailureWithExplanation e;
128 if (!AddWithArtificialVariable(*pexpr, e))
130 #ifdef CL_DEBUG_FAILURES
131 cerr << "Failed solve! Could not Add constraint.\n"
134 RemoveConstraintInternal(pcn);
138 throw ExCLRequiredFailure();
142 _fNeedsSolving = true;
144 if (pcn->IsEditConstraint())
146 ClEditConstraint *pcnEdit = dynamic_cast<ClEditConstraint *>(pcn);
147 ClVariable clv = pcnEdit->variable();
148 _editInfoList.push_back(new ClEditInfo(clv, pcnEdit, clvEplus, clvEminus,
154 Optimize(_objective);
155 SetExternalVariables();
162 // Add weak stays to the x and y parts of each point. These have
163 // increasing weights so that the solver will try to satisfy the x
164 // and y stays on the same point, rather than the x stay on one and
165 // the y stay on another.
167 ClSimplexSolver::AddPointStays(const vector<const ClPoint *> &listOfPoints,
168 const ClStrength &strength)
171 Tracer TRACER(__FUNCTION__);
174 vector<const ClPoint *>::const_iterator it = listOfPoints.begin();
176 static const double multiplier = 2.0;
177 for ( ; it != listOfPoints.end(); ++it )
179 AddPointStay((*it)->X(),(*it)->Y(),strength,weight);
180 weight *= multiplier;
186 ClSimplexSolver::AddPointStay(const ClPoint &clp, const ClStrength &strength, double weight)
188 AddPointStay(clp.X(),clp.Y(),strength,weight);
194 ClSimplexSolver::RemoveEditVarsTo(unsigned int n)
196 queue<ClVariable> qclv;
197 ClVarSet sclvStillEditing; // Set of edit variables that we need to *not* remove
198 #ifdef DEBUG_NESTED_EDITS
199 cerr << __FUNCTION__ << " " << n << endl;
202 for ( ClEditInfoList::const_iterator it = _editInfoList.begin();
203 (it != _editInfoList.end() && _editInfoList.size() != static_cast<unsigned int>(n));
206 const ClEditInfo *pcei = (*it);
208 #ifdef DEBUG_NESTED_EDITS
209 cerr << __FUNCTION__ << "Checking " << pcei->_clv
210 << ", index = " << i << endl;
213 qclv.push(pcei->_clv);
215 sclvStillEditing.insert(pcei->_clv);
217 while (!qclv.empty())
219 ClVariable clv = qclv.front();
220 // only remove the variable if it's not in the set of variable
221 // from a previous nested outer edit
226 // The end edit needs to only get rid of the edits on w,h
227 // not the ones on x,y
228 if (sclvStillEditing.find(clv) == sclvStillEditing.end())
230 #ifdef DEBUG_NESTED_EDITS
231 cerr << __FUNCTION__ << ": Removing " << clv << endl;
235 #ifdef DEBUG_NESTED_EDITS
238 cerr << __FUNCTION__ << ": Not removing " << clv << endl;
243 while (_editInfoList.size() > n) {
244 _editInfoList.pop_back();
251 /* A predicate used for remove_if */
252 class VarInVarSet : public unary_function<ClVariable,bool> {
254 VarInVarSet(const ClVarSet &clvset) :
256 _setEnd(clvset.end())
259 bool operator ()(ClVariable clv) const {
260 return (_set.find(clv) != _setEnd);
264 const ClVarSet &_set;
265 const ClVarSet::iterator _setEnd;
270 // Remove the constraint cn from the tableau
271 // Also remove any error variable associated with cn
273 ClSimplexSolver::RemoveConstraintInternal(const ClConstraint *const pcn)
276 Tracer TRACER(__FUNCTION__);
277 cerr << "(" << *pcn << ")" << endl;
280 // We are about to remove a constraint. There may be some stay
281 // constraints that were unsatisfied previously -- if we just
282 // removed the constraint these could come into play. Instead,
283 // Reset all of the stays so that things should stay where they are
285 _fNeedsSolving = true;
287 ResetStayConstants();
289 // remove any error variables from the objective function
290 ClLinearExpression *pzRow = RowExpression(_objective);
293 cerr << _errorVars << endl << endl;
296 ClConstraintToVarSetMap::iterator
297 it_eVars = _errorVars.find(pcn);
298 bool fFoundErrorVar = (it_eVars != _errorVars.end());
302 ClVarSet &eVars = (*it_eVars).second;
303 ClVarSet::iterator it = eVars.begin();
304 for ( ; it != eVars.end(); ++it )
306 const ClLinearExpression *pexpr = RowExpression(*it);
309 pzRow->AddVariable(*it,-pcn->weight() * pcn->strength().symbolicWeight().AsDouble(),
313 { // the error variable was in the basis
314 pzRow->AddExpression(*pexpr,-pcn->weight() * pcn->strength().symbolicWeight().AsDouble(),
320 ClConstraintToVarMap::iterator
321 it_marker = _markerVars.find(pcn);
322 if (it_marker == _markerVars.end())
323 { // could not find the constraint
324 throw ExCLConstraintNotFound();
326 // try to make the marker variable basic if it isn't already
327 const ClVariable marker = (*it_marker).second;
328 _markerVars.erase(it_marker);
329 _constraintsMarked.erase(marker);
331 cerr << "Looking to remove var " << marker << endl;
333 if (!FIsBasicVar(marker))
334 { // not in the basis, so need to do some work
335 // first choose which variable to move out of the basis
336 // only consider restricted basic variables
337 ClVarSet &col = _columns[marker];
338 ClVarSet::iterator it_col = col.begin();
340 cerr << "Must Pivot -- columns are " << col << endl;
343 ClVariable exitVar = clvNil;
344 bool fExitVarSet = false;
345 double minRatio = 0.0;
346 for ( ; it_col != col.end(); ++it_col)
348 const ClVariable &v = *it_col;
349 if (v.IsRestricted() )
351 const ClLinearExpression *pexpr = RowExpression(v);
352 assert(pexpr != NULL );
353 Number coeff = pexpr->CoefficientFor(marker);
355 cerr << "Marker " << marker << "'s coefficient in " << *pexpr << " is "
358 // only consider negative coefficients
361 Number r = - pexpr->Constant() / coeff;
362 if (!fExitVarSet || r < minRatio)
371 // if we didn't set exitvar above, then either the marker
372 // variable has a positive coefficient in all equations, or it
373 // only occurs in equations for unrestricted variables. If it
374 // does occur in an equation for a restricted variable, pick the
375 // equation that gives the smallest ratio. (The row with the
376 // marker variable will become infeasible, but all the other rows
377 // will still be feasible; and we will be dropping the row with
378 // the marker variable. In effect we are removing the
379 // non-negativity restriction on the marker variable.)
383 cerr << "exitVar did not get set" << endl;
385 it_col = col.begin();
386 for ( ; it_col != col.end(); ++it_col)
388 ClVariable v = *it_col;
389 if (v.IsRestricted() )
391 const ClLinearExpression *pexpr = RowExpression(v);
392 assert(pexpr != NULL);
393 Number coeff = pexpr->CoefficientFor(marker);
394 Number r = pexpr->Constant() / coeff;
395 if (!fExitVarSet || r < minRatio)
406 { // exitVar is still nil
407 // If col is empty, then exitVar doesn't occur in any equations,
408 // so just remove it. Otherwise pick an exit var from among the
409 // unrestricted variables whose equation involves the marker var
412 RemoveColumn(marker);
416 // A. Beurive' Tue Sep 14 18:26:05 CEST 1999
417 // Don't pick the objective, or it will be removed!
418 it_col = col.begin();
419 for ( ; it_col != col.end(); ++it_col)
421 ClVariable v = *it_col;
429 assert(fExitVarSet == true);
435 Pivot(marker,exitVar);
439 if (FIsBasicVar(marker))
441 ClLinearExpression *pexpr = RemoveRow(marker);
443 cerr << "delete@ " << pexpr << endl;
448 // Delete any error variables. If cn is an inequality, it also
449 // contains a slack variable; but we use that as the marker variable
450 // and so it has been deleted when we removed its row.
453 ClVarSet &eVars = (*it_eVars).second;
454 ClVarSet::iterator it = eVars.begin();
455 for ( ; it != eVars.end(); ++it )
457 ClVariable v = (*it);
465 if (pcn->isStayConstraint())
467 // iterate over the stay{Plus,Minus}ErrorVars and remove those
468 // variables v in those vectors that are also in set eVars
471 ClVarSet &eVars = (*it_eVars).second;
473 .erase(remove_if(_stayPlusErrorVars.begin(),_stayPlusErrorVars.end(),
475 _stayPlusErrorVars.end());
477 .erase(remove_if(_stayMinusErrorVars.begin(),_stayMinusErrorVars.end(),
479 _stayMinusErrorVars.end());
482 else if (pcn->IsEditConstraint())
484 const ClEditConstraint *pcnEdit = dynamic_cast<const ClEditConstraint *>(pcn);
485 const ClVariable clv = pcnEdit->variable();
486 ClEditInfo *pcei = PEditInfoFromClv(clv);
488 ClVariable clvEditMinus = pcei->_clvEditMinus;
489 RemoveColumn(clvEditMinus); // clvEditPlus is a marker var and gets removed later
491 _editInfoList.remove(pcei);
496 // This code is not needed since the variables are deleted
497 // when they are removed from the row --
498 // leaving it in results in double deletions
499 // delete the constraint's error variables
500 // ClVarSet &evars_set = (*it_eVars).second;
501 // ClVarSet::const_iterator it_set = evars_set.begin();
502 // for ( ; it_set != evars_set.end(); ++it_set)
506 _errorVars.erase((*it_eVars).first);
511 Optimize(_objective);
512 SetExternalVariables();
519 // Re-initialize this solver from the original constraints, thus
520 // getting rid of any accumulated numerical problems. (Actually,
521 // Alan hasn't observed any such problems yet, but here's the method
524 ClSimplexSolver::Reset()
527 Tracer TRACER(__FUNCTION__);
528 cerr << "()" << endl;
530 // FIXGJB -- can postpone writing this for a while
531 // gotta be careful, though, as it's a likely place for
532 // a memory leak to sneak in
537 // Re-solve the cuurent collection of constraints, given the new
538 // values for the edit variables that have already been
539 // suggested (see SuggestValue() method)
541 ClSimplexSolver::Resolve()
542 { // CODE DUPLICATED ABOVE
544 Tracer TRACER(__FUNCTION__);
547 SetExternalVariables();
548 _infeasibleRows.clear();
549 if (_fResetStayConstantsAutomatically)
550 ResetStayConstants();
554 ClSimplexSolver::SuggestValue(ClVariable v, Number x)
557 Tracer TRACER(__FUNCTION__);
559 ClEditInfo *pcei = PEditInfoFromClv(v);
563 std::stringstream ss;
564 ss << "SuggestValue for variable " << v << ", but var is not an edit variable" << ends;
565 throw ExCLEditMisuse(ss.str().c_str());
567 throw ExCLEditMisuse(v.Name().c_str());
570 ClVariable clvEditPlus = pcei->_clvEditPlus;
571 ClVariable clvEditMinus = pcei->_clvEditMinus;
572 Number delta = x - pcei->_prevEditConstant;
573 pcei->_prevEditConstant = x;
574 DeltaEditConstant(delta,clvEditPlus,clvEditMinus);
578 // Re-solve the cuurent collection of constraints, given the new
579 // values for the edit variables that have already been
580 // suggested (see SuggestValue() method)
581 // This is not guaranteed to work if you remove an edit constraint
582 // from the middle of the edit constraints you added
583 // (e.g., edit A, edit B, edit C, remove B -> this will fail!)
586 ClSimplexSolver::Resolve(const vector<Number> &newEditConstants)
588 ClEditInfoList::iterator it = _editInfoList.begin();
590 for (; i < newEditConstants.size() && it != _editInfoList.end(); ++it, ++i)
592 ClEditInfo *pcei = (*it);
593 SuggestValue(pcei->_clv,newEditConstants[i]);
601 // Add the constraint expr=0 to the inequality tableau using an
602 // artificial variable. To do this, create an artificial variable
603 // av and Add av=expr to the inequality tableau, then make av be 0.
604 // (Raise an exception if we can't attain av=0 -- and prepare explanation)
606 ClSimplexSolver::AddWithArtificialVariable(ClLinearExpression &expr,
607 ExCLRequiredFailureWithExplanation &e)
610 Tracer TRACER(__FUNCTION__);
611 cerr << "(" << expr << ")" << endl;
614 // Allocate the objects on the heap because the objects
615 // will remain in the tableau if we throw an exception,
616 // and that will result in the destructor cleaning up
618 ClSlackVariable *pav = new ClSlackVariable(++_artificialCounter,"a");
619 ClObjectiveVariable *paz = new ClObjectiveVariable("az");
620 ClLinearExpression *pazRow = new ClLinearExpression(expr);
621 // the artificial objective is av, which we know is equal to expr
622 // (which contains only parametric variables)
625 cerr << "aC = " << _artificialCounter
626 << "\nDeletes = " << _cArtificialVarsDeleted << endl;
629 cerr << __FUNCTION__ << " before addRow-s:\n"
633 // the artificial objective is av, which we know is equal to expr
634 // (which contains only parametric variables)
636 // objective is treated as a row in the tableau,
637 // so do the substitution for its value (we are minimizing
638 // the artificial variable)
639 // this row will be removed from the tableau after optimizing
640 addRow(*paz,*pazRow);
642 // now Add the normal row to the tableau -- when artifical
643 // variable is minimized to 0 (if possible)
644 // this row remains in the tableau to maintain the constraint
645 // we are trying to Add
649 cerr << __FUNCTION__ << " after addRow-s:\n"
653 // try to Optimize az to 0
654 // note we are *not* optimizing the real objective, but optimizing
655 // the artificial objective to see if the error in the constraint
656 // we are adding can be set to 0
659 // Careful, we want to get the Expression that is in
660 // the tableau, not the one we initialized it with!
661 ClLinearExpression *pazTableauRow = RowExpression(*paz);
663 cerr << "pazTableauRow->Constant() == " << pazTableauRow->Constant() << endl;
666 // Check that we were able to make the objective value 0
667 // If not, the original constraint was not satisfiable
668 if (!ClApprox(pazTableauRow->Constant(),0.0))
670 BuildExplanation(e, paz, pazTableauRow);
671 // remove the artificial objective row that we just
673 delete RemoveRow(*paz);
674 // and delete the artificial objective variable that we also added above
679 // see if av is a basic variable
680 const ClLinearExpression *pe = RowExpression(*pav);
683 // FIXGJB: do we ever even get here?
684 // Find another variable in this row and Pivot, so that av becomes parametric
685 // If there isn't another variable in the row then
686 // the tableau contains the equation av = 0 -- just delete av's row
687 if (pe->IsConstant())
689 // FIXGJB: do we ever get here?
690 assert(ClApprox(pe->Constant(),0.0));
691 delete RemoveRow(*pav);
692 // remove the temporary objective function
693 // FIXGJB may need this too: delete RemoveRow(*paz);
696 ++_cArtificialVarsDeleted;
700 ClVariable entryVar = pe->AnyPivotableVariable();
701 if (entryVar.IsNil())
703 BuildExplanation(e, *pav, pe);
704 return false; /* required failure */
706 Pivot(entryVar, *pav);
708 // now av should be parametric
709 assert(RowExpression(*pav) == NULL);
713 ++_cArtificialVarsDeleted;
715 // remove the temporary objective function
716 delete RemoveRow(*paz);
722 // Using the given equation (av = cle) build an explanation which
723 // implicates all constraints used to construct the equation. That
724 // is, everything for which the variables in the equation are markers.
725 void ClSimplexSolver::BuildExplanation(ExCLRequiredFailureWithExplanation &e,
727 const ClLinearExpression *pcle)
729 ClVarToConstraintMap::iterator it_cn;
730 it_cn = _constraintsMarked.find(av);
731 if (it_cn != _constraintsMarked.end())
733 e.AddConstraint((*it_cn).second);
736 assert(pcle != NULL);
738 const ClVarToNumberMap & terms = pcle->Terms();
739 ClVarToNumberMap::const_iterator it_term;
740 for (it_term = terms.begin(); it_term != terms.end(); it_term++)
742 it_cn = _constraintsMarked.find((*it_term).first);
743 if (it_cn != _constraintsMarked.end())
745 e.AddConstraint((*it_cn).second);
752 // We are trying to Add the constraint expr=0 to the appropriate
753 // tableau. Try to Add expr directly to the tableaus without
754 // creating an artificial variable. Return true if successful and
757 ClSimplexSolver::TryAddingDirectly(ClLinearExpression &expr)
760 Tracer TRACER(__FUNCTION__);
761 cerr << "(" << expr << ")" << endl;
763 ClVariable subject = ChooseSubject(expr);
764 if (subject.get_pclv() == NULL )
767 cerr << "- returning false" << endl;
771 expr.NewSubject(subject);
772 if (ColumnsHasKey(subject))
774 SubstituteOut(subject,expr);
776 addRow(subject,expr);
778 cerr << "- returning true" << endl;
780 return true; // successfully added directly
784 // We are trying to Add the constraint expr=0 to the tableaux. Try
785 // to choose a subject (a variable to become basic) from among the
786 // current variables in expr. If expr contains any unrestricted
787 // variables, then we must choose an unrestricted variable as the
788 // subject. Also, if the subject is new to the solver we won't have
789 // to do any substitutions, so we prefer new variables to ones that
790 // are currently noted as parametric. If expr contains only
791 // restricted variables, if there is a restricted variable with a
792 // negative coefficient that is new to the solver we can make that
793 // the subject. Otherwise we can't find a subject, so return nil.
794 // (In this last case we have to Add an artificial variable and use
795 // that variable as the subject -- this is done outside this method
798 // Note: in checking for variables that are new to the solver, we
799 // ignore whether a variable occurs in the objective function, since
800 // new slack variables are added to the objective function by
801 // 'NewExpression:', which is called before this method.
803 ClSimplexSolver::ChooseSubject(ClLinearExpression &expr)
806 Tracer TRACER(__FUNCTION__);
807 cerr << "(" << expr << ")" << endl;
809 ClVariable subject(clvNil); // the current best subject, if any
811 // true iff we have found a subject that is an unrestricted variable
812 bool foundUnrestricted = false;
814 // true iff we have found a restricted variable that is new to the
815 // solver (except for being in the obj. function) and that has a
816 // negative coefficient
817 bool foundNewRestricted = false;
819 const ClVarToNumberMap &terms = expr.Terms();
820 ClVarToNumberMap::const_iterator it = terms.begin();
821 for ( ; it != terms.end(); ++it )
823 ClVariable v = (*it).first;
824 Number c = (*it).second;
826 if (foundUnrestricted)
828 // We have already found an unrestricted variable. The only
829 // time we will want to use v instead of the current choice
830 // 'subject' is if v is unrestricted and new to the solver and
831 // 'subject' isn't new. If this is the case just pick v
832 // immediately and return.
833 if (!v.IsRestricted())
835 if (!ColumnsHasKey(v))
840 { // we haven't found an restricted variable yet
841 if (v.IsRestricted())
843 // v is restricted. If we have already found a suitable
844 // restricted variable just stick with that. Otherwise, if v
845 // is new to the solver and has a negative coefficient pick
846 // it. Regarding being new to the solver -- if the variable
847 // occurs only in the objective function we regard it as being
848 // new to the solver, since error variables are added to the
849 // objective function when we make the Expression. We also
850 // never pick a dummy variable here.
851 if (!foundNewRestricted && !v.IsDummy() && c < 0.0)
853 const ClTableauColumnsMap &col = Columns();
854 ClTableauColumnsMap::const_iterator it_col = col.find(v);
855 if (it_col == col.end() ||
856 ( col.size() == 1 && ColumnsHasKey(_objective) ) )
859 foundNewRestricted = true;
865 // v is unrestricted.
866 // If v is also new to the solver just pick it now
868 foundUnrestricted = true;
872 if (!subject.IsNil())
876 // Make one last check -- if all of the variables in expr are dummy
877 // variables, then we can pick a dummy variable as the subject
880 for ( ; it != terms.end(); ++it )
882 ClVariable v = (*it).first;
883 Number c = (*it).second;
885 return clvNil; // nope, no luck
886 // if v is new to the solver, tentatively make it the subject
887 if (!ColumnsHasKey(v))
894 // If we get this far, all of the variables in the Expression should
895 // be dummy variables. If the Constant is nonzero we are trying to
896 // Add an unsatisfiable required constraint. (Remember that dummy
897 // variables must take on a value of 0.) Otherwise, if the Constant
898 // is Zero, multiply by -1 if necessary to make the coefficient for
899 // the subject negative."
900 if (!ClApprox(expr.Constant(),0.0))
902 #ifdef CL_DEBUG_FAILURES
903 cerr << "required failure in choose subject:\n"
908 ExCLRequiredFailureWithExplanation e;
909 BuildExplanation(e, clvNil, &expr);
913 throw ExCLRequiredFailure();
922 // Each of the non-required edits will be represented by an equation
924 // v = c + eplus - eminus
925 // where v is the variable with the edit, c is the previous edit
926 // value, and eplus and eminus are slack variables that hold the
927 // error in satisfying the edit constraint. We are about to change
928 // something, and we want to fix the constants in the equations
929 // representing the edit constraints. If one of eplus and eminus is
930 // basic, the other must occur only in the Expression for that basic
931 // error variable. (They can't both be basic.) Fix the Constant in
932 // this Expression. Otherwise they are both nonbasic. Find all of
933 // the expressions in which they occur, and fix the constants in
934 // those. See the UIST paper for details.
935 // (This comment was for resetEditConstants(), but that is now
936 // gone since it was part of the screwey vector-based interface
937 // to resolveing. --02/15/99 gjb)
939 ClSimplexSolver::DeltaEditConstant(Number delta,
940 ClVariable plusErrorVar,
941 ClVariable minusErrorVar)
944 Tracer TRACER(__FUNCTION__);
945 cerr << "(" << delta << ", " << plusErrorVar << ", " << minusErrorVar << ")" << endl;
947 // first check if the plusErrorVar is basic
948 ClLinearExpression *pexprPlus = RowExpression(plusErrorVar);
949 if (pexprPlus != NULL )
951 pexprPlus->IncrementConstant(delta);
952 // error variables are always restricted
953 // so the row is infeasible if the Constant is negative
954 if (pexprPlus->Constant() < 0.0)
956 _infeasibleRows.insert(plusErrorVar);
960 // check if minusErrorVar is basic
961 ClLinearExpression *pexprMinus = RowExpression(minusErrorVar);
962 if (pexprMinus != NULL)
964 pexprMinus->IncrementConstant(-delta);
965 if (pexprMinus->Constant() < 0.0)
967 _infeasibleRows.insert(minusErrorVar);
971 // Neither is basic. So they must both be nonbasic, and will both
972 // occur in exactly the same expressions. Find all the expressions
973 // in which they occur by finding the column for the minusErrorVar
974 // (it doesn't matter whether we look for that one or for
975 // plusErrorVar). Fix the constants in these expressions.
976 ClVarSet &columnVars = _columns[minusErrorVar];
977 ClVarSet::iterator it = columnVars.begin();
978 for (; it != columnVars.end(); ++it)
980 ClVariable basicVar = *it;
981 ClLinearExpression *pexpr = RowExpression(basicVar);
982 assert(pexpr != NULL );
983 double c = pexpr->CoefficientFor(minusErrorVar);
984 pexpr->IncrementConstant(c*delta);
985 if (basicVar.IsRestricted() && pexpr->Constant() < 0.0)
987 _infeasibleRows.insert(basicVar);
992 // We have set new values for the constants in the edit constraints.
993 // Re-Optimize using the dual simplex algorithm.
995 ClSimplexSolver::DualOptimize()
998 Tracer TRACER(__FUNCTION__);
999 cerr << "()" << endl;
1001 const ClLinearExpression *pzRow = RowExpression(_objective);
1002 // need to handle infeasible rows
1003 while (!_infeasibleRows.empty())
1005 ClVarSet::iterator it_exitVar = _infeasibleRows.begin();
1006 ClVariable exitVar = *it_exitVar;
1007 _infeasibleRows.erase(it_exitVar);
1008 ClVariable entryVar;
1009 // exitVar might have become basic after some other pivoting
1010 // so allow for the case of its not being there any longer
1011 ClLinearExpression *pexpr = RowExpression(exitVar);
1014 // make sure the row is still not feasible
1015 if (pexpr->Constant() < 0.0)
1017 double ratio = DBL_MAX;
1019 ClVarToNumberMap &terms = pexpr->Terms();
1020 ClVarToNumberMap::iterator it = terms.begin();
1021 for ( ; it != terms.end(); ++it )
1023 ClVariable v = (*it).first;
1024 Number c = (*it).second;
1025 if (c > 0.0 && v.IsPivotable())
1027 Number zc = pzRow->CoefficientFor(v);
1028 r = zc/c; // FIXGJB r:= zc/c or Zero, as ClSymbolicWeight-s
1036 if (ratio == DBL_MAX)
1039 ss << "ratio == nil (DBL_MAX)" << ends;
1040 throw ExCLInternalError(ss.str().c_str());
1042 Pivot(entryVar,exitVar);
1048 // Make a new linear Expression representing the constraint cn,
1049 // replacing any basic variables with their defining expressions.
1050 // Normalize if necessary so that the Constant is non-negative. If
1051 // the constraint is non-required give its error variables an
1052 // appropriate weight in the objective function.
1053 ClLinearExpression *
1054 ClSimplexSolver::NewExpression(const ClConstraint *pcn,
1056 ClVariable &clvEplus,
1057 ClVariable &clvEminus,
1058 Number &prevEConstant)
1061 Tracer TRACER(__FUNCTION__);
1062 cerr << "(" << *pcn << ")" << endl;
1063 cerr << "cn.IsInequality() == " << pcn->IsInequality() << endl;
1064 cerr << "cn.IsRequired() == " << pcn->IsRequired() << endl;
1066 const ClLinearExpression &cnExpr = pcn->Expression();
1067 cl_auto_ptr<ClLinearExpression> pexpr ( new ClLinearExpression(cnExpr.Constant()) );
1068 cl_auto_ptr<ClSlackVariable> pslackVar;
1069 cl_auto_ptr<ClDummyVariable> pdummyVar;
1070 cl_auto_ptr<ClSlackVariable> peminus(0);
1071 cl_auto_ptr<ClSlackVariable> peplus(0);
1072 const ClVarToNumberMap &cnTerms = cnExpr.Terms();
1073 ClVarToNumberMap::const_iterator it = cnTerms.begin();
1074 for ( ; it != cnTerms.end(); ++it)
1076 ClVariable v = (*it).first;
1077 Number c = (*it).second;
1078 const ClLinearExpression *pe = RowExpression(v);
1081 pexpr->AddVariable(v,c);
1085 pexpr->AddExpression(*pe,c);
1089 // Add slack and error variables as needed
1090 if (pcn->IsInequality())
1092 // cn is an inequality, so Add a slack variable. The original
1093 // constraint is expr>=0, so that the resulting equality is
1094 // expr-slackVar=0. If cn is also non-required Add a negative
1095 // error variable, giving
1096 // expr-slackVar = -errorVar, in other words
1097 // expr-slackVar+errorVar=0.
1098 // Since both of these variables are newly created we can just Add
1099 // them to the Expression (they can't be basic).
1101 ReinitializeAutoPtr(pslackVar,new ClSlackVariable (_slackCounter, "s"));
1102 pexpr->setVariable(*pslackVar,-1);
1103 // index the constraint under its slack variable and vice-versa
1104 _markerVars[pcn] = pslackVar.get();
1105 _constraintsMarked[pslackVar.get()] = pcn;
1107 if (!pcn->IsRequired())
1110 ReinitializeAutoPtr(peminus,new ClSlackVariable (_slackCounter, "em"));
1111 pexpr->setVariable(*peminus,1.0);
1112 // Add emnius to the objective function with the appropriate weight
1113 ClLinearExpression *pzRow = RowExpression(_objective);
1114 // FIXGJB: pzRow->AddVariable(eminus,pcn->strength().symbolicWeight() * pcn->weight());
1115 ClSymbolicWeight sw = pcn->strength().symbolicWeight().Times(pcn->weight());
1116 pzRow->setVariable(*peminus,sw.AsDouble());
1117 _errorVars[pcn].insert(peminus.get());
1118 NoteAddedVariable(*peminus,_objective);
1122 { // cn is an equality
1123 if (pcn->IsRequired())
1125 // Add a dummy variable to the Expression to serve as a marker
1126 // for this constraint. The dummy variable is never allowed to
1127 // enter the basis when pivoting.
1129 ReinitializeAutoPtr(pdummyVar,new ClDummyVariable (_dummyCounter, "d"));
1130 pexpr->setVariable(*pdummyVar,1.0);
1131 _markerVars[pcn] = pdummyVar.get();
1132 _constraintsMarked[pdummyVar.get()] = pcn;
1134 cerr << "Adding dummyVar == d" << _dummyCounter << endl;
1139 // cn is a non-required equality. Add a positive and a negative
1140 // error variable, making the resulting constraint
1141 // expr = eplus - eminus,
1142 // in other words: expr-eplus+eminus=0
1144 ReinitializeAutoPtr(peplus,new ClSlackVariable (_slackCounter, "ep"));
1145 ReinitializeAutoPtr(peminus,new ClSlackVariable (_slackCounter, "em"));
1147 pexpr->setVariable(*peplus,-1.0);
1148 pexpr->setVariable(*peminus,1.0);
1149 // index the constraint under one of the error variables
1150 _markerVars[pcn] = peplus.get();
1151 _constraintsMarked[peplus.get()] = pcn;
1153 ClLinearExpression *pzRow = RowExpression(_objective);
1154 // FIXGJB: pzRow->AddVariable(eplus,pcn->strength().symbolicWeight() * pcn->weight());
1155 ClSymbolicWeight sw = pcn->strength().symbolicWeight().Times(pcn->weight());
1156 double swCoeff = sw.AsDouble();
1160 cerr << "sw == " << sw << endl
1161 << "cn == " << *pcn << endl;
1162 cerr << "adding " << *peplus << " and " << *peminus
1163 << " with swCoeff == " << swCoeff << endl;
1166 pzRow->setVariable(*peplus,swCoeff);
1167 NoteAddedVariable(*peplus,_objective);
1168 // FIXGJB: pzRow->AddVariable(eminus,pcn->strength().symbolicWeight() * pcn->weight());
1169 pzRow->setVariable(*peminus,swCoeff);
1170 NoteAddedVariable(*peminus,_objective);
1171 _errorVars[pcn].insert(peminus.get());
1172 _errorVars[pcn].insert(peplus.get());
1173 if (pcn->isStayConstraint())
1175 _stayPlusErrorVars.push_back(peplus.get());
1176 _stayMinusErrorVars.push_back(peminus.get());
1178 else if (pcn->IsEditConstraint())
1180 clvEplus = peplus.get();
1181 clvEminus = peminus.get();
1182 prevEConstant = cnExpr.Constant();
1187 // the Constant in the Expression should be non-negative.
1188 // If necessary normalize the Expression by multiplying by -1
1189 if (pexpr->Constant() < 0)
1192 cerr << "NewExpression's Constant is " << pexpr->Constant() << ", < 0, so flipping" << endl;
1194 pexpr->MultiplyMe(-1);
1197 cerr << "- returning " << *pexpr << endl;
1199 // Terrible Name -- release() does *not* delete the object,
1200 // only makes sure that the destructor won't delete the object
1201 // (it releases the cl_auto_ptr from the responsibility of deleting the object)
1202 pslackVar.release();
1203 pdummyVar.release();
1206 return pexpr.release();
1209 // Minimize the value of the objective. (The tableau should already
1212 ClSimplexSolver::Optimize(ClVariable zVar)
1215 Tracer TRACER(__FUNCTION__);
1216 cerr << "(" << zVar << ")\n"
1219 ClLinearExpression *pzRow = RowExpression(zVar);
1220 assert(pzRow != NULL);
1221 ClVariable entryVar = clvNil;
1222 ClVariable exitVar = clvNil;
1225 Number objectiveCoeff = 0;
1226 // Find the most negative coefficient in the objective function
1227 // (ignoring the non-pivotable dummy variables). If all
1228 // coefficients are positive we're done
1229 ClVarToNumberMap &terms = pzRow->Terms();
1230 ClVarToNumberMap::iterator it = terms.begin();
1231 for (; it != terms.end(); ++it)
1233 ClVariable v = (*it).first;
1234 Number c = (*it).second;
1235 if (v.IsPivotable() && c < objectiveCoeff)
1239 // A. Beurive' Tue Jul 13 23:03:05 CEST 1999 Why the most
1240 // negative? I encountered unending cycles of pivots!
1244 // if all coefficients were positive (or if the objective
1245 // function has no pivotable variables)
1246 // we are at an optimum
1247 if (objectiveCoeff >= -_epsilon)
1250 cerr << "entryVar == " << entryVar << ", "
1251 << "objectiveCoeff == " << objectiveCoeff
1255 // choose which variable to move out of the basis
1256 // Only consider pivotable basic variables
1257 // (i.e. restricted, non-dummy variables)
1258 double minRatio = DBL_MAX;
1259 ClVarSet &columnVars = _columns[entryVar];
1260 ClVarSet::iterator it_rowvars = columnVars.begin();
1262 for (; it_rowvars != columnVars.end(); ++it_rowvars)
1264 ClVariable v = *it_rowvars;
1266 cerr << "Checking " << v << endl;
1268 if (v.IsPivotable())
1270 const ClLinearExpression *pexpr = RowExpression(v);
1271 Number coeff = pexpr->CoefficientFor(entryVar);
1272 // only consider negative coefficients
1275 r = - pexpr->Constant() / coeff;
1279 cerr << "New minRatio == " << r << endl;
1287 // If minRatio is still nil at this point, it means that the
1288 // objective function is unbounded, i.e. it can become
1289 // arbitrarily negative. This should never happen in this
1291 if (minRatio == DBL_MAX)
1294 ss << "objective function is unbounded!" << ends;
1295 throw ExCLInternalError(ss.str().c_str());
1297 Pivot(entryVar, exitVar);
1299 cerr << "After Optimize:\n"
1305 // Do a Pivot. Move entryVar into the basis (i.e. make it a basic variable),
1306 // and move exitVar out of the basis (i.e., make it a parametric variable)
1308 ClSimplexSolver::Pivot(ClVariable entryVar, ClVariable exitVar)
1311 Tracer TRACER(__FUNCTION__);
1312 cerr << "(" << entryVar << ", " << exitVar << ")" << endl;
1315 // the entryVar might be non-pivotable if we're doing a RemoveConstraint --
1316 // otherwise it should be a pivotable variable -- enforced at call sites,
1319 // expr is the Expression for the exit variable (about to leave the basis) --
1320 // so that the old tableau includes the equation:
1322 ClLinearExpression *pexpr = RemoveRow(exitVar);
1324 // Compute an Expression for the entry variable. Since expr has
1325 // been deleted from the tableau we can destructively modify it to
1326 // build this Expression.
1327 pexpr->ChangeSubject(exitVar,entryVar);
1328 SubstituteOut(entryVar,*pexpr);
1330 if (entryVar.IsExternal())
1332 // entry var is no longer a parametric variable since we're moving
1333 // it into the basis
1334 _externalParametricVars.erase(entryVar);
1336 addRow(entryVar,*pexpr);
1341 // Each of the non-required stays will be represented by an equation
1343 // v = c + eplus - eminus
1344 // where v is the variable with the stay, c is the previous value of
1345 // v, and eplus and eminus are slack variables that hold the error
1346 // in satisfying the stay constraint. We are about to change
1347 // something, and we want to fix the constants in the equations
1348 // representing the stays. If both eplus and eminus are nonbasic
1349 // they have value 0 in the current solution, meaning the previous
1350 // stay was exactly satisfied. In this case nothing needs to be
1351 // changed. Otherwise one of them is basic, and the other must
1352 // occur only in the Expression for that basic error variable.
1353 // Reset the Constant in this Expression to 0.
1355 ClSimplexSolver::ResetStayConstants()
1358 Tracer TRACER(__FUNCTION__);
1359 cerr << "()" << endl;
1361 ClVarVector::const_iterator
1362 itStayPlusErrorVars = _stayPlusErrorVars.begin();
1363 ClVarVector::const_iterator
1364 itStayMinusErrorVars = _stayMinusErrorVars.begin();
1366 for ( ; itStayPlusErrorVars != _stayPlusErrorVars.end();
1367 ++itStayPlusErrorVars, ++itStayMinusErrorVars )
1369 ClLinearExpression *pexpr = RowExpression(*itStayPlusErrorVars);
1372 pexpr = RowExpression(*itStayMinusErrorVars);
1376 pexpr->Set_constant(0.0);
1381 // Set the external variables known to this solver to their appropriate values.
1382 // Set each external basic variable to its value, and set each
1383 // external parametric variable to 0. (It isn't clear that we will
1384 // ever have external parametric variables -- every external
1385 // variable should either have a stay on it, or have an equation
1386 // that defines it in terms of other external variables that do have
1387 // stays. For the moment I'll put this in though.) Variables that
1388 // are internal to the solver don't actually store values -- their
1389 // values are just implicit in the tableu -- so we don't need to set
1392 ClSimplexSolver::SetExternalVariables()
1395 Tracer TRACER(__FUNCTION__);
1400 // FIXGJB -- oughta check some invariants here
1402 // Set external parametric variables first
1403 // in case I've screwed up
1404 ClVarSet::iterator itParVars = _externalParametricVars.begin();
1405 for ( ; itParVars != _externalParametricVars.end(); ++itParVars )
1407 ClVariable v = *itParVars;
1409 // defensively skip it if it is basic -- ChangeValue is virtual
1410 // so don't want to call it twice; this should never
1416 cerr << __FUNCTION__ << "Error: variable " << v
1417 << " in _externalParametricVars is basic" << endl;
1418 cerr << "Row is: " << *RowExpression(v) << endl;
1426 // Only iterate over the rows w/ external variables
1427 ClVarSet::iterator itRowVars = _externalRows.begin();
1428 for ( ; itRowVars != _externalRows.end() ; ++itRowVars )
1430 ClVariable v = *itRowVars;
1431 ClLinearExpression *pexpr = RowExpression(v);
1432 ChangeClv(v,pexpr->Constant());
1435 _fNeedsSolving = false;
1436 if (_pfnResolveCallback)
1437 _pfnResolveCallback(this);
1442 PrintTo(ostream &xo, const ClVarVector &varlist)
1444 ClVarVector::const_iterator it = varlist.begin();
1445 xo << varlist.size() << ":" << "[ ";
1446 if (it != varlist.end())
1451 for (; it != varlist.end(); ++it)
1459 ostream &operator<<(ostream &xo, const ClVarVector &varlist)
1460 { return PrintTo(xo,varlist); }
1464 PrintTo(ostream &xo, const ClConstraintToVarSetMap &mapCnToVarSet)
1466 ClConstraintToVarSetMap::const_iterator it = mapCnToVarSet.begin();
1467 for ( ; it != mapCnToVarSet.end(); ++it) {
1468 const ClConstraint *pcn = (*it).first;
1469 const ClVarSet &set = (*it).second;
1470 xo << "CN: " << pcn << *pcn << ":: " << set << endl;
1475 ostream &operator <<(ostream &xo, const ClConstraintToVarSetMap &mapCnToVarSet)
1476 { return PrintTo(xo,mapCnToVarSet); }
1481 ClSimplexSolver::PrintOn(ostream &xo) const
1483 ClTableau::PrintOn(xo);
1485 xo << "_stayPlusErrorVars: "
1486 << _stayPlusErrorVars << endl;
1487 xo << "_stayMinusErrorVars: "
1488 << _stayMinusErrorVars << endl;
1489 xo << "_editInfoList:\n"
1490 << _editInfoList << endl;
1496 ClSimplexSolver::PrintInternalInfo(ostream &xo) const
1498 ClTableau::PrintInternalInfo(xo);
1499 xo << "; edvars: " << _editInfoList.size();
1501 printExternalVariablesTo(xo);
1505 ostream &operator<<(ostream &xo, const ClSimplexSolver &clss)
1507 return clss.PrintOn(xo);
1513 ClSimplexSolver::FIsConstraintSatisfied(const ClConstraint *const pcn) const
1515 ClConstraintToVarMap::const_iterator it_marker = _markerVars.find(pcn);
1516 if (it_marker == _markerVars.end())
1517 { // could not find the constraint
1518 throw ExCLConstraintNotFound();
1522 bool fCnsays = pcn->FIsSatisfied();
1525 ClConstraintToVarSetMap::const_iterator it_eVars = _errorVars.find(pcn);
1527 if (it_eVars != _errorVars.end())
1529 const ClVarSet &eVars = (*it_eVars).second;
1530 ClVarSet::const_iterator it = eVars.begin();
1531 for ( ; it != eVars.end(); ++it )
1533 const ClLinearExpression *pexpr = RowExpression(*it);
1534 if (pexpr != NULL && !ClApprox(pexpr->Constant(),0.0))
1538 cerr << __FUNCTION__ << ": constraint says satisfiable, but solver does not" << endl;
1547 cerr << __FUNCTION__ << ": solver says satisfiable, but constraint does not" << endl;
1556 ostream &PrintTo(ostream &xo, const ClSimplexSolver::ClEditInfoList &listPEditInfo)
1558 ClSimplexSolver::ClEditInfoList::const_iterator it = listPEditInfo.begin();
1559 for ( ; it != listPEditInfo.end(); ++it) {
1560 const ClSimplexSolver::ClEditInfo *pcei = (*it);
1561 xo << *pcei << endl;
1567 ostream &operator<<(ostream &xo, const ClSimplexSolver::ClEditInfoList &listPEditInfo)
1568 { return PrintTo(xo,listPEditInfo); }
1572 // A. Beurive' Tue Jul 6 17:03:32 CEST 1999
1574 ClSimplexSolver::ChangeStrengthAndWeight(ClConstraint *pcn, const ClStrength &strength, double weight)
1576 ClConstraintToVarSetMap::iterator it_eVars = _errorVars.find(pcn);
1577 // Only for constraints that already have error variables (i.e. non-required constraints)
1578 assert(it_eVars != _errorVars.end());
1580 ClLinearExpression *pzRow = RowExpression(_objective);
1582 Number old_coeff = pcn->weight() * pcn->strength().symbolicWeight().AsDouble();
1583 pcn->setStrength(strength);
1584 pcn->setWeight(weight);
1585 Number new_coeff = pcn->weight() * pcn->strength().symbolicWeight().AsDouble();
1587 if (new_coeff != old_coeff)
1590 cerr << "Changing strength and/or weight for constraint: " << endl << *pcn << endl;
1591 cerr << "Updating objective row from:" << endl << *pzRow << endl;
1593 ClVarSet &eVars = (*it_eVars).second;
1594 ClVarSet::iterator it = eVars.begin();
1595 for ( ; it != eVars.end(); ++it )
1597 const ClLinearExpression *pexpr = RowExpression(*it);
1600 pzRow->AddVariable(*it,-old_coeff,_objective,*this);
1601 pzRow->AddVariable(*it,new_coeff,_objective,*this);
1605 pzRow->AddExpression(*pexpr,-old_coeff,_objective,*this);
1606 pzRow->AddExpression(*pexpr,new_coeff,_objective,*this);
1610 cerr << "to: " << endl << *pzRow << endl;
1615 Optimize(_objective);
1616 SetExternalVariables();
1621 // A. Beurive' Tue Jul 6 17:03:42 CEST 1999
1623 ClSimplexSolver::ChangeStrength(ClConstraint *pcn, const ClStrength &strength)
1625 ChangeStrengthAndWeight(pcn,strength,pcn->weight());
1628 // A. Beurive' Tue Jul 6 17:03:42 CEST 1999
1630 ClSimplexSolver::ChangeWeight(ClConstraint *pcn, double weight)
1632 ChangeStrengthAndWeight(pcn,pcn->strength(),weight);