xref: /haiku/src/libs/linprog/LinearSpec.cpp (revision fc691d7de2182d23659b86d87c9c36b0feaa6b40)
1 /*
2  * Copyright 2007-2008, Christof Lutteroth, lutteroth@cs.auckland.ac.nz
3  * Copyright 2007-2008, James Kim, jkim202@ec.auckland.ac.nz
4  * Copyright 2010, Clemens Zeidler <haiku@clemens-zeidler.de>
5  * Distributed under the terms of the MIT License.
6  */
7 
8 
9 #include "LinearSpec.h"
10 
11 
12 /**
13  * Constructor.
14  * Creates a new specification for a linear programming problem.
15  */
16 LinearSpec::LinearSpec()
17 	:
18 	fLpPresolved(NULL),
19 	fOptimization(MINIMIZE),
20 	fObjFunction(new SummandList()),
21 	fResult(ERROR),
22 	fObjectiveValue(NAN),
23 	fSolvingTime(NAN)
24 {
25 	fLP = make_lp(0, 0);
26 	if (fLP == NULL)
27 		printf("Couldn't construct a new model.");
28 	set_verbose(fLP, 1);
29 }
30 
31 
32 /**
33  * Destructor.
34  * Removes the specification and deletes all constraints,
35  * objective function summands and variables.
36  */
37 LinearSpec::~LinearSpec()
38 {
39 	RemovePresolved();
40 	for (int32 i = 0; i < fConstraints.CountItems(); i++)
41 		delete (Constraint*)fConstraints.ItemAt(i);
42 	for (int32 i = 0; i < fObjFunction->CountItems(); i++)
43 		delete (Summand*)fObjFunction->ItemAt(i);
44 	while (fVariables.CountItems() > 0)
45 		RemoveVariable(fVariables.ItemAt(0));
46 
47 	delete_lp(fLP);
48 
49 	delete fObjFunction;
50 }
51 
52 
53 /**
54  * Adds a new variable to the specification.
55  *
56  * @return the new variable
57  */
58 Variable*
59 LinearSpec::AddVariable()
60 {
61 	Variable* variable = new Variable(this);
62 	if (!variable)
63 		return NULL;
64 	if (!AddVariable(variable)) {
65 		delete variable;
66 		return NULL;
67 	}
68 
69 	return variable;
70 }
71 
72 
73 bool
74 LinearSpec::AddVariable(Variable* variable)
75 {
76 	double d = 0;
77 	int i = 0;
78 
79 	if (!fVariables.AddItem(variable))
80 		return false;
81 	if (add_columnex(fLP, 0, &d, &i) == 0) {
82 		fVariables.RemoveItem(variable);
83 		return false;
84 	}
85 
86 	if (!SetRange(variable, -20000, 20000)) {
87 		RemoveVariable(variable, false);
88 		return false;
89 	}
90 
91 	variable->fIsValid = true;
92 	return true;
93 }
94 
95 
96 bool
97 LinearSpec::RemoveVariable(Variable* variable, bool deleteVariable)
98 {
99 	int32 index = IndexOf(variable);
100 	if (index < 0)
101 		return false;
102 
103 	if (!del_column(fLP, index))
104 		return false;
105 	fVariables.RemoveItemAt(index - 1);
106 	variable->Invalidate();
107 
108 	if (deleteVariable)
109 		delete variable;
110 	return true;
111 }
112 
113 
114 int32
115 LinearSpec::IndexOf(const Variable* variable) const
116 {
117 	int32 i = fVariables.IndexOf(variable);
118 	if (i == -1) {
119 		printf("Variable 0x%p not part of fLS->Variables().\n", variable);
120 		return -1;
121 	}
122 	return i + 1;
123 }
124 
125 
126 bool
127 LinearSpec::SetRange(Variable* variable, double min, double max)
128 {
129 	return set_bounds(fLP, IndexOf(variable), min, max);
130 }
131 
132 
133 /**
134  * Adds a new hard linear constraint to the specification.
135  *
136  * @param coeffs		the constraint's coefficients
137  * @param vars			the constraint's variables
138  * @param op			the constraint's operand
139  * @param rightSide	the constant value on the constraint's right side
140  * @return the new constraint
141  */
142 Constraint*
143 LinearSpec::AddConstraint(SummandList* summands, OperatorType op,
144 	double rightSide)
145 {
146 	Constraint* c = new Constraint(this, summands, op, rightSide,
147 		INFINITY, INFINITY);
148 	RemovePresolved();
149 	return c;
150 }
151 
152 
153 /**
154  * Adds a new hard linear constraint to the specification with a single summand.
155  *
156  * @param coeff1		the constraint's first coefficient
157  * @param var1			the constraint's first variable
158  * @param op			the constraint's operand
159  * @param rightSide	the constant value on the constraint's right side
160  * @return the new constraint
161  */
162 Constraint*
163 LinearSpec::AddConstraint(double coeff1, Variable* var1,
164 	OperatorType op, double rightSide)
165 {
166 	SummandList* summands = new SummandList(1);
167 	summands->AddItem(new Summand(coeff1, var1));
168 	Constraint* c = new Constraint(this, summands, op, rightSide,
169 		INFINITY, INFINITY);
170 	RemovePresolved();
171 	return c;
172 }
173 
174 
175 /**
176  * Adds a new hard linear constraint to the specification with two summands.
177  *
178  * @param coeff1		the constraint's first coefficient
179  * @param var1			the constraint's first variable
180  * @param coeff2		the constraint's second coefficient
181  * @param var2			the constraint's second variable
182  * @param op			the constraint's operand
183  * @param rightSide	the constant value on the constraint's right side
184  * @return the new constraint
185  */
186 Constraint*
187 LinearSpec::AddConstraint(double coeff1, Variable* var1,
188 	double coeff2, Variable* var2, OperatorType op, double rightSide)
189 {
190 	SummandList* summands = new SummandList(2);
191 	summands->AddItem(new Summand(coeff1, var1));
192 	summands->AddItem(new Summand(coeff2, var2));
193 	Constraint* c = new Constraint(this, summands, op, rightSide,
194 		INFINITY, INFINITY);
195 	RemovePresolved();
196 	return c;
197 }
198 
199 
200 /**
201  * Adds a new hard linear constraint to the specification with three summands.
202  *
203  * @param coeff1		the constraint's first coefficient
204  * @param var1			the constraint's first variable
205  * @param coeff2		the constraint's second coefficient
206  * @param var2			the constraint's second variable
207  * @param coeff3		the constraint's third coefficient
208  * @param var3			the constraint's third variable
209  * @param op			the constraint's operand
210  * @param rightSide	the constant value on the constraint's right side
211  * @return the new constraint
212  */
213 Constraint*
214 LinearSpec::AddConstraint(double coeff1, Variable* var1,
215 		double coeff2, Variable* var2, double coeff3, Variable* var3,
216 		OperatorType op, double rightSide)
217 {
218 	SummandList* summands = new SummandList(3);
219 	summands->AddItem(new Summand(coeff1, var1));
220 	summands->AddItem(new Summand(coeff2, var2));
221 	summands->AddItem(new Summand(coeff3, var3));
222 	Constraint* c = new Constraint(this, summands, op, rightSide,
223 		INFINITY, INFINITY);
224 	RemovePresolved();
225 	return c;
226 }
227 
228 
229 /**
230  * Adds a new hard linear constraint to the specification with four summands.
231  *
232  * @param coeff1		the constraint's first coefficient
233  * @param var1			the constraint's first variable
234  * @param coeff2		the constraint's second coefficient
235  * @param var2			the constraint's second variable
236  * @param coeff3		the constraint's third coefficient
237  * @param var3			the constraint's third variable
238  * @param coeff4		the constraint's fourth coefficient
239  * @param var4			the constraint's fourth variable
240  * @param op			the constraint's operand
241  * @param rightSide	the constant value on the constraint's right side
242  * @return the new constraint
243  */
244 Constraint*
245 LinearSpec::AddConstraint(double coeff1, Variable* var1,
246 		double coeff2, Variable* var2, double coeff3, Variable* var3,
247 		double coeff4, Variable* var4, OperatorType op, double rightSide)
248 {
249 	SummandList* summands = new SummandList(3);
250 	summands->AddItem(new Summand(coeff1, var1));
251 	summands->AddItem(new Summand(coeff2, var2));
252 	summands->AddItem(new Summand(coeff3, var3));
253 	summands->AddItem(new Summand(coeff4, var4));
254 	Constraint* c = new Constraint(this, summands, op, rightSide,
255 		INFINITY, INFINITY);
256 	RemovePresolved();
257 	return c;
258 }
259 
260 
261 /**
262  * Adds a new soft linear constraint to the specification.
263  * i.e. a constraint that does not always have to be satisfied.
264  *
265  * @param coeffs		the constraint's coefficients
266  * @param vars			the constraint's variables
267  * @param op			the constraint's operand
268  * @param rightSide		the constant value on the constraint's right side
269  * @param penaltyNeg	the coefficient penalizing negative deviations from the exact solution
270  * @param penaltyPos	the coefficient penalizing positive deviations from the exact solution
271  */
272 Constraint*
273 LinearSpec::AddConstraint(SummandList* summands, OperatorType op,
274 		double rightSide, double penaltyNeg, double penaltyPos)
275 {
276 	Constraint* c = new Constraint(this, summands, op, rightSide,
277 			penaltyNeg, penaltyPos);
278 	RemovePresolved();
279 	return c;
280 }
281 
282 
283 /**
284  * Adds a new soft linear constraint to the specification with a single summand.
285  *
286  * @param coeff1		the constraint's first coefficient
287  * @param var1			the constraint's first variable
288  * @param op			the constraint's operand
289  * @param rightSide		the constant value on the constraint's right side
290  * @param penaltyNeg	the coefficient penalizing negative deviations from the exact solution
291  * @param penaltyPos	the coefficient penalizing positive deviations from the exact solution
292  */
293 Constraint*
294 LinearSpec::AddConstraint(double coeff1, Variable* var1,
295 		OperatorType op, double rightSide, double penaltyNeg, double penaltyPos)
296 {
297 	SummandList* summands = new SummandList(1);
298 	summands->AddItem(new Summand(coeff1, var1));
299 	Constraint* c = new Constraint(this, summands, op, rightSide,
300 		penaltyNeg, penaltyPos);
301 	RemovePresolved();
302 	return c;
303 }
304 
305 
306 /**
307  * Adds a new soft linear constraint to the specification with two summands.
308  *
309  * @param coeff1		the constraint's first coefficient
310  * @param var1			the constraint's first variable
311  * @param coeff2		the constraint's second coefficient
312  * @param var2			the constraint's second variable
313  * @param op			the constraint's operand
314  * @param rightSide		the constant value on the constraint's right side
315  * @param penaltyNeg	the coefficient penalizing negative deviations from the exact solution
316  * @param penaltyPos	the coefficient penalizing positive deviations from the exact solution
317  */
318 Constraint*
319 LinearSpec::AddConstraint(double coeff1, Variable* var1,
320 	double coeff2, Variable* var2, OperatorType op, double rightSide,
321 	double penaltyNeg, double penaltyPos)
322 {
323 	SummandList* summands = new SummandList(2);
324 	summands->AddItem(new Summand(coeff1, var1));
325 	summands->AddItem(new Summand(coeff2, var2));
326 	Constraint* c = new Constraint(this, summands, op, rightSide,
327 		penaltyNeg, penaltyPos);
328 	RemovePresolved();
329 	return c;
330 }
331 
332 
333 /**
334  * Adds a new soft linear constraint to the specification with three summands.
335  *
336  * @param coeff1		the constraint's first coefficient
337  * @param var1			the constraint's first variable
338  * @param coeff2		the constraint's second coefficient
339  * @param var2			the constraint's second variable
340  * @param coeff3		the constraint's third coefficient
341  * @param var3			the constraint's third variable
342  * @param op			the constraint's operand
343  * @param rightSide		the constant value on the constraint's right side
344  * @param penaltyNeg	the coefficient penalizing negative deviations from the exact solution
345  * @param penaltyPos	the coefficient penalizing positive deviations from the exact solution
346  */
347 Constraint*
348 LinearSpec::AddConstraint(double coeff1, Variable* var1,
349 	double coeff2, Variable* var2, double coeff3, Variable* var3,
350 	OperatorType op, double rightSide, double penaltyNeg, double penaltyPos)
351 {
352 	SummandList* summands = new SummandList(2);
353 	summands->AddItem(new Summand(coeff1, var1));
354 	summands->AddItem(new Summand(coeff2, var2));
355 	summands->AddItem(new Summand(coeff3, var3));
356 	Constraint* c = new Constraint(this, summands, op, rightSide,
357 		penaltyNeg, penaltyPos);
358 	RemovePresolved();
359 	return c;
360 }
361 
362 
363 /**
364  * Adds a new soft linear constraint to the specification with four summands.
365  *
366  * @param coeff1		the constraint's first coefficient
367  * @param var1			the constraint's first variable
368  * @param coeff2		the constraint's second coefficient
369  * @param var2			the constraint's second variable
370  * @param coeff3		the constraint's third coefficient
371  * @param var3			the constraint's third variable
372  * @param coeff4		the constraint's fourth coefficient
373  * @param var4			the constraint's fourth variable
374  * @param op			the constraint's operand
375  * @param rightSide		the constant value on the constraint's right side
376  * @param penaltyNeg	the coefficient penalizing negative deviations from the exact solution
377  * @param penaltyPos	the coefficient penalizing positive deviations from the exact solution
378  */
379 Constraint*
380 LinearSpec::AddConstraint(double coeff1, Variable* var1,
381 	double coeff2, Variable* var2, double coeff3, Variable* var3,
382 	double coeff4, Variable* var4, OperatorType op, double rightSide,
383 	double penaltyNeg, double penaltyPos)
384 {
385 	SummandList* summands = new SummandList(2);
386 	summands->AddItem(new Summand(coeff1, var1));
387 	summands->AddItem(new Summand(coeff2, var2));
388 	summands->AddItem(new Summand(coeff3, var3));
389 	summands->AddItem(new Summand(coeff4, var4));
390 	Constraint* c = new Constraint(this, summands, op, rightSide,
391 		penaltyNeg, penaltyPos);
392 	RemovePresolved();
393 	return c;
394 }
395 
396 
397 /**
398  * Adds a new penalty function to the specification.
399  *
400  * @param var	the penalty function's variable
401  * @param xs	the penalty function's sampling points
402  * @param gs	the penalty function's gradients
403  * @return the new penalty function
404  */
405 PenaltyFunction*
406 LinearSpec::AddPenaltyFunction(Variable* var, BList* xs, BList* gs)
407 {
408 	return new PenaltyFunction(this, var, xs, gs);
409 }
410 
411 
412 /**
413  * Gets the objective function.
414  *
415  * @return SummandList containing the objective function's summands
416  */
417 SummandList*
418 LinearSpec::ObjectiveFunction()
419 {
420 	return fObjFunction;
421 }
422 
423 
424 SummandList*
425 LinearSpec::ReplaceObjectiveFunction(SummandList* objFunction)
426 {
427 	SummandList* list = fObjFunction;
428 	fObjFunction = objFunction;
429 	UpdateObjectiveFunction();
430 	return list;
431 }
432 
433 
434 /**
435  * Sets a new objective function.
436  *
437  * @param summands	SummandList containing the objective function's summands
438  */
439 void
440 LinearSpec::SetObjectiveFunction(SummandList* objFunction)
441 {
442 	for (int32 i = 0; i < fObjFunction->CountItems(); i++)
443 		delete (Summand*)fObjFunction->ItemAt(i);
444 	delete fObjFunction;
445 
446 	fObjFunction = objFunction;
447 	UpdateObjectiveFunction();
448 }
449 
450 
451 /**
452  * Updates the internal representation of the objective function.
453  * Must be called whenever the summands of the objective function are changed.
454  */
455 void
456 LinearSpec::UpdateObjectiveFunction()
457 {
458 	int32 size = fObjFunction->CountItems();
459 	double coeffs[size];
460 	int varIndexes[size];
461 	Summand* current;
462 	for (int32 i = 0; i < size; i++) {
463 		current = (Summand*)fObjFunction->ItemAt(i);
464 		coeffs[i] = current->Coeff();
465 		varIndexes[i] = current->Var()->Index();
466 	}
467 
468 	if (!set_obj_fnex(fLP, size, &coeffs[0], &varIndexes[0]))
469 		printf("Error in set_obj_fnex.\n");
470 
471 	RemovePresolved();
472 }
473 
474 
475 /**
476  * Remove a cached presolved model, if existent.
477  * This is automatically done each time after the model has been changed,
478  * to avoid an old cached presolved model getting out of sync.
479  */
480 void
481 LinearSpec::RemovePresolved()
482 {
483 	if (fLpPresolved == NULL)
484 		return;
485 	delete_lp(fLpPresolved);
486 	fLpPresolved = NULL;
487 }
488 
489 
490 /**
491  * Creates and caches a simplified version of the linear programming problem,
492  * where redundant rows, columns and constraints are removed,
493  * if it has not been created before.
494  * Then, the simplified problem is solved.
495  *
496  * @return the result of the solving attempt
497  */
498 ResultType
499 LinearSpec::Presolve()
500 {
501 	bigtime_t start, end;
502 	start = system_time();
503 
504 	if (fLpPresolved == NULL) {
505 		fLpPresolved = copy_lp(fLP);
506 		set_presolve(fLpPresolved, PRESOLVE_ROWS | PRESOLVE_COLS
507 			| PRESOLVE_LINDEP, get_presolveloops(fLpPresolved));
508 	}
509 
510 	fResult = (ResultType)solve(fLpPresolved);
511 	fObjectiveValue = get_objective(fLpPresolved);
512 
513 	if (fResult == OPTIMAL) {
514 		int32 size = fVariables.CountItems();
515 		for (int32 i = 0; i < size; i++) {
516 			Variable* current = (Variable*)fVariables.ItemAt(i);
517 			current->SetValue(get_var_primalresult(fLpPresolved,
518 					get_Norig_rows(fLpPresolved) + current->Index()));
519 		}
520 	}
521 
522 	end = system_time();
523 	fSolvingTime = (end - start) / 1000.0;
524 
525 	return fResult;
526 }
527 
528 
529 /**
530  * Tries to solve the linear programming problem.
531  * If a cached simplified version of the problem exists, it is used instead.
532  *
533  * @return the result of the solving attempt
534  */
535 ResultType
536 LinearSpec::Solve()
537 {
538 	if (fLpPresolved != NULL)
539 		return Presolve();
540 
541 	bigtime_t start, end;
542 	start = system_time();
543 
544 	fResult = (ResultType)solve(fLP);
545 	fObjectiveValue = get_objective(fLP);
546 
547 	if (fResult == OPTIMAL) {
548 		int32 size = fVariables.CountItems();
549 		double x[size];
550 		if (!get_variables(fLP, &x[0]))
551 			printf("Error in get_variables.\n");
552 
553 		int32 i = 0;
554 		while (i < size) {
555 			((Variable*)fVariables.ItemAt(i))->SetValue(x[i]);
556 			i++;
557 		}
558 	}
559 
560 	end = system_time();
561 	fSolvingTime = (end - start) / 1000.0;
562 
563 	return fResult;
564 }
565 
566 
567 /**
568  * Writes the specification into a text file.
569  * The file will be overwritten if it exists.
570  *
571  * @param fname	the file name
572  */
573 void
574 LinearSpec::Save(const char* fileName)
575 {
576 	// TODO: Constness should be fixed in liblpsolve API.
577 	write_lp(fLP, const_cast<char*>(fileName));
578 }
579 
580 
581 /**
582  * Gets the current optimization.
583  * The default is minimization.
584  *
585  * @return the current optimization
586  */
587 OptimizationType
588 LinearSpec::Optimization() const
589 {
590 	return fOptimization;
591 }
592 
593 
594 /**
595  * Sets whether the solver should minimize or maximize the objective function.
596  * The default is minimization.
597  *
598  * @param optimization	the optimization type
599  */
600 void
601 LinearSpec::SetOptimization(OptimizationType value)
602 {
603 	fOptimization = value;
604 	if (fOptimization == MINIMIZE)
605 		set_minim(fLP);
606 	else
607 		set_maxim(fLP);
608 }
609 
610 
611 /**
612  * Gets the constraints.
613  *
614  * @return the constraints
615  */
616 const ConstraintList&
617 LinearSpec::Constraints() const
618 {
619 	return fConstraints;
620 }
621 
622 
623 /**
624  * Gets the result type.
625  *
626  * @return the result type
627  */
628 ResultType
629 LinearSpec::Result() const
630 {
631 	return fResult;
632 }
633 
634 
635 /**
636  * Gets the objective value.
637  *
638  * @return the objective value
639  */
640 double
641 LinearSpec::ObjectiveValue() const
642 {
643 	return fObjectiveValue;
644 }
645 
646 
647 /**
648  * Gets the solving time.
649  *
650  * @return the solving time
651  */
652 double
653 LinearSpec::SolvingTime() const
654 {
655 	return fSolvingTime;
656 }
657 
658 
659 LinearSpec::operator BString() const
660 {
661 	BString string;
662 	GetString(string);
663 	return string;
664 }
665 
666 
667 void
668 LinearSpec::GetString(BString& string) const
669 {
670 	string << "LinearSpec " << (int32)this << ":\n";
671 	for (int i = 0; i < fVariables.CountItems(); i++) {
672 		Variable* variable = static_cast<Variable*>(fVariables.ItemAt(i));
673 		variable->GetString(string);
674 		string << "=" << (float)variable->Value() << " ";
675 	}
676 	string << "\n";
677 	for (int i = 0; i < fConstraints.CountItems(); i++) {
678 		Constraint* c = static_cast<Constraint*>(fConstraints.ItemAt(i));
679 		string << i << ": ";
680 		c->GetString(string);
681 		string << "\n";
682 	}
683 	string << "Result=";
684 	if (fResult==-1)
685 		string << "ERROR";
686 	else if (fResult==0)
687 		string << "OPTIMAL";
688 	else if (fResult==1)
689 		string << "SUBOPTIMAL";
690 	else if (fResult==2)
691 		string << "INFEASIBLE";
692 	else if (fResult==3)
693 		string << "UNBOUNDED";
694 	else if (fResult==4)
695 		string << "DEGENERATE";
696 	else if (fResult==5)
697 		string << "NUMFAILURE";
698 	else
699 		string << fResult;
700 	string << " SolvingTime=" << (float)fSolvingTime << "ms";
701 }
702 
703