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