Actual source code: ssls.c

  1: #include "ssls.h"


  4: /*------------------------------------------------------------*/
  7: PetscErrorCode TaoSetFromOptions_SSLS(TaoSolver tao)
  8: {
  9:   TAO_SSLS *ssls = (TAO_SSLS *)tao->data;
 11:   PetscBool flg;

 14:   PetscOptionsHead("Semismooth method with a linesearch for "
 15:                           "complementarity problems"); 
 16:   PetscOptionsReal("-ssls_delta", "descent test fraction", "",
 17:                          ssls->delta, &(ssls->delta), &flg);
 18:   PetscOptionsReal("-ssls_rho", "descent test power", "",
 19:                          ssls->rho, &(ssls->rho), &flg);
 20:   TaoLineSearchSetFromOptions(tao->linesearch);
 21:   PetscOptionsTail(); 
 22:   return(0);
 23: }

 25: /*------------------------------------------------------------*/
 28: PetscErrorCode TaoView_SSLS(TaoSolver tao, PetscViewer pv)
 29: {

 33:   return(0);
 34: }

 36: /*------------------------------------------------------------*/
 39: PetscErrorCode Tao_SSLS_Function(TaoLineSearch ls, Vec X, PetscReal *fcn, void *ptr) 
 40: {
 41:   TaoSolver tao = (TaoSolver)ptr;
 42:   TAO_SSLS *ssls = (TAO_SSLS *)tao->data;

 46:   
 47:   TaoComputeConstraints(tao, X, tao->constraints); 
 48:   VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff); 
 49:   VecNorm(ssls->ff,NORM_2,&ssls->merit); 
 50:   *fcn = 0.5*ssls->merit*ssls->merit;
 51:   return(0);
 52: }

 54: /*------------------------------------------------------------*/
 57: PetscErrorCode Tao_SSLS_FunctionGradient(TaoLineSearch ls, Vec X, PetscReal *fcn,  Vec G, void *ptr)
 58: {
 59:   TaoSolver tao = (TaoSolver)ptr;
 60:   TAO_SSLS *ssls = (TAO_SSLS *)tao->data;


 65:   TaoComputeConstraints(tao, X, tao->constraints); 
 66:   VecFischer(X,tao->constraints,tao->XL,tao->XU,ssls->ff); 
 67:   VecNorm(ssls->ff,NORM_2,&ssls->merit); 
 68:   *fcn = 0.5*ssls->merit*ssls->merit;

 70:   TaoComputeJacobian(tao, tao->solution, &tao->jacobian, &tao->jacobian_pre, &ssls->matflag); 
 71:   
 72:   D_Fischer(tao->jacobian, tao->solution, tao->constraints, 
 73:                    tao->XL, tao->XU, ssls->t1, ssls->t2, 
 74:                    ssls->da, ssls->db); 
 75:   MatDiagonalScale(tao->jacobian,ssls->db,PETSC_NULL); 
 76:   MatDiagonalSet(tao->jacobian,ssls->da,ADD_VALUES); 
 77:   MatMultTranspose(tao->jacobian,ssls->ff,G); 

 79:   return(0);
 80: }