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