ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@4188 f3b2605a-c512-4ea7-a41b-209d697bcdaa

Browse Source
This commit is contained in:
sjplimp
2010-05-21 20:05:47 +00:00
parent 278d6cb95c
commit 7a5cf1480a
7 changed files with 94 additions and 94 deletions
Download Patch File Download Diff File Expand all files Collapse all files

6
src/min.cpp
View File

@ -640,7 +640,7 @@ double Min::fnorm_sqr()
double *fatom;
double local_norm2_sqr = 0.0;
for (i = 0; i < n3; i++) local_norm2_sqr += f[i]*f[i];
for (i = 0; i < nvec; i++) local_norm2_sqr += fvec[i]*fvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];
@ -670,8 +670,8 @@ double Min::fnorm_inf()
double *fatom;
double local_norm_inf = 0.0;
for (i = 0; i < n3; i++)
local_norm_inf = MAX(fabs(f[i]),local_norm_inf);
for (i = 0; i < nvec; i++)
local_norm_inf = MAX(fabs(fvec[i]),local_norm_inf);
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];

6
src/min.h
View File

@ -66,9 +66,9 @@ class Min : protected Pointers {
double ndoftotal; // total dof for entire problem
int n3; // local atomic dof
double *x; // variables for atomic dof, as 1d vector
double *f; // force vector for atomic dof, as 1d vector
int nvec; // local atomic dof = length of xvec
double *xvec; // variables for atomic dof, as 1d vector
double *fvec; // force vector for atomic dof, as 1d vector
int nextra_global; // # of extra global dof due to fixes
double *fextra; // force vector for extra global dof

16
src/min_cg.cpp
View File

@ -56,7 +56,7 @@ int MinCG::iterate(int niter_max)
// initialize working vectors
for (i = 0; i < n3; i++) h[i] = g[i] = f[i];
for (i = 0; i < nvec; i++) h[i] = g[i] = fvec[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];
@ -94,9 +94,9 @@ int MinCG::iterate(int niter_max)
// force tolerance criterion
dot[0] = dot[1] = 0.0;
for (i = 0; i < n3; i++) {
dot[0] += f[i]*f[i];
dot[1] += f[i]*g[i];
for (i = 0; i < nvec; i++) {
dot[0] += fvec[i]*fvec[i];
dot[1] += fvec[i]*g[i];
}
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
@ -126,8 +126,8 @@ int MinCG::iterate(int niter_max)
if ((niter+1) % nlimit == 0) beta = 0.0;
gg = dotall[0];
for (i = 0; i < n3; i++) {
g[i] = f[i];
for (i = 0; i < nvec; i++) {
g[i] = fvec[i];
h[i] = g[i] + beta*h[i];
}
if (nextra_atom)
@ -150,7 +150,7 @@ int MinCG::iterate(int niter_max)
// reinitialize CG if new search direction h is not downhill
dot[0] = 0.0;
for (i = 0; i < n3; i++) dot[0] += g[i]*h[i];
for (i = 0; i < nvec; i++) dot[0] += g[i]*h[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
gatom = gextra_atom[m];
@ -164,7 +164,7 @@ int MinCG::iterate(int niter_max)
dotall[0] += gextra[i]*hextra[i];
if (dotall[0] <= 0.0) {
for (i = 0; i < n3; i++) h[i] = g[i];
for (i = 0; i < nvec; i++) h[i] = g[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
gatom = gextra_atom[m];

125
src/min_hftn.cpp
View File

@ -165,12 +165,12 @@ void MinHFTN::setup_style()
void MinHFTN::reset_vectors()
{
n3 = 3 * atom->nlocal;
nvec = 3 * atom->nlocal;
//---- ATOMIC DEGREES OF FREEDOM.
if (n3 > 0) {
x = atom->x[0];
f = atom->f[0];
if (nvec > 0) {
xvec = atom->x[0];
fvec = atom->f[0];
}
for (int i = 0; i < NUM_HFTN_ATOM_BASED_VECTORS; i++)
_daAVectors[i] = fix_minimize->request_vector (i);
@ -251,8 +251,8 @@ int MinHFTN::execute_hftn_(const bool bPrintProgress,
//---- SAVE ATOM POSITIONS BEFORE AN ITERATION.
fix_minimize->store_box();
for (int i = 0; i < n3; i++)
_daAVectors[VEC_XK][i] = atom->x[0][i];
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_XK][i] = xvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -268,7 +268,7 @@ int MinHFTN::execute_hftn_(const bool bPrintProgress,
double dXInf = calc_xinf_using_mpi_();
//---- FIND THE NUMBER OF UNKNOWNS.
int nLocalNumUnknowns = n3 + nextra_atom;
int nLocalNumUnknowns = nvec + nextra_atom;
MPI_Allreduce (&nLocalNumUnknowns, &_nNumUnknowns,
1, MPI_INT, MPI_SUM, world);
@ -316,8 +316,8 @@ int MinHFTN::execute_hftn_(const bool bPrintProgress,
//---- THERE WAS AN ERROR. RESTORE TO LAST ACCEPTED STEP.
if (nextra_global)
modify->min_step (0.0, _daExtraGlobal[VEC_CG_P]);
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_XK][i];
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_XK][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -358,8 +358,8 @@ int MinHFTN::execute_hftn_(const bool bPrintProgress,
//---- MOVE TO THE NEW POINT AND EVALUATE ENERGY AND FORCES.
//---- THIS IS THE PLACE WHERE energy_force IS ALLOWED TO RESET.
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_XK][i] + _daAVectors[VEC_CG_P][i];
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_XK][i] + _daAVectors[VEC_CG_P][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -429,8 +429,8 @@ int MinHFTN::execute_hftn_(const bool bPrintProgress,
fix_minimize->store_box();
modify->min_clearstore();
for (int i = 0; i < n3; i++)
_daAVectors[VEC_XK][i] = atom->x[0][i];
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_XK][i] = xvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -478,8 +478,8 @@ int MinHFTN::execute_hftn_(const bool bPrintProgress,
//---- RESTORE THE LAST X_K POSITION.
if (nextra_global)
modify->min_step (0.0, _daExtraGlobal[VEC_CG_P]);
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_XK][i];
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_XK][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -605,7 +605,7 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
for (int i = 0; i < nextra_global; i++)
_daExtraGlobal[VEC_CG_P][i] = 0.0;
}
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_CG_P][i] = 0.0;
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -647,9 +647,9 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
_daExtraGlobal[VEC_CG_D][i] = fextra[i];
}
}
for (int i = 0; i < n3; i++) {
_daAVectors[VEC_CG_R][i] = atom->f[0][i];
_daAVectors[VEC_CG_D][i] = atom->f[0][i];
for (int i = 0; i < nvec; i++) {
_daAVectors[VEC_CG_R][i] = fvec[i];
_daAVectors[VEC_CG_D][i] = fvec[i];
}
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -709,8 +709,8 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
//---- MOVE TO X_K AND COMPUTE ENERGY AND FORCES.
if (nextra_global)
modify->min_step (0.0, _daExtraGlobal[VEC_CG_P]);
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_XK][i];
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_XK][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -736,7 +736,7 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
pGlobal[i] += tau * dGlobal[i];
}
}
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_CG_P][i] += tau * _daAVectors[VEC_CG_D][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -770,7 +770,7 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
pGlobal[i] += dAlpha * dGlobal[i];
}
}
for (int i = 0; i < n3; i++) {
for (int i = 0; i < nvec; i++) {
_daAVectors[VEC_DIF1][i] = _daAVectors[VEC_CG_P][i];
_daAVectors[VEC_CG_P][i] += dAlpha * _daAVectors[VEC_CG_D][i];
}
@ -818,7 +818,7 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
for (int i = 0; i < nextra_global; i++)
rGlobal[i] -= dAlpha * hdGlobal[i];
}
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_CG_R][i] -= dAlpha * _daAVectors[VEC_CG_HD][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -847,7 +847,7 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
for (int i = 0; i < nextra_global; i++)
dGlobal[i] = rGlobal[i] + dBeta * dGlobal[i];
}
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_CG_D][i] = _daAVectors[VEC_CG_R][i]
+ dBeta * _daAVectors[VEC_CG_D][i];
if (nextra_atom) {
@ -877,8 +877,8 @@ bool MinHFTN::compute_inner_cg_step_(const double dTrustRadius,
double MinHFTN::calc_xinf_using_mpi_(void) const
{
double dXInfLocal = 0.0;
for (int i = 0; i < n3; i++)
dXInfLocal = MAX (dXInfLocal, fabs (atom->x[0][i]));
for (int i = 0; i < nvec; i++)
dXInfLocal = MAX(dXInfLocal,fabs(xvec[i]));
double dXInf;
MPI_Allreduce (&dXInfLocal, &dXInf, 1, MPI_DOUBLE, MPI_MAX, world);
@ -909,7 +909,7 @@ double MinHFTN::calc_dot_prod_using_mpi_(const int nIx1,
const int nIx2) const
{
double dDotLocal = 0.0;
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
dDotLocal += _daAVectors[nIx1][i] * _daAVectors[nIx2][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -945,13 +945,13 @@ double MinHFTN::calc_grad_dot_v_using_mpi_(const int nIx) const
//---- REMEMBER THAT FORCES = -GRADIENT.
double dGradDotVLocal = 0.0;
for (int i = 0; i < n3; i++)
dGradDotVLocal += - _daAVectors[nIx][i] * atom->f[0][i];
for (int i = 0; i < nvec; i++)
dGradDotVLocal += - _daAVectors[nIx][i] * fvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
for (int m = 0; m < nextra_atom; m++) {
double * iAtom = _daExtraAtom[nIx][m];
int n = extra_nlen[m];
for (int i = 0; i < n; i++)
for (int i = 0; i < n; i++)
dGradDotVLocal += - iAtom[i] * fextra_atom[m][i];
}
}
@ -978,7 +978,7 @@ void MinHFTN::calc_dhd_dd_using_mpi_(double & dDHD,
{
double dDHDLocal = 0.0;
double dDDLocal = 0.0;
for (int i = 0; i < n3; i++) {
for (int i = 0; i < nvec; i++) {
dDHDLocal += _daAVectors[VEC_CG_D][i] * _daAVectors[VEC_CG_HD][i];
dDDLocal += _daAVectors[VEC_CG_D][i] * _daAVectors[VEC_CG_D][i];
}
@ -1023,7 +1023,7 @@ void MinHFTN::calc_ppnew_pdold_using_mpi_(double & dPnewDotPnew,
{
double dPnewDotPnewLocal = 0.0;
double dPoldDotDLocal = 0.0;
for (int i = 0; i < n3; i++) {
for (int i = 0; i < nvec; i++) {
dPnewDotPnewLocal
+= _daAVectors[VEC_CG_P][i] * _daAVectors[VEC_CG_P][i];
dPoldDotDLocal
@ -1073,7 +1073,7 @@ void MinHFTN::calc_plengths_using_mpi_(double & dStepLength2,
{
double dPPLocal = 0.0;
double dPInfLocal = 0.0;
for (int i = 0; i < n3; i++) {
for (int i = 0; i < nvec; i++) {
dPPLocal += _daAVectors[VEC_CG_P][i] * _daAVectors[VEC_CG_P][i];
dPInfLocal = MAX (dPInfLocal, fabs (_daAVectors[VEC_CG_P][i]));
}
@ -1147,7 +1147,7 @@ bool MinHFTN::step_exceeds_DMAX_(void) const
double dAlpha = dmax * sqrt (_nNumUnknowns);
double dPInfLocal = 0.0;
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
dPInfLocal = MAX (dPInfLocal, fabs (_daAVectors[VEC_CG_P][i]));
double dPInf;
MPI_Allreduce (&dPInfLocal, &dPInf, 1, MPI_DOUBLE, MPI_MAX, world);
@ -1198,7 +1198,7 @@ void MinHFTN::adjust_step_to_tau_(const double tau)
for (int i = 0; i < nextra_global; i++)
pGlobal[i] = d1Global[i] + (tau * dGlobal[i]);
}
for (int i = 0; i < n3; i++) {
for (int i = 0; i < nvec; i++) {
_daAVectors[VEC_CG_P][i] = _daAVectors[VEC_DIF1][i]
+ (tau * _daAVectors[VEC_CG_D][i]);
}
@ -1313,7 +1313,7 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
{
//---- COMPUTE THE MAGNITUDE OF THE DIRECTION VECTOR: |p|_2.
double dDirNorm2SqrdLocal = 0.0;
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
dDirNorm2SqrdLocal
+= _daAVectors[nIxDir][i] * _daAVectors[nIxDir][i];
if (nextra_atom) {
@ -1337,7 +1337,7 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- IF THE STEP IS TOO SMALL, RETURN ZERO FOR THE DERIVATIVE.
if (dDirNorm2 == 0.0) {
for (int i = 0; i < n3; i++)
for (int i = 0; i < nvec; i++)
_daAVectors[nIxResult][i] = 0.0;
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -1370,8 +1370,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- SAVE A COPY OF x.
fix_minimize->store_box();
for (int i = 0; i < n3; i++)
_daAVectors[VEC_DIF1][i] = atom->x[0][i];
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_DIF1][i] = xvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -1389,8 +1389,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- EVALUATE FORCES AT x + eps*p.
if (nextra_global)
modify->min_step (dEps, _daExtraGlobal[nIxDir]);
for (int i = 0; i < n3; i++)
atom->x[0][i] += dEps * _daAVectors[nIxDir][i];
for (int i = 0; i < nvec; i++)
xvec[i] += dEps * _daAVectors[nIxDir][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -1408,8 +1408,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
for (int i = 0; i < nextra_global; i++)
_daExtraGlobal[VEC_DIF2][i] = fextra[i];
}
for (int i = 0; i < n3; i++)
_daAVectors[VEC_DIF2][i] = atom->f[0][i];
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_DIF2][i] = fvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * d2Atom = _daExtraAtom[VEC_DIF2][m];
@ -1424,8 +1424,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
modify->min_step (0.0, _daExtraGlobal[VEC_DIF1]);
modify->min_popstore();
}
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_DIF1][i];
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_DIF1][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -1440,9 +1440,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- COMPUTE THE DIFFERENCE VECTOR: [grad(x + eps*p) - grad(x)] / eps.
//---- REMEMBER THAT FORCES = -GRADIENT.
for (int i = 0; i < n3; i++)
_daAVectors[nIxResult][i]
= (atom->f[0][i] - _daAVectors[VEC_DIF2][i]) / dEps;
for (int i = 0; i < nvec; i++)
_daAVectors[nIxResult][i] = (fvec[i] - _daAVectors[VEC_DIF2][i]) / dEps;
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * iAtom = _daExtraAtom[nIxResult][m];
@ -1465,8 +1464,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- SAVE A COPY OF x.
fix_minimize->store_box();
for (int i = 0; i < n3; i++)
_daAVectors[VEC_DIF1][i] = atom->x[0][i];
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_DIF1][i] = xvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -1484,8 +1483,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- EVALUATE FORCES AT x + eps*p.
if (nextra_global)
modify->min_step (dEps, _daExtraGlobal[nIxDir]);
for (int i = 0; i < n3; i++)
atom->x[0][i] += dEps * _daAVectors[nIxDir][i];
for (int i = 0; i < nvec; i++)
xvec[i] += dEps * _daAVectors[nIxDir][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];
@ -1503,8 +1502,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
for (int i = 0; i < nextra_global; i++)
_daExtraGlobal[VEC_DIF2][i] = fextra[i];
}
for (int i = 0; i < n3; i++)
_daAVectors[VEC_DIF2][i] = atom->f[0][i];
for (int i = 0; i < nvec; i++)
_daAVectors[VEC_DIF2][i] = fvec[i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * d2Atom = _daExtraAtom[VEC_DIF2][m];
@ -1517,8 +1516,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
//---- EVALUATE FORCES AT x - eps*p.
if (nextra_global)
modify->min_step (-dEps, _daExtraGlobal[nIxDir]);
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_DIF1][i]
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_DIF1][i]
- dEps * _daAVectors[nIxDir][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
@ -1542,9 +1541,9 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
for (int i = 0; i < nextra_global; i++)
iGlobal[i] = (fextra[i] - d2Global[i]) / (2.0 + dEps);
}
for (int i = 0; i < n3; i++)
_daAVectors[nIxResult][i]
= (atom->f[0][i] - _daAVectors[VEC_DIF2][i]) / (2.0 * dEps);
for (int i = 0; i < nvec; i++)
_daAVectors[nIxResult][i] =
(fvec[i] - _daAVectors[VEC_DIF2][i]) / (2.0 * dEps);
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * iAtom = _daExtraAtom[nIxResult][m];
@ -1561,8 +1560,8 @@ void MinHFTN::evaluate_dir_der_(const bool bUseForwardDiffs,
modify->min_step (0.0, _daExtraGlobal[VEC_DIF1]);
modify->min_popstore();
}
for (int i = 0; i < n3; i++)
atom->x[0][i] = _daAVectors[VEC_DIF1][i];
for (int i = 0; i < nvec; i++)
xvec[i] = _daAVectors[VEC_DIF1][i];
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
double * xatom = xextra_atom[m];

28
src/min_linesearch.cpp
View File

@ -133,9 +133,9 @@ void MinLineSearch::reset_vectors()
{
// atomic dof
n3 = 3 * atom->nlocal;
if (n3) x = atom->x[0];
if (n3) f = atom->f[0];
nvec = 3 * atom->nlocal;
if (nvec) xvec = atom->x[0];
if (nvec) fvec = atom->f[0];
x0 = fix_minimize->request_vector(0);
g = fix_minimize->request_vector(1);
h = fix_minimize->request_vector(2);
@ -188,7 +188,7 @@ int MinLineSearch::linemin_backtrack(double eoriginal, double &alpha,
// if search direction is not downhill, exit with error
fdothme = 0.0;
for (i = 0; i < n3; i++) fdothme += f[i]*h[i];
for (i = 0; i < nvec; i++) fdothme += fvec[i]*h[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];
@ -211,7 +211,7 @@ int MinLineSearch::linemin_backtrack(double eoriginal, double &alpha,
// if all search dir components are already 0.0, exit with error
hme = 0.0;
for (i = 0; i < n3; i++) hme = MAX(hme,fabs(h[i]));
for (i = 0; i < nvec; i++) hme = MAX(hme,fabs(h[i]));
MPI_Allreduce(&hme,&hmaxall,1,MPI_DOUBLE,MPI_MAX,world);
alpha = MIN(ALPHA_MAX,dmax/hmaxall);
if (nextra_atom)
@ -235,7 +235,7 @@ int MinLineSearch::linemin_backtrack(double eoriginal, double &alpha,
// store box and values of all dof at start of linesearch
fix_minimize->store_box();
for (i = 0; i < n3; i++) x0[i] = x[i];
for (i = 0; i < nvec; i++) x0[i] = xvec[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
xatom = xextra_atom[m];
@ -338,7 +338,7 @@ int MinLineSearch::linemin_quadratic(double eoriginal, double &alpha,
// if search direction is not downhill, exit with error
fdothme = 0.0;
for (i = 0; i < n3; i++) fdothme += f[i]*h[i];
for (i = 0; i < nvec; i++) fdothme += fvec[i]*h[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];
@ -361,7 +361,7 @@ int MinLineSearch::linemin_quadratic(double eoriginal, double &alpha,
// if all search dir components are already 0.0, exit with error
hme = 0.0;
for (i = 0; i < n3; i++) hme = MAX(hme,fabs(h[i]));
for (i = 0; i < nvec; i++) hme = MAX(hme,fabs(h[i]));
MPI_Allreduce(&hme,&hmaxall,1,MPI_DOUBLE,MPI_MAX,world);
alphamax = MIN(ALPHA_MAX,dmax/hmaxall);
if (nextra_atom)
@ -386,7 +386,7 @@ int MinLineSearch::linemin_quadratic(double eoriginal, double &alpha,
// store box and values of all dof at start of linesearch
fix_minimize->store_box();
for (i = 0; i < n3; i++) x0[i] = x[i];
for (i = 0; i < nvec; i++) x0[i] = xvec[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
xatom = xextra_atom[m];
@ -419,9 +419,9 @@ int MinLineSearch::linemin_quadratic(double eoriginal, double &alpha,
// compute new fh, alpha, delfh
dot[0] = dot[1] = 0.0;
for (i = 0; i < n3; i++) {
dot[0] += f[i]*f[i];
dot[1] += f[i]*h[i];
for (i = 0; i < nvec; i++) {
dot[0] += fvec[i]*fvec[i];
dot[1] += fvec[i]*h[i];
}
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
@ -516,7 +516,7 @@ double MinLineSearch::alpha_step(double alpha, int resetflag, int &nfunc)
// reset to starting point
if (nextra_global) modify->min_step(0.0,hextra);
for (i = 0; i < n3; i++) x[i] = x0[i];
for (i = 0; i < nvec; i++) xvec[i] = x0[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
xatom = xextra_atom[m];
@ -529,7 +529,7 @@ double MinLineSearch::alpha_step(double alpha, int resetflag, int &nfunc)
if (alpha > 0.0) {
if (nextra_global) modify->min_step(alpha,hextra);
for (i = 0; i < n3; i++) x[i] += alpha*h[i];
for (i = 0; i < nvec; i++) xvec[i] += alpha*h[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
xatom = xextra_atom[m];

4
src/min_sd.cpp
View File

@ -45,7 +45,7 @@ int MinSD::iterate(int niter_max)
// initialize working vectors
for (i = 0; i < n3; i++) h[i] = f[i];
for (i = 0; i < nvec; i++) h[i] = fvec[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];
@ -85,7 +85,7 @@ int MinSD::iterate(int niter_max)
// set new search direction h to f = -Grad(x)
for (i = 0; i < n3; i++) h[i] = f[i];
for (i = 0; i < nvec; i++) h[i] = fvec[i];
if (nextra_atom)
for (m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];

3
src/run.cpp
View File

@ -55,7 +55,8 @@ void Run::command(int narg, char **arg)
int preflag = 1;
int postflag = 1;
int nevery = 0;
int first,last,ncommands;
int ncommands = 0;
int first,last;
int iarg = 1;
while (iarg < narg) {