Actual source code: lmvm.c
1: /*$Id: lmvm.c 1.99 05/05/10 15:21:38-05:00 sarich@zorak.(none) $*/
3: #include "lmvm.h" /*I "tao_solver.h" I*/
4: #include "src/matrix/lmvmmat.h" /*I "tao_solver.h" I*/
6: static int TaoDestroy_LMVM(TAO_SOLVER, void*);
8: /* --- TESTING ------------------------------------------------------- */
11: int TaoInitializeLMVMmatrix(TAO_SOLVER tao, TaoVec *HV){
12: int info;
13: TAO_LMVM *lmvm;
14:
15: TaoFunctionBegin;
16: info = TaoGetSolverContext(tao,"tao_lmvm",(void**)&lmvm); CHKERRQ(info);
17: if (lmvm && lmvm->M){
18: info = lmvm->M->InitialApproximation( HV );CHKERRQ(info);
19: }
20: TaoFunctionReturn(0);
21: }
25: int TaoSetUp_LMVM(TAO_SOLVER tao, void* solver)
26: {
27: TAO_LMVM *lmvm = (TAO_LMVM *)solver;
28: TaoVec *xx;
29: int info, lm;
31: TaoFunctionBegin;
33: info = TaoGetSolution(tao,&xx);CHKERRQ(info);
35: info = TaoLMVMGetSize(tao,&lm); CHKERRQ(info);
36: lmvm->M = new TaoLMVMMat(xx,lm);
38: info = xx->Clone(&lmvm->G); CHKERRQ(info);
39: info = xx->Clone(&lmvm->Work); CHKERRQ(info);
40: info = xx->Clone(&lmvm->DX); CHKERRQ(info);
42: info = TaoSetLagrangianGradientVector(tao,lmvm->G);CHKERRQ(info);
43: info = TaoSetStepDirectionVector(tao,lmvm->DX);CHKERRQ(info);
45: TaoFunctionReturn(0);
46: }
48: /* ---------------------------------------------------------- */
51: static int TaoDestroy_LMVM(TAO_SOLVER tao, void*solver)
52: {
53: int info;
54: TAO_LMVM *lmvm = (TAO_LMVM *)solver;
56: TaoFunctionBegin;
58: info=TaoVecDestroy(lmvm->G);CHKERRQ(info);
59: info=TaoVecDestroy(lmvm->DX);CHKERRQ(info);
60: info=TaoVecDestroy(lmvm->Work);CHKERRQ(info);
61: info=TaoMatDestroy(lmvm->M);CHKERRQ(info);
63: info = TaoSetLagrangianGradientVector(tao,0);CHKERRQ(info);
64: info = TaoSetStepDirectionVector(tao,0);CHKERRQ(info);
65: lmvm->DX=0;lmvm->Work=0;lmvm->M=0;
67: TaoFunctionReturn(0);
68: }
70: /*------------------------------------------------------------*/
73: static int TaoSolve_LMVM(TAO_SOLVER tao, void *solver)
74: {
75: TAO_LMVM *lmvm = (TAO_LMVM *)solver;
76: TaoVec *xx, *gg=lmvm->G; /* gradient vector */
77: TaoVec *dxdx=lmvm->DX, *ww=lmvm->Work;
78: int info,lsflag=0,iter=0;
79: double f, gdx, gnorm, dxnorm, r0, step=0;
80: TaoTerminateReason reason;
81: TaoTruth success;
83: TaoFunctionBegin;
85: info = TaoGetSolution(tao,&xx);CHKERRQ(info);
87: info = TaoComputeMeritFunctionGradient(tao,xx,&f,gg);CHKERRQ(info);
88: info = gg->Norm2(&gnorm); CHKERRQ(info);
90: /* Enter loop */
91: while (1){
93: /* Test for convergence */
94: info=TaoMonitor(tao,iter++,f,gnorm,0,step,&reason);CHKERRQ(info);
95: if (reason!=TAO_CONTINUE_ITERATING) break;
97: /* Compute step direction */
98: info = lmvm->M->update(xx,gg); CHKERRQ(info);
99: info = lmvm->M->Solve(gg,dxdx,&success); CHKERRQ(info);
100: info = dxdx->Negate(); CHKERRQ(info);
101: info = dxdx->Dot(gg,&gdx); CHKERRQ(info);
102: if (gdx>=0){
103: info = dxdx->ScaleCopyFrom(-1.0,gg); CHKERRQ(info);
104: gdx = -gnorm*gnorm;
105: }
107: /* Line Search */
108: info = dxdx->Norm2(&dxnorm); CHKERRQ(info);
109: info = TaoGetInitialTrustRegionRadius(tao,&r0); CHKERRQ(info);
110: step = TaoMin(1.0,r0/dxnorm);
111: if (step < 1.0){
112: info = PetscLogInfo((tao,"Restrict step size to %g so the change in solution <= %g\n",step,r0)); CHKERRQ(info);
113: }
114: info = TaoLineSearchApply(tao,xx,gg,dxdx,ww,&f,&step,&lsflag);
115: CHKERRQ(info);
116: info = gg->Norm2(&gnorm);CHKERRQ(info);
118: }
120: TaoFunctionReturn(0);
122: }
125: /*------------------------------------------------------------*/
128: static int TaoSetOptions_LMVM(TAO_SOLVER tao, void *solver)
129: {
130: int info;
132: TaoFunctionBegin;
133: info = TaoOptionsHead("Limited Memory Variable Metric method for unconstrained optimization");CHKERRQ(info);
134: info = TaoLineSearchSetFromOptions(tao);CHKERRQ(info);
135: info = TaoOptionsTail();CHKERRQ(info);
137: TaoFunctionReturn(0);
138: }
140: /* ---------------------------------------------------------- */
143: int TaoLMVMGetX0(TAO_SOLVER tao,TaoVec *x0)
144: {
145: TAO_LMVM *lmvm;
146: int info;
147: TaoFunctionBegin;
148: TaoValidHeaderSpecific(tao,TAO_COOKIE,1);
149: info=TaoGetSolverContext(tao,"tao_lmvm",(void**)&lmvm);CHKERRQ(info);
150: if (lmvm && lmvm->M){
151: info=lmvm->M->GetX0(x0);CHKERRQ(info);
152: }
153: TaoFunctionReturn(0);
154: }
156: /*------------------------------------------------------------*/
159: static int TaoView_LMVM(TAO_SOLVER tao,void *solver)
160: {
161: int info;
163: TaoFunctionBegin;
164: info = TaoLineSearchView(tao);CHKERRQ(info);
165: TaoFunctionReturn(0);
166: }
168: /* ---------------------------------------------------------- */
173: int TaoCreate_LMVM(TAO_SOLVER tao)
174: {
175: TAO_LMVM *lmvm;
176: int info;
178: TaoFunctionBegin;
180: info = TaoNew(TAO_LMVM,&lmvm); CHKERRQ(info);
181: info = PetscLogObjectMemory(tao,sizeof(TAO_LMVM)); CHKERRQ(info);
183: info=TaoSetTaoSolveRoutine(tao,TaoSolve_LMVM,(void*)lmvm); CHKERRQ(info);
184: info=TaoSetTaoSetUpDownRoutines(tao,TaoSetUp_LMVM,TaoDestroy_LMVM); CHKERRQ(info);
185: info=TaoSetTaoOptionsRoutine(tao,TaoSetOptions_LMVM); CHKERRQ(info);
186: info=TaoSetTaoViewRoutine(tao,TaoView_LMVM); CHKERRQ(info);
188: info = TaoSetMaximumIterates(tao,2000); CHKERRQ(info);
189: info = TaoSetMaximumFunctionEvaluations(tao,4000); CHKERRQ(info);
190: info = TaoSetTolerances(tao,1e-4,1e-4,0,0); CHKERRQ(info);
191: info = TaoLMVMSetSize(tao,5); CHKERRQ(info);
192: info = TaoSetTrustRegionRadius(tao,1.0e16); CHKERRQ(info);
193: info = TaoCreateMoreThuenteLineSearch(tao,0,0.99); CHKERRQ(info);
194: lmvm->DX=0;lmvm->Work=0;lmvm->M=0;
196: TaoFunctionReturn(0);
197: }