Actual source code: ssfls.c
1: #include "src/complementarity/impls/ssls/ssls.h"
3: /*------------------------------------------------------------*/
6: static int TaoSolve_SSFLS(TAO_SOLVER tao, void *solver)
7: {
8: TAO_SSLS *ssls = (TAO_SSLS *)solver;
9: // TaoLinearSolver *lsolver;
10: TaoVec *x, *l, *u, *ff, *dpsi, *d, *w;
11: TaoMat *J;
12: double psi, ndpsi, normd, innerd, t=0;
13: double delta, rho;
14: int iter=0, info;
15: TaoTerminateReason reason;
16: TaoTruth flag;
18: TaoFunctionBegin;
20: // Assume that Setup has been called!
21: // Set the structure for the Jacobian and create a linear solver.
22:
23: delta = ssls->delta;
24: rho = ssls->rho;
26: info = TaoGetSolution(tao, &x); CHKERRQ(info);
27: l=ssls->xl;
28: u=ssls->xu;
29: info = TaoEvaluateVariableBounds(tao,l,u); CHKERRQ(info);
30: info = x->Median(l,x,u); CHKERRQ(info);
31: info = TaoGetJacobian(tao, &J); CHKERRQ(info);
33: ff = ssls->ff;
34: dpsi = ssls->dpsi;
35: d = ssls->d;
36: w = ssls->w;
38: info = x->PointwiseMaximum(x, l); CHKERRQ(info);
39: info = x->PointwiseMinimum(x, u); CHKERRQ(info);
40: info = TaoSetMeritFunction(tao, Tao_SSLS_Function, Tao_SSLS_FunctionGradient,
41: TAO_NULL, TAO_NULL, ssls); CHKERRQ(info);
43: info = J->SetStructure(TAOMAT_UNSYMMETRIC); CHKERRQ(info);
45: // Calculate the function value and fischer function value at the
46: // current iterate
47: info = TaoComputeMeritFunctionGradient(tao, x, &psi, dpsi); CHKERRQ(info);
48: info = dpsi->Norm2(&ndpsi);
50: while (1) {
51: info=PetscLogInfo((tao, "TaoSolve_SSFLS: %d: merit: %5.4e, ndpsi: %5.4e\n",
52: iter, ssls->merit, ndpsi));CHKERRQ(info);
54: // Check the termination criteria
55: info = TaoMonitor(tao,iter++,ssls->merit,ndpsi,0.0,t,&reason);
56: CHKERRQ(info);
57: if (reason!=TAO_CONTINUE_ITERATING) break;
59: // Calculate direction. (Really negative of newton direction. Therefore,
60: // rest of the code uses -d.)
61: info = TaoPreLinearSolve(tao, J); CHKERRQ(info);
62: info = TaoLinearSolve(tao, J, ff, d, &flag); CHKERRQ(info);
63:
64: info = w->CopyFrom(d); CHKERRQ(info);
65: info = w->Negate(); CHKERRQ(info);
66: info = w->BoundGradientProjection(w,l, x, u);
68: info = w->Norm2(&normd); CHKERRQ(info);
69: info = w->Dot(dpsi, &innerd); CHKERRQ(info);
71: // Make sure that we have a descent direction
72: if (innerd >= -delta*pow(normd, rho)) {
73: info = PetscLogInfo((tao, "TaoSolve_SSFLS: %d: newton direction not descent\n", iter)); CHKERRQ(info);
74: info = d->CopyFrom(dpsi); CHKERRQ(info);
75: info = w->Dot(dpsi, &innerd); CHKERRQ(info);
76: }
77: info = d->Negate(); CHKERRQ(info);
78: innerd = -innerd;
80: t = 1;
81: info = TaoLineSearchApply(tao, x, dpsi, d, w,
82: &psi, &t, &tao->lsflag); CHKERRQ(info);
83: info = dpsi->Norm2(&ndpsi);
84: }
85: TaoFunctionReturn(0);
86: }
88: /* ---------------------------------------------------------- */
92: int TaoCreate_SSFLS(TAO_SOLVER tao)
93: {
94: TAO_SSLS *ssls;
95: int info;
97: TaoFunctionBegin;
99: info = TaoNew(TAO_SSLS,&ssls); CHKERRQ(info);
100: info = PetscLogObjectMemory(tao, sizeof(TAO_SSLS)); CHKERRQ(info);
102: ssls->delta = 1e-10;
103: ssls->rho = 2.1;
105: info=TaoSetTaoSolveRoutine(tao,TaoSolve_SSFLS,(void*)ssls); CHKERRQ(info);
106: info=TaoSetTaoSetUpDownRoutines(tao,TaoSetUp_SSLS,TaoSetDown_SSLS); CHKERRQ(info);
107: info=TaoSetTaoOptionsRoutine(tao,TaoSetOptions_SSLS); CHKERRQ(info);
108: info=TaoSetTaoViewRoutine(tao,TaoView_SSLS); CHKERRQ(info);
110: info = TaoCreateProjectedArmijoLineSearch(tao); CHKERRQ(info);
112: info = TaoSetMaximumIterates(tao,2000); CHKERRQ(info);
113: info = TaoSetMaximumFunctionEvaluations(tao,4000); CHKERRQ(info);
115: info = TaoSetTolerances(tao,0,0,0,0); CHKERRQ(info);
116: info = TaoSetGradientTolerances(tao,1.0e-16,0.0,0.0); CHKERRQ(info);
117: info = TaoSetFunctionLowerBound(tao,1.0e-8); CHKERRQ(info);
119: TaoFunctionReturn(0);
120: }