From 82a3df7d1d75eee339b28ad878b5a29b5ab8d338 Mon Sep 17 00:00:00 2001 From: Bertram Kopf <bertram@ep1.rub.de> Date: Wed, 18 Sep 2024 20:29:35 +0200 Subject: [PATCH] added option to minimize w/ Adam algorithm; added setenv-script for using MINUIT2 standalone version from ROOT --- AppUtils/AppBase.cc | 30 +++-- .../benchmarkCoupledChannel_runscriptAdam | 26 ++++ JamrootAL9_Min2Standalone | 78 +++++++++++ MinFunctions/AdamMinimizer.cc | 124 ++++++++++++++++++ MinFunctions/AdamMinimizer.hh | 74 +++++++++++ PwaApps/coupledChannelApp.cc | 2 +- PwaApps/singleChannelApp.cc | 2 +- SetEnv_rubMin2Standalone | 12 ++ .../pbarpTopi0pi0eta/benchmark_runscriptAdam | 19 +++ 9 files changed, 352 insertions(+), 15 deletions(-) create mode 100755 Examples/benchmark/coupledChannel/benchmarkCoupledChannel_runscriptAdam create mode 100644 JamrootAL9_Min2Standalone create mode 100644 MinFunctions/AdamMinimizer.cc create mode 100644 MinFunctions/AdamMinimizer.hh create mode 100755 SetEnv_rubMin2Standalone create mode 100755 benchmark/pbarpTopi0pi0eta/benchmark_runscriptAdam diff --git a/AppUtils/AppBase.cc b/AppUtils/AppBase.cc index f8d9ec41..1c500418 100644 --- a/AppUtils/AppBase.cc +++ b/AppUtils/AppBase.cc @@ -76,6 +76,7 @@ #include "MinFunctions/PwaFcnServerGradNumSlow.hh" #include "MinFunctions/AbsPawianMinimizer.hh" #include "MinFunctions/EvoMinimizer.hh" +#include "MinFunctions/AdamMinimizer.hh" #include "MinFunctions/MinuitMinimizer.hh" #include "FitParams/AbsPawianParameters.hh" @@ -882,31 +883,34 @@ void AppBase::fitServerMode(std::shared_ptr<AbsPawianParameters> upar){ InfoMsg << "Closing server." << endmsg; } - else if(GlobalEnv::instance()->parser()->mode()=="serverGradientNum" || GlobalEnv::instance()->parser()->mode()=="serverGradientNumSlow"){ + else if(GlobalEnv::instance()->parser()->mode()=="serverGradientNum" || GlobalEnv::instance()->parser()->mode()=="serverGradientNumSlow" || GlobalEnv::instance()->parser()->mode()=="serverAdamNum"){ std::shared_ptr<AbsFcn<FCNGradientBase>> absFcn; std::shared_ptr<NetworkServer> theServer(new NetworkServer(GlobalEnv::instance()->parser()->serverPort(), GlobalEnv::instance()->parser()->noOfClients(), numEventMap, GlobalEnv::instance()->parser()-> clientNumberWeights())); - if (GlobalEnv::instance()->parser()->mode()=="serverGradientNum") absFcn=std::shared_ptr<AbsFcn<FCNGradientBase>>(new PwaFcnServer<FCNGradientBase>(theServer)); + if (GlobalEnv::instance()->parser()->mode()=="serverGradientNum" || GlobalEnv::instance()->parser()->mode()=="serverAdamNum") absFcn=std::shared_ptr<AbsFcn<FCNGradientBase>>(new PwaFcnServer<FCNGradientBase>(theServer)); else absFcn=std::shared_ptr<AbsFcn<FCNGradientBase>>(new PwaFcnServerGradNumSlow(theServer)); theServer->WaitForFirstClientLogin(); - std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>> absMinimizerPtr; - if (GlobalEnv::instance()->useCovMatrix()){ - // std::shared_ptr<PwaCovMatrix> theCovMatrix=GlobalEnv::instance()->pwaCovMatrix(); - // const std::vector<double> dataFromMnCovMat=theCovMatrix->dataFromMnCovMatrix(); - // //int nrow=std::sqrt(2.*dataFromMnCovMat.size()+1./4.)-0.5; - // unsigned int nrow=theCovMatrix->nRow(); - // InfoMsg << "dataFromMnCovMat.size(): " << dataFromMnCovMat.size() << " nrow: " << nrow << endmsg; - // std::shared_ptr<MnUserCovariance> theMnUserCov(new MnUserCovariance(dataFromMnCovMat, nrow)); - absMinimizerPtr = std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>>(new MinuitMinimizer<FCNGradientBase>(absFcn, upar, theMnUserCov)); + if (GlobalEnv::instance()->parser()->mode()=="serverGradientNum" || GlobalEnv::instance()->parser()->mode()=="serverGradientNumSlow"){ + if (GlobalEnv::instance()->useCovMatrix()){ + // std::shared_ptr<PwaCovMatrix> theCovMatrix=GlobalEnv::instance()->pwaCovMatrix(); + // const std::vector<double> dataFromMnCovMat=theCovMatrix->dataFromMnCovMatrix(); + // //int nrow=std::sqrt(2.*dataFromMnCovMat.size()+1./4.)-0.5; + // unsigned int nrow=theCovMatrix->nRow(); + // InfoMsg << "dataFromMnCovMat.size(): " << dataFromMnCovMat.size() << " nrow: " << nrow << endmsg; + // std::shared_ptr<MnUserCovariance> theMnUserCov(new MnUserCovariance(dataFromMnCovMat, nrow)); + absMinimizerPtr = std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>>(new MinuitMinimizer<FCNGradientBase>(absFcn, upar, theMnUserCov)); + } + else absMinimizerPtr = std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>>(new MinuitMinimizer<FCNGradientBase>(absFcn, upar)); + } + else if (GlobalEnv::instance()->parser()->mode()=="serverAdamNum"){ + absMinimizerPtr = std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>>(new AdamMinimizer(absFcn, upar)); } - else absMinimizerPtr = std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>>(new MinuitMinimizer<FCNGradientBase>(absFcn, upar)); - // std::shared_ptr<AbsPawianMinimizer<FCNGradientBase>> absMinimizerPtr(new MinuitMinimizer<FCNGradientBase>(absFcn, upar)); absMinimizerPtr->minimize(); diff --git a/Examples/benchmark/coupledChannel/benchmarkCoupledChannel_runscriptAdam b/Examples/benchmark/coupledChannel/benchmarkCoupledChannel_runscriptAdam new file mode 100755 index 00000000..d6246832 --- /dev/null +++ b/Examples/benchmark/coupledChannel/benchmarkCoupledChannel_runscriptAdam @@ -0,0 +1,26 @@ +#!/bin/sh + +echo Starting server for coupled channel pwa ... +$TOP_DIR/bin/coupledChannelApp -c global.cfg --pbarpFiles pbarpToEtaEtaPiReaction.cfg --pbarpFiles pbarpToKKPiReaction.cfg --epemFiles psiToGamKKpi.cfg --resFiles D0ToKmpippi0.cfg --mode serverAdamNum !>& logAdamNum.server & + +sleep 15 + +echo "Starting first client for channel pbar p -> eta eta pi0 ..." +$TOP_DIR/bin/coupledChannelApp -c global.cfg --pbarpFiles pbarpToEtaEtaPiReaction.cfg --pbarpFiles pbarpToKKPiReaction.cfg --epemFiles psiToGamKKpi.cfg --resFiles D0ToKmpippi0.cfg --mode client !>& log.client1 & + +sleep 5 + +echo "Starting first client for channel pbar p -> K+ K- pi0 ..." +$TOP_DIR/bin/coupledChannelApp -c global.cfg --pbarpFiles pbarpToEtaEtaPiReaction.cfg --pbarpFiles pbarpToKKPiReaction.cfg --epemFiles psiToGamKKpi.cfg --resFiles D0ToKmpippi0.cfg --mode client !>& log.client2 & + +sleep 5 + +echo "Starting first client for channel e+ e- -> J/psi -> chi_c1 gamma -> (K+ K- pi0) gamma ..." +$TOP_DIR/bin/coupledChannelApp -c global.cfg --pbarpFiles pbarpToEtaEtaPiReaction.cfg --pbarpFiles pbarpToKKPiReaction.cfg --epemFiles psiToGamKKpi.cfg --resFiles D0ToKmpippi0.cfg --mode client !>& log.client3 & + +sleep 5 + +echo "Starting first client for D0 -> K- pi+ pi0 ..." +$TOP_DIR/bin/coupledChannelApp -c global.cfg --pbarpFiles pbarpToEtaEtaPiReaction.cfg --pbarpFiles pbarpToKKPiReaction.cfg --epemFiles psiToGamKKpi.cfg --resFiles D0ToKmpippi0.cfg --mode client !>& log.client4 & + +echo done diff --git a/JamrootAL9_Min2Standalone b/JamrootAL9_Min2Standalone new file mode 100644 index 00000000..332ce8cf --- /dev/null +++ b/JamrootAL9_Min2Standalone @@ -0,0 +1,78 @@ +import os ; +import testing ; + +path-constant TOP : . ; +local extern = [ os.environ extern ] ; + +local ROOTSYS = [ os.environ ROOTSYS ] ; +local rlibs = [ SHELL "$(ROOTSYS)/bin/root-config --libs" ] ; +ROOTLIBS = [ MATCH "(.*)[\n]" : $(rlibs) ] ; + +BOOSTLIBS = + -lboost_chrono + -lboost_date_time + -lboost_filesystem + -lboost_program_options + -lboost_regex + -lboost_serialization + -lboost_system + -lboost_timer + -lboost_unit_test_framework + -lrt + ; + +lib boost_test : : <name>boost_unit_test_framework ; + +project : + requirements <include>./ + <include>$(TOP) + <include>$(BOOSTINCLUDE) + <include>$(extern)/Minuit2-6.32.04-Standalone/include/ +# <include>$(extern)/include + <include>$(ROOTSYS)/include + <cxxflags>--std=c++17 + <cxxflags>-ftemplate-depth=256 + <cxxflags>-DBOOST_BIND_GLOBAL_PLACEHOLDERS + <link>shared + <linkflags>-L$(TOP)/lib + <linkflags>$(ROOTLIBS) + <linkflags>-L$(BOOSTLIBPATH) + <linkflags>$(BOOSTLIBS) + <linkflags>-lgsl + <linkflags>-lgslcblas + <cxxflags>-fPIC + <cxxflags>-Wall + ; + +actions rootlibs +{ + $(ROOTSYS)/bin/root-config --libs +} + +# lib Minuit2 : : <file>$(extern)/libgcc920/libMinuit2.so : : ; +lib Minuit2 : : <file>$(extern)/Minuit2-6.32.04-Standalone/lib/libMinuit2.so : : ; + +build-project ErrLogger ; +build-project qft++ ; +build-project qft++Extension ; +build-project ConfigParser ; +build-project Particle ; +build-project Utils ; +build-project Event ; +build-project FitParams ; +build-project PwaDynamics ; +build-project PwaUtils ; +build-project MinFunctions ; +build-project AppUtils ; +build-project pbarpUtils ; +build-project epemUtils ; +build-project resUtils ; +build-project ggUtils ; +# build-project gammapUtils ; +build-project pipiScatteringUtils ; +build-project gslUtils ; +build-project qaErrorExtract ; +build-project PspGen ; +build-project Examples ; +build-project KMatrixExtract ; +build-project PwaApps ; diff --git a/MinFunctions/AdamMinimizer.cc b/MinFunctions/AdamMinimizer.cc new file mode 100644 index 00000000..4f1b531e --- /dev/null +++ b/MinFunctions/AdamMinimizer.cc @@ -0,0 +1,124 @@ +//************************************************************************// +// // +// Copyright 2014 Maya Elspass (maya.elspass@ruhr-uni-bochum.de) // +// Bertram Kopf (bertram@ep1.rub.de) // +// - Ruhr-Universität Bochum // +// // +// This file is part of Pawian. // +// // +// Pawian is free software: you can redistribute it and/or modify // +// it under the terms of the GNU General Public License as published by // +// the Free Software Foundation, either version 3 of the License, or // +// (at your option) any later version. // +// // +// Pawian is distributed in the hope that it will be useful, // +// but WITHOUT ANY WARRANTY; without even the implied warranty of // +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the // +// GNU General Public License for more details. // +// // +// You should have received a copy of the GNU General Public License // +// along with Pawian. If not, see <http://www.gnu.org/licenses/>. // +// // +//************************************************************************// + +#include <fstream> +#include <sstream> +#include <string> +#include <math.h> +#include <iomanip> + + +#include "MinFunctions/AdamMinimizer.hh" +#include "PwaUtils/GlobalEnv.hh" +#include "ErrLogger/ErrLogger.hh" +#include "PwaUtils/GlobalEnv.hh" +#include "ConfigParser/ParserBase.hh" + +AdamMinimizer::AdamMinimizer(std::shared_ptr<AbsFcn<FCNGradientBase>> theAbsFcnPtr, std::shared_ptr<AbsPawianParameters> upar) : + AbsPawianMinimizer<FCNGradientBase>(theAbsFcnPtr, upar) + ,_max_iterations(1000) + ,_s(upar->Params().size(),0.) + ,_v(upar->Params().size(),0.) + ,_learning_rate(0.2) + ,_iterations(0) + ,_currentPawianParams(std::shared_ptr<AbsPawianParameters>(upar->Clone())) + ,_bestLH(1.e20) + ,_bestPawianParams(std::shared_ptr<AbsPawianParameters>(upar->Clone())) + ,_noItersWoImprovement(0) +{ +} + +AdamMinimizer::~AdamMinimizer() +{ +} + +void AdamMinimizer::minimize(){ + if(_iterations==0){ + double currentLH=_absFcn->operator()(_bestPawianParams->Params()); + if(currentLH<_bestLH) _bestLH=currentLH; + InfoMsg << "best LH: " << _bestLH << "\tnoIters: " << _iterations << "\tnoItersWoImprovement: " << _noItersWoImprovement << endmsg; + } + while(_iterations <= _max_iterations && _noItersWoImprovement<100){ + if(_noItersWoImprovement>0 && _noItersWoImprovement%10 == 0){ //continue with parameters of the best fit + _currentPawianParams=std::shared_ptr<AbsPawianParameters>(_bestPawianParams->Clone()); + } + std::vector<double> derivatives = _absFcn->Gradient(_currentPawianParams->Params()); + + updateParameters(_currentPawianParams, derivatives, _s, _v, _iterations); + double currentLH=_absFcn->operator()(_currentPawianParams->Params()); + if(currentLH<_bestLH){ + _bestLH=currentLH; + _bestPawianParams=std::shared_ptr<AbsPawianParameters>(_currentPawianParams->Clone()); + _noItersWoImprovement=0; + } + else ++_noItersWoImprovement; + ++_iterations; + _learning_rate *= 0.999; + InfoMsg << "best LH: " << _bestLH << "\tnoIters: " << _iterations << "\tnoItersWoImprovement: " << _noItersWoImprovement << endmsg; + } +} + +void AdamMinimizer::updateParameters(std::shared_ptr<AbsPawianParameters> pawianParams, std::vector<double>& gradients, std::vector<double>& s, std::vector<double>& v, int t){ + double beta1=0.9; + double beta2=0.99; + double epsilon=1.*pow(10., -8); + + for(unsigned int i = 0; i < pawianParams->Params().size(); ++i){ + if (pawianParams->IsFixed(i)) continue; + _s.at(i) = beta1 * _s.at(i) + (1.0 - beta1) * gradients.at(i); + _v.at(i) = beta2 * _v.at(i) + (1.0 - beta2) * gradients.at(i) * gradients.at(i); + + double s_hat = _s.at(i) / (1.0 - pow(beta1 , (t + 1))); + double v_hat = _v.at(i) / (1.0 - pow(beta2 , (t + 1))); + + double newVal = pawianParams->Value(i) - _learning_rate * s_hat / (std::sqrt(v_hat) + epsilon); + if(pawianParams->HasLimits(i)){ + if(newVal>pawianParams->UpperLimit(i)) newVal=pawianParams->UpperLimit(i); + else if(newVal<pawianParams->LowerLimit(i)) newVal=pawianParams->LowerLimit(i); + } + pawianParams->SetValue(i,newVal); + } +} + + +void AdamMinimizer::printFitResult(double evtWeightSumData){ + InfoMsg << "\n\n********************** Final fit parameters *************************" << endmsg; + _bestPawianParams->print(std::cout, true); + InfoMsg << "\nFinal LH: " << _bestLH; + + /////////////////////////////////////////////////////////////////////////////////////////////// + // calculate AIC, BIC criteria and output selected wave contrib + /////////////////////////////////////////////////////////////////////////////////////////////// + unsigned int noOfFreeFitParams=_bestPawianParams->VariableParameters(); + + double BICcriterion=2.*_bestLH+noOfFreeFitParams*log(evtWeightSumData); + double AICcriterion=2.*_bestLH+2.*noOfFreeFitParams; + double AICccriterion=AICcriterion+2.*noOfFreeFitParams*(noOfFreeFitParams+1) + / (evtWeightSumData-noOfFreeFitParams-1); + InfoMsg << "\nnoOfFreeFitParams:\t" <<noOfFreeFitParams << endmsg; + InfoMsg << "evtWeightSumData:\t" <<evtWeightSumData << endmsg; + InfoMsg << "BIC:\t" << BICcriterion << endmsg; + InfoMsg << "AIC:\t" << AICcriterion << endmsg; + InfoMsg << "AICc:\t" << AICccriterion << endmsg; +} + diff --git a/MinFunctions/AdamMinimizer.hh b/MinFunctions/AdamMinimizer.hh new file mode 100644 index 00000000..ea72023f --- /dev/null +++ b/MinFunctions/AdamMinimizer.hh @@ -0,0 +1,74 @@ +//************************************************************************// +// // +// Copyright 2014 Maya Elspass (maya.elspass@ruhr-uni-bochum.de) // +// Bertram Kopf (bertram@ep1.rub.de) // +// - Ruhr-Universität Bochum // +// // +// This file is part of Pawian. // +// // +// Pawian is free software: you can redistribute it and/or modify // +// it under the terms of the GNU General Public License as published by // +// the Free Software Foundation, either version 3 of the License, or // +// (at your option) any later version. // +// // +// Pawian is distributed in the hope that it will be useful, // +// but WITHOUT ANY WARRANTY; without even the implied warranty of // +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the // +// GNU General Public License for more details. // +// // +// You should have received a copy of the GNU General Public License // +// along with Pawian. If not, see <http://www.gnu.org/licenses/>. // +// // +//************************************************************************// + + +#pragma once + +#include <iostream> +#include <fstream> +#include <string> +#include <vector> +#include <memory> + +#include "MinFunctions/AbsFcn.hh" +#include "MinFunctions/AbsPawianMinimizer.hh" +#include "Minuit2/FCNGradientBase.h" + +using namespace ROOT::Minuit2; + +class AdamMinimizer : public AbsPawianMinimizer<FCNGradientBase>{ + +public: + + // create/copy/destroy: + + ///Constructor + AdamMinimizer(std::shared_ptr<AbsFcn<FCNGradientBase>> theAbsFcnPtr, std::shared_ptr<AbsPawianParameters> upar); + + + ///Destructor + virtual ~AdamMinimizer(); + + //Funktionen + virtual std::string type() {return "AdamMinimizer";}; + virtual void minimize(); + virtual void printFitResult(double evtWeightSumData); + + // Getters: + +protected: + const unsigned int _max_iterations; + std::vector<double> _s; + std::vector<double> _v; + double _learning_rate; + unsigned int _iterations; + std::shared_ptr<AbsPawianParameters> _currentPawianParams; + + double _bestLH; + std::shared_ptr<AbsPawianParameters> _bestPawianParams; + unsigned int _noItersWoImprovement; + + void updateParameters(std::shared_ptr<AbsPawianParameters> pawianParams, std::vector<double>& gradients, std::vector<double>& s, std::vector<double>& v, int t); + + private: +}; diff --git a/PwaApps/coupledChannelApp.cc b/PwaApps/coupledChannelApp.cc index feff405c..c7ae292f 100644 --- a/PwaApps/coupledChannelApp.cc +++ b/PwaApps/coupledChannelApp.cc @@ -160,7 +160,7 @@ int main(int __argc,char *__argv[]){ return 1; } - if(mode == "server" || mode == "serverGradientNum" || mode == "serverGradientNumSlow" || mode == "evoserver" || mode=="serverQA"){ + if(mode == "server" || mode == "serverGradientNum" || mode == "serverGradientNumSlow" || mode == "serverAdamNum" || mode == "evoserver" || mode=="serverQA"){ theAppBase.fitServerMode(startPawianParams); return 1; } diff --git a/PwaApps/singleChannelApp.cc b/PwaApps/singleChannelApp.cc index a945cdb7..074b4567 100644 --- a/PwaApps/singleChannelApp.cc +++ b/PwaApps/singleChannelApp.cc @@ -201,7 +201,7 @@ int main(int __argc,char *__argv[]){ // Disable output buffering setvbuf(stdout, NULL, _IONBF, 0); - if(mode == "server" || mode == "serverGradientNum" || mode == "serverGradientNumSlow" || mode == "evoserver"){ + if(mode == "server" || mode == "serverGradientNum" || mode == "serverGradientNumSlow" || mode == "serverAdamNum" || mode == "evoserver"){ theAppBase.fitServerMode(startPawianParams); return 1; } diff --git a/SetEnv_rubMin2Standalone b/SetEnv_rubMin2Standalone new file mode 100755 index 00000000..0ee54295 --- /dev/null +++ b/SetEnv_rubMin2Standalone @@ -0,0 +1,12 @@ +# source /opt/root/6-28.04-AL9.2-gcc12.2.0/bin/thisroot.csh +source /opt/root/6-32.04-AL9.4-gcc12.2.0/bin/thisroot.csh +setenv TOP_DIR `pwd | sed -e 's/\/nfs//'` +setenv extern /data/sleipnir/PANDA/PWA +setenv LD_LIBRARY_PATH ${TOP_DIR}/lib:${extern}/libgcc920-opt +setenv BOOST_BUILD_PATH /usr/share/boost-build +setenv PATH ${TOP_DIR}/bin:${PATH} +setenv KMAT_DIR /data/duldul/bertram/KMatStore/ +setenv EVT_DIR /data/duldul/bertram/EvtStore/ + +rm Jamroot +ln -s JamrootAL9_Min2Standalone Jamroot diff --git a/benchmark/pbarpTopi0pi0eta/benchmark_runscriptAdam b/benchmark/pbarpTopi0pi0eta/benchmark_runscriptAdam new file mode 100755 index 00000000..e94acdcb --- /dev/null +++ b/benchmark/pbarpTopi0pi0eta/benchmark_runscriptAdam @@ -0,0 +1,19 @@ +#!/bin/sh + +echo Starting server ... + +$TOP_DIR/bin/singleChannelApp --pbarpFiles pbarpReactionDefault.cfg --mode serverAdamNum --paramFile $TOP_DIR/benchmark/pbarpTopi0pi0eta/defaultparams.dat !>& serverAdamNum.log & + +sleep 5 + +echo Starting clients ... + +$TOP_DIR/bin/singleChannelApp --pbarpFiles pbarpReactionDefault.cfg --mode client --paramFile $TOP_DIR/benchmark/pbarpTopi0pi0eta/defaultparams.dat !>& client1.log & + +$TOP_DIR/bin/singleChannelApp --pbarpFiles pbarpReactionDefault.cfg --mode client --paramFile $TOP_DIR/benchmark/pbarpTopi0pi0eta/defaultparams.dat !>& client2.log & + +$TOP_DIR/bin/singleChannelApp --pbarpFiles pbarpReactionDefault.cfg --mode client --paramFile $TOP_DIR/benchmark/pbarpTopi0pi0eta/defaultparams.dat !>& client3.log & + +$TOP_DIR/bin/singleChannelApp --pbarpFiles pbarpReactionDefault.cfg --mode client --paramFile $TOP_DIR/benchmark/pbarpTopi0pi0eta/defaultparams.dat !>& client4.log & + +echo done -- GitLab