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