Revision | 57956c5db62c87f7ce2a19342e66cae88289ba1d (tree) |
---|---|
Time | 2012-11-28 13:40:12 |
Author | Katsuhiko Nishimra <ktns.87@gmai...> |
Commiter | Katsuhiko Nishimra |
tmp
@@ -66,7 +66,7 @@ GDIIS::~GDIIS(){ | ||
66 | 66 | } |
67 | 67 | } |
68 | 68 | |
69 | -bool GDIIS::DoGDIIS(double* vectorError, | |
69 | +bool GDIIS::TryGDIIS(double* vectorError, | |
70 | 70 | double* vectorPosition, |
71 | 71 | double const* vectorRefStep){ |
72 | 72 | this->Update(vectorError, vectorPosition); |
@@ -128,13 +128,19 @@ bool GDIIS::CalcGDIIS(double* vectorError, | ||
128 | 128 | } |
129 | 129 | matrixGDIIS[numErrors][numErrors] = 0; |
130 | 130 | |
131 | + // If only one error vector is given, following routine is meaningless. | |
132 | + if(numErrors <= 1){ | |
133 | + return false; | |
134 | + } | |
135 | + | |
131 | 136 | double* vectorCoefs = NULL; |
132 | 137 | try{ |
133 | 138 | // Solve DIIS equation |
134 | 139 | MallocerFreer::GetInstance()->Malloc(&vectorCoefs, numErrors+1); |
135 | - for(int i=0;i<numErrors;i++){ | |
136 | - vectorCoefs[i]=0.0; | |
137 | - } | |
140 | + // Initialized in Malloc | |
141 | + // for(int i=0;i<numErrors;i++){ | |
142 | + // vectorCoefs[i]=0.0; | |
143 | + // } | |
138 | 144 | vectorCoefs[numErrors]=1.0; |
139 | 145 | try{ |
140 | 146 | MolDS_wrappers::Lapack::GetInstance()->Dsysv(matrixGDIIS, |
@@ -150,12 +156,6 @@ bool GDIIS::CalcGDIIS(double* vectorError, | ||
150 | 156 | return false; |
151 | 157 | } |
152 | 158 | |
153 | - // If only one error vector is given, following routine is meaningless. | |
154 | - if(numErrors <= 1){ | |
155 | - MallocerFreer::GetInstance()->Free(&vectorCoefs, numErrors+1); | |
156 | - return false; | |
157 | - } | |
158 | - | |
159 | 159 | // If Lagrange multiplier is too small; |
160 | 160 | if(-vectorCoefs[numErrors] < 1e-8){ |
161 | 161 | this->OutputLog((formatTooSmallLagrangeMultiplier % -vectorCoefs[numErrors]).str()); |
@@ -177,23 +177,23 @@ bool GDIIS::CalcGDIIS(double* vectorError, | ||
177 | 177 | |
178 | 178 | // Calculate cosine of the angle between GDIIS step and vectorRefStep |
179 | 179 | // and lengths of the vectors. |
180 | - double innerprod = 0, norm2gdiis = 0, norm2ref = 0; | |
180 | + double innerprod = 0, normSquaregdiis = 0, normSquareref = 0; | |
181 | 181 | for(int i=0;i<this->sizeErrorVector;i++){ |
182 | 182 | double diff = vectorPosition[i] - listPositions.back()[i]; |
183 | 183 | innerprod += diff * vectorRefStep[i]; |
184 | - norm2gdiis += diff * diff; | |
185 | - norm2ref += vectorRefStep[i] * vectorRefStep[i]; | |
184 | + normSquaregdiis += diff * diff; | |
185 | + normSquareref += vectorRefStep[i] * vectorRefStep[i]; | |
186 | 186 | } |
187 | - double cosine = innerprod/sqrt(norm2gdiis*norm2ref); | |
187 | + double cosine = innerprod/sqrt(normSquaregdiis*normSquareref); | |
188 | 188 | |
189 | 189 | // If length of the GDIIS step is larger than reference step * 10 |
190 | - if(norm2gdiis >= norm2ref * 100){ | |
190 | + if(normSquaregdiis >= normSquareref * 100){ | |
191 | 191 | // Rollback vectorPosition and vectorError to original value |
192 | 192 | for(int i=0; i<this->sizeErrorVector; i++){ |
193 | 193 | vectorError[i] = listErrors.back()[i]; |
194 | 194 | vectorPosition[i] = listPositions.back()[i]; |
195 | 195 | } |
196 | - this->OutputLog((formatTooLargeGDIISStep % sqrt(norm2gdiis) % sqrt(norm2ref)).str()); | |
196 | + this->OutputLog((formatTooLargeGDIISStep % sqrt(normSquaregdiis) % sqrt(normSquareref)).str()); | |
197 | 197 | MallocerFreer::GetInstance()->Free(&vectorCoefs, numErrors+1); |
198 | 198 | // and recalculate GDIIS step without the oldest data |
199 | 199 | this->DiscardOldest(); |
@@ -223,7 +223,7 @@ bool GDIIS::CalcGDIIS(double* vectorError, | ||
223 | 223 | return true; |
224 | 224 | } |
225 | 225 | |
226 | -bool GDIIS::DoGDIIS(double *vectorError, Molecule& molecule, double const* vectorRefStep){ | |
226 | +bool GDIIS::TryGDIIS(double *vectorError, Molecule& molecule, double const* vectorRefStep){ | |
227 | 227 | double** matrixPosition = NULL; |
228 | 228 | try{ |
229 | 229 | MallocerFreer::GetInstance()->Malloc(&matrixPosition, molecule.GetNumberAtoms(), CartesianType_end); |
@@ -233,7 +233,7 @@ bool GDIIS::DoGDIIS(double *vectorError, Molecule& molecule, double const* vecto | ||
233 | 233 | matrixPosition[i][j] = atom->GetXyz()[j]; |
234 | 234 | } |
235 | 235 | } |
236 | - if(this->DoGDIIS(vectorError,&matrixPosition[0][0],vectorRefStep)){ | |
236 | + if(this->TryGDIIS(vectorError,&matrixPosition[0][0],vectorRefStep)==SUCCEED){ | |
237 | 237 | for(int i=0;i<molecule.GetNumberAtoms();i++){ |
238 | 238 | Atom* atom = molecule.GetAtom(i); |
239 | 239 | for(int j=0;j<CartesianType_end;j++){ |
@@ -255,6 +255,7 @@ bool GDIIS::DoGDIIS(double *vectorError, Molecule& molecule, double const* vecto | ||
255 | 255 | } |
256 | 256 | |
257 | 257 | double GDIIS::MinCosine(){ |
258 | + //reference | |
258 | 259 | static const double inf = std::numeric_limits<double>::infinity(); |
259 | 260 | static const double mincos[] = {inf, inf, 0.97, 0.84, 0.71, 0.67, 0.62, 0.56, 0.49, 0.41}; |
260 | 261 | static const int nummincos = sizeof(mincos)/sizeof(mincos[0]); |