Added new dftParameters for choosing ion relaxation solver and other …

…cleanups
ferroteam · Sep 11, 2019 · 11a382f · 11a382f
1 parent 367dddb
commit 11a382f
Show file tree

Hide file tree

Showing 6 changed files with 606 additions and 580 deletions.
diff --git a/include/cgPRPNonLinearSolver.h b/include/cgPRPNonLinearSolver.h
@@ -211,7 +211,7 @@ namespace dftfe {
 
  /// damping parameter (0,1] to be multiplied with the steepest descent direction,
  /// which controls the initial guess to the line search iteration.
- const double d_lineSearchDampingParameter;
+ double d_lineSearchDampingParameter;
 
  //parallel objects
  MPI_Comm mpi_communicator;

diff --git a/include/dftParameters.h b/include/dftParameters.h
@@ -47,7 +47,7 @@ namespace dftfe {
  extern bool autoUserMeshParams;
  extern double meshSizeInnerBall, meshSizeOuterBall;
  extern double chebyshevTolerance, topfrac;
- extern std::string mixingMethod ;
+ extern std::string mixingMethod,ionOptSolver;
 
 
  extern bool isIonOpt, isCellOpt, isIonForce, isCellStress;

diff --git a/src/dft/moveAtoms.cc b/src/dft/moveAtoms.cc
@@ -373,7 +373,7 @@ void dftClass<FEOrder>::updateAtomPositionsAndMoveMesh(const std::vector<Tensor<
  const std::pair<bool,double> meshQualityMetrics=d_gaussianMovePar.moveMeshTwoStep(controlPointLocationsInitialMove,d_controlPointLocationsCurrentMove,controlPointDisplacementsInitialMove,controlPointDisplacementsCurrentMove,d_gaussianConstantAutoMove,forcePtr->getGaussianGeneratorParameter());
  double factor;
  if(maximumForceToBeRelaxed >= 1e-03)
- factor = 1.35;
+ factor = 2.00;
  else if(maximumForceToBeRelaxed < 1e-03 && maximumForceToBeRelaxed >= 1e-04)
  factor = 1.25;
  else if(maximumForceToBeRelaxed < 1e-04)

diff --git a/src/geoOpt/geoOptIon.cc b/src/geoOpt/geoOptIon.cc
@@ -133,25 +133,43 @@ namespace dftfe {
  if(getNumberUnknowns()>0)
  {
  nonLinearSolver::ReturnValueType cgReturn=nonLinearSolver::FAILURE;
+ bool cgSuccess;
 
- /*if (dftParameters::chkType>=1 && dftParameters::restartFromChk)
+ if (dftParameters::chkType>=1 && dftParameters::restartFromChk && dftParameters::ionOptSolver == "CGPRP")
  cgReturn=cgSolver.solve(*this,std::string("ionRelaxCG.chk"),true);
- else if (dftParameters::chkType>=1 && !dftParameters::restartFromChk)
+ else if (dftParameters::chkType>=1 && !dftParameters::restartFromChk && dftParameters::ionOptSolver == "CGPRP")
  cgReturn=cgSolver.solve(*this,std::string("ionRelaxCG.chk"));
+ else if (dftParameters::ionOptSolver == "CGPRP")
+ cgReturn=cgSolver.solve(*this);
+ else if(dftParameters::ionOptSolver == "LBFGS")
+ {
+ cg_descent.set_step(0.8);
+ cg_descent.set_lbfgs(true);
+ cgSuccess = cg_descent.run(*this);
+ if(this_mpi_process == 0)
+ cg_descent.set_PrintLevel(2);
+ }
  else
- cgReturn=cgSolver.solve(*this);*/
- //if(getNumberUnknowns() < 12)
- // cg_descent.set_memory(0);
+ {
+ cg_descent.set_step(0.8);
+ cgSuccess = cg_descent.run(*this);
+ if(this_mpi_process == 0)
+ cg_descent.set_PrintLevel(2);
+ }
 
- cg_descent.set_step(0.8);
 
- if(this_mpi_process == 0)
- cg_descent.set_PrintLevel(2);
 
- cg_descent.set_AWolfe(true);
- bool cgDescentSucess = cg_descent.run(*this);
+ //cg_descent.set_memory(57);
+ // cg_descent.set_lbfgs(true);
+ //cg_descent.set_step(0.8);
 
- if (cgReturn == nonLinearSolver::SUCCESS || cgDescentSucess)
+ //if(this_mpi_process == 0)
+ //cg_descent.set_PrintLevel(2);
+
+ //cg_descent.set_AWolfe(true);
+ //bool cgDescentSucess = cg_descent.run(*this);
+
+ if (cgReturn == nonLinearSolver::SUCCESS || cgSuccess)
  {
  pcout<< " ...Ion force relaxation completed as maximum force magnitude is less than FORCE TOL: "<< dftParameters::forceRelaxTol<<", total number of ion position updates: "<<d_totalUpdateCalls<<std::endl;
 
@@ -185,7 +203,7 @@ namespace dftfe {
  }
  pcout<<"-----------------------------------------------------------------------------------"<<std::endl;
  }
- else if (cgReturn == nonLinearSolver::FAILURE || !cgDescentSucess)
+ else if (cgReturn == nonLinearSolver::FAILURE || !cgSuccess)
  {
  pcout<< " ...Ion force relaxation failed "<<std::endl;
 
@@ -249,9 +267,6 @@ namespace dftfe {
 
  }
 
-
-
-
  template<unsigned int FEOrder>
  void geoOptIon<FEOrder>::precondition(std::vector<double> & s,
  const std::vector<double> & gradient) const

diff --git a/src/solvers/cgSolvers/cgPRPNonLinearSolver.cc b/src/solvers/cgSolvers/cgPRPNonLinearSolver.cc
@@ -498,6 +498,8 @@ namespace dftfe {
 
  }
 
+ d_lineSearchDampingParameter = alpha;
+
  //
  //
  //