1 /* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 /*
3 * This file is part of the LibreOffice project.
4 *
5 * This Source Code Form is subject to the terms of the Mozilla Public
6 * License, v. 2.0. If a copy of the MPL was not distributed with this
7 * file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 *
9 * This file incorporates work covered by the following license notice:
10 *
11 * Licensed to the Apache Software Foundation (ASF) under one or more
12 * contributor license agreements. See the NOTICE file distributed
13 * with this work for additional information regarding copyright
14 * ownership. The ASF licenses this file to you under the Apache
15 * License, Version 2.0 (the "License"); you may not use this file
16 * except in compliance with the License. You may obtain a copy of
17 * the License at http://www.apache.org/licenses/LICENSE-2.0 .
18 */
19
20 #include <CoinMP.h>
21 #include <CoinError.hpp>
22
23 #include "SolverComponent.hxx"
24 #include <strings.hrc>
25
26 #include <com/sun/star/frame/XModel.hpp>
27 #include <com/sun/star/table/CellAddress.hpp>
28
29 #include <rtl/math.hxx>
30 #include <stdexcept>
31 #include <vector>
32 #include <float.h>
33
34 namespace com::sun::star::uno { class XComponentContext; }
35
36 using namespace com::sun::star;
37
38 namespace {
39
40 class CoinMPSolver : public SolverComponent
41 {
42 public:
CoinMPSolver()43 CoinMPSolver() {}
44
45 private:
46 virtual void SAL_CALL solve() override;
getImplementationName()47 virtual OUString SAL_CALL getImplementationName() override
48 {
49 return "com.sun.star.comp.Calc.CoinMPSolver";
50 }
getComponentDescription()51 virtual OUString SAL_CALL getComponentDescription() override
52 {
53 return SolverComponent::GetResourceString( RID_COINMP_SOLVER_COMPONENT );
54 }
55 };
56
57 }
58
solve()59 void SAL_CALL CoinMPSolver::solve()
60 {
61 uno::Reference<frame::XModel> xModel( mxDoc, uno::UNO_QUERY_THROW );
62
63 maStatus.clear();
64 mbSuccess = false;
65
66 xModel->lockControllers();
67
68 // collect variables in vector (?)
69
70 auto aVariableCells = comphelper::sequenceToContainer<std::vector<table::CellAddress>>(maVariables);
71 size_t nVariables = aVariableCells.size();
72 size_t nVar = 0;
73
74 // collect all dependent cells
75
76 ScSolverCellHashMap aCellsHash;
77 aCellsHash[maObjective].reserve( nVariables + 1 ); // objective function
78
79 for (const auto& rConstr : std::as_const(maConstraints))
80 {
81 table::CellAddress aCellAddr = rConstr.Left;
82 aCellsHash[aCellAddr].reserve( nVariables + 1 ); // constraints: left hand side
83
84 if ( rConstr.Right >>= aCellAddr )
85 aCellsHash[aCellAddr].reserve( nVariables + 1 ); // constraints: right hand side
86 }
87
88 // set all variables to zero
89 //! store old values?
90 //! use old values as initial values?
91 for ( const auto& rVarCell : aVariableCells )
92 {
93 SolverComponent::SetValue( mxDoc, rVarCell, 0.0 );
94 }
95
96 // read initial values from all dependent cells
97 for ( auto& rEntry : aCellsHash )
98 {
99 double fValue = SolverComponent::GetValue( mxDoc, rEntry.first );
100 rEntry.second.push_back( fValue ); // store as first element, as-is
101 }
102
103 // loop through variables
104 for ( const auto& rVarCell : aVariableCells )
105 {
106 SolverComponent::SetValue( mxDoc, rVarCell, 1.0 ); // set to 1 to examine influence
107
108 // read value change from all dependent cells
109 for ( auto& rEntry : aCellsHash )
110 {
111 double fChanged = SolverComponent::GetValue( mxDoc, rEntry.first );
112 double fInitial = rEntry.second.front();
113 rEntry.second.push_back( fChanged - fInitial );
114 }
115
116 SolverComponent::SetValue( mxDoc, rVarCell, 2.0 ); // minimal test for linearity
117
118 for ( const auto& rEntry : aCellsHash )
119 {
120 double fInitial = rEntry.second.front();
121 double fCoeff = rEntry.second.back(); // last appended: coefficient for this variable
122 double fTwo = SolverComponent::GetValue( mxDoc, rEntry.first );
123
124 bool bLinear = rtl::math::approxEqual( fTwo, fInitial + 2.0 * fCoeff ) ||
125 rtl::math::approxEqual( fInitial, fTwo - 2.0 * fCoeff );
126 // second comparison is needed in case fTwo is zero
127 if ( !bLinear )
128 maStatus = SolverComponent::GetResourceString( RID_ERROR_NONLINEAR );
129 }
130
131 SolverComponent::SetValue( mxDoc, rVarCell, 0.0 ); // set back to zero for examining next variable
132 }
133
134 xModel->unlockControllers();
135
136 if ( !maStatus.isEmpty() )
137 return;
138
139 //
140 // build parameter arrays for CoinMP
141 //
142
143 // set objective function
144
145 const std::vector<double>& rObjCoeff = aCellsHash[maObjective];
146 std::unique_ptr<double[]> pObjectCoeffs(new double[nVariables]);
147 for (nVar=0; nVar<nVariables; nVar++)
148 pObjectCoeffs[nVar] = rObjCoeff[nVar+1];
149 double nObjectConst = rObjCoeff[0]; // constant term of objective
150
151 // add rows
152
153 size_t nRows = maConstraints.getLength();
154 size_t nCompSize = nVariables * nRows;
155 std::unique_ptr<double[]> pCompMatrix(new double[nCompSize]); // first collect all coefficients, row-wise
156 for (size_t i=0; i<nCompSize; i++)
157 pCompMatrix[i] = 0.0;
158
159 std::unique_ptr<double[]> pRHS(new double[nRows]);
160 std::unique_ptr<char[]> pRowType(new char[nRows]);
161 for (size_t i=0; i<nRows; i++)
162 {
163 pRHS[i] = 0.0;
164 pRowType[i] = 'N';
165 }
166
167 for (sal_Int32 nConstrPos = 0; nConstrPos < maConstraints.getLength(); ++nConstrPos)
168 {
169 // integer constraints are set later
170 sheet::SolverConstraintOperator eOp = maConstraints[nConstrPos].Operator;
171 if ( eOp == sheet::SolverConstraintOperator_LESS_EQUAL ||
172 eOp == sheet::SolverConstraintOperator_GREATER_EQUAL ||
173 eOp == sheet::SolverConstraintOperator_EQUAL )
174 {
175 double fDirectValue = 0.0;
176 bool bRightCell = false;
177 table::CellAddress aRightAddr;
178 const uno::Any& rRightAny = maConstraints[nConstrPos].Right;
179 if ( rRightAny >>= aRightAddr )
180 bRightCell = true; // cell specified as right-hand side
181 else
182 rRightAny >>= fDirectValue; // constant value
183
184 table::CellAddress aLeftAddr = maConstraints[nConstrPos].Left;
185
186 const std::vector<double>& rLeftCoeff = aCellsHash[aLeftAddr];
187 double* pValues = &pCompMatrix[nConstrPos * nVariables];
188 for (nVar=0; nVar<nVariables; nVar++)
189 pValues[nVar] = rLeftCoeff[nVar+1];
190
191 // if left hand cell has a constant term, put into rhs value
192 double fRightValue = -rLeftCoeff[0];
193
194 if ( bRightCell )
195 {
196 const std::vector<double>& rRightCoeff = aCellsHash[aRightAddr];
197 // modify pValues with rhs coefficients
198 for (nVar=0; nVar<nVariables; nVar++)
199 pValues[nVar] -= rRightCoeff[nVar+1];
200
201 fRightValue += rRightCoeff[0]; // constant term
202 }
203 else
204 fRightValue += fDirectValue;
205
206 switch ( eOp )
207 {
208 case sheet::SolverConstraintOperator_LESS_EQUAL: pRowType[nConstrPos] = 'L'; break;
209 case sheet::SolverConstraintOperator_GREATER_EQUAL: pRowType[nConstrPos] = 'G'; break;
210 case sheet::SolverConstraintOperator_EQUAL: pRowType[nConstrPos] = 'E'; break;
211 default:
212 OSL_ENSURE( false, "unexpected enum type" );
213 }
214 pRHS[nConstrPos] = fRightValue;
215 }
216 }
217
218 // Find non-zero coefficients, column-wise
219
220 std::unique_ptr<int[]> pMatrixBegin(new int[nVariables+1]);
221 std::unique_ptr<int[]> pMatrixCount(new int[nVariables]);
222 std::unique_ptr<double[]> pMatrix(new double[nCompSize]); // not always completely used
223 std::unique_ptr<int[]> pMatrixIndex(new int[nCompSize]);
224 int nMatrixPos = 0;
225 for (nVar=0; nVar<nVariables; nVar++)
226 {
227 int nBegin = nMatrixPos;
228 for (size_t nRow=0; nRow<nRows; nRow++)
229 {
230 double fCoeff = pCompMatrix[ nRow * nVariables + nVar ]; // row-wise
231 if ( fCoeff != 0.0 )
232 {
233 pMatrix[nMatrixPos] = fCoeff;
234 pMatrixIndex[nMatrixPos] = nRow;
235 ++nMatrixPos;
236 }
237 }
238 pMatrixBegin[nVar] = nBegin;
239 pMatrixCount[nVar] = nMatrixPos - nBegin;
240 }
241 pMatrixBegin[nVariables] = nMatrixPos;
242 pCompMatrix.reset();
243
244 // apply settings to all variables
245
246 std::unique_ptr<double[]> pLowerBounds(new double[nVariables]);
247 std::unique_ptr<double[]> pUpperBounds(new double[nVariables]);
248 for (nVar=0; nVar<nVariables; nVar++)
249 {
250 pLowerBounds[nVar] = mbNonNegative ? 0.0 : -DBL_MAX;
251 pUpperBounds[nVar] = DBL_MAX;
252
253 // bounds could possibly be further restricted from single-cell constraints
254 }
255
256 std::unique_ptr<char[]> pColType(new char[nVariables]);
257 for (nVar=0; nVar<nVariables; nVar++)
258 pColType[nVar] = mbInteger ? 'I' : 'C';
259
260 // apply single-var integer constraints
261
262 for (const auto& rConstr : std::as_const(maConstraints))
263 {
264 sheet::SolverConstraintOperator eOp = rConstr.Operator;
265 if ( eOp == sheet::SolverConstraintOperator_INTEGER ||
266 eOp == sheet::SolverConstraintOperator_BINARY )
267 {
268 table::CellAddress aLeftAddr = rConstr.Left;
269 // find variable index for cell
270 for (nVar=0; nVar<nVariables; nVar++)
271 if ( AddressEqual( aVariableCells[nVar], aLeftAddr ) )
272 {
273 if ( eOp == sheet::SolverConstraintOperator_INTEGER )
274 pColType[nVar] = 'I';
275 else
276 {
277 pColType[nVar] = 'B';
278 pLowerBounds[nVar] = 0.0;
279 pUpperBounds[nVar] = 1.0;
280 }
281 }
282 }
283 }
284
285 int nObjectSense = mbMaximize ? SOLV_OBJSENS_MAX : SOLV_OBJSENS_MIN;
286
287 HPROB hProb = CoinCreateProblem("");
288 int nResult = CoinLoadProblem( hProb, nVariables, nRows, nMatrixPos, 0,
289 nObjectSense, nObjectConst, pObjectCoeffs.get(),
290 pLowerBounds.get(), pUpperBounds.get(), pRowType.get(), pRHS.get(), nullptr,
291 pMatrixBegin.get(), pMatrixCount.get(), pMatrixIndex.get(), pMatrix.get(),
292 nullptr, nullptr, nullptr );
293 if (nResult == SOLV_CALL_SUCCESS)
294 {
295 nResult = CoinLoadInteger( hProb, pColType.get() );
296 }
297
298 pColType.reset();
299 pMatrixIndex.reset();
300 pMatrix.reset();
301 pMatrixCount.reset();
302 pMatrixBegin.reset();
303 pUpperBounds.reset();
304 pLowerBounds.reset();
305 pRowType.reset();
306 pRHS.reset();
307 pObjectCoeffs.reset();
308
309 CoinSetRealOption( hProb, COIN_REAL_MAXSECONDS, mnTimeout );
310 CoinSetRealOption( hProb, COIN_REAL_MIPMAXSEC, mnTimeout );
311
312 // TODO: handle (or remove) settings: epsilon, B&B depth
313
314 // solve model
315
316 if (nResult == SOLV_CALL_SUCCESS)
317 {
318 nResult = CoinCheckProblem( hProb );
319 }
320
321 if (nResult == SOLV_CALL_SUCCESS)
322 {
323 try
324 {
325 nResult = CoinOptimizeProblem( hProb, 0 );
326 }
327 catch (const CoinError& e)
328 {
329 CoinUnloadProblem(hProb);
330 throw std::runtime_error(e.message());
331 }
332 }
333
334 mbSuccess = ( nResult == SOLV_CALL_SUCCESS );
335 if ( mbSuccess )
336 {
337 // get solution
338
339 maSolution.realloc( nVariables );
340 CoinGetSolutionValues( hProb, maSolution.getArray(), nullptr, nullptr, nullptr );
341 mfResultValue = CoinGetObjectValue( hProb );
342 }
343 else
344 {
345 int nSolutionStatus = CoinGetSolutionStatus( hProb );
346 if ( nSolutionStatus == 1 )
347 maStatus = SolverComponent::GetResourceString( RID_ERROR_INFEASIBLE );
348 else if ( nSolutionStatus == 2 )
349 maStatus = SolverComponent::GetResourceString( RID_ERROR_UNBOUNDED );
350 // TODO: detect timeout condition and report as RID_ERROR_TIMEOUT
351 // (currently reported as infeasible)
352 }
353
354 CoinUnloadProblem( hProb );
355 }
356
357 extern "C" SAL_DLLPUBLIC_EXPORT css::uno::XInterface *
com_sun_star_comp_Calc_CoinMPSolver_get_implementation(css::uno::XComponentContext *,css::uno::Sequence<css::uno::Any> const &)358 com_sun_star_comp_Calc_CoinMPSolver_get_implementation(
359 css::uno::XComponentContext *,
360 css::uno::Sequence<css::uno::Any> const &)
361 {
362 return cppu::acquire(new CoinMPSolver());
363 }
364
365 /* vim:set shiftwidth=4 softtabstop=4 expandtab: */
366