1 2 #include <../src/snes/impls/vi/ss/vissimpl.h> /*I "petscsnes.h" I*/ 3 4 /* 5 SNESVIComputeMeritFunction - Evaluates the merit function for the mixed complementarity problem. 6 7 Input Parameter: 8 . phi - the semismooth function 9 10 Output Parameter: 11 . merit - the merit function 12 . phinorm - ||phi|| 13 14 Notes: 15 The merit function for the mixed complementarity problem is defined as 16 merit = 0.5*phi^T*phi 17 */ 18 #undef __FUNCT__ 19 #define __FUNCT__ "SNESVIComputeMeritFunction" 20 static PetscErrorCode SNESVIComputeMeritFunction(Vec phi, PetscReal *merit,PetscReal *phinorm) 21 { 22 PetscErrorCode ierr; 23 24 PetscFunctionBegin; 25 ierr = VecNormBegin(phi,NORM_2,phinorm);CHKERRQ(ierr); 26 ierr = VecNormEnd(phi,NORM_2,phinorm);CHKERRQ(ierr); 27 28 *merit = 0.5*(*phinorm)*(*phinorm); 29 PetscFunctionReturn(0); 30 } 31 32 PETSC_STATIC_INLINE PetscScalar Phi(PetscScalar a,PetscScalar b) 33 { 34 return a + b - PetscSqrtScalar(a*a + b*b); 35 } 36 37 PETSC_STATIC_INLINE PetscScalar DPhi(PetscScalar a,PetscScalar b) 38 { 39 if ((PetscAbsScalar(a) >= 1.e-6) || (PetscAbsScalar(b) >= 1.e-6)) return 1.0 - a/ PetscSqrtScalar(a*a + b*b); 40 else return .5; 41 } 42 43 /* 44 SNESVIComputeFunction - Reformulates a system of nonlinear equations in mixed complementarity form to a system of nonlinear equations in semismooth form. 45 46 Input Parameters: 47 . snes - the SNES context 48 . X - current iterate 49 . functx - user defined function context 50 51 Output Parameters: 52 . phi - Semismooth function 53 54 */ 55 #undef __FUNCT__ 56 #define __FUNCT__ "SNESVIComputeFunction" 57 static PetscErrorCode SNESVIComputeFunction(SNES snes,Vec X,Vec phi,void *functx) 58 { 59 PetscErrorCode ierr; 60 SNES_VINEWTONSSLS *vi = (SNES_VINEWTONSSLS*)snes->data; 61 Vec Xl = snes->xl,Xu = snes->xu,F = snes->vec_func; 62 PetscScalar *phi_arr,*f_arr,*l,*u; 63 const PetscScalar *x_arr; 64 PetscInt i,nlocal; 65 66 PetscFunctionBegin; 67 ierr = (*vi->computeuserfunction)(snes,X,F,functx);CHKERRQ(ierr); 68 ierr = VecGetLocalSize(X,&nlocal);CHKERRQ(ierr); 69 ierr = VecGetArrayRead(X,&x_arr);CHKERRQ(ierr); 70 ierr = VecGetArray(F,&f_arr);CHKERRQ(ierr); 71 ierr = VecGetArray(Xl,&l);CHKERRQ(ierr); 72 ierr = VecGetArray(Xu,&u);CHKERRQ(ierr); 73 ierr = VecGetArray(phi,&phi_arr);CHKERRQ(ierr); 74 75 for (i=0; i < nlocal; i++) { 76 if ((PetscRealPart(l[i]) <= PETSC_NINFINITY) && (PetscRealPart(u[i]) >= PETSC_INFINITY)) { /* no constraints on variable */ 77 phi_arr[i] = f_arr[i]; 78 } else if (PetscRealPart(l[i]) <= PETSC_NINFINITY) { /* upper bound on variable only */ 79 phi_arr[i] = -Phi(u[i] - x_arr[i],-f_arr[i]); 80 } else if (PetscRealPart(u[i]) >= PETSC_INFINITY) { /* lower bound on variable only */ 81 phi_arr[i] = Phi(x_arr[i] - l[i],f_arr[i]); 82 } else if (l[i] == u[i]) { 83 phi_arr[i] = l[i] - x_arr[i]; 84 } else { /* both bounds on variable */ 85 phi_arr[i] = Phi(x_arr[i] - l[i],-Phi(u[i] - x_arr[i],-f_arr[i])); 86 } 87 } 88 89 ierr = VecRestoreArrayRead(X,&x_arr);CHKERRQ(ierr); 90 ierr = VecRestoreArray(F,&f_arr);CHKERRQ(ierr); 91 ierr = VecRestoreArray(Xl,&l);CHKERRQ(ierr); 92 ierr = VecRestoreArray(Xu,&u);CHKERRQ(ierr); 93 ierr = VecRestoreArray(phi,&phi_arr);CHKERRQ(ierr); 94 PetscFunctionReturn(0); 95 } 96 97 /* 98 SNESVIComputeBsubdifferentialVectors - Computes the diagonal shift (Da) and row scaling (Db) vectors needed for the 99 the semismooth jacobian. 100 */ 101 #undef __FUNCT__ 102 #define __FUNCT__ "SNESVIComputeBsubdifferentialVectors" 103 PetscErrorCode SNESVIComputeBsubdifferentialVectors(SNES snes,Vec X,Vec F,Mat jac,Vec Da,Vec Db) 104 { 105 PetscErrorCode ierr; 106 PetscScalar *l,*u,*x,*f,*da,*db,da1,da2,db1,db2; 107 PetscInt i,nlocal; 108 109 PetscFunctionBegin; 110 ierr = VecGetArray(X,&x);CHKERRQ(ierr); 111 ierr = VecGetArray(F,&f);CHKERRQ(ierr); 112 ierr = VecGetArray(snes->xl,&l);CHKERRQ(ierr); 113 ierr = VecGetArray(snes->xu,&u);CHKERRQ(ierr); 114 ierr = VecGetArray(Da,&da);CHKERRQ(ierr); 115 ierr = VecGetArray(Db,&db);CHKERRQ(ierr); 116 ierr = VecGetLocalSize(X,&nlocal);CHKERRQ(ierr); 117 118 for (i=0; i< nlocal; i++) { 119 if ((PetscRealPart(l[i]) <= PETSC_NINFINITY) && (PetscRealPart(u[i]) >= PETSC_INFINITY)) { /* no constraints on variable */ 120 da[i] = 0; 121 db[i] = 1; 122 } else if (PetscRealPart(l[i]) <= PETSC_NINFINITY) { /* upper bound on variable only */ 123 da[i] = DPhi(u[i] - x[i], -f[i]); 124 db[i] = DPhi(-f[i],u[i] - x[i]); 125 } else if (PetscRealPart(u[i]) >= PETSC_INFINITY) { /* lower bound on variable only */ 126 da[i] = DPhi(x[i] - l[i], f[i]); 127 db[i] = DPhi(f[i],x[i] - l[i]); 128 } else if (l[i] == u[i]) { /* fixed variable */ 129 da[i] = 1; 130 db[i] = 0; 131 } else { /* upper and lower bounds on variable */ 132 da1 = DPhi(x[i] - l[i], -Phi(u[i] - x[i], -f[i])); 133 db1 = DPhi(-Phi(u[i] - x[i], -f[i]),x[i] - l[i]); 134 da2 = DPhi(u[i] - x[i], -f[i]); 135 db2 = DPhi(-f[i],u[i] - x[i]); 136 da[i] = da1 + db1*da2; 137 db[i] = db1*db2; 138 } 139 } 140 141 ierr = VecRestoreArray(X,&x);CHKERRQ(ierr); 142 ierr = VecRestoreArray(F,&f);CHKERRQ(ierr); 143 ierr = VecRestoreArray(snes->xl,&l);CHKERRQ(ierr); 144 ierr = VecRestoreArray(snes->xu,&u);CHKERRQ(ierr); 145 ierr = VecRestoreArray(Da,&da);CHKERRQ(ierr); 146 ierr = VecRestoreArray(Db,&db);CHKERRQ(ierr); 147 PetscFunctionReturn(0); 148 } 149 150 /* 151 SNESVIComputeJacobian - Computes the jacobian of the semismooth function.The Jacobian for the semismooth function is an element of the B-subdifferential of the Fischer-Burmeister function for complementarity problems. 152 153 Input Parameters: 154 . Da - Diagonal shift vector for the semismooth jacobian. 155 . Db - Row scaling vector for the semismooth jacobian. 156 157 Output Parameters: 158 . jac - semismooth jacobian 159 . jac_pre - optional preconditioning matrix 160 161 Notes: 162 The semismooth jacobian matrix is given by 163 jac = Da + Db*jacfun 164 where Db is the row scaling matrix stored as a vector, 165 Da is the diagonal perturbation matrix stored as a vector 166 and jacfun is the jacobian of the original nonlinear function. 167 */ 168 #undef __FUNCT__ 169 #define __FUNCT__ "SNESVIComputeJacobian" 170 PetscErrorCode SNESVIComputeJacobian(Mat jac, Mat jac_pre,Vec Da, Vec Db) 171 { 172 PetscErrorCode ierr; 173 174 /* Do row scaling and add diagonal perturbation */ 175 ierr = MatDiagonalScale(jac,Db,NULL);CHKERRQ(ierr); 176 ierr = MatDiagonalSet(jac,Da,ADD_VALUES);CHKERRQ(ierr); 177 if (jac != jac_pre) { /* If jac and jac_pre are different */ 178 ierr = MatDiagonalScale(jac_pre,Db,NULL);CHKERRQ(ierr); 179 ierr = MatDiagonalSet(jac_pre,Da,ADD_VALUES);CHKERRQ(ierr); 180 } 181 PetscFunctionReturn(0); 182 } 183 184 /* 185 SNESVIComputeMeritFunctionGradient - Computes the gradient of the merit function psi. 186 187 Input Parameters: 188 phi - semismooth function. 189 H - semismooth jacobian 190 191 Output Parameters: 192 dpsi - merit function gradient 193 194 Notes: 195 The merit function gradient is computed as follows 196 dpsi = H^T*phi 197 */ 198 #undef __FUNCT__ 199 #define __FUNCT__ "SNESVIComputeMeritFunctionGradient" 200 PetscErrorCode SNESVIComputeMeritFunctionGradient(Mat H, Vec phi, Vec dpsi) 201 { 202 PetscErrorCode ierr; 203 204 PetscFunctionBegin; 205 ierr = MatMultTranspose(H,phi,dpsi);CHKERRQ(ierr); 206 PetscFunctionReturn(0); 207 } 208 209 210 211 /* 212 SNESSolve_VINEWTONSSLS - Solves the complementarity problem with a semismooth Newton 213 method using a line search. 214 215 Input Parameters: 216 . snes - the SNES context 217 218 Application Interface Routine: SNESSolve() 219 220 Notes: 221 This implements essentially a semismooth Newton method with a 222 line search. The default line search does not do any line search 223 but rather takes a full Newton step. 224 225 Developer Note: the code in this file should be slightly modified so that this routine need not exist and the SNESSolve_NEWTONLS() routine is called directly with the appropriate wrapped function and Jacobian evaluations 226 227 */ 228 #undef __FUNCT__ 229 #define __FUNCT__ "SNESSolve_VINEWTONSSLS" 230 PetscErrorCode SNESSolve_VINEWTONSSLS(SNES snes) 231 { 232 SNES_VINEWTONSSLS *vi = (SNES_VINEWTONSSLS*)snes->data; 233 PetscErrorCode ierr; 234 PetscInt maxits,i,lits; 235 SNESLineSearchReason lssucceed; 236 PetscReal gnorm,xnorm=0,ynorm; 237 Vec Y,X,F; 238 KSPConvergedReason kspreason; 239 DM dm; 240 DMSNES sdm; 241 242 PetscFunctionBegin; 243 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 244 ierr = DMGetDMSNES(dm,&sdm);CHKERRQ(ierr); 245 246 vi->computeuserfunction = sdm->ops->computefunction; 247 sdm->ops->computefunction = SNESVIComputeFunction; 248 249 snes->numFailures = 0; 250 snes->numLinearSolveFailures = 0; 251 snes->reason = SNES_CONVERGED_ITERATING; 252 253 maxits = snes->max_its; /* maximum number of iterations */ 254 X = snes->vec_sol; /* solution vector */ 255 F = snes->vec_func; /* residual vector */ 256 Y = snes->work[0]; /* work vectors */ 257 258 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 259 snes->iter = 0; 260 snes->norm = 0.0; 261 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 262 263 ierr = SNESVIProjectOntoBounds(snes,X);CHKERRQ(ierr); 264 ierr = SNESComputeFunction(snes,X,vi->phi);CHKERRQ(ierr); 265 if (snes->domainerror) { 266 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 267 sdm->ops->computefunction = vi->computeuserfunction; 268 PetscFunctionReturn(0); 269 } 270 /* Compute Merit function */ 271 ierr = SNESVIComputeMeritFunction(vi->phi,&vi->merit,&vi->phinorm);CHKERRQ(ierr); 272 273 ierr = VecNormBegin(X,NORM_2,&xnorm);CHKERRQ(ierr); /* xnorm <- ||x|| */ 274 ierr = VecNormEnd(X,NORM_2,&xnorm);CHKERRQ(ierr); 275 SNESCheckFunctionNorm(snes,vi->merit); 276 277 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 278 snes->norm = vi->phinorm; 279 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 280 ierr = SNESLogConvergenceHistory(snes,vi->phinorm,0);CHKERRQ(ierr); 281 ierr = SNESMonitor(snes,0,vi->phinorm);CHKERRQ(ierr); 282 283 /* test convergence */ 284 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,vi->phinorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 285 if (snes->reason) { 286 sdm->ops->computefunction = vi->computeuserfunction; 287 PetscFunctionReturn(0); 288 } 289 290 for (i=0; i<maxits; i++) { 291 292 /* Call general purpose update function */ 293 if (snes->ops->update) { 294 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 295 } 296 297 /* Solve J Y = Phi, where J is the semismooth jacobian */ 298 299 /* Get the jacobian -- note that the function must be the original function for snes_fd and snes_fd_color to work for this*/ 300 sdm->ops->computefunction = vi->computeuserfunction; 301 ierr = SNESComputeJacobian(snes,X,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 302 sdm->ops->computefunction = SNESVIComputeFunction; 303 304 /* Get the diagonal shift and row scaling vectors */ 305 ierr = SNESVIComputeBsubdifferentialVectors(snes,X,F,snes->jacobian,vi->Da,vi->Db);CHKERRQ(ierr); 306 /* Compute the semismooth jacobian */ 307 ierr = SNESVIComputeJacobian(snes->jacobian,snes->jacobian_pre,vi->Da,vi->Db);CHKERRQ(ierr); 308 /* Compute the merit function gradient */ 309 ierr = SNESVIComputeMeritFunctionGradient(snes->jacobian,vi->phi,vi->dpsi);CHKERRQ(ierr); 310 ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 311 ierr = KSPSolve(snes->ksp,vi->phi,Y);CHKERRQ(ierr); 312 ierr = KSPGetConvergedReason(snes->ksp,&kspreason);CHKERRQ(ierr); 313 314 if (kspreason < 0) { 315 if (++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) { 316 ierr = PetscInfo2(snes,"iter=%D, number linear solve failures %D greater than current SNES allowed, stopping solve\n",snes->iter,snes->numLinearSolveFailures);CHKERRQ(ierr); 317 snes->reason = SNES_DIVERGED_LINEAR_SOLVE; 318 break; 319 } 320 } 321 ierr = KSPGetIterationNumber(snes->ksp,&lits);CHKERRQ(ierr); 322 snes->linear_its += lits; 323 ierr = PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);CHKERRQ(ierr); 324 /* 325 if (snes->ops->precheck) { 326 PetscBool changed_y = PETSC_FALSE; 327 ierr = (*snes->ops->precheck)(snes,X,Y,snes->precheck,&changed_y);CHKERRQ(ierr); 328 } 329 330 if (PetscLogPrintInfo) { 331 ierr = SNESVICheckResidual_Private(snes,snes->jacobian,F,Y,G,W);CHKERRQ(ierr); 332 } 333 */ 334 /* Compute a (scaled) negative update in the line search routine: 335 Y <- X - lambda*Y 336 and evaluate G = function(Y) (depends on the line search). 337 */ 338 ierr = VecCopy(Y,snes->vec_sol_update);CHKERRQ(ierr); 339 ynorm = 1; gnorm = vi->phinorm; 340 ierr = SNESLineSearchApply(snes->linesearch, X, vi->phi, &gnorm, Y);CHKERRQ(ierr); 341 ierr = SNESLineSearchGetReason(snes->linesearch, &lssucceed);CHKERRQ(ierr); 342 ierr = SNESLineSearchGetNorms(snes->linesearch, &xnorm, &gnorm, &ynorm);CHKERRQ(ierr); 343 ierr = PetscInfo4(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n",(double)vi->phinorm,(double)gnorm,(double)ynorm,(int)lssucceed);CHKERRQ(ierr); 344 if (snes->reason == SNES_DIVERGED_FUNCTION_COUNT) break; 345 if (snes->domainerror) { 346 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 347 sdm->ops->computefunction = vi->computeuserfunction; 348 PetscFunctionReturn(0); 349 } 350 if (lssucceed) { 351 if (++snes->numFailures >= snes->maxFailures) { 352 PetscBool ismin; 353 snes->reason = SNES_DIVERGED_LINE_SEARCH; 354 ierr = SNESVICheckLocalMin_Private(snes,snes->jacobian,vi->phi,X,gnorm,&ismin);CHKERRQ(ierr); 355 if (ismin) snes->reason = SNES_DIVERGED_LOCAL_MIN; 356 break; 357 } 358 } 359 /* Update function and solution vectors */ 360 vi->phinorm = gnorm; 361 vi->merit = 0.5*vi->phinorm*vi->phinorm; 362 /* Monitor convergence */ 363 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 364 snes->iter = i+1; 365 snes->norm = vi->phinorm; 366 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 367 ierr = SNESLogConvergenceHistory(snes,snes->norm,lits);CHKERRQ(ierr); 368 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 369 /* Test for convergence, xnorm = || X || */ 370 if (snes->ops->converged != SNESConvergedSkip) { ierr = VecNorm(X,NORM_2,&xnorm);CHKERRQ(ierr); } 371 ierr = (*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,vi->phinorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 372 if (snes->reason) break; 373 } 374 if (i == maxits) { 375 ierr = PetscInfo1(snes,"Maximum number of iterations has been reached: %D\n",maxits);CHKERRQ(ierr); 376 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 377 } 378 sdm->ops->computefunction = vi->computeuserfunction; 379 PetscFunctionReturn(0); 380 } 381 382 /* -------------------------------------------------------------------------- */ 383 /* 384 SNESSetUp_VINEWTONSSLS - Sets up the internal data structures for the later use 385 of the SNES nonlinear solver. 386 387 Input Parameter: 388 . snes - the SNES context 389 390 Application Interface Routine: SNESSetUp() 391 392 Notes: 393 For basic use of the SNES solvers, the user need not explicitly call 394 SNESSetUp(), since these actions will automatically occur during 395 the call to SNESSolve(). 396 */ 397 #undef __FUNCT__ 398 #define __FUNCT__ "SNESSetUp_VINEWTONSSLS" 399 PetscErrorCode SNESSetUp_VINEWTONSSLS(SNES snes) 400 { 401 PetscErrorCode ierr; 402 SNES_VINEWTONSSLS *vi = (SNES_VINEWTONSSLS*) snes->data; 403 404 PetscFunctionBegin; 405 ierr = SNESSetUp_VI(snes);CHKERRQ(ierr); 406 ierr = VecDuplicate(snes->vec_sol, &vi->dpsi);CHKERRQ(ierr); 407 ierr = VecDuplicate(snes->vec_sol, &vi->phi);CHKERRQ(ierr); 408 ierr = VecDuplicate(snes->vec_sol, &vi->Da);CHKERRQ(ierr); 409 ierr = VecDuplicate(snes->vec_sol, &vi->Db);CHKERRQ(ierr); 410 ierr = VecDuplicate(snes->vec_sol, &vi->z);CHKERRQ(ierr); 411 ierr = VecDuplicate(snes->vec_sol, &vi->t);CHKERRQ(ierr); 412 PetscFunctionReturn(0); 413 } 414 /* -------------------------------------------------------------------------- */ 415 #undef __FUNCT__ 416 #define __FUNCT__ "SNESReset_VINEWTONSSLS" 417 PetscErrorCode SNESReset_VINEWTONSSLS(SNES snes) 418 { 419 SNES_VINEWTONSSLS *vi = (SNES_VINEWTONSSLS*) snes->data; 420 PetscErrorCode ierr; 421 422 PetscFunctionBegin; 423 ierr = SNESReset_VI(snes);CHKERRQ(ierr); 424 ierr = VecDestroy(&vi->dpsi);CHKERRQ(ierr); 425 ierr = VecDestroy(&vi->phi);CHKERRQ(ierr); 426 ierr = VecDestroy(&vi->Da);CHKERRQ(ierr); 427 ierr = VecDestroy(&vi->Db);CHKERRQ(ierr); 428 ierr = VecDestroy(&vi->z);CHKERRQ(ierr); 429 ierr = VecDestroy(&vi->t);CHKERRQ(ierr); 430 PetscFunctionReturn(0); 431 } 432 433 /* -------------------------------------------------------------------------- */ 434 /* 435 SNESSetFromOptions_VINEWTONSSLS - Sets various parameters for the SNESVI method. 436 437 Input Parameter: 438 . snes - the SNES context 439 440 Application Interface Routine: SNESSetFromOptions() 441 */ 442 #undef __FUNCT__ 443 #define __FUNCT__ "SNESSetFromOptions_VINEWTONSSLS" 444 static PetscErrorCode SNESSetFromOptions_VINEWTONSSLS(PetscOptionItems *PetscOptionsObject,SNES snes) 445 { 446 PetscErrorCode ierr; 447 SNESLineSearch linesearch; 448 449 PetscFunctionBegin; 450 ierr = SNESSetFromOptions_VI(PetscOptionsObject,snes);CHKERRQ(ierr); 451 ierr = PetscOptionsHead(PetscOptionsObject,"SNES semismooth method options");CHKERRQ(ierr); 452 ierr = PetscOptionsTail();CHKERRQ(ierr); 453 /* set up the default line search */ 454 if (!snes->linesearch) { 455 ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr); 456 ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHBT);CHKERRQ(ierr); 457 ierr = SNESLineSearchBTSetAlpha(linesearch, 0.0);CHKERRQ(ierr); 458 } 459 PetscFunctionReturn(0); 460 } 461 462 463 /* -------------------------------------------------------------------------- */ 464 /*MC 465 SNESVINEWTONSSLS - Semi-smooth solver for variational inequalities based on Newton's method 466 467 Options Database: 468 + -snes_type <vinewtonssls,vinewtonrsls> a semi-smooth solver, a reduced space active set method 469 - -snes_vi_monitor - prints the number of active constraints at each iteration. 470 471 Level: beginner 472 473 References: 474 + 1. - T. S. Munson, F. Facchinei, M. C. Ferris, A. Fischer, and C. Kanzow. The semismooth 475 algorithm for large scale complementarity problems. INFORMS Journal on Computing, 13 (2001). 476 - 2. - T. S. Munson, and S. Benson. Flexible Complementarity Solvers for Large Scale 477 Applications, Optimization Methods and Software, 21 (2006). 478 479 .seealso: SNESVISetVariableBounds(), SNESVISetComputeVariableBounds(), SNESCreate(), SNES, SNESSetType(), SNESVINEWTONRSLS, SNESNEWTONTR, SNESLineSearchSetType(),SNESLineSearchSetPostCheck(), SNESLineSearchSetPreCheck() 480 481 M*/ 482 #undef __FUNCT__ 483 #define __FUNCT__ "SNESCreate_VINEWTONSSLS" 484 PETSC_EXTERN PetscErrorCode SNESCreate_VINEWTONSSLS(SNES snes) 485 { 486 PetscErrorCode ierr; 487 SNES_VINEWTONSSLS *vi; 488 489 PetscFunctionBegin; 490 snes->ops->reset = SNESReset_VINEWTONSSLS; 491 snes->ops->setup = SNESSetUp_VINEWTONSSLS; 492 snes->ops->solve = SNESSolve_VINEWTONSSLS; 493 snes->ops->destroy = SNESDestroy_VI; 494 snes->ops->setfromoptions = SNESSetFromOptions_VINEWTONSSLS; 495 snes->ops->view = NULL; 496 497 snes->usesksp = PETSC_TRUE; 498 snes->usespc = PETSC_FALSE; 499 500 snes->alwayscomputesfinalresidual = PETSC_FALSE; 501 502 ierr = PetscNewLog(snes,&vi);CHKERRQ(ierr); 503 snes->data = (void*)vi; 504 505 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESVISetVariableBounds_C",SNESVISetVariableBounds_VI);CHKERRQ(ierr); 506 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESVISetComputeVariableBounds_C",SNESVISetComputeVariableBounds_VI);CHKERRQ(ierr); 507 PetscFunctionReturn(0); 508 } 509 510