sccomp/source/solver/SwarmSolver.cxx

   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  */
  10
  11 #include <sal/config.h>
  12
  13 #include <com/sun/star/frame/XModel.hpp>
  14 #include <com/sun/star/container/XIndexAccess.hpp>
  15 #include <com/sun/star/sheet/XSpreadsheetDocument.hpp>
  16 #include <com/sun/star/sheet/XSpreadsheet.hpp>
  17 #include <com/sun/star/sheet/XSolver.hpp>
  18 #include <com/sun/star/sheet/XSolverDescription.hpp>
  19 #include <com/sun/star/table/CellAddress.hpp>
  20 #include <com/sun/star/table/CellContentType.hpp>
  21 #include <com/sun/star/table/XCell.hpp>
  22 #include <com/sun/star/lang/XServiceInfo.hpp>
  23
  24 #include <rtl/math.hxx>
  25 #include <cppuhelper/implbase.hxx>
  26 #include <cppuhelper/supportsservice.hxx>
  27
  28 #include <comphelper/compbase.hxx>
  29 #include <comphelper/propertycontainer2.hxx>
  30 #include <comphelper/proparrhlp.hxx>
  31
  32 #include <cmath>
  33 #include <vector>
  34 #include <limits>
  35 #include <chrono>
  36 #include <random>
  37
  38 #include <unotools/resmgr.hxx>
  39
  40 #include "DifferentialEvolution.hxx"
  41 #include "ParticelSwarmOptimization.hxx"
  42
  43 #include <strings.hrc>
  44
  45 namespace com::sun::star::uno
  46 {
  47 class XComponentContext;
  48 }
  49
  50 using namespace css;
  51
  52 namespace
  53 {
  54 struct Bound
  55 {
  56     double lower;
  57     double upper;
  58
  59     Bound()
  60         // float bounds should be low/high enough for all practical uses
  61         // otherwise we are too far away from the solution
  62         : lower(std::numeric_limits<float>::lowest())
  63         , upper(std::numeric_limits<float>::max())
  64     {
  65     }
  66
  67     void updateBound(sheet::SolverConstraintOperator eOp, double fValue)
  68     {
  69         if (eOp == sheet::SolverConstraintOperator_LESS_EQUAL)
  70         {
  71             // if we set the bound multiple times use the one which includes both values
  72             // for example bound values 100, 120, 150 -> use 100 -> the lowest one
  73             if (fValue < upper)
  74                 upper = fValue;
  75         }
  76         else if (eOp == sheet::SolverConstraintOperator_GREATER_EQUAL)
  77         {
  78             if (fValue > lower)
  79                 lower = fValue;
  80         }
  81         else if (eOp == sheet::SolverConstraintOperator_EQUAL)
  82         {
  83             lower = fValue;
  84             upper = fValue;
  85         }
  86     }
  87 };
  88
  89 enum
  90 {
  91     PROP_NONNEGATIVE,
  92     PROP_INTEGER,
  93     PROP_TIMEOUT,
  94     PROP_ALGORITHM,
  95 };
  96
  97 } // end anonymous namespace
  98
  99 typedef comphelper::WeakImplHelper<sheet::XSolver, sheet::XSolverDescription, lang::XServiceInfo>
 100     SwarmSolver_Base;
 101
 102 namespace
 103 {
 104 class SwarmSolver : public comphelper::OPropertyContainer2,
 105                     public comphelper::OPropertyArrayUsageHelper<SwarmSolver>,
 106                     public SwarmSolver_Base
 107 {
 108 private:
 109     uno::Reference<sheet::XSpreadsheetDocument> mxDocument;
 110     table::CellAddress maObjective;
 111     uno::Sequence<table::CellAddress> maVariables;
 112     uno::Sequence<sheet::SolverConstraint> maConstraints;
 113     bool mbMaximize;
 114
 115     // set via XPropertySet
 116     bool mbNonNegative;
 117     bool mbInteger;
 118     sal_Int32 mnTimeout;
 119     sal_Int32 mnAlgorithm;
 120
 121     // results
 122     bool mbSuccess;
 123     double mfResultValue;
 124
 125     uno::Sequence<double> maSolution;
 126     OUString maStatus;
 127
 128     std::vector<Bound> maBounds;
 129     std::vector<sheet::SolverConstraint> maNonBoundedConstraints;
 130
 131 private:
 132     static OUString getResourceString(TranslateId aId);
 133
 134     uno::Reference<table::XCell> getCell(const table::CellAddress& rPosition);
 135     void setValue(const table::CellAddress& rPosition, double fValue);
 136     double getValue(const table::CellAddress& rPosition);
 137
 138 public:
 139     SwarmSolver()
 140         : mbMaximize(true)
 141         , mbNonNegative(false)
 142         , mbInteger(false)
 143         , mnTimeout(60000)
 144         , mnAlgorithm(0)
 145         , mbSuccess(false)
 146         , mfResultValue(0.0)
 147     {
 148         registerProperty(u"NonNegative"_ustr, PROP_NONNEGATIVE, 0, &mbNonNegative,
 149                          cppu::UnoType<decltype(mbNonNegative)>::get());
 150         registerProperty(u"Integer"_ustr, PROP_INTEGER, 0, &mbInteger,
 151                          cppu::UnoType<decltype(mbInteger)>::get());
 152         registerProperty(u"Timeout"_ustr, PROP_TIMEOUT, 0, &mnTimeout,
 153                          cppu::UnoType<decltype(mnTimeout)>::get());
 154         registerProperty(u"Algorithm"_ustr, PROP_ALGORITHM, 0, &mnAlgorithm,
 155                          cppu::UnoType<decltype(mnAlgorithm)>::get());
 156     }
 157
 158     DECLARE_XINTERFACE()
 159     DECLARE_XTYPEPROVIDER()
 160
 161     virtual uno::Reference<beans::XPropertySetInfo> SAL_CALL getPropertySetInfo() override
 162     {
 163         return createPropertySetInfo(getInfoHelper());
 164     }
 165     // OPropertySetHelper
 166     virtual cppu::IPropertyArrayHelper& getInfoHelper() override { return *getArrayHelper(); }
 167     // OPropertyArrayUsageHelper
 168     virtual cppu::IPropertyArrayHelper* createArrayHelper() const override
 169     {
 170         uno::Sequence<beans::Property> aProperties;
 171         describeProperties(aProperties);
 172         return new cppu::OPropertyArrayHelper(aProperties);
 173     }
 174
 175     // XSolver
 176     virtual uno::Reference<sheet::XSpreadsheetDocument> SAL_CALL getDocument() override
 177     {
 178         return mxDocument;
 179     }
 180     virtual void SAL_CALL
 181     setDocument(const uno::Reference<sheet::XSpreadsheetDocument>& rDocument) override
 182     {
 183         mxDocument = rDocument;
 184     }
 185
 186     virtual table::CellAddress SAL_CALL getObjective() override { return maObjective; }
 187     virtual void SAL_CALL setObjective(const table::CellAddress& rObjective) override
 188     {
 189         maObjective = rObjective;
 190     }
 191
 192     virtual uno::Sequence<table::CellAddress> SAL_CALL getVariables() override
 193     {
 194         return maVariables;
 195     }
 196     virtual void SAL_CALL setVariables(const uno::Sequence<table::CellAddress>& rVariables) override
 197     {
 198         maVariables = rVariables;
 199     }
 200
 201     virtual uno::Sequence<sheet::SolverConstraint> SAL_CALL getConstraints() override
 202     {
 203         return maConstraints;
 204     }
 205     virtual void SAL_CALL
 206     setConstraints(const uno::Sequence<sheet::SolverConstraint>& rConstraints) override
 207     {
 208         maConstraints = rConstraints;
 209     }
 210
 211     virtual sal_Bool SAL_CALL getMaximize() override { return mbMaximize; }
 212     virtual void SAL_CALL setMaximize(sal_Bool bMaximize) override { mbMaximize = bMaximize; }
 213
 214     virtual sal_Bool SAL_CALL getSuccess() override { return mbSuccess; }
 215     virtual double SAL_CALL getResultValue() override { return mfResultValue; }
 216
 217     virtual uno::Sequence<double> SAL_CALL getSolution() override { return maSolution; }
 218
 219     virtual void SAL_CALL solve() override;
 220
 221     // XSolverDescription
 222     virtual OUString SAL_CALL getComponentDescription() override
 223     {
 224         return SwarmSolver::getResourceString(RID_SWARM_SOLVER_COMPONENT);
 225     }
 226
 227     virtual OUString SAL_CALL getStatusDescription() override { return maStatus; }
 228
 229     virtual OUString SAL_CALL getPropertyDescription(const OUString& rPropertyName) override
 230     {
 231         TranslateId pResId;
 232         switch (getInfoHelper().getHandleByName(rPropertyName))
 233         {
 234             case PROP_NONNEGATIVE:
 235                 pResId = RID_PROPERTY_NONNEGATIVE;
 236                 break;
 237             case PROP_INTEGER:
 238                 pResId = RID_PROPERTY_INTEGER;
 239                 break;
 240             case PROP_TIMEOUT:
 241                 pResId = RID_PROPERTY_TIMEOUT;
 242                 break;
 243             case PROP_ALGORITHM:
 244                 pResId = RID_PROPERTY_ALGORITHM;
 245                 break;
 246             default:
 247                 break;
 248         }
 249         return SwarmSolver::getResourceString(pResId);
 250     }
 251
 252     // XServiceInfo
 253     virtual OUString SAL_CALL getImplementationName() override
 254     {
 255         return u"com.sun.star.comp.Calc.SwarmSolver"_ustr;
 256     }
 257
 258     sal_Bool SAL_CALL supportsService(const OUString& rServiceName) override
 259     {
 260         return cppu::supportsService(this, rServiceName);
 261     }
 262
 263     uno::Sequence<OUString> SAL_CALL getSupportedServiceNames() override
 264     {
 265         return { u"com.sun.star.sheet.Solver"_ustr };
 266     }
 267
 268 private:
 269     void applyVariables(std::vector<double> const& rVariables);
 270     bool doesViolateConstraints();
 271
 272 public:
 273     double calculateFitness(std::vector<double> const& rVariables);
 274     size_t getDimensionality() const;
 275     void initializeVariables(std::vector<double>& rVariables, std::mt19937& rGenerator);
 276     double clampVariable(size_t nVarIndex, double fValue);
 277     double boundVariable(size_t nVarIndex, double fValue);
 278 };
 279 }
 280
 281 OUString SwarmSolver::getResourceString(TranslateId aId)
 282 {
 283     if (!aId)
 284         return OUString();
 285
 286     return Translate::get(aId, Translate::Create("scc"));
 287 }
 288
 289 uno::Reference<table::XCell> SwarmSolver::getCell(const table::CellAddress& rPosition)
 290 {
 291     uno::Reference<container::XIndexAccess> xSheets(mxDocument->getSheets(), uno::UNO_QUERY);
 292     uno::Reference<sheet::XSpreadsheet> xSheet(xSheets->getByIndex(rPosition.Sheet),
 293                                                uno::UNO_QUERY);
 294     return xSheet->getCellByPosition(rPosition.Column, rPosition.Row);
 295 }
 296
 297 void SwarmSolver::setValue(const table::CellAddress& rPosition, double fValue)
 298 {
 299     getCell(rPosition)->setValue(fValue);
 300 }
 301
 302 double SwarmSolver::getValue(const table::CellAddress& rPosition)
 303 {
 304     return getCell(rPosition)->getValue();
 305 }
 306
 307 IMPLEMENT_FORWARD_XINTERFACE2(SwarmSolver, SwarmSolver_Base, comphelper::OPropertyContainer2)
 308 IMPLEMENT_FORWARD_XTYPEPROVIDER2(SwarmSolver, SwarmSolver_Base, comphelper::OPropertyContainer2)
 309
 310 void SwarmSolver::applyVariables(std::vector<double> const& rVariables)
 311 {
 312     for (sal_Int32 i = 0; i < maVariables.getLength(); ++i)
 313     {
 314         setValue(maVariables[i], rVariables[i]);
 315     }
 316 }
 317
 318 double SwarmSolver::calculateFitness(std::vector<double> const& rVariables)
 319 {
 320     applyVariables(rVariables);
 321
 322     if (doesViolateConstraints())
 323         return std::numeric_limits<float>::lowest();
 324
 325     double x = getValue(maObjective);
 326
 327     if (mbMaximize)
 328         return x;
 329     else
 330         return -x;
 331 }
 332
 333 void SwarmSolver::initializeVariables(std::vector<double>& rVariables, std::mt19937& rGenerator)
 334 {
 335     int nTry = 1;
 336     bool bConstraintsOK = false;
 337
 338     while (!bConstraintsOK && nTry < 10)
 339     {
 340         size_t noVariables(maVariables.getLength());
 341
 342         rVariables.resize(noVariables);
 343
 344         for (size_t i = 0; i < noVariables; ++i)
 345         {
 346             Bound const& rBound = maBounds[i];
 347             if (mbInteger)
 348             {
 349                 sal_Int64 intLower(rBound.lower);
 350                 sal_Int64 intUpper(rBound.upper);
 351                 std::uniform_int_distribution<sal_Int64> random(intLower, intUpper);
 352                 rVariables[i] = double(random(rGenerator));
 353             }
 354             else
 355             {
 356                 std::uniform_real_distribution<double> random(rBound.lower, rBound.upper);
 357                 rVariables[i] = random(rGenerator);
 358             }
 359         }
 360
 361         applyVariables(rVariables);
 362
 363         bConstraintsOK = !doesViolateConstraints();
 364         nTry++;
 365     }
 366 }
 367
 368 double SwarmSolver::clampVariable(size_t nVarIndex, double fValue)
 369 {
 370     Bound const& rBound = maBounds[nVarIndex];
 371     double fResult = std::clamp(fValue, rBound.lower, rBound.upper);
 372
 373     if (mbInteger)
 374         return std::trunc(fResult);
 375
 376     return fResult;
 377 }
 378
 379 double SwarmSolver::boundVariable(size_t nVarIndex, double fValue)
 380 {
 381     Bound const& rBound = maBounds[nVarIndex];
 382     // double fResult = std::max(std::min(fValue, rBound.upper), rBound.lower);
 383     double fResult = fValue;
 384     while (fResult < rBound.lower || fResult > rBound.upper)
 385     {
 386         if (fResult < rBound.lower)
 387             fResult = rBound.upper - (rBound.lower - fResult);
 388         if (fResult > rBound.upper)
 389             fResult = (fResult - rBound.upper) + rBound.lower;
 390     }
 391
 392     if (mbInteger)
 393         return std::trunc(fResult);
 394
 395     return fResult;
 396 }
 397
 398 size_t SwarmSolver::getDimensionality() const { return maVariables.getLength(); }
 399
 400 bool SwarmSolver::doesViolateConstraints()
 401 {
 402     for (const sheet::SolverConstraint& rConstraint : maNonBoundedConstraints)
 403     {
 404         double fLeftValue = getValue(rConstraint.Left);
 405         double fRightValue = 0.0;
 406
 407         table::CellAddress aCellAddress;
 408
 409         if (rConstraint.Right >>= aCellAddress)
 410         {
 411             fRightValue = getValue(aCellAddress);
 412         }
 413         else if (rConstraint.Right >>= fRightValue)
 414         {
 415             // empty
 416         }
 417         else
 418         {
 419             return false;
 420         }
 421
 422         sheet::SolverConstraintOperator eOp = rConstraint.Operator;
 423         switch (eOp)
 424         {
 425             case sheet::SolverConstraintOperator_LESS_EQUAL:
 426             {
 427                 if (fLeftValue > fRightValue)
 428                     return true;
 429             }
 430             break;
 431             case sheet::SolverConstraintOperator_GREATER_EQUAL:
 432             {
 433                 if (fLeftValue < fRightValue)
 434                     return true;
 435             }
 436             break;
 437             case sheet::SolverConstraintOperator_EQUAL:
 438             {
 439                 if (!rtl::math::approxEqual(fLeftValue, fRightValue))
 440                     return true;
 441             }
 442             break;
 443             default:
 444                 break;
 445         }
 446     }
 447     return false;
 448 }
 449
 450 namespace
 451 {
 452 template <typename SwarmAlgorithm> class SwarmRunner
 453 {
 454 private:
 455     SwarmAlgorithm& mrAlgorithm;
 456     double mfTimeout;
 457
 458     static constexpr size_t mnPopulationSize = 40;
 459     static constexpr int constNumberOfGenerationsWithoutChange = 50;
 460
 461     std::chrono::high_resolution_clock::time_point maStart;
 462     std::chrono::high_resolution_clock::time_point maEnd;
 463
 464 public:
 465     SwarmRunner(SwarmAlgorithm& rAlgorithm)
 466         : mrAlgorithm(rAlgorithm)
 467         , mfTimeout(5000)
 468     {
 469     }
 470
 471     void setTimeout(double fTimeout) { mfTimeout = fTimeout; }
 472
 473     std::vector<double> const& solve()
 474     {
 475         using std::chrono::duration_cast;
 476         using std::chrono::high_resolution_clock;
 477         using std::chrono::milliseconds;
 478
 479         mrAlgorithm.initialize();
 480
 481         maEnd = maStart = high_resolution_clock::now();
 482
 483         int nLastChange = 0;
 484
 485         while ((mrAlgorithm.getGeneration() - nLastChange) < constNumberOfGenerationsWithoutChange
 486                && duration_cast<milliseconds>(maEnd - maStart).count() < mfTimeout)
 487         {
 488             bool bChange = mrAlgorithm.next();
 489
 490             if (bChange)
 491                 nLastChange = mrAlgorithm.getGeneration();
 492
 493             maEnd = high_resolution_clock::now();
 494         }
 495         return mrAlgorithm.getResult();
 496     }
 497 };
 498 }
 499
 500 void SAL_CALL SwarmSolver::solve()
 501 {
 502     uno::Reference<frame::XModel> xModel(mxDocument, uno::UNO_QUERY_THROW);
 503
 504     maStatus.clear();
 505     mbSuccess = false;
 506     if (!maVariables.getLength())
 507         return;
 508
 509     maBounds.resize(maVariables.getLength());
 510
 511     xModel->lockControllers();
 512
 513     if (mbNonNegative)
 514     {
 515         for (Bound& rBound : maBounds)
 516             rBound.lower = 0;
 517     }
 518
 519     // Determine variable bounds
 520     for (sheet::SolverConstraint const& rConstraint : maConstraints)
 521     {
 522         table::CellAddress aLeftCellAddress = rConstraint.Left;
 523         sheet::SolverConstraintOperator eOp = rConstraint.Operator;
 524
 525         size_t index = 0;
 526         bool bFoundVariable = false;
 527         for (const table::CellAddress& rVariableCell : maVariables)
 528         {
 529             if (aLeftCellAddress == rVariableCell)
 530             {
 531                 bFoundVariable = true;
 532                 table::CellAddress aCellAddress;
 533                 double fValue;
 534
 535                 if (rConstraint.Right >>= aCellAddress)
 536                 {
 537                     uno::Reference<table::XCell> xCell = getCell(aCellAddress);
 538                     if (xCell->getType() == table::CellContentType_VALUE)
 539                     {
 540                         maBounds[index].updateBound(eOp, xCell->getValue());
 541                     }
 542                     else
 543                     {
 544                         maNonBoundedConstraints.push_back(rConstraint);
 545                     }
 546                 }
 547                 else if (rConstraint.Right >>= fValue)
 548                 {
 549                     maBounds[index].updateBound(eOp, fValue);
 550                 }
 551             }
 552             index++;
 553         }
 554         if (!bFoundVariable)
 555             maNonBoundedConstraints.push_back(rConstraint);
 556     }
 557
 558     std::vector<double> aSolution;
 559
 560     if (mnAlgorithm == 0)
 561     {
 562         DifferentialEvolutionAlgorithm<SwarmSolver> aDE(*this, 50);
 563         SwarmRunner<DifferentialEvolutionAlgorithm<SwarmSolver>> aEvolution(aDE);
 564         aEvolution.setTimeout(mnTimeout);
 565         aSolution = aEvolution.solve();
 566     }
 567     else
 568     {
 569         ParticleSwarmOptimizationAlgorithm<SwarmSolver> aPSO(*this, 100);
 570         SwarmRunner<ParticleSwarmOptimizationAlgorithm<SwarmSolver>> aSwarmSolver(aPSO);
 571         aSwarmSolver.setTimeout(mnTimeout);
 572         aSolution = aSwarmSolver.solve();
 573     }
 574
 575     xModel->unlockControllers();
 576
 577     mbSuccess = true;
 578
 579     maSolution.realloc(aSolution.size());
 580     std::copy(aSolution.begin(), aSolution.end(), maSolution.getArray());
 581 }
 582
 583 extern "C" SAL_DLLPUBLIC_EXPORT uno::XInterface*
 584 com_sun_star_comp_Calc_SwarmSolver_get_implementation(uno::XComponentContext*,
 585                                                       uno::Sequence<uno::Any> const&)
 586 {
 587     return cppu::acquire(new SwarmSolver());
 588 }
 589
 590 /* vim:set shiftwidth=4 softtabstop=4 expandtab: */