BBA.cxx

/* Copyright (C) 2021-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
   SPDX-License-Identifier: GPL-3.0-only
   Authors: Sergey Gorbunov [committer], Martin Beyer */

/// @file BBA.cxx
/// @author Sergey Gorbunov
/// @date 2021-08-06

#include "bba/BBA.h"

#include <iomanip>
#include <iostream>
#include <random>

#include "bba/Parameter.h"

namespace bba
{
  void BBA::printInfo()
  {
    std::cout << "\n\n ------------------------------ \n"
              << "     bba package version 1.0 \n"
              << " ------------------------------\n\n"
              << std::endl;
  }

  void BBA::align()
  {
    double Fbest           = fF(fP);
    double Foriginal       = Fbest;
    bool isAnythingChanged = 0;
    do {
      isAnythingChanged = 0;
      for (unsigned int ip = 0; ip < fPar.size(); ip++) {
        Parameter& par = fPar[ip];
        double& p      = fP[ip];
        if (!par.IsActive()) { continue; }
        double step             = par.GetStep();
        bool isParameterChanged = 0;
        double FbestOld         = Fbest;
        if (p + step <= par.GetBoundaryMax()) {
          double pStore = p;
          p += step;
          double Fnew = fF(fP);
          if (Fnew < Fbest) {
            Fbest              = Fnew;
            isParameterChanged = 1;
          }
          else {
            p = pStore;
          }
        }
        if (!isParameterChanged && (p - step >= par.GetBoundaryMin())) {
          double pStore = p;
          p -= step;
          double Fnew = fF(fP);
          if (Fnew < Fbest) {
            Fbest              = Fnew;
            isParameterChanged = 1;
          }
          else {
            p = pStore;
          }
        }
        if (isParameterChanged) {
          std::cout << "BBA: cost " << std::setprecision(15) << std::setw(20) << Fbest
                    << ";     diff to the initial cost " << std::setprecision(15) << std::setw(20) << Fbest - Foriginal
                    << ";     dF= " << std::setw(30) << std::setprecision(25) << std::fixed << Fbest - FbestOld
                    << std::endl;
          /*
          std::cout << " parameters: ";
          for (unsigned int i = 0; i < fP.size(); i++) {
            std::cout << fP[i] << " ";
          }
          std::cout << std::endl;
          */

          step = 2 * par.GetStep();
          if (step > par.GetStepMax()) { step = par.GetStepMax(); }
          par.SetStep(step);
          isAnythingChanged = 1;
        }
        else {
          step = 0.5 * par.GetStep();
          if (step < par.GetStepMin()) { step = par.GetStepMin(); }
          if (step < par.GetStep()) { isAnythingChanged = 1; }
          par.SetStep(step);
        }
      }  // ip

      std::cout << " BBA: next round over parameters.. " << std::endl;

    } while (isAnythingChanged);
  }

  void BBA::randomSearch(int n)
  {
    std::random_device rd;
    std::mt19937 gen(rd());
    double Fbest = fF(fP);
    std::vector<double> Pbest(fP.size());
    int i = 0;
    do {
      for (unsigned int ip = 0; ip < fPar.size(); ip++) {
        Parameter& par = fPar[ip];
        double& p      = fP[ip];
        if (!par.IsActive()) { continue; }
        std::uniform_real_distribution<double> dist(par.GetBoundaryMin(), par.GetBoundaryMax());
        p = dist(gen);
      }
      double F = fF(fP);
      if (F < Fbest) {
        Fbest = F;
        for (unsigned int ip = 0; ip < fP.size(); ip++) {
          Pbest[ip] = fP[ip];
        }
        std::cout << "============ New best ============" << std::endl;
      }
      std::cout << " F= " << F << " parameters: ";
      for (unsigned int ip = 0; ip < fP.size(); ip++) {
        std::cout << fP[ip] << " ";
      }
      std::cout << std::endl;
      i++;
    } while (i < n);

    std::cout << "Best F= " << Fbest << " parameters: ";
    for (unsigned int ip = 0; ip < Pbest.size(); ip++) {
      std::cout << Pbest[ip] << " ";
    }
    std::cout << std::endl;
  }

}  // namespace bba