Commit: 79a98bb1a167945549973f6a2124708cac6e45ec - MolDS (git)

Commit MetaInfo

• Revision: 79a98bb1a167945549973f6a2124708cac6e45ec (tree)
• Time: 2014-01-11 21:21:41
• Author: Katsuhiko Nishimra <ktns.87@gmai...>
• Commiter: Katsuhiko Nishimra

Log Message

Factor Out Optimizer::PrepareState. #32881

git-svn-id: https://svn.sourceforge.jp/svnroot/molds/branches/refactor_opt@1643 1136aad2-a195-0410-b898-f5ea1d11b9d8

Change Summary

• modified: src/optimization/BFGS.cpp
• modified: src/optimization/BFGS.h
• modified: src/optimization/ConjugateGradient.cpp
• modified: src/optimization/ConjugateGradient.h
• modified: src/optimization/Optimizer.h
• modified: src/optimization/SteepestDescent.h

Incremental Difference

--- a/src/optimization/BFGS.cpp
+++ b/src/optimization/BFGS.cpp
@@ -135,9 +135,9 @@ void BFGS::SearchMinimum(boost::shared_ptr<ElectronicStructure> electronicStruct
    double maxGradientThreshold = Parameters::GetInstance()->GetMaxGradientOptimization();
    double rmsGradientThreshold = Parameters::GetInstance()->GetRmsGradientOptimization();
    const int dimension = molecule.GetAtomVect().size()*CartesianType_end;
-   double       trustRadius      = Parameters::GetInstance()->GetInitialTrustRadiusOptimization();
-   const double maxNormStep      = Parameters::GetInstance()->GetMaxNormStepOptimization();
    BFGSState state(molecule);
+   state.SetTrustRadius(Parameters::GetInstance()->GetInitialTrustRadiusOptimization());
+   state.SetMaxNormStep(Parameters::GetInstance()->GetMaxNormStepOptimization());
 
    // initial calculation
    bool requireGuess = true;
@@ -152,23 +152,7 @@ void BFGS::SearchMinimum(boost::shared_ptr<ElectronicStructure> electronicStruct
    for(int s=0; s<totalSteps; s++){
       this->OutputOptimizationStepMessage(s);
 
-      // Store old Force data
-      MolDS_wrappers::Blas::GetInstance()->Dcopy(dimension,
-                                                 static_cast<double const*>(state.GetVectorForce()),
-                                                 state.GetVectorOldForce());
-
-      this->StoreMolecularGeometry(state.GetMatrixOldCoordinatesRef(), molecule);
-
-      // Level shift Hessian redundant modes
-      this->ShiftHessianRedundantMode(state.GetMatrixHessian(), molecule);
-
-      // Limit the trustRadius to maxNormStep
-      trustRadius=min(trustRadius,maxNormStep);
-
-      //Calculate RFO step
-      this->CalcRFOStep(state.GetVectorStep(), state.GetMatrixHessian(), state.GetVectorForce(), trustRadius, dimension);
-
-      double approximateChange = this->ApproximateEnergyChange(dimension, state.GetMatrixHessian(), state.GetVectorForce(), state.GetVectorStep());
+      this->PrepareState(state, molecule, electronicStructure, elecState);
 
       // Take a RFO step
       bool doLineSearch = false;
@@ -189,7 +173,7 @@ void BFGS::SearchMinimum(boost::shared_ptr<ElectronicStructure> electronicStruct
       }
       this->OutputMoleculeElectronicStructure(electronicStructure, molecule, this->CanOutputLogs());
 
-      this->UpdateTrustRadius(trustRadius, approximateChange, state.GetInitialEnergy(), state.GetCurrentEnergy());
+      this->UpdateTrustRadius(state.GetTrustRadiusRef(), state.GetApproximateChange(), state.GetInitialEnergy(), state.GetCurrentEnergy());
 
       // check convergence
       if(this->SatisfiesConvergenceCriterion(state.GetMatrixForce(),
@@ -234,6 +218,35 @@ void BFGS::InitializeState(OptimizerState &stateOrig, const Molecule& molecule)
    MolDS_wrappers::Blas::GetInstance()->Dcopy(dimension, &one, 0, &state.GetMatrixHessian()[0][0], dimension+1);
 }
 
+void BFGS::PrepareState(OptimizerState& stateOrig,
+                        const MolDS_base::Molecule& molecule,
+                        const boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
+                        const int elecState) const{
+   BFGSState& state = stateOrig.CastRef<BFGSState>();
+   const MolDS_wrappers::molds_blas_int dimension = molecule.GetAtomVect().size()*CartesianType_end;
+   // Store old Force data
+   MolDS_wrappers::Blas::GetInstance()->Dcopy(dimension,
+         static_cast<double const*>(state.GetVectorForce()),
+         state.GetVectorOldForce());
+
+   this->StoreMolecularGeometry(state.GetMatrixOldCoordinatesRef(), molecule);
+
+   // Level shift Hessian redundant modes
+   this->ShiftHessianRedundantMode(state.GetMatrixHessian(), molecule);
+
+   // Limit the trustRadius to maxNormStep
+   state.SetTrustRadius(min(state.GetTrustRadius(),state.GetMaxNormStep()));
+
+   //Calculate RFO step
+   this->CalcRFOStep(state.GetVectorStep(), state.GetMatrixHessian(), state.GetVectorForce(), state.GetTrustRadius(), dimension);
+
+   state.SetApproximateChange(this->ApproximateEnergyChange(dimension,
+                                                            state.GetMatrixHessian(),
+                                                            state.GetVectorForce(),
+                                                            state.GetVectorStep())
+                             );
+}
+
 void BFGS::CalcRFOStep(double* vectorStep,
                        double const* const* matrixHessian,
                        double const* vectorForce,

--- a/src/optimization/BFGS.h
+++ b/src/optimization/BFGS.h
@@ -31,6 +31,9 @@ private:
       double** matrixOldCoordinates;
       double*  vectorOldCoordinates;
       double** matrixDisplacement;
+      double   approximateChange;
+      double   trustRadius;
+      double   maxNormStep;
       size_t numAtoms;
    private:
       template<class vector>
@@ -48,6 +51,13 @@ private:
       double**  GetMatrixOldCoordinates()   {return this->matrixOldCoordinates;}
       double**& GetMatrixOldCoordinatesRef(){return this->matrixOldCoordinates;}
       double**  GetMatrixDisplacement()     {return this->matrixDisplacement;}
+      double    GetApproximateChange()      {return this->approximateChange;}
+      double    GetTrustRadius()            {return this->trustRadius;}
+      double&   GetTrustRadiusRef()         {return this->trustRadius;}
+      double    GetMaxNormStep()            {return this->maxNormStep;}
+      void      SetApproximateChange(double approximateChange){this->approximateChange = approximateChange;}
+      void      SetTrustRadius(double trustRadius)            {this->trustRadius       = trustRadius;}
+      void      SetMaxNormStep(double maxNormStep)            {this->maxNormStep       = maxNormStep;}
    };
 public:
    BFGS();
@@ -82,6 +92,10 @@ private:
 
 protected:
    void InitializeState(OptimizerState &state, const MolDS_base::Molecule& molecule) const;
+   void PrepareState(OptimizerState& state,
+                             const MolDS_base::Molecule& molecule,
+                             const boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
+                             const int elecState) const;
    void CalcRFOStep(double* vectorStep,
                     double const* const* matrixHessian,
                     double const* vectorForce,

--- a/src/optimization/ConjugateGradient.cpp
+++ b/src/optimization/ConjugateGradient.cpp
@@ -101,6 +101,9 @@ void ConjugateGradient::SearchMinimum(boost::shared_ptr<ElectronicStructure> ele
    // conugate gradient loop
    for(int s=0; s<totalSteps; s++){
       this->OutputOptimizationStepMessage(s);
+
+      this->PrepareState(state, molecule, electronicStructure, elecState);
+
       state.SetInitialEnergy(state.GetCurrentEnergy());
 
       // do line search

--- a/src/optimization/ConjugateGradient.h
+++ b/src/optimization/ConjugateGradient.h
@@ -49,6 +49,10 @@ private:
                       double* lineSearchedEnergy,
                       bool* obainesOptimizedStructure) const;
    void InitializeState(OptimizerState &state, const MolDS_base::Molecule& molecule) const;
+   virtual void PrepareState(OptimizerState& state,
+                             const MolDS_base::Molecule& molecule,
+                             const boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
+                             const int elecState) const{};
    void UpdateSearchDirection(OptimizerState& state,
                               boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
                               const MolDS_base::Molecule& molecule,

--- a/src/optimization/Optimizer.h
+++ b/src/optimization/Optimizer.h
@@ -104,6 +104,10 @@ private:
                               double* lineSearchedEnergy,
                               bool* obainesOptimizedStructure) const = 0;
    virtual void InitializeState(OptimizerState &state, const MolDS_base::Molecule& molecule) const = 0;
+   virtual void PrepareState(OptimizerState& state,
+                             const MolDS_base::Molecule& molecule,
+                             const boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
+                             const int elecState) const = 0;
    virtual const std::string& OptimizationStepMessage() const = 0;
 };
 

--- a/src/optimization/SteepestDescent.h
+++ b/src/optimization/SteepestDescent.h
@@ -37,6 +37,10 @@ private:
                       double* lineSearchedEnergy,
                       bool* obainesOptimizedStructure) const;
    void InitializeState(OptimizerState&, const MolDS_base::Molecule&) const{}
+   virtual void PrepareState(OptimizerState& state,
+                             const MolDS_base::Molecule& molecule,
+                             const boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
+                             const int elecState) const{};
    void UpdateSearchDirection(OptimizerState& state,
                               boost::shared_ptr<MolDS_base::ElectronicStructure> electronicStructure,
                               const MolDS_base::Molecule& molecule,