Actual source code: ssls.c
petsc-dev 2014-02-02
1: #include <../src/tao/complementarity/impls/ssls/ssls.h>
3: /*------------------------------------------------------------*/
6: PetscErrorCode TaoSetFromOptions_SSLS(Tao tao)
7: {
8: TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
10: PetscBool flg;
13: PetscOptionsHead("Semismooth method with a linesearch for complementarity problems");
14: PetscOptionsReal("-ssls_delta", "descent test fraction", "",ssls->delta, &(ssls->delta), &flg);
15: PetscOptionsReal("-ssls_rho", "descent test power", "",ssls->rho, &(ssls->rho), &flg);
16: TaoLineSearchSetFromOptions(tao->linesearch);
17: KSPSetFromOptions(tao->ksp);
18: PetscOptionsTail();
19: return(0);
20: }
22: /*------------------------------------------------------------*/
25: PetscErrorCode TaoView_SSLS(Tao tao, PetscViewer pv)
26: {
30: return(0);
31: }
33: /*------------------------------------------------------------*/
36: PetscErrorCode Tao_SSLS_Function(TaoLineSearch ls, Vec X, PetscReal *fcn, void *ptr)
37: {
38: Tao tao = (Tao)ptr;
39: TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
44: TaoComputeConstraints(tao, X, tao->constraints);
45: VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff);
46: VecNorm(ssls->ff,NORM_2,&ssls->merit);
47: *fcn = 0.5*ssls->merit*ssls->merit;
48: return(0);
49: }
51: /*------------------------------------------------------------*/
54: PetscErrorCode Tao_SSLS_FunctionGradient(TaoLineSearch ls, Vec X, PetscReal *fcn, Vec G, void *ptr)
55: {
56: Tao tao = (Tao)ptr;
57: TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
62: TaoComputeConstraints(tao, X, tao->constraints);
63: VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff);
64: VecNorm(ssls->ff,NORM_2,&ssls->merit);
65: *fcn = 0.5*ssls->merit*ssls->merit;
67: TaoComputeJacobian(tao, tao->solution, &tao->jacobian, &tao->jacobian_pre, &ssls->matflag);
69: D_Fischer(tao->jacobian, tao->solution, tao->constraints,tao->XL, tao->XU, ssls->t1, ssls->t2,ssls->da, ssls->db);
70: MatDiagonalScale(tao->jacobian,ssls->db,NULL);
71: MatDiagonalSet(tao->jacobian,ssls->da,ADD_VALUES);
72: MatMultTranspose(tao->jacobian,ssls->ff,G);
74: return(0);
75: }