1 /* Program usage: mpiexec -n 1 minsurf1 [-help] [all TAO options] */ 2 3 /* Include "petsctao.h" so we can use TAO solvers. */ 4 #include <petsctao.h> 5 6 static char help[] = 7 "This example demonstrates use of the TAO package to \n\ 8 solve an unconstrained minimization problem. This example is based on a \n\ 9 problem from the MINPACK-2 test suite. Given a rectangular 2-D domain and \n\ 10 boundary values along the edges of the domain, the objective is to find the\n\ 11 surface with the minimal area that satisfies the boundary conditions.\n\ 12 The command line options are:\n\ 13 -mx <xg>, where <xg> = number of grid points in the 1st coordinate direction\n\ 14 -my <yg>, where <yg> = number of grid points in the 2nd coordinate direction\n\ 15 -start <st>, where <st> =0 for zero vector, and an average of the boundary conditions otherwise \n\n"; 16 17 /* 18 User-defined application context - contains data needed by the 19 application-provided call-back routines, FormFunctionGradient() 20 and FormHessian(). 21 */ 22 typedef struct { 23 PetscInt mx, my; /* discretization in x, y directions */ 24 PetscReal *bottom, *top, *left, *right; /* boundary values */ 25 Mat H; 26 } AppCtx; 27 28 /* -------- User-defined Routines --------- */ 29 30 static PetscErrorCode MSA_BoundaryConditions(AppCtx*); 31 static PetscErrorCode MSA_InitialPoint(AppCtx*,Vec); 32 static PetscErrorCode QuadraticH(AppCtx*,Vec,Mat); 33 PetscErrorCode FormFunctionGradient(Tao,Vec,PetscReal*,Vec,void*); 34 PetscErrorCode FormHessian(Tao,Vec,Mat,Mat,void*); 35 36 int main(int argc, char **argv) 37 { 38 PetscInt N; /* Size of vector */ 39 PetscMPIInt size; /* Number of processors */ 40 Vec x; /* solution, gradient vectors */ 41 KSP ksp; /* PETSc Krylov subspace method */ 42 PetscBool flg; /* A return value when checking for user options */ 43 Tao tao; /* Tao solver context */ 44 AppCtx user; /* user-defined work context */ 45 46 /* Initialize TAO,PETSc */ 47 PetscCall(PetscInitialize(&argc, &argv,(char *)0,help)); 48 49 PetscCallMPI(MPI_Comm_size(MPI_COMM_WORLD,&size)); 50 PetscCheck(size == 1,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"Incorrect number of processors"); 51 52 /* Specify default dimension of the problem */ 53 user.mx = 4; user.my = 4; 54 55 /* Check for any command line arguments that override defaults */ 56 PetscCall(PetscOptionsGetInt(NULL,NULL,"-mx",&user.mx,&flg)); 57 PetscCall(PetscOptionsGetInt(NULL,NULL,"-my",&user.my,&flg)); 58 59 PetscCall(PetscPrintf(PETSC_COMM_SELF,"\n---- Minimum Surface Area Problem -----\n")); 60 PetscCall(PetscPrintf(PETSC_COMM_SELF,"mx: %" PetscInt_FMT " my: %" PetscInt_FMT " \n\n",user.mx,user.my)); 61 62 /* Calculate any derived values from parameters */ 63 N = user.mx*user.my; 64 65 /* Create TAO solver and set desired solution method */ 66 PetscCall(TaoCreate(PETSC_COMM_SELF,&tao)); 67 PetscCall(TaoSetType(tao,TAOLMVM)); 68 69 /* Initialize minsurf application data structure for use in the function evaluations */ 70 PetscCall(MSA_BoundaryConditions(&user)); /* Application specific routine */ 71 72 /* 73 Create a vector to hold the variables. Compute an initial solution. 74 Set this vector, which will be used by TAO. 75 */ 76 PetscCall(VecCreateSeq(PETSC_COMM_SELF,N,&x)); 77 PetscCall(MSA_InitialPoint(&user,x)); /* Application specific routine */ 78 PetscCall(TaoSetSolution(tao,x)); /* A TAO routine */ 79 80 /* Provide TAO routines for function, gradient, and Hessian evaluation */ 81 PetscCall(TaoSetObjectiveAndGradient(tao,NULL,FormFunctionGradient,(void *)&user)); 82 83 /* Create a matrix data structure to store the Hessian. This structure will be used by TAO */ 84 PetscCall(MatCreateSeqAIJ(PETSC_COMM_SELF,N,N,7,NULL,&(user.H))); 85 PetscCall(MatSetOption(user.H,MAT_SYMMETRIC,PETSC_TRUE)); 86 PetscCall(TaoSetHessian(tao,user.H,user.H,FormHessian,(void *)&user)); 87 88 /* Check for any TAO command line options */ 89 PetscCall(TaoSetFromOptions(tao)); 90 91 /* Limit the number of iterations in the KSP linear solver */ 92 PetscCall(TaoGetKSP(tao,&ksp)); 93 if (ksp) PetscCall(KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,user.mx*user.my)); 94 95 /* SOLVE THE APPLICATION */ 96 PetscCall(TaoSolve(tao)); 97 98 PetscCall(TaoDestroy(&tao)); 99 PetscCall(VecDestroy(&x)); 100 PetscCall(MatDestroy(&user.H)); 101 PetscCall(PetscFree(user.bottom)); 102 PetscCall(PetscFree(user.top)); 103 PetscCall(PetscFree(user.left)); 104 PetscCall(PetscFree(user.right)); 105 106 PetscCall(PetscFinalize()); 107 return 0; 108 } 109 110 /* -------------------------------------------------------------------- */ 111 112 /* FormFunctionGradient - Evaluates function and corresponding gradient. 113 114 Input Parameters: 115 . tao - the Tao context 116 . X - input vector 117 . userCtx - optional user-defined context, as set by TaoSetFunctionGradient() 118 119 Output Parameters: 120 . fcn - the newly evaluated function 121 . G - vector containing the newly evaluated gradient 122 */ 123 PetscErrorCode FormFunctionGradient(Tao tao,Vec X,PetscReal *fcn,Vec G,void *userCtx) 124 { 125 AppCtx *user = (AppCtx *) userCtx; 126 PetscInt i,j,row; 127 PetscInt mx=user->mx, my=user->my; 128 PetscReal rhx=mx+1, rhy=my+1; 129 PetscReal hx=1.0/(mx+1),hy=1.0/(my+1), hydhx=hy/hx, hxdhy=hx/hy, area=0.5*hx*hy, ft=0; 130 PetscReal f1,f2,f3,f4,f5,f6,d1,d2,d3,d4,d5,d6,d7,d8,xc,xl,xr,xt,xb,xlt,xrb; 131 PetscReal df1dxc,df2dxc,df3dxc,df4dxc,df5dxc,df6dxc; 132 PetscReal zero=0.0; 133 PetscScalar *g; 134 const PetscScalar *x; 135 136 PetscFunctionBeginUser; 137 PetscCall(VecSet(G, zero)); 138 139 PetscCall(VecGetArrayRead(X,&x)); 140 PetscCall(VecGetArray(G,&g)); 141 142 /* Compute function over the locally owned part of the mesh */ 143 for (j=0; j<my; j++) { 144 for (i=0; i< mx; i++) { 145 row=(j)*mx + (i); 146 xc = x[row]; 147 xlt=xrb=xl=xr=xb=xt=xc; 148 if (i==0) { /* left side */ 149 xl = user->left[j+1]; 150 xlt = user->left[j+2]; 151 } else { 152 xl = x[row-1]; 153 } 154 155 if (j==0) { /* bottom side */ 156 xb = user->bottom[i+1]; 157 xrb = user->bottom[i+2]; 158 } else { 159 xb = x[row-mx]; 160 } 161 162 if (i+1 == mx) { /* right side */ 163 xr = user->right[j+1]; 164 xrb = user->right[j]; 165 } else { 166 xr = x[row+1]; 167 } 168 169 if (j+1==0+my) { /* top side */ 170 xt = user->top[i+1]; 171 xlt = user->top[i]; 172 }else { 173 xt = x[row+mx]; 174 } 175 176 if (i>0 && j+1<my) { 177 xlt = x[row-1+mx]; 178 } 179 if (j>0 && i+1<mx) { 180 xrb = x[row+1-mx]; 181 } 182 183 d1 = (xc-xl); 184 d2 = (xc-xr); 185 d3 = (xc-xt); 186 d4 = (xc-xb); 187 d5 = (xr-xrb); 188 d6 = (xrb-xb); 189 d7 = (xlt-xl); 190 d8 = (xt-xlt); 191 192 df1dxc = d1*hydhx; 193 df2dxc = (d1*hydhx + d4*hxdhy); 194 df3dxc = d3*hxdhy; 195 df4dxc = (d2*hydhx + d3*hxdhy); 196 df5dxc = d2*hydhx; 197 df6dxc = d4*hxdhy; 198 199 d1 *= rhx; 200 d2 *= rhx; 201 d3 *= rhy; 202 d4 *= rhy; 203 d5 *= rhy; 204 d6 *= rhx; 205 d7 *= rhy; 206 d8 *= rhx; 207 208 f1 = PetscSqrtScalar(1.0 + d1*d1 + d7*d7); 209 f2 = PetscSqrtScalar(1.0 + d1*d1 + d4*d4); 210 f3 = PetscSqrtScalar(1.0 + d3*d3 + d8*d8); 211 f4 = PetscSqrtScalar(1.0 + d3*d3 + d2*d2); 212 f5 = PetscSqrtScalar(1.0 + d2*d2 + d5*d5); 213 f6 = PetscSqrtScalar(1.0 + d4*d4 + d6*d6); 214 215 ft = ft + (f2 + f4); 216 217 df1dxc /= f1; 218 df2dxc /= f2; 219 df3dxc /= f3; 220 df4dxc /= f4; 221 df5dxc /= f5; 222 df6dxc /= f6; 223 224 g[row] = (df1dxc+df2dxc+df3dxc+df4dxc+df5dxc+df6dxc)/2.0; 225 } 226 } 227 228 for (j=0; j<my; j++) { /* left side */ 229 d3 = (user->left[j+1] - user->left[j+2])*rhy; 230 d2 = (user->left[j+1] - x[j*mx])*rhx; 231 ft = ft+PetscSqrtScalar(1.0 + d3*d3 + d2*d2); 232 } 233 234 for (i=0; i<mx; i++) { /* bottom */ 235 d2 = (user->bottom[i+1]-user->bottom[i+2])*rhx; 236 d3 = (user->bottom[i+1]-x[i])*rhy; 237 ft = ft+PetscSqrtScalar(1.0 + d3*d3 + d2*d2); 238 } 239 240 for (j=0; j< my; j++) { /* right side */ 241 d1 = (x[(j+1)*mx-1]-user->right[j+1])*rhx; 242 d4 = (user->right[j]-user->right[j+1])*rhy; 243 ft = ft+PetscSqrtScalar(1.0 + d1*d1 + d4*d4); 244 } 245 246 for (i=0; i<mx; i++) { /* top side */ 247 d1 = (x[(my-1)*mx + i] - user->top[i+1])*rhy; 248 d4 = (user->top[i+1] - user->top[i])*rhx; 249 ft = ft+PetscSqrtScalar(1.0 + d1*d1 + d4*d4); 250 } 251 252 /* Bottom left corner */ 253 d1 = (user->left[0]-user->left[1])*rhy; 254 d2 = (user->bottom[0]-user->bottom[1])*rhx; 255 ft += PetscSqrtScalar(1.0 + d1*d1 + d2*d2); 256 257 /* Top right corner */ 258 d1 = (user->right[my+1] - user->right[my])*rhy; 259 d2 = (user->top[mx+1] - user->top[mx])*rhx; 260 ft += PetscSqrtScalar(1.0 + d1*d1 + d2*d2); 261 262 (*fcn)=ft*area; 263 264 /* Restore vectors */ 265 PetscCall(VecRestoreArrayRead(X,&x)); 266 PetscCall(VecRestoreArray(G,&g)); 267 PetscCall(PetscLogFlops(67.0*mx*my)); 268 PetscFunctionReturn(0); 269 } 270 271 /* ------------------------------------------------------------------- */ 272 /* 273 FormHessian - Evaluates the Hessian matrix. 274 275 Input Parameters: 276 . tao - the Tao context 277 . x - input vector 278 . ptr - optional user-defined context, as set by TaoSetHessian() 279 280 Output Parameters: 281 . H - Hessian matrix 282 . Hpre - optionally different preconditioning matrix 283 . flg - flag indicating matrix structure 284 285 */ 286 PetscErrorCode FormHessian(Tao tao,Vec X,Mat H, Mat Hpre, void *ptr) 287 { 288 AppCtx *user = (AppCtx *) ptr; 289 290 PetscFunctionBeginUser; 291 /* Evaluate the Hessian entries*/ 292 PetscCall(QuadraticH(user,X,H)); 293 PetscFunctionReturn(0); 294 } 295 296 /* ------------------------------------------------------------------- */ 297 /* 298 QuadraticH - Evaluates the Hessian matrix. 299 300 Input Parameters: 301 . user - user-defined context, as set by TaoSetHessian() 302 . X - input vector 303 304 Output Parameter: 305 . H - Hessian matrix 306 */ 307 PetscErrorCode QuadraticH(AppCtx *user, Vec X, Mat Hessian) 308 { 309 PetscInt i,j,k,row; 310 PetscInt mx=user->mx, my=user->my; 311 PetscInt col[7]; 312 PetscReal hx=1.0/(mx+1), hy=1.0/(my+1), hydhx=hy/hx, hxdhy=hx/hy; 313 PetscReal rhx=mx+1, rhy=my+1; 314 PetscReal f1,f2,f3,f4,f5,f6,d1,d2,d3,d4,d5,d6,d7,d8,xc,xl,xr,xt,xb,xlt,xrb; 315 PetscReal hl,hr,ht,hb,hc,htl,hbr; 316 const PetscScalar *x; 317 PetscReal v[7]; 318 319 PetscFunctionBeginUser; 320 /* Get pointers to vector data */ 321 PetscCall(VecGetArrayRead(X,&x)); 322 323 /* Initialize matrix entries to zero */ 324 PetscCall(MatZeroEntries(Hessian)); 325 326 /* Set various matrix options */ 327 PetscCall(MatSetOption(Hessian,MAT_IGNORE_OFF_PROC_ENTRIES,PETSC_TRUE)); 328 329 /* Compute Hessian over the locally owned part of the mesh */ 330 for (i=0; i< mx; i++) { 331 for (j=0; j<my; j++) { 332 333 row=(j)*mx + (i); 334 335 xc = x[row]; 336 xlt=xrb=xl=xr=xb=xt=xc; 337 338 /* Left side */ 339 if (i==0) { 340 xl= user->left[j+1]; 341 xlt = user->left[j+2]; 342 } else { 343 xl = x[row-1]; 344 } 345 346 if (j==0) { 347 xb=user->bottom[i+1]; 348 xrb = user->bottom[i+2]; 349 } else { 350 xb = x[row-mx]; 351 } 352 353 if (i+1 == mx) { 354 xr=user->right[j+1]; 355 xrb = user->right[j]; 356 } else { 357 xr = x[row+1]; 358 } 359 360 if (j+1==my) { 361 xt=user->top[i+1]; 362 xlt = user->top[i]; 363 }else { 364 xt = x[row+mx]; 365 } 366 367 if (i>0 && j+1<my) { 368 xlt = x[row-1+mx]; 369 } 370 if (j>0 && i+1<mx) { 371 xrb = x[row+1-mx]; 372 } 373 374 d1 = (xc-xl)*rhx; 375 d2 = (xc-xr)*rhx; 376 d3 = (xc-xt)*rhy; 377 d4 = (xc-xb)*rhy; 378 d5 = (xrb-xr)*rhy; 379 d6 = (xrb-xb)*rhx; 380 d7 = (xlt-xl)*rhy; 381 d8 = (xlt-xt)*rhx; 382 383 f1 = PetscSqrtScalar(1.0 + d1*d1 + d7*d7); 384 f2 = PetscSqrtScalar(1.0 + d1*d1 + d4*d4); 385 f3 = PetscSqrtScalar(1.0 + d3*d3 + d8*d8); 386 f4 = PetscSqrtScalar(1.0 + d3*d3 + d2*d2); 387 f5 = PetscSqrtScalar(1.0 + d2*d2 + d5*d5); 388 f6 = PetscSqrtScalar(1.0 + d4*d4 + d6*d6); 389 390 hl = (-hydhx*(1.0+d7*d7)+d1*d7)/(f1*f1*f1)+(-hydhx*(1.0+d4*d4)+d1*d4)/(f2*f2*f2); 391 hr = (-hydhx*(1.0+d5*d5)+d2*d5)/(f5*f5*f5)+(-hydhx*(1.0+d3*d3)+d2*d3)/(f4*f4*f4); 392 ht = (-hxdhy*(1.0+d8*d8)+d3*d8)/(f3*f3*f3)+(-hxdhy*(1.0+d2*d2)+d2*d3)/(f4*f4*f4); 393 hb = (-hxdhy*(1.0+d6*d6)+d4*d6)/(f6*f6*f6)+(-hxdhy*(1.0+d1*d1)+d1*d4)/(f2*f2*f2); 394 395 hbr = -d2*d5/(f5*f5*f5) - d4*d6/(f6*f6*f6); 396 htl = -d1*d7/(f1*f1*f1) - d3*d8/(f3*f3*f3); 397 398 hc = hydhx*(1.0+d7*d7)/(f1*f1*f1) + hxdhy*(1.0+d8*d8)/(f3*f3*f3) + hydhx*(1.0+d5*d5)/(f5*f5*f5) + hxdhy*(1.0+d6*d6)/(f6*f6*f6) + 399 (hxdhy*(1.0+d1*d1)+hydhx*(1.0+d4*d4)-2*d1*d4)/(f2*f2*f2) + (hxdhy*(1.0+d2*d2)+hydhx*(1.0+d3*d3)-2*d2*d3)/(f4*f4*f4); 400 401 hl*=0.5; hr*=0.5; ht*=0.5; hb*=0.5; hbr*=0.5; htl*=0.5; hc*=0.5; 402 403 k=0; 404 if (j>0) { 405 v[k]=hb; col[k]=row - mx; k++; 406 } 407 408 if (j>0 && i < mx -1) { 409 v[k]=hbr; col[k]=row - mx+1; k++; 410 } 411 412 if (i>0) { 413 v[k]= hl; col[k]=row - 1; k++; 414 } 415 416 v[k]= hc; col[k]=row; k++; 417 418 if (i < mx-1) { 419 v[k]= hr; col[k]=row+1; k++; 420 } 421 422 if (i>0 && j < my-1) { 423 v[k]= htl; col[k] = row+mx-1; k++; 424 } 425 426 if (j < my-1) { 427 v[k]= ht; col[k] = row+mx; k++; 428 } 429 430 /* 431 Set matrix values using local numbering, which was defined 432 earlier, in the main routine. 433 */ 434 PetscCall(MatSetValues(Hessian,1,&row,k,col,v,INSERT_VALUES)); 435 } 436 } 437 438 /* Restore vectors */ 439 PetscCall(VecRestoreArrayRead(X,&x)); 440 441 /* Assemble the matrix */ 442 PetscCall(MatAssemblyBegin(Hessian,MAT_FINAL_ASSEMBLY)); 443 PetscCall(MatAssemblyEnd(Hessian,MAT_FINAL_ASSEMBLY)); 444 445 PetscCall(PetscLogFlops(199.0*mx*my)); 446 PetscFunctionReturn(0); 447 } 448 449 /* ------------------------------------------------------------------- */ 450 /* 451 MSA_BoundaryConditions - Calculates the boundary conditions for 452 the region. 453 454 Input Parameter: 455 . user - user-defined application context 456 457 Output Parameter: 458 . user - user-defined application context 459 */ 460 static PetscErrorCode MSA_BoundaryConditions(AppCtx * user) 461 { 462 PetscInt i,j,k,limit=0; 463 PetscInt maxits=5; 464 PetscInt mx=user->mx,my=user->my; 465 PetscInt bsize=0, lsize=0, tsize=0, rsize=0; 466 PetscReal one=1.0, two=2.0, three=3.0, tol=1e-10; 467 PetscReal fnorm,det,hx,hy,xt=0,yt=0; 468 PetscReal u1,u2,nf1,nf2,njac11,njac12,njac21,njac22; 469 PetscReal b=-0.5, t=0.5, l=-0.5, r=0.5; 470 PetscReal *boundary; 471 472 PetscFunctionBeginUser; 473 bsize=mx+2; lsize=my+2; rsize=my+2; tsize=mx+2; 474 475 PetscCall(PetscMalloc1(bsize,&user->bottom)); 476 PetscCall(PetscMalloc1(tsize,&user->top)); 477 PetscCall(PetscMalloc1(lsize,&user->left)); 478 PetscCall(PetscMalloc1(rsize,&user->right)); 479 480 hx= (r-l)/(mx+1); hy=(t-b)/(my+1); 481 482 for (j=0; j<4; j++) { 483 if (j==0) { 484 yt=b; 485 xt=l; 486 limit=bsize; 487 boundary=user->bottom; 488 } else if (j==1) { 489 yt=t; 490 xt=l; 491 limit=tsize; 492 boundary=user->top; 493 } else if (j==2) { 494 yt=b; 495 xt=l; 496 limit=lsize; 497 boundary=user->left; 498 } else { /* if (j==3) */ 499 yt=b; 500 xt=r; 501 limit=rsize; 502 boundary=user->right; 503 } 504 505 for (i=0; i<limit; i++) { 506 u1=xt; 507 u2=-yt; 508 for (k=0; k<maxits; k++) { 509 nf1=u1 + u1*u2*u2 - u1*u1*u1/three-xt; 510 nf2=-u2 - u1*u1*u2 + u2*u2*u2/three-yt; 511 fnorm=PetscSqrtScalar(nf1*nf1+nf2*nf2); 512 if (fnorm <= tol) break; 513 njac11=one+u2*u2-u1*u1; 514 njac12=two*u1*u2; 515 njac21=-two*u1*u2; 516 njac22=-one - u1*u1 + u2*u2; 517 det = njac11*njac22-njac21*njac12; 518 u1 = u1-(njac22*nf1-njac12*nf2)/det; 519 u2 = u2-(njac11*nf2-njac21*nf1)/det; 520 } 521 522 boundary[i]=u1*u1-u2*u2; 523 if (j==0 || j==1) { 524 xt=xt+hx; 525 } else { /* if (j==2 || j==3) */ 526 yt=yt+hy; 527 } 528 } 529 } 530 PetscFunctionReturn(0); 531 } 532 533 /* ------------------------------------------------------------------- */ 534 /* 535 MSA_InitialPoint - Calculates the initial guess in one of three ways. 536 537 Input Parameters: 538 . user - user-defined application context 539 . X - vector for initial guess 540 541 Output Parameters: 542 . X - newly computed initial guess 543 */ 544 static PetscErrorCode MSA_InitialPoint(AppCtx * user, Vec X) 545 { 546 PetscInt start=-1,i,j; 547 PetscReal zero=0.0; 548 PetscBool flg; 549 550 PetscFunctionBeginUser; 551 PetscCall(VecSet(X, zero)); 552 PetscCall(PetscOptionsGetInt(NULL,NULL,"-start",&start,&flg)); 553 554 if (flg && start==0) { /* The zero vector is reasonable */ 555 PetscCall(VecSet(X, zero)); 556 } else { /* Take an average of the boundary conditions */ 557 PetscInt row; 558 PetscInt mx=user->mx,my=user->my; 559 PetscReal *x; 560 561 /* Get pointers to vector data */ 562 PetscCall(VecGetArray(X,&x)); 563 /* Perform local computations */ 564 for (j=0; j<my; j++) { 565 for (i=0; i< mx; i++) { 566 row=(j)*mx + (i); 567 x[row] = (((j+1)*user->bottom[i+1]+(my-j+1)*user->top[i+1])/(my+2)+ ((i+1)*user->left[j+1]+(mx-i+1)*user->right[j+1])/(mx+2))/2.0; 568 } 569 } 570 /* Restore vectors */ 571 PetscCall(VecRestoreArray(X,&x)); 572 } 573 PetscFunctionReturn(0); 574 } 575 576 /*TEST 577 578 build: 579 requires: !complex 580 581 test: 582 args: -tao_smonitor -tao_type nls -mx 10 -my 8 -tao_gatol 1.e-4 583 requires: !single 584 585 test: 586 suffix: 2 587 args: -tao_smonitor -tao_type bmrm -mx 10 -my 8 -tao_gatol 1.e-3 588 requires: !single 589 590 test: 591 suffix: 3 592 args: -tao_smonitor -tao_type lmvm -mx 10 -my 8 -tao_gatol 1.e-3 593 requires: !single 594 595 test: 596 suffix: 4 597 args: -tao_smonitor -tao_type bntr -mx 10 -my 8 -tao_gatol 1.e-4 598 requires: !single 599 600 test: 601 suffix: 4_ew 602 args: -tao_smonitor -tao_type bntr -tao_ksp_ew -mx 10 -my 8 -tao_gatol 1.e-4 603 requires: !single 604 605 test: 606 suffix: 5 607 args: -tao_smonitor -tao_type bntl -mx 10 -my 8 -tao_gatol 1.e-4 608 requires: !single 609 610 test: 611 suffix: 5_ew 612 args: -tao_smonitor -tao_type bntl -tao_ksp_ew -mx 10 -my 8 -tao_gatol 1.e-4 613 requires: !single 614 615 test: 616 suffix: 6 617 args: -tao_smonitor -tao_type bnls -mx 10 -my 8 -tao_gatol 1.e-4 618 requires: !single 619 620 test: 621 suffix: 6_ew 622 args: -tao_smonitor -tao_type bnls -tao_ksp_ew -tao_bnk_ksp_ew_version 3 -mx 10 -my 8 -tao_gatol 1.e-4 623 requires: !single 624 625 test: 626 suffix: 7 627 args: -tao_smonitor -tao_type bntr -mx 10 -my 8 -tao_bnk_max_cg_its 3 -tao_gatol 1.e-4 628 requires: !single 629 630 test: 631 suffix: 8 632 args: -tao_smonitor -tao_type bntl -mx 10 -my 8 -tao_bnk_max_cg_its 3 -tao_gatol 1.e-4 633 requires: !single 634 635 test: 636 suffix: 9 637 args: -tao_smonitor -tao_type bnls -mx 10 -my 8 -tao_bnk_max_cg_its 3 -tao_gatol 1.e-4 638 requires: !single 639 640 test: 641 suffix: 10 642 args: -tao_smonitor -tao_type bnls -mx 10 -my 8 -tao_bnk_max_cg_its 3 -tao_gatol 1.e-4 -tao_mf_hessian 643 644 test: 645 suffix: 11 646 args: -tao_smonitor -tao_type bntr -mx 10 -my 8 -tao_gatol 1.e-4 -tao_mf_hessian 647 requires: !single 648 649 test: 650 suffix: 12 651 args: -tao_smonitor -tao_type bntl -mx 10 -my 8 -tao_gatol 1.e-4 -tao_mf_hessian 652 requires: !single 653 654 TEST*/ 655