Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-2023 The plumed team
3 : (see the PEOPLE file at the root of the distribution for a list of names)
4 :
5 : See http://www.plumed.org for more information.
6 :
7 : This file is part of plumed, version 2.
8 :
9 : plumed is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : plumed is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 : #ifndef __PLUMED_gridtools_GridSearch_h
23 : #define __PLUMED_gridtools_GridSearch_h
24 :
25 : #include "tools/MinimiseBase.h"
26 : #include "GridVessel.h"
27 : #include <iostream>
28 : #include <memory>
29 :
30 : namespace PLMD {
31 : namespace gridtools {
32 :
33 : template <class FCLASS>
34 1 : class GridSearch {
35 : private:
36 : /// This is the pointer to the member function in the energy
37 : /// calculating class that calculates the energy
38 : typedef double(FCLASS::*engf_pointer)( const std::vector<double>& p, std::vector<double>& der );
39 : FCLASS* myclass_func;
40 : std::unique_ptr<GridVessel> mygrid;
41 : std::unique_ptr<GridVessel> myfgrid;
42 : public:
43 1 : GridSearch( const std::vector<double>& mmin, const std::vector<double>& mmax, const std::vector<unsigned>& ng, const std::vector<unsigned>& nfg, FCLASS* funcc ) :
44 1 : myclass_func( funcc )
45 : {
46 : // Create the grid objects
47 1 : std::string nstr, vstring="COMPONENTS=func COORDINATES=x1";
48 2 : for(unsigned i=1; i<mmin.size(); ++i) { Tools::convert(i+1,nstr); vstring += ",x" + nstr; }
49 2 : vstring += " PBC=F"; for(unsigned i=1; i<mmin.size(); ++i) vstring += ",F";
50 2 : vesselbase::VesselOptions da("mygrid","",-1,vstring,NULL);
51 1 : Keywords keys; gridtools::GridVessel::registerKeywords( keys );
52 1 : vesselbase::VesselOptions dar( da, keys );
53 2 : mygrid=Tools::make_unique<GridVessel>(dar);
54 2 : if( nfg[0]>0 ) myfgrid=Tools::make_unique<GridVessel>(dar);
55 :
56 : // Now setup the min and max values for the grid
57 1 : std::vector<std::string> gmin( nfg.size() ), gmax( nfg.size() ); std::vector<double> dummy_spacing;
58 3 : for(unsigned i=0; i<nfg.size(); ++i) { Tools::convert(mmin[i],gmin[i]); Tools::convert(mmax[i],gmax[i]); }
59 1 : mygrid->setBounds( gmin, gmax, ng, dummy_spacing ); mygrid->resize();
60 1 : if( myfgrid ) myfgrid->setBounds( gmin, gmax, nfg, dummy_spacing );
61 2 : }
62 : bool minimise( std::vector<double>& p, engf_pointer myfunc );
63 : };
64 :
65 : template <class FCLASS>
66 200 : bool GridSearch<FCLASS>::minimise( std::vector<double>& p, engf_pointer myfunc ) {
67 200 : std::vector<double> der( p.size() ); std::vector<double> coords( p.size() );
68 200 : double initial_eng = (myclass_func->*myfunc)( p, der );
69 200 : mygrid->getGridPointCoordinates( 0, coords );
70 200 : double emin=(myclass_func->*myfunc)( coords, der );
71 200 : mygrid->setValueAndDerivatives( 0, 0, emin, der ); unsigned pmin=0;
72 24200 : for(unsigned i=1; i<mygrid->getNumberOfPoints(); ++i) {
73 24000 : mygrid->getGridPointCoordinates( i, coords );
74 24000 : double eng = (myclass_func->*myfunc)( coords, der );
75 24000 : mygrid->setValueAndDerivatives( i, 0, eng, der );
76 24000 : if( eng<emin ) { emin=eng; pmin=i; }
77 : }
78 : // This prevents division by zero
79 : mygrid->setNorm( 1.0 );
80 :
81 200 : if( myfgrid ) {
82 200 : myfgrid->getGridPointCoordinates( 0, coords ); pmin=0;
83 200 : double emin=mygrid->getValueAndDerivatives( coords, 0, der );
84 520200 : for(unsigned i=1; i<myfgrid->getNumberOfPoints(); ++i) {
85 520000 : myfgrid->getGridPointCoordinates( i, coords );
86 520000 : double eng = mygrid->getValueAndDerivatives( coords, 0, der );
87 520000 : if( eng<emin ) { emin=eng; pmin=i; }
88 : }
89 200 : myfgrid->getGridPointCoordinates( pmin, coords );
90 200 : double checkEng = (myclass_func->*myfunc)( coords, der );
91 200 : if( checkEng<initial_eng ) {
92 22 : myfgrid->getGridPointCoordinates( pmin, p );
93 : return true;
94 : } else {
95 : return false;
96 : }
97 : }
98 :
99 0 : if( emin<initial_eng ) {
100 0 : mygrid->getGridPointCoordinates( pmin, p );
101 : return true;
102 : } else {
103 : return false;
104 : }
105 : }
106 :
107 : }
108 : }
109 : #endif
110 :
|
