1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifdef MINIMIZE_CLASS
MinimizeStyle(hftn,MinHFTN)
#else
#ifndef LMP_MIN_HFTN_H
#define LMP_MIN_HFTN_H
#include "min.h"
namespace LAMMPS_NS
{
class MinHFTN : public Min
{
public:
MinHFTN (LAMMPS *);
~MinHFTN (void);
void init();
void setup_style();
void reset_vectors();
int iterate (int);
private:
//---- THE ALGORITHM NEEDS TO STORE THIS MANY ATOM-BASED VECTORS,
//---- IN ADDITION TO ATOM POSITIONS AND THE FORCE VECTOR.
enum {
VEC_XK=0, //-- ATOM POSITIONS AT SUBITER START
VEC_CG_P, //-- STEP p IN CG SUBITER
VEC_CG_D, //-- DIRECTION d IN CG SUBITER
VEC_CG_HD, //-- HESSIAN-VECTOR PRODUCT Hd
VEC_CG_R, //-- RESIDUAL r IN CG SUBITER
VEC_DIF1, //-- FOR FINITE DIFFERENCING
VEC_DIF2, //-- FOR FINITE DIFFERENCING
NUM_HFTN_ATOM_BASED_VECTORS
};
//---- ATOM-BASED STORAGE VECTORS.
double * _daAVectors[NUM_HFTN_ATOM_BASED_VECTORS];
double ** _daExtraAtom[NUM_HFTN_ATOM_BASED_VECTORS];
//---- GLOBAL DOF STORAGE. ELEMENT [0] IS X0 (XK), NOT USED IN THIS ARRAY.
double * _daExtraGlobal[NUM_HFTN_ATOM_BASED_VECTORS];
int _nNumUnknowns;
FILE * _fpPrint;
int execute_hftn_ (const bool bPrintProgress,
const double dInitialEnergy,
const double dInitialForce2,
double & dFinalEnergy,
double & dFinalForce2);
bool compute_inner_cg_step_ (const double dTrustRadius,
const double dForceTol,
const int nMaxEvals,
const bool bHaveEvalAtXin,
const double dEnergyAtXin,
const double dForce2AtXin,
double & dEnergyAtXout,
double & dForce2AtXout,
int & nStepType,
double & dStepLength2,
double & dStepLengthInf);
double calc_xinf_using_mpi_ (void) const;
double calc_dot_prod_using_mpi_ (const int nIx1,
const int nIx2) const;
double calc_grad_dot_v_using_mpi_ (const int nIx) const;
void calc_dhd_dd_using_mpi_ (double & dDHD,
double & dDD) const;
void calc_ppnew_pdold_using_mpi_ (double & dPnewDotPnew,
double & dPoldDotD) const;
void calc_plengths_using_mpi_ (double & dStepLength2,
double & dStepLengthInf) const;
bool step_exceeds_TR_ (const double dTrustRadius,
const double dPP,
const double dPD,
const double dDD,
double & dTau) const;
bool step_exceeds_DMAX_ (void) const;
void adjust_step_to_tau_ (const double tau);
double compute_to_tr_ (const double dPP,
const double dPD,
const double dDD,
const double dTrustRadius,
const bool bConsiderBothRoots,
const double dDHD,
const double dPdotHD,
const double dGradDotD) const;
void evaluate_dir_der_ (const bool bUseForwardDiffs,
const int nIxDir,
const int nIxResult,
const bool bEvaluateAtX,
double & dNewEnergy);
void open_hftn_print_file_ (void);
void hftn_print_line_ (const bool bIsStepAccepted,
const int nIteration,
const int nTotalEvals,
const double dEnergy,
const double dForce2,
const int nStepType,
const double dTrustRadius,
const double dStepLength2,
const double dActualRed,
const double dPredictedRed) const;
void close_hftn_print_file_ (void);
};
}
#endif
#endif