1 #include <../src/tao/complementarity/impls/ssls/ssls.h>
2
3 /*------------------------------------------------------------*/
TaoSetFromOptions_SSLS(PetscOptionItems * PetscOptionsObject,Tao tao)4 PetscErrorCode TaoSetFromOptions_SSLS(PetscOptionItems *PetscOptionsObject,Tao tao)
5 {
6 TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
7 PetscErrorCode ierr;
8
9 PetscFunctionBegin;
10 ierr = PetscOptionsHead(PetscOptionsObject,"Semismooth method with a linesearch for complementarity problems");CHKERRQ(ierr);
11 ierr = PetscOptionsReal("-ssls_delta", "descent test fraction", "",ssls->delta, &ssls->delta, NULL);CHKERRQ(ierr);
12 ierr = PetscOptionsReal("-ssls_rho", "descent test power", "",ssls->rho, &ssls->rho, NULL);CHKERRQ(ierr);
13 ierr = TaoLineSearchSetFromOptions(tao->linesearch);CHKERRQ(ierr);
14 ierr = KSPSetFromOptions(tao->ksp);CHKERRQ(ierr);
15 ierr = PetscOptionsTail();CHKERRQ(ierr);
16 PetscFunctionReturn(0);
17 }
18
19 /*------------------------------------------------------------*/
TaoView_SSLS(Tao tao,PetscViewer pv)20 PetscErrorCode TaoView_SSLS(Tao tao, PetscViewer pv)
21 {
22 PetscFunctionBegin;
23 PetscFunctionReturn(0);
24 }
25
26 /*------------------------------------------------------------*/
Tao_SSLS_Function(TaoLineSearch ls,Vec X,PetscReal * fcn,void * ptr)27 PetscErrorCode Tao_SSLS_Function(TaoLineSearch ls, Vec X, PetscReal *fcn, void *ptr)
28 {
29 Tao tao = (Tao)ptr;
30 TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
31 PetscErrorCode ierr;
32
33 PetscFunctionBegin;
34 ierr = TaoComputeConstraints(tao, X, tao->constraints);CHKERRQ(ierr);
35 ierr = VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff);CHKERRQ(ierr);
36 ierr = VecNorm(ssls->ff,NORM_2,&ssls->merit);CHKERRQ(ierr);
37 *fcn = 0.5*ssls->merit*ssls->merit;
38 PetscFunctionReturn(0);
39 }
40
41 /*------------------------------------------------------------*/
Tao_SSLS_FunctionGradient(TaoLineSearch ls,Vec X,PetscReal * fcn,Vec G,void * ptr)42 PetscErrorCode Tao_SSLS_FunctionGradient(TaoLineSearch ls, Vec X, PetscReal *fcn, Vec G, void *ptr)
43 {
44 Tao tao = (Tao)ptr;
45 TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
46 PetscErrorCode ierr;
47
48 PetscFunctionBegin;
49 ierr = TaoComputeConstraints(tao, X, tao->constraints);CHKERRQ(ierr);
50 ierr = VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff);CHKERRQ(ierr);
51 ierr = VecNorm(ssls->ff,NORM_2,&ssls->merit);CHKERRQ(ierr);
52 *fcn = 0.5*ssls->merit*ssls->merit;
53
54 ierr = TaoComputeJacobian(tao,tao->solution,tao->jacobian,tao->jacobian_pre);CHKERRQ(ierr);
55
56 ierr = MatDFischer(tao->jacobian, tao->solution, tao->constraints,tao->XL, tao->XU, ssls->t1, ssls->t2,ssls->da, ssls->db);CHKERRQ(ierr);
57 ierr = MatDiagonalScale(tao->jacobian,ssls->db,NULL);CHKERRQ(ierr);
58 ierr = MatDiagonalSet(tao->jacobian,ssls->da,ADD_VALUES);CHKERRQ(ierr);
59 ierr = MatMultTranspose(tao->jacobian,ssls->ff,G);CHKERRQ(ierr);
60 PetscFunctionReturn(0);
61 }
62
63