modify line for spin_cubic, spin_none. edit docs a bit.

doc/src/min_modify.txt

@@ -84,12 +84,12 @@ The choice of a line search algorithm for the {spin_oso_cg} and
 {spin_oso_lbfgs} styles can be specified via the {line} keyword.
 The {spin_cubic} and  {spin_none} only make sense when one of those 
 two minimization styles is declared.
-
+The {spin_cubic} performs the line search based on a cubic interpolation
-The {spin_cubic} keyword activates the line search procedure when 
+of the energy along the search direction. The {spin_none} keyword
-the {spin_oso_cg} algorithm is used.
+deactivates the line search procedure.
-
+The {spin_none} is a default value for {line} keyword apart from the case when
-The {spin_none} keyword deactivates the line search procedure when 
+single-replica calculations are performed with {spin_oso_lbfgs} that
-the {spin_oso_lbfgs} algorithm is used.
+uses {spin_cubic} line search.
 
 [Restrictions:] The line search procedure of styles
 {spin_oso_cg} and {spin_oso_lbfgs} cannot be used for magnetic

doc/src/min_spin.txt

@@ -18,7 +18,7 @@ min_style spin_oso_lbfgs :pre
 [Examples:]
 
 min_style 	spin_oso_lbfgs 
-min_modify      discrete_factor 10.0 line spin_none :pre 
+min_modify      line spin_none discrete_factor 10.0 :pre
 
 [Description:]
 
@@ -56,11 +56,20 @@ Style {spin_oso_cg} defines an orthogonal spin optimization
 The "min_modify"_min_modify.html command can be used to
 couple the {spin_oso_cg} to a line search procedure, and to modify the 
 discretization factor {discrete_factor}.
+By defualt, the style {spin_oso_cg} does not employ line search procedure and
+and uses the adaptive time-step technique in the same way as style {spin}.
 
 Style {spin_oso_lbfgs} defines an orthogonal spin optimization
 (OSO) combined to a limited-memory Broyden-Fletcher-Goldfarb-Shanno 
-(LBFGS) algorithm. 
+(L-BFGS) algorithm.
-By default, style {spin_oso_lbfgs} uses a line search procedure.
+By default, style {spin_oso_lbfgs} uses a line search procedure
+based on cubic interpolation for
+a single-replica calculation, and it does not use line search procedure
+for a multireplica calculation (such as in case of GNEB calculation).
+If the line search procedure is not used then the discrete factor defines
+the maximum root mean squared rotation angle of spins by equation {pi/(5*Kappa)}.
+The default value for Kappa is 10.
+
 The "min_modify"_min_modify.html command can be used to
 deactivate the line search procedure, and to modify the
 discretization factor {discrete_factor}. 

doc/src/minimize.txt

@@ -106,10 +106,10 @@ the number of total force evaluations exceeds {maxeval} :ul
 NOTE: the "minimization style"_min_style.html {spin},
 {spin_oso_cg}, and {spin_oso_lbfgs} replace
 the force tolerance {ftol} by a torque tolerance.
-The minimization procedure stops if the 2-norm (length) of the 
+The minimization procedure stops if the 2-norm (length) of the torque vector on atom
-global torque vector (defined as the cross product between the 
+(defined as the cross product between the
-spins and their precession vectors omega) is less than {ftol}, 
+atomic spin and its precession vectors omega) is less than {ftol},
-or if any of the other criteria are met.
+or if any of the other criteria are met. Torque have the same units as the energy.
 
 NOTE: You can also use the "fix halt"_fix_halt.html command to specify
 a general criterion for exiting a minimization, that is a calculation

doc/src/neb_spin.txt

@@ -59,7 +59,7 @@ performance speed-up you would see with one or more physical
 processors per replica.  See the "Howto replica"_Howto_replica.html
 doc page for further discussion.
 
-NOTE: As explained below, a GNEB calculation performs a damped dynamics
+NOTE: As explained below, a GNEB calculation performs a
 minimization across all the replicas. One of the "spin"_min_spin.html
 style minimizers has to be defined in your input script. 
 
@@ -170,8 +170,8 @@ command is issued.
 :line
 
 A NEB calculation proceeds in two stages, each of which is a
-minimization procedure, performed via damped dynamics.  To enable
+minimization procedure.  To enable
-this, you must first define a damped spin dynamics
+this, you must first define a
 "min_style"_min_style.html, using either the {spin},
 {spin_oso_cg}, or {spin_oso_lbfgs} style (see
 "min_spin"_min_spin.html for more information).  

examples/SPIN/spinmin/in.spinmin_cg.bfo

@@ -42,13 +42,13 @@ variable 	magnorm   equal c_out_mag[4]
 variable 	emag      equal c_out_mag[5]
 variable 	tmag      equal c_out_mag[6]
 
-thermo          50
+thermo          100
 thermo_style    custom step time v_magnorm v_emag v_tmag temp etotal
 thermo_modify   format float %20.15g
 
 compute 	outsp all property/atom spx spy spz sp fmx fmy fmz
 dump 		1 all custom 50 dump.lammpstrj type x y z c_outsp[1] c_outsp[2] c_outsp[3] c_outsp[4] c_outsp[5] c_outsp[6] c_outsp[7]
 
-min_style	spin/oso_cg
+min_style	spin_oso_cg
-min_modify 	discrete_factor 10.0 line spin_cubic
+# min_modify 	line spin_none discrete_factor 10.0
-minimize        1.0e-10 1.0e-10 10000 1000
+minimize        1.0e-10 1.0e-7 1000 1000

examples/SPIN/spinmin/in.spinmin_lbfgs.bfo

@@ -49,6 +49,6 @@ thermo_modify   format float %20.15g
 compute 	outsp all property/atom spx spy spz sp fmx fmy fmz
 dump 		1 all custom 50 dump.lammpstrj type x y z c_outsp[1] c_outsp[2] c_outsp[3] c_outsp[4] c_outsp[5] c_outsp[6] c_outsp[7]
 
-min_style	spin/oso_lbfgs
+min_style	spin_oso_lbfgs
-min_modify 	discrete_factor 10.0 line spin_none
+min_modify 	line spin_cubic discrete_factor 10.0
-minimize        1.0e-15 1.0e-10 10000 1000
+minimize        1.0e-15 1.0e-7 10000 1000

src/MAKE/Makefile.serial

@@ -7,7 +7,7 @@ SHELL = /bin/sh
 # specify flags and libraries needed for your compiler
 
 CC =		g++
-CCFLAGS =	-g -O3
+CCFLAGS =	-g -O3 -Wall
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M
 

src/SPIN/min_spin_oso_cg.cpp

@@ -63,7 +63,7 @@ static const char cite_minstyle_spin_oso_cg[] =
 /* ---------------------------------------------------------------------- */
 
 MinSpinOSO_CG::MinSpinOSO_CG(LAMMPS *lmp) :
-  Min(lmp), g_old(NULL), g_cur(NULL), p_s(NULL)
+  Min(lmp), g_old(NULL), g_cur(NULL), p_s(NULL), sp_copy(NULL)
 {
   if (lmp->citeme) lmp->citeme->add(cite_minstyle_spin_oso_cg);
   nlocal_max = 0;
@@ -99,6 +99,13 @@ void MinSpinOSO_CG::init()
 
   Min::init();
 
+  if (linestyle == 3 && nreplica == 1){
+    use_line_search = 1;
+  }
+  else{
+    use_line_search = 0;
+  }
+
   dts = dt = update->dt;
   last_negative = update->ntimestep;
 
@@ -132,13 +139,6 @@ void MinSpinOSO_CG::setup_style()
 
 int MinSpinOSO_CG::modify_param(int narg, char **arg)
 {
-  if (strcmp(arg[0],"line_search") == 0) {
-    if (narg < 2) error->all(FLERR,"Illegal min_modify command");
-    use_line_search = force->numeric(FLERR,arg[1]);
-    if (nreplica > 1 && use_line_search)
-      error->all(FLERR,"Illegal fix_modify command, cannot use NEB and line search together");
-    return 2;
-  }
   if (strcmp(arg[0],"discrete_factor") == 0) {
     if (narg < 2) error->all(FLERR,"Illegal fix_modify command");
     discrete_factor = force->numeric(FLERR,arg[1]);
@@ -181,7 +181,6 @@ int MinSpinOSO_CG::iterate(int maxiter)
   double der_e_cur_tmp = 0.0;
 
   if (nlocal_max < nlocal) {
-    nlocal_max = nlocal;
     local_iter = 0;
     nlocal_max = nlocal;
     memory->grow(g_old,3*nlocal_max,"min/spin/oso/cg:g_old");
@@ -205,8 +204,9 @@ int MinSpinOSO_CG::iterate(int maxiter)
     if (use_line_search) {
 
       // here we need to do line search
-      if (local_iter == 0)
+      if (local_iter == 0){
         calc_gradient();
+      }
 
       calc_search_direction();
       der_e_cur = 0.0;
@@ -228,24 +228,15 @@ int MinSpinOSO_CG::iterate(int maxiter)
     else{
 
       // here we don't do line search
-      // but use cutoff rotation angle
       // if gneb calc., nreplica > 1
       // then calculate gradients and advance spins
       // of intermediate replicas only
-
-      if (nreplica > 1) {
-      if(ireplica != 0 && ireplica != nreplica-1)
       calc_gradient();
       calc_search_direction();
       advance_spins();
-      } else{
+      neval++;
-      calc_gradient();
-      calc_search_direction();
-      advance_spins();
-      }
       eprevious = ecurrent;
       ecurrent = energy_force(0);
-      neval++;
     }
 
     // energy tolerance criterion
@@ -336,10 +327,18 @@ void MinSpinOSO_CG::calc_search_direction()
   double g2_global = 0.0;
   double g2old_global = 0.0;
 
+  double factor = 1.0;
+
+  // for multiple replica do not move end points
+  if (nreplica > 1)
+    if (ireplica == 0 || ireplica == nreplica - 1)
+      factor = 0.0;
+
+
   if (local_iter == 0 || local_iter % 5 == 0){ 	// steepest descent direction
     for (int i = 0; i < 3 * nlocal; i++) {
-      p_s[i] = -g_cur[i];
+      p_s[i] = -g_cur[i] * factor;
-      g_old[i] = g_cur[i];
+      g_old[i] = g_cur[i] * factor;
     }
   } else { 				// conjugate direction
     for (int i = 0; i < 3 * nlocal; i++) {
@@ -354,9 +353,9 @@ void MinSpinOSO_CG::calc_search_direction()
     MPI_Allreduce(&g2old,&g2old_global,1,MPI_DOUBLE,MPI_SUM,world);
 
     // Sum over all replicas. Good for GNEB.
-    if (update->multireplica == 1) {
+    if (nreplica > 1) {
-      g2 = g2_global;
+      g2 = g2_global * factor;
-      g2old = g2old_global;
+      g2old = g2old_global * factor;
       MPI_Allreduce(&g2,&g2_global,1,MPI_DOUBLE,MPI_SUM,universe->uworld);
       MPI_Allreduce(&g2old,&g2old_global,1,MPI_DOUBLE,MPI_SUM,universe->uworld);
     }
@@ -364,8 +363,8 @@ void MinSpinOSO_CG::calc_search_direction()
     else beta = g2_global / g2old_global;
     // calculate conjugate direction
     for (int i = 0; i < 3 * nlocal; i++) {
-      p_s[i] = (beta * p_s[i] - g_cur[i]);
+      p_s[i] = (beta * p_s[i] - g_cur[i]) * factor;
-      g_old[i] = g_cur[i];
+      g_old[i] = g_cur[i] * factor;
     }
   }
 
@@ -380,8 +379,6 @@ void MinSpinOSO_CG::advance_spins()
 {
   int nlocal = atom->nlocal;
   double **sp = atom->sp;
-  double **fm = atom->fm;
-  double tdampx, tdampy, tdampz;
   double rot_mat[9]; // exponential of matrix made of search direction
   double s_new[3];
 
@@ -477,7 +474,7 @@ void MinSpinOSO_CG::rodrigues_rotation(const double *upp_tr, double *out)
 
   A = cos(theta);
   B = sin(theta);
-  D = 1 - A;
+  D = 1.0 - A;
   x = upp_tr[0]/theta;
   y = upp_tr[1]/theta;
   z = upp_tr[2]/theta;
@@ -529,7 +526,7 @@ void MinSpinOSO_CG::make_step(double c, double *energy_and_der)
   double rot_mat[9]; // exponential of matrix made of search direction
   double s_new[3];
   double **sp = atom->sp;
-  double der_e_cur_tmp = 0.0;;
+  double der_e_cur_tmp = 0.0;
 
   for (int i = 0; i < nlocal; i++) {
 
@@ -629,7 +626,8 @@ int MinSpinOSO_CG::awc(double der_phi_0, double phi_0, double der_phi_j, double
   double delta = 0.1;
   double sigma = 0.9;
 
-  if ((phi_j<=phi_0+eps*fabs(phi_0)) && ((2.0*delta-1.0) * der_phi_0>=der_phi_j>=sigma*der_phi_0))
+  if ((phi_j<=phi_0+eps*fabs(phi_0)) &&
+      ((2.0*delta-1.0) * der_phi_0>=der_phi_j>=sigma*der_phi_0))
     return 1;
   else
     return 0;

src/SPIN/min_spin_oso_cg.h

@@ -13,7 +13,7 @@
 
 #ifdef MINIMIZE_CLASS
 
-MinimizeStyle(spin/oso_cg, MinSpinOSO_CG)
+MinimizeStyle(spin_oso_cg, MinSpinOSO_CG)
 
 #else
 
@@ -39,8 +39,8 @@ class MinSpinOSO_CG: public Min {
     int ireplica,nreplica; // for neb
     double *spvec;		// variables for atomic dof, as 1d vector
     double *fmvec;		// variables for atomic dof, as 1d vector
-    double *g_cur;  	// current gradient vector
     double *g_old;  	// gradient vector at previous step
+    double *g_cur;  	// current gradient vector
     double *p_s;  		// search direction vector
     double **sp_copy;   // copy of the spins
     int local_iter;     // for neb

src/SPIN/min_spin_oso_lbfgs.cpp

@@ -107,6 +107,13 @@ void MinSpinOSO_LBFGS::init()
 
   Min::init();
 
+  if (linestyle != 4 && nreplica == 1){
+    use_line_search = 1;
+  }
+  else{
+    use_line_search = 0;
+  }
+
   last_negative = update->ntimestep;
 
   // allocate tables
@@ -143,14 +150,6 @@ void MinSpinOSO_LBFGS::setup_style()
 
 int MinSpinOSO_LBFGS::modify_param(int narg, char **arg)
 {
-
-  if (strcmp(arg[0],"line_search") == 0) {
-    if (narg < 2) error->all(FLERR,"Illegal min_modify command");
-    use_line_search = force->numeric(FLERR,arg[1]);
-    if (nreplica > 1 && use_line_search)
-      error->all(FLERR,"Illegal min_modify command, cannot use NEB and line search together");
-    return 2;
-  }
   if (strcmp(arg[0],"discrete_factor") == 0) {
     if (narg < 2) error->all(FLERR,"Illegal min_modify command");
     double discrete_factor;
@@ -221,8 +220,11 @@ int MinSpinOSO_LBFGS::iterate(int maxiter)
     if (use_line_search) {
 
       // here we need to do line search
-      if (local_iter == 0)
+      if (local_iter == 0){
+        eprevious = ecurrent;
+        ecurrent = energy_force(0);
         calc_gradient();
+      }
 
       calc_search_direction();
       der_e_cur = 0.0;
@@ -248,19 +250,11 @@ int MinSpinOSO_LBFGS::iterate(int maxiter)
       // if gneb calc., nreplica > 1
       // then calculate gradients and advance spins
       // of intermediate replicas only
-
-      if (nreplica > 1) {
-      if(ireplica != 0 && ireplica != nreplica-1)
-      calc_gradient();
-      calc_search_direction();
-      advance_spins();
-      } else{
-      calc_gradient();
-      calc_search_direction();
-      advance_spins();
-      }
       eprevious = ecurrent;
       ecurrent = energy_force(0);
+      calc_gradient();
+      calc_search_direction();
+      advance_spins();
       neval++;
     }
 
@@ -398,7 +392,7 @@ void MinSpinOSO_LBFGS::calc_search_direction()
     }
     MPI_Allreduce(&dyds, &dyds_global, 1, MPI_DOUBLE, MPI_SUM, world);
 
-    if (update->multireplica == 1) {
+    if (nreplica > 1) {
       dyds_global *= factor;
       dyds = dyds_global;
       MPI_Allreduce(&dyds, &dyds_global, 1,MPI_DOUBLE,MPI_SUM,universe->uworld);
@@ -437,7 +431,7 @@ void MinSpinOSO_LBFGS::calc_search_direction()
         sq += ds[c_ind][i] * q[i];
       }
       MPI_Allreduce(&sq,&sq_global,1,MPI_DOUBLE,MPI_SUM,world);
-      if (update->multireplica == 1) {
+      if (nreplica > 1) {
         sq_global *= factor;
         sq = sq_global;
         MPI_Allreduce(&sq,&sq_global,1,MPI_DOUBLE,MPI_SUM,universe->uworld);
@@ -460,7 +454,7 @@ void MinSpinOSO_LBFGS::calc_search_direction()
       yy += dy[m_index][i] * dy[m_index][i];
     }
     MPI_Allreduce(&yy,&yy_global,1,MPI_DOUBLE,MPI_SUM,world);
-    if (update->multireplica == 1) {
+    if (nreplica > 1) {
       yy_global *= factor;
       yy = yy_global;
       MPI_Allreduce(&yy,&yy_global,1,MPI_DOUBLE,MPI_SUM,universe->uworld);
@@ -493,7 +487,7 @@ void MinSpinOSO_LBFGS::calc_search_direction()
       }
 
       MPI_Allreduce(&yr,&yr_global,1,MPI_DOUBLE,MPI_SUM,world);
-      if (update->multireplica == 1) {
+      if (nreplica > 1) {
         yr_global *= factor;
         yr = yr_global;
         MPI_Allreduce(&yr,&yr_global,1,MPI_DOUBLE,MPI_SUM,universe->uworld);
@@ -668,7 +662,7 @@ void MinSpinOSO_LBFGS::make_step(double c, double *energy_and_der)
   double rot_mat[9]; // exponential of matrix made of search direction
   double s_new[3];
   double **sp = atom->sp;
-  double der_e_cur_tmp = 0.0;;
+  double der_e_cur_tmp = 0.0;
 
   for (int i = 0; i < nlocal; i++) {
 
@@ -784,12 +778,12 @@ double MinSpinOSO_LBFGS::maximum_rotation(double *p)
   for (int i = 0; i < 3 * nlocal; i++) norm2 += p[i] * p[i];
 
   MPI_Allreduce(&norm2,&norm2_global,1,MPI_DOUBLE,MPI_SUM,world);
-  if (update->multireplica == 1) {
+  if (nreplica > 1) {
     norm2 = norm2_global;
     MPI_Allreduce(&norm2,&norm2_global,1,MPI_DOUBLE,MPI_SUM,universe->uworld);
   }
   MPI_Allreduce(&nlocal,&ntotal,1,MPI_INT,MPI_SUM,world);
-  if (update->multireplica == 1) {
+  if (nreplica > 1) {
     nlocal = ntotal;
     MPI_Allreduce(&nlocal,&ntotal,1,MPI_INT,MPI_SUM,universe->uworld);
   }

src/SPIN/min_spin_oso_lbfgs.h

@@ -13,7 +13,7 @@
 
 #ifdef MINIMIZE_CLASS
 
-MinimizeStyle(spin/oso_lbfgs, MinSpinOSO_LBFGS)
+MinimizeStyle(spin_oso_lbfgs, MinSpinOSO_LBFGS)
 
 #else