1 2 static char help[] = "Basic equation for generator stability analysis.\n"; 3 4 /*F 5 6 \begin{eqnarray} 7 \frac{d \theta}{dt} = \omega_b (\omega - \omega_s) 8 \frac{2 H}{\omega_s}\frac{d \omega}{dt} & = & P_m - P_max \sin(\theta) -D(\omega - \omega_s)\\ 9 \end{eqnarray} 10 11 Ensemble of initial conditions 12 ./ex2 -ensemble -ts_monitor_draw_solution_phase -1,-3,3,3 -ts_adapt_dt_max .01 -ts_monitor -ts_type rosw -pc_type lu -ksp_type preonly 13 14 Fault at .1 seconds 15 ./ex2 -ts_monitor_draw_solution_phase .42,.95,.6,1.05 -ts_adapt_dt_max .01 -ts_monitor -ts_type rosw -pc_type lu -ksp_type preonly 16 17 Initial conditions same as when fault is ended 18 ./ex2 -u 0.496792,1.00932 -ts_monitor_draw_solution_phase .42,.95,.6,1.05 -ts_adapt_dt_max .01 -ts_monitor -ts_type rosw -pc_type lu -ksp_type preonly 19 20 F*/ 21 22 /* 23 Include "petscts.h" so that we can use TS solvers. Note that this 24 file automatically includes: 25 petscsys.h - base PETSc routines petscvec.h - vectors 26 petscmat.h - matrices 27 petscis.h - index sets petscksp.h - Krylov subspace methods 28 petscviewer.h - viewers petscpc.h - preconditioners 29 petscksp.h - linear solvers 30 */ 31 32 #include <petsctao.h> 33 #include <petscts.h> 34 35 typedef struct { 36 TS ts; 37 PetscScalar H, D, omega_b, omega_s, Pmax, Pm, E, V, X, u_s, c; 38 PetscInt beta; 39 PetscReal tf, tcl, dt; 40 } AppCtx; 41 42 PetscErrorCode FormFunction(Tao, Vec, PetscReal *, void *); 43 PetscErrorCode FormGradient(Tao, Vec, Vec, void *); 44 45 /* 46 Defines the ODE passed to the ODE solver 47 */ 48 static PetscErrorCode RHSFunction(TS ts, PetscReal t, Vec U, Vec F, AppCtx *ctx) { 49 PetscScalar *f, Pmax; 50 const PetscScalar *u; 51 52 PetscFunctionBegin; 53 /* The next three lines allow us to access the entries of the vectors directly */ 54 PetscCall(VecGetArrayRead(U, &u)); 55 PetscCall(VecGetArray(F, &f)); 56 if ((t > ctx->tf) && (t < ctx->tcl)) Pmax = 0.0; /* A short-circuit on the generator terminal that drives the electrical power output (Pmax*sin(delta)) to 0 */ 57 else Pmax = ctx->Pmax; 58 59 f[0] = ctx->omega_b * (u[1] - ctx->omega_s); 60 f[1] = (-Pmax * PetscSinScalar(u[0]) - ctx->D * (u[1] - ctx->omega_s) + ctx->Pm) * ctx->omega_s / (2.0 * ctx->H); 61 62 PetscCall(VecRestoreArrayRead(U, &u)); 63 PetscCall(VecRestoreArray(F, &f)); 64 PetscFunctionReturn(0); 65 } 66 67 /* 68 Defines the Jacobian of the ODE passed to the ODE solver. See TSSetIJacobian() for the meaning of a and the Jacobian. 69 */ 70 static PetscErrorCode RHSJacobian(TS ts, PetscReal t, Vec U, Mat A, Mat B, AppCtx *ctx) { 71 PetscInt rowcol[] = {0, 1}; 72 PetscScalar J[2][2], Pmax; 73 const PetscScalar *u; 74 75 PetscFunctionBegin; 76 PetscCall(VecGetArrayRead(U, &u)); 77 if ((t > ctx->tf) && (t < ctx->tcl)) Pmax = 0.0; /* A short-circuit on the generator terminal that drives the electrical power output (Pmax*sin(delta)) to 0 */ 78 else Pmax = ctx->Pmax; 79 80 J[0][0] = 0; 81 J[0][1] = ctx->omega_b; 82 J[1][1] = -ctx->D * ctx->omega_s / (2.0 * ctx->H); 83 J[1][0] = -Pmax * PetscCosScalar(u[0]) * ctx->omega_s / (2.0 * ctx->H); 84 85 PetscCall(MatSetValues(A, 2, rowcol, 2, rowcol, &J[0][0], INSERT_VALUES)); 86 PetscCall(VecRestoreArrayRead(U, &u)); 87 88 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 89 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 90 if (A != B) { 91 PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY)); 92 PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY)); 93 } 94 PetscFunctionReturn(0); 95 } 96 97 static PetscErrorCode RHSJacobianP(TS ts, PetscReal t, Vec X, Mat A, void *ctx0) { 98 PetscInt row[] = {0, 1}, col[] = {0}; 99 PetscScalar J[2][1]; 100 AppCtx *ctx = (AppCtx *)ctx0; 101 102 PetscFunctionBeginUser; 103 J[0][0] = 0; 104 J[1][0] = ctx->omega_s / (2.0 * ctx->H); 105 PetscCall(MatSetValues(A, 2, row, 1, col, &J[0][0], INSERT_VALUES)); 106 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 107 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 108 PetscFunctionReturn(0); 109 } 110 111 static PetscErrorCode CostIntegrand(TS ts, PetscReal t, Vec U, Vec R, AppCtx *ctx) { 112 PetscScalar *r; 113 const PetscScalar *u; 114 115 PetscFunctionBegin; 116 PetscCall(VecGetArrayRead(U, &u)); 117 PetscCall(VecGetArray(R, &r)); 118 r[0] = ctx->c * PetscPowScalarInt(PetscMax(0., u[0] - ctx->u_s), ctx->beta); 119 PetscCall(VecRestoreArray(R, &r)); 120 PetscCall(VecRestoreArrayRead(U, &u)); 121 PetscFunctionReturn(0); 122 } 123 124 static PetscErrorCode DRDUJacobianTranspose(TS ts, PetscReal t, Vec U, Mat DRDU, Mat B, AppCtx *ctx) { 125 PetscScalar ru[1]; 126 const PetscScalar *u; 127 PetscInt row[] = {0}, col[] = {0}; 128 129 PetscFunctionBegin; 130 PetscCall(VecGetArrayRead(U, &u)); 131 ru[0] = ctx->c * ctx->beta * PetscPowScalarInt(PetscMax(0., u[0] - ctx->u_s), ctx->beta - 1); 132 PetscCall(VecRestoreArrayRead(U, &u)); 133 PetscCall(MatSetValues(DRDU, 1, row, 1, col, ru, INSERT_VALUES)); 134 PetscCall(MatAssemblyBegin(DRDU, MAT_FINAL_ASSEMBLY)); 135 PetscCall(MatAssemblyEnd(DRDU, MAT_FINAL_ASSEMBLY)); 136 PetscFunctionReturn(0); 137 } 138 139 static PetscErrorCode DRDPJacobianTranspose(TS ts, PetscReal t, Vec U, Mat DRDP, AppCtx *ctx) { 140 PetscFunctionBegin; 141 PetscCall(MatZeroEntries(DRDP)); 142 PetscCall(MatAssemblyBegin(DRDP, MAT_FINAL_ASSEMBLY)); 143 PetscCall(MatAssemblyEnd(DRDP, MAT_FINAL_ASSEMBLY)); 144 PetscFunctionReturn(0); 145 } 146 147 PetscErrorCode ComputeSensiP(Vec lambda, Vec mu, AppCtx *ctx) { 148 PetscScalar *y, sensip; 149 const PetscScalar *x; 150 151 PetscFunctionBegin; 152 PetscCall(VecGetArrayRead(lambda, &x)); 153 PetscCall(VecGetArray(mu, &y)); 154 sensip = 1. / PetscSqrtScalar(1. - (ctx->Pm / ctx->Pmax) * (ctx->Pm / ctx->Pmax)) / ctx->Pmax * x[0] + y[0]; 155 y[0] = sensip; 156 PetscCall(VecRestoreArray(mu, &y)); 157 PetscCall(VecRestoreArrayRead(lambda, &x)); 158 PetscFunctionReturn(0); 159 } 160 161 int main(int argc, char **argv) { 162 Vec p; 163 PetscScalar *x_ptr; 164 PetscMPIInt size; 165 AppCtx ctx; 166 Vec lowerb, upperb; 167 Tao tao; 168 KSP ksp; 169 PC pc; 170 Vec U, lambda[1], mu[1]; 171 Mat A; /* Jacobian matrix */ 172 Mat Jacp; /* Jacobian matrix */ 173 Mat DRDU, DRDP; 174 PetscInt n = 2; 175 TS quadts; 176 177 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 178 Initialize program 179 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 180 PetscFunctionBeginUser; 181 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 182 PetscFunctionBeginUser; 183 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 184 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!"); 185 186 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 187 Set runtime options 188 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 189 PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Swing equation options", ""); 190 { 191 ctx.beta = 2; 192 ctx.c = PetscRealConstant(10000.0); 193 ctx.u_s = PetscRealConstant(1.0); 194 ctx.omega_s = PetscRealConstant(1.0); 195 ctx.omega_b = PetscRealConstant(120.0) * PETSC_PI; 196 ctx.H = PetscRealConstant(5.0); 197 PetscCall(PetscOptionsScalar("-Inertia", "", "", ctx.H, &ctx.H, NULL)); 198 ctx.D = PetscRealConstant(5.0); 199 PetscCall(PetscOptionsScalar("-D", "", "", ctx.D, &ctx.D, NULL)); 200 ctx.E = PetscRealConstant(1.1378); 201 ctx.V = PetscRealConstant(1.0); 202 ctx.X = PetscRealConstant(0.545); 203 ctx.Pmax = ctx.E * ctx.V / ctx.X; 204 PetscCall(PetscOptionsScalar("-Pmax", "", "", ctx.Pmax, &ctx.Pmax, NULL)); 205 ctx.Pm = PetscRealConstant(1.0194); 206 PetscCall(PetscOptionsScalar("-Pm", "", "", ctx.Pm, &ctx.Pm, NULL)); 207 ctx.tf = PetscRealConstant(0.1); 208 ctx.tcl = PetscRealConstant(0.2); 209 PetscCall(PetscOptionsReal("-tf", "Time to start fault", "", ctx.tf, &ctx.tf, NULL)); 210 PetscCall(PetscOptionsReal("-tcl", "Time to end fault", "", ctx.tcl, &ctx.tcl, NULL)); 211 } 212 PetscOptionsEnd(); 213 214 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 215 Create necessary matrix and vectors 216 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 217 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 218 PetscCall(MatSetSizes(A, n, n, PETSC_DETERMINE, PETSC_DETERMINE)); 219 PetscCall(MatSetType(A, MATDENSE)); 220 PetscCall(MatSetFromOptions(A)); 221 PetscCall(MatSetUp(A)); 222 223 PetscCall(MatCreateVecs(A, &U, NULL)); 224 225 PetscCall(MatCreate(PETSC_COMM_WORLD, &Jacp)); 226 PetscCall(MatSetSizes(Jacp, PETSC_DECIDE, PETSC_DECIDE, 2, 1)); 227 PetscCall(MatSetFromOptions(Jacp)); 228 PetscCall(MatSetUp(Jacp)); 229 PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, &DRDP)); 230 PetscCall(MatSetUp(DRDP)); 231 PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, 1, 2, NULL, &DRDU)); 232 PetscCall(MatSetUp(DRDU)); 233 234 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 235 Create timestepping solver context 236 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 237 PetscCall(TSCreate(PETSC_COMM_WORLD, &ctx.ts)); 238 PetscCall(TSSetProblemType(ctx.ts, TS_NONLINEAR)); 239 PetscCall(TSSetEquationType(ctx.ts, TS_EQ_ODE_EXPLICIT)); /* less Jacobian evaluations when adjoint BEuler is used, otherwise no effect */ 240 PetscCall(TSSetType(ctx.ts, TSRK)); 241 PetscCall(TSSetRHSFunction(ctx.ts, NULL, (TSRHSFunction)RHSFunction, &ctx)); 242 PetscCall(TSSetRHSJacobian(ctx.ts, A, A, (TSRHSJacobian)RHSJacobian, &ctx)); 243 PetscCall(TSSetExactFinalTime(ctx.ts, TS_EXACTFINALTIME_MATCHSTEP)); 244 245 PetscCall(MatCreateVecs(A, &lambda[0], NULL)); 246 PetscCall(MatCreateVecs(Jacp, &mu[0], NULL)); 247 PetscCall(TSSetCostGradients(ctx.ts, 1, lambda, mu)); 248 PetscCall(TSSetRHSJacobianP(ctx.ts, Jacp, RHSJacobianP, &ctx)); 249 250 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 251 Set solver options 252 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 253 PetscCall(TSSetMaxTime(ctx.ts, PetscRealConstant(1.0))); 254 PetscCall(TSSetTimeStep(ctx.ts, PetscRealConstant(0.01))); 255 PetscCall(TSSetFromOptions(ctx.ts)); 256 257 PetscCall(TSGetTimeStep(ctx.ts, &ctx.dt)); /* save the stepsize */ 258 259 PetscCall(TSCreateQuadratureTS(ctx.ts, PETSC_TRUE, &quadts)); 260 PetscCall(TSSetRHSFunction(quadts, NULL, (TSRHSFunction)CostIntegrand, &ctx)); 261 PetscCall(TSSetRHSJacobian(quadts, DRDU, DRDU, (TSRHSJacobian)DRDUJacobianTranspose, &ctx)); 262 PetscCall(TSSetRHSJacobianP(quadts, DRDP, (TSRHSJacobianP)DRDPJacobianTranspose, &ctx)); 263 PetscCall(TSSetSolution(ctx.ts, U)); 264 265 /* Create TAO solver and set desired solution method */ 266 PetscCall(TaoCreate(PETSC_COMM_WORLD, &tao)); 267 PetscCall(TaoSetType(tao, TAOBLMVM)); 268 269 /* 270 Optimization starts 271 */ 272 /* Set initial solution guess */ 273 PetscCall(VecCreateSeq(PETSC_COMM_WORLD, 1, &p)); 274 PetscCall(VecGetArray(p, &x_ptr)); 275 x_ptr[0] = ctx.Pm; 276 PetscCall(VecRestoreArray(p, &x_ptr)); 277 278 PetscCall(TaoSetSolution(tao, p)); 279 /* Set routine for function and gradient evaluation */ 280 PetscCall(TaoSetObjective(tao, FormFunction, (void *)&ctx)); 281 PetscCall(TaoSetGradient(tao, NULL, FormGradient, (void *)&ctx)); 282 283 /* Set bounds for the optimization */ 284 PetscCall(VecDuplicate(p, &lowerb)); 285 PetscCall(VecDuplicate(p, &upperb)); 286 PetscCall(VecGetArray(lowerb, &x_ptr)); 287 x_ptr[0] = 0.; 288 PetscCall(VecRestoreArray(lowerb, &x_ptr)); 289 PetscCall(VecGetArray(upperb, &x_ptr)); 290 x_ptr[0] = PetscRealConstant(1.1); 291 PetscCall(VecRestoreArray(upperb, &x_ptr)); 292 PetscCall(TaoSetVariableBounds(tao, lowerb, upperb)); 293 294 /* Check for any TAO command line options */ 295 PetscCall(TaoSetFromOptions(tao)); 296 PetscCall(TaoGetKSP(tao, &ksp)); 297 if (ksp) { 298 PetscCall(KSPGetPC(ksp, &pc)); 299 PetscCall(PCSetType(pc, PCNONE)); 300 } 301 302 /* SOLVE THE APPLICATION */ 303 PetscCall(TaoSolve(tao)); 304 305 PetscCall(VecView(p, PETSC_VIEWER_STDOUT_WORLD)); 306 /* Free TAO data structures */ 307 PetscCall(TaoDestroy(&tao)); 308 PetscCall(VecDestroy(&p)); 309 PetscCall(VecDestroy(&lowerb)); 310 PetscCall(VecDestroy(&upperb)); 311 312 PetscCall(TSDestroy(&ctx.ts)); 313 PetscCall(VecDestroy(&U)); 314 PetscCall(MatDestroy(&A)); 315 PetscCall(MatDestroy(&Jacp)); 316 PetscCall(MatDestroy(&DRDU)); 317 PetscCall(MatDestroy(&DRDP)); 318 PetscCall(VecDestroy(&lambda[0])); 319 PetscCall(VecDestroy(&mu[0])); 320 PetscCall(PetscFinalize()); 321 return 0; 322 } 323 324 /* ------------------------------------------------------------------ */ 325 /* 326 FormFunction - Evaluates the function 327 328 Input Parameters: 329 tao - the Tao context 330 X - the input vector 331 ptr - optional user-defined context, as set by TaoSetObjectiveAndGradient() 332 333 Output Parameters: 334 f - the newly evaluated function 335 */ 336 PetscErrorCode FormFunction(Tao tao, Vec P, PetscReal *f, void *ctx0) { 337 AppCtx *ctx = (AppCtx *)ctx0; 338 TS ts = ctx->ts; 339 Vec U; /* solution will be stored here */ 340 PetscScalar *u; 341 PetscScalar *x_ptr; 342 Vec q; 343 344 PetscCall(VecGetArrayRead(P, (const PetscScalar **)&x_ptr)); 345 ctx->Pm = x_ptr[0]; 346 PetscCall(VecRestoreArrayRead(P, (const PetscScalar **)&x_ptr)); 347 348 /* reset time */ 349 PetscCall(TSSetTime(ts, 0.0)); 350 /* reset step counter, this is critical for adjoint solver */ 351 PetscCall(TSSetStepNumber(ts, 0)); 352 /* reset step size, the step size becomes negative after TSAdjointSolve */ 353 PetscCall(TSSetTimeStep(ts, ctx->dt)); 354 /* reinitialize the integral value */ 355 PetscCall(TSGetCostIntegral(ts, &q)); 356 PetscCall(VecSet(q, 0.0)); 357 358 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 359 Set initial conditions 360 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 361 PetscCall(TSGetSolution(ts, &U)); 362 PetscCall(VecGetArray(U, &u)); 363 u[0] = PetscAsinScalar(ctx->Pm / ctx->Pmax); 364 u[1] = PetscRealConstant(1.0); 365 PetscCall(VecRestoreArray(U, &u)); 366 367 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 368 Solve nonlinear system 369 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 370 PetscCall(TSSolve(ts, U)); 371 PetscCall(TSGetCostIntegral(ts, &q)); 372 PetscCall(VecGetArray(q, &x_ptr)); 373 *f = -ctx->Pm + x_ptr[0]; 374 PetscCall(VecRestoreArray(q, &x_ptr)); 375 return 0; 376 } 377 378 PetscErrorCode FormGradient(Tao tao, Vec P, Vec G, void *ctx0) { 379 AppCtx *ctx = (AppCtx *)ctx0; 380 TS ts = ctx->ts; 381 Vec U; /* solution will be stored here */ 382 PetscReal ftime; 383 PetscInt steps; 384 PetscScalar *u; 385 PetscScalar *x_ptr, *y_ptr; 386 Vec *lambda, q, *mu; 387 388 PetscCall(VecGetArrayRead(P, (const PetscScalar **)&x_ptr)); 389 ctx->Pm = x_ptr[0]; 390 PetscCall(VecRestoreArrayRead(P, (const PetscScalar **)&x_ptr)); 391 392 /* reset time */ 393 PetscCall(TSSetTime(ts, 0.0)); 394 /* reset step counter, this is critical for adjoint solver */ 395 PetscCall(TSSetStepNumber(ts, 0)); 396 /* reset step size, the step size becomes negative after TSAdjointSolve */ 397 PetscCall(TSSetTimeStep(ts, ctx->dt)); 398 /* reinitialize the integral value */ 399 PetscCall(TSGetCostIntegral(ts, &q)); 400 PetscCall(VecSet(q, 0.0)); 401 402 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 403 Set initial conditions 404 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 405 PetscCall(TSGetSolution(ts, &U)); 406 PetscCall(VecGetArray(U, &u)); 407 u[0] = PetscAsinScalar(ctx->Pm / ctx->Pmax); 408 u[1] = PetscRealConstant(1.0); 409 PetscCall(VecRestoreArray(U, &u)); 410 411 /* Set up to save trajectory before TSSetFromOptions() so that TSTrajectory options can be captured */ 412 PetscCall(TSSetSaveTrajectory(ts)); 413 PetscCall(TSSetFromOptions(ts)); 414 415 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 416 Solve nonlinear system 417 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 418 PetscCall(TSSolve(ts, U)); 419 420 PetscCall(TSGetSolveTime(ts, &ftime)); 421 PetscCall(TSGetStepNumber(ts, &steps)); 422 423 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 424 Adjoint model starts here 425 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 426 PetscCall(TSGetCostGradients(ts, NULL, &lambda, &mu)); 427 /* Set initial conditions for the adjoint integration */ 428 PetscCall(VecGetArray(lambda[0], &y_ptr)); 429 y_ptr[0] = 0.0; 430 y_ptr[1] = 0.0; 431 PetscCall(VecRestoreArray(lambda[0], &y_ptr)); 432 PetscCall(VecGetArray(mu[0], &x_ptr)); 433 x_ptr[0] = PetscRealConstant(-1.0); 434 PetscCall(VecRestoreArray(mu[0], &x_ptr)); 435 436 PetscCall(TSAdjointSolve(ts)); 437 PetscCall(TSGetCostIntegral(ts, &q)); 438 PetscCall(ComputeSensiP(lambda[0], mu[0], ctx)); 439 PetscCall(VecCopy(mu[0], G)); 440 return 0; 441 } 442 443 /*TEST 444 445 build: 446 requires: !complex 447 448 test: 449 args: -viewer_binary_skip_info -ts_adapt_type none -tao_monitor -tao_gatol 0.0 -tao_grtol 1.e-3 -tao_converged_reason 450 451 test: 452 suffix: 2 453 args: -viewer_binary_skip_info -ts_adapt_type none -tao_monitor -tao_gatol 0.0 -tao_grtol 1.e-3 -tao_converged_reason -tao_test_gradient 454 455 TEST*/ 456