Actual source code: ssls.c

petsc-3.5.4 2015-05-23
Report Typos and Errors
  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: {
 28:   return(0);
 29: }

 31: /*------------------------------------------------------------*/
 34: PetscErrorCode Tao_SSLS_Function(TaoLineSearch ls, Vec X, PetscReal *fcn, void *ptr)
 35: {
 36:   Tao            tao = (Tao)ptr;
 37:   TAO_SSLS       *ssls = (TAO_SSLS *)tao->data;

 41:   TaoComputeConstraints(tao, X, tao->constraints);
 42:   VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff);
 43:   VecNorm(ssls->ff,NORM_2,&ssls->merit);
 44:   *fcn = 0.5*ssls->merit*ssls->merit;
 45:   return(0);
 46: }

 48: /*------------------------------------------------------------*/
 51: PetscErrorCode Tao_SSLS_FunctionGradient(TaoLineSearch ls, Vec X, PetscReal *fcn,  Vec G, void *ptr)
 52: {
 53:   Tao            tao = (Tao)ptr;
 54:   TAO_SSLS       *ssls = (TAO_SSLS *)tao->data;

 58:   TaoComputeConstraints(tao, X, tao->constraints);
 59:   VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff);
 60:   VecNorm(ssls->ff,NORM_2,&ssls->merit);
 61:   *fcn = 0.5*ssls->merit*ssls->merit;

 63:   TaoComputeJacobian(tao,tao->solution,tao->jacobian,tao->jacobian_pre);

 65:   MatDFischer(tao->jacobian, tao->solution, tao->constraints,tao->XL, tao->XU, ssls->t1, ssls->t2,ssls->da, ssls->db);
 66:   MatDiagonalScale(tao->jacobian,ssls->db,NULL);
 67:   MatDiagonalSet(tao->jacobian,ssls->da,ADD_VALUES);
 68:   MatMultTranspose(tao->jacobian,ssls->ff,G);
 69:   return(0);
 70: }