sccomp/source/solver/ParticelSwarmOptimization.hxx

   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 #ifndef INCLUDED_SCCOMP_SOURCE_PARTICLESWARM_HXX
  12 #define INCLUDED_SCCOMP_SOURCE_PARTICLESWARM_HXX
  13
  14 #include <vector>
  15 #include <random>
  16 #include <limits>
  17
  18 struct Particle
  19 {
  20     Particle(size_t nDimensionality)
  21         : mVelocity(nDimensionality)
  22         , mPosition(nDimensionality)
  23         , mCurrentFitness(std::numeric_limits<float>::lowest())
  24         , mBestPosition(nDimensionality)
  25         , mBestFitness(std::numeric_limits<float>::lowest())
  26     {
  27     }
  28
  29     std::vector<double> mVelocity;
  30
  31     std::vector<double> mPosition;
  32     double mCurrentFitness;
  33
  34     std::vector<double> mBestPosition;
  35     double mBestFitness;
  36 };
  37
  38 template <typename DataProvider> class ParticleSwarmOptimizationAlgorithm
  39 {
  40 private:
  41     // inertia
  42     static constexpr double constWeight = 0.729;
  43     // cognitive coefficient
  44     static constexpr double c1 = 1.49445;
  45     // social  coefficient
  46     static constexpr double c2 = 1.49445;
  47
  48     static constexpr double constAcceptedPrecision = 0.000000001;
  49
  50     DataProvider& mrDataProvider;
  51
  52     size_t const mnNumOfParticles;
  53
  54     std::vector<Particle> maSwarm;
  55
  56     std::random_device maRandomDevice;
  57     std::mt19937 maGenerator;
  58     size_t mnDimensionality;
  59
  60     std::uniform_real_distribution<> maRandom01;
  61
  62     std::vector<double> maBestPosition;
  63     double mfBestFitness;
  64     int mnGeneration;
  65     int mnLastChange;
  66
  67 public:
  68     ParticleSwarmOptimizationAlgorithm(DataProvider& rDataProvider, size_t nNumOfParticles)
  69         : mrDataProvider(rDataProvider)
  70         , mnNumOfParticles(nNumOfParticles)
  71         , maGenerator(maRandomDevice())
  72         , mnDimensionality(mrDataProvider.getDimensionality())
  73         , maRandom01(0.0, 1.0)
  74         , maBestPosition(mnDimensionality)
  75         , mfBestFitness(std::numeric_limits<float>::lowest())
  76         , mnGeneration(0)
  77         , mnLastChange(0)
  78     {
  79     }
  80
  81     std::vector<double> const& getResult() { return maBestPosition; }
  82
  83     int getGeneration() { return mnGeneration; }
  84
  85     int getLastChange() { return mnLastChange; }
  86
  87     void initialize()
  88     {
  89         mnGeneration = 0;
  90         mnLastChange = 0;
  91         maSwarm.clear();
  92
  93         mfBestFitness = std::numeric_limits<float>::lowest();
  94
  95         maSwarm.reserve(mnNumOfParticles);
  96         for (size_t i = 0; i < mnNumOfParticles; i++)
  97         {
  98             maSwarm.emplace_back(mnDimensionality);
  99             Particle& rParticle = maSwarm.back();
 100
 101             mrDataProvider.initializeVariables(rParticle.mPosition, maGenerator);
 102             mrDataProvider.initializeVariables(rParticle.mVelocity, maGenerator);
 103
 104             for (size_t k = 0; k < mnDimensionality; k++)
 105             {
 106                 rParticle.mPosition[k] = mrDataProvider.clampVariable(k, rParticle.mPosition[k]);
 107             }
 108
 109             rParticle.mCurrentFitness = mrDataProvider.calculateFitness(rParticle.mPosition);
 110
 111             for (size_t k = 0; k < mnDimensionality; k++)
 112             {
 113                 rParticle.mPosition[k] = mrDataProvider.clampVariable(k, rParticle.mPosition[k]);
 114             }
 115
 116             std::copy(rParticle.mPosition.begin(), rParticle.mPosition.end(),
 117                       rParticle.mBestPosition.begin());
 118             rParticle.mBestFitness = rParticle.mCurrentFitness;
 119
 120             if (rParticle.mCurrentFitness > mfBestFitness)
 121             {
 122                 mfBestFitness = rParticle.mCurrentFitness;
 123                 std::copy(rParticle.mPosition.begin(), rParticle.mPosition.end(),
 124                           maBestPosition.begin());
 125             }
 126         }
 127     }
 128
 129     bool next()
 130     {
 131         bool bBestChanged = false;
 132
 133         for (Particle& rParticle : maSwarm)
 134         {
 135             double fRandom1 = maRandom01(maGenerator);
 136             double fRandom2 = maRandom01(maGenerator);
 137
 138             for (size_t k = 0; k < mnDimensionality; k++)
 139             {
 140                 rParticle.mVelocity[k]
 141                     = (constWeight * rParticle.mVelocity[k])
 142                       + (c1 * fRandom1 * (rParticle.mBestPosition[k] - rParticle.mPosition[k]))
 143                       + (c2 * fRandom2 * (maBestPosition[k] - rParticle.mPosition[k]));
 144
 145                 mrDataProvider.clampVariable(k, rParticle.mVelocity[k]);
 146
 147                 rParticle.mPosition[k] += rParticle.mVelocity[k];
 148                 rParticle.mPosition[k] = mrDataProvider.clampVariable(k, rParticle.mPosition[k]);
 149             }
 150
 151             rParticle.mCurrentFitness = mrDataProvider.calculateFitness(rParticle.mPosition);
 152
 153             if (rParticle.mCurrentFitness > rParticle.mBestFitness)
 154             {
 155                 rParticle.mBestFitness = rParticle.mCurrentFitness;
 156                 std::copy(rParticle.mPosition.begin(), rParticle.mPosition.end(),
 157                           rParticle.mBestPosition.begin());
 158             }
 159
 160             if (rParticle.mCurrentFitness > mfBestFitness)
 161             {
 162                 if (std::abs(rParticle.mCurrentFitness - mfBestFitness) > constAcceptedPrecision)
 163                 {
 164                     bBestChanged = true;
 165                     mnLastChange = mnGeneration;
 166                 }
 167                 std::copy(rParticle.mPosition.begin(), rParticle.mPosition.end(),
 168                           maBestPosition.begin());
 169                 mfBestFitness = rParticle.mCurrentFitness;
 170             }
 171         }
 172         mnGeneration++;
 173         return bBestChanged;
 174     }
 175 };
 176
 177 #endif
 178
 179 /* vim:set shiftwidth=4 softtabstop=4 expandtab: */