ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/ConjugateMinimizer.cpp
(Generate patch)

Comparing trunk/OOPSE/libmdtools/ConjugateMinimizer.cpp (file contents):
Revision 1011 by tim, Tue Feb 3 20:47:10 2004 UTC vs.
Revision 1064 by tim, Tue Feb 24 15:44:45 2004 UTC

# Line 1 | Line 1
1   #include "ConjugateMinimizer.hpp"
2 + #include "Utility.hpp"
3 + ConjugateMinimizerBase::ConjugateMinimizerBase(NLModel1* nlmodel, MinimizerParameterSet* param)
4 +                                : MinimizerUsingLineSearch(param){
5 +  int dim;
6 +  
7 +  model = nlmodel;
8 +  //set the dimension
9 +  
10 + #ifndef IS_MPI
11 +  dim = model->getDim();
12 + #else
13 +  dim = model->getDim();
14 + #endif
15 +  prevGrad.resize(dim);    
16 +  gradient.resize(dim);
17 +  prevDirection.resize(dim);
18 +  direction.resize(dim);
19 + }
20 +
21   bool ConjugateMinimizerBase::isSolvable(){
22  
23    //conjuage gradient can only solve unconstrained nonlinear model
24  
25 <  if (!model->hasConstraint())
25 >  if (!model->hasConstraints())
26      return true;
27    else
28      return false;
29   }
30  
12 void ConjugateMinimizerBase::Init(){
13
14 }
15
31   void ConjugateMinimizerBase::printMinizerInfo(){
32  
33   }
34  
35 < void ConjugateMinimizerBase::Minimize(){
35 > void ConjugateMinimizerBase::minimize(){
36    int maxIteration;
37 <  int nextRestIter;
37 >  int nextResetIter;
38    int resetFrq;
39    int nextWriteIter;
40    int writeFrq;
41 +  int lsStatus;
42 +  double gamma;
43 +  double lamda;
44    
45 +  
46    if (!isSolvable()){
47 <    cout << "ConjugateMinimizerBase Error: This nonlinear model can not be solved by " << methodName <<endl;
47 >    cout << "ConjugateMinimizerBase Error: This nonlinear model can not be solved by " << minimizerName <<endl;
48  
49      exit(1);
50    }
# Line 33 | Line 52 | void ConjugateMinimizerBase::Minimize(){
52    printMinizerInfo();
53  
54    resetFrq = paramSet->getResetFrq();
55 <  nextRestIter = resetFrq;
55 >  nextResetIter = resetFrq;
56  
57    writeFrq = paramSet->getWriteFrq();
58    nextWriteIter = writeFrq;
40  
41  prevGrad = model->calcGrad();
59  
60 <  direction = preGrad;
60 >  minX = model->getX();
61 >  gradient = model->calcGrad();
62 >
63 >  for(int i = 0; i < direction.size(); i++)
64 >    direction[i] = -gradient[i];
65    
66    maxIteration = paramSet->getMaxIteration();
67  
68 <  for(currentIter = 0;currentIter < maxIteration; currentIter++){
68 >  for(currentIter = 1;currentIter <= maxIteration; currentIter++){
69  
70 <    // perform line search to minimize f(x + stepSize * direction) where stepSize > 0
71 <    lsMinimizer->minimize(direction, 0.0, 1.0);
70 >    // perform line search to minimize f(x + lamda * direction) where stepSize > 0
71 >    lsMinimizer->minimize(direction, 0.0, 0.01);
72      
73      lsStatus = lsMinimizer->getMinimizationStatus();
74 +
75 +    lamda = lsMinimizer->getMinVar();
76      
77      if(lsStatus ==MINSTATUS_ERROR){
78 <      minStatus = MINSTATUS_ERROR;
78 >      if (lamda == 0){
79 >
80 >          for(int i = 0; i < direction.size(); i++)
81 >          direction[i] = -prevGrad[i];
82 >          
83 >        continue;
84 >      }
85 >        minStatus = MINSTATUS_ERROR;
86        return;
87      }
88 <    
89 <    prevMinX = minX;
90 <    minX = minX + lsMinimizer->getMinVar() * direction;
88 >    else{
89 >      prevMinX = minX;
90 >    }
91  
92 +    for(int i = 0; i < direction.size(); i++)
93 +      minX[i] = minX[i] + lamda * direction[i];
94 +
95      //calculate the gradient
96      prevGrad = gradient;
97 <    
97 >
98 >    model->setX(minX);
99      gradient = model->calcGrad();
100  
101 +    minX = model->getX();
102      // stop if converge
103 <    convStatus = checkConvergence();
69 <    if (convStatus == ){
103 >    if (checkConvergence() > 0){
104        writeOut(minX, currentIter);
105  
106        minStatus = MINSTATUS_CONVERGE;
# Line 75 | Line 109 | void ConjugateMinimizerBase::Minimize(){
109            
110  
111      //calculate the
112 <    gamma = calcGamma(grad, preGrad);
112 >    gamma = calcGamma(gradient, prevGrad);
113  
114      // update new direction
115      prevDirection = direction;
82    direction += gamma * direction;
116  
117 +    for(int i = 0; i < direction.size(); i++)  
118 +      direction[i] = -gradient[i] + gamma * direction[i];
119 +
120      //
121      if (currentIter == nextWriteIter){
122        nextWriteIter += writeFrq;
# Line 96 | Line 132 | void ConjugateMinimizerBase::Minimize(){
132  
133    // if writeFrq is not a multipiler of maxIteration, we need to write the final result
134    // otherwise, we already write it inside the loop, just skip it
135 <  if(currentIter != (nextWriteIter - writeFrq))
135 >  if(currentIter - 1 != (nextWriteIter - writeFrq))
136      writeOut(minX, currentIter);
137  
138    minStatus = MINSTATUS_MAXITER;
# Line 107 | Line 143 | int ConjugateMinimizerBase::checkConvergence(){
143  
144    //test absolute gradient tolerance
145    
146 <  if (norm2(gradient) < paramSet->absGradTol)
146 >  if (sqrt(dotProduct(gradient, gradient)) < paramSet->getGTol())
147      return 1;
148    else
149      return -1;
# Line 118 | Line 154 | double FRCGMinimizer::calcGamma(vector<double>& newGra
154   }
155  
156   double FRCGMinimizer::calcGamma(vector<double>& newGrad, vector<double>& oldGrad){
157 <  return norm2(newGrad) / norm2(oldGrad);
157 >  return dotProduct(newGrad, newGrad) / dotProduct(oldGrad, newGrad);
158   }
159  
160   double PRCGMinimizer::calcGamma(vector<double>& newGrad, vector<double>& oldGrad){
161    double gamma;
162    vector<double> deltaGrad;
163 +  
164 +  for(int i = 0; i < newGrad.size(); i++)
165 +    deltaGrad.push_back(newGrad[i] - oldGrad[i]);
166  
167 <  deltaGrad = newGrad - oldGrad;
129 <
130 <  return norm(deltaGrad, newGrad) / norm2(oldGrad);
167 >  return dotProduct(deltaGrad, newGrad) / dotProduct(oldGrad, oldGrad);
168    
169   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines