1 2 static char help[] = "Solves a simple time-dependent linear PDE (the heat equation).\n\ 3 Input parameters include:\n\ 4 -m <points>, where <points> = number of grid points\n\ 5 -time_dependent_rhs : Treat the problem as having a time-dependent right-hand side\n\ 6 -use_ifunc : Use IFunction/IJacobian interface\n\ 7 -debug : Activate debugging printouts\n\ 8 -nox : Deactivate x-window graphics\n\n"; 9 10 /* ------------------------------------------------------------------------ 11 12 This program solves the one-dimensional heat equation (also called the 13 diffusion equation), 14 u_t = u_xx, 15 on the domain 0 <= x <= 1, with the boundary conditions 16 u(t,0) = 0, u(t,1) = 0, 17 and the initial condition 18 u(0,x) = sin(6*pi*x) + 3*sin(2*pi*x). 19 This is a linear, second-order, parabolic equation. 20 21 We discretize the right-hand side using finite differences with 22 uniform grid spacing h: 23 u_xx = (u_{i+1} - 2u_{i} + u_{i-1})/(h^2) 24 We then demonstrate time evolution using the various TS methods by 25 running the program via 26 ex3 -ts_type <timestepping solver> 27 28 We compare the approximate solution with the exact solution, given by 29 u_exact(x,t) = exp(-36*pi*pi*t) * sin(6*pi*x) + 30 3*exp(-4*pi*pi*t) * sin(2*pi*x) 31 32 Notes: 33 This code demonstrates the TS solver interface to two variants of 34 linear problems, u_t = f(u,t), namely 35 - time-dependent f: f(u,t) is a function of t 36 - time-independent f: f(u,t) is simply f(u) 37 38 The parallel version of this code is ts/tutorials/ex4.c 39 40 ------------------------------------------------------------------------- */ 41 42 /* 43 Include "petscts.h" so that we can use TS solvers. Note that this file 44 automatically includes: 45 petscsys.h - base PETSc routines petscvec.h - vectors 46 petscmat.h - matrices 47 petscis.h - index sets petscksp.h - Krylov subspace methods 48 petscviewer.h - viewers petscpc.h - preconditioners 49 petscksp.h - linear solvers petscsnes.h - nonlinear solvers 50 */ 51 52 #include <petscts.h> 53 #include <petscdraw.h> 54 55 /* 56 User-defined application context - contains data needed by the 57 application-provided call-back routines. 58 */ 59 typedef struct { 60 Vec solution; /* global exact solution vector */ 61 PetscInt m; /* total number of grid points */ 62 PetscReal h; /* mesh width h = 1/(m-1) */ 63 PetscBool debug; /* flag (1 indicates activation of debugging printouts) */ 64 PetscViewer viewer1, viewer2; /* viewers for the solution and error */ 65 PetscReal norm_2, norm_max; /* error norms */ 66 Mat A; /* RHS mat, used with IFunction interface */ 67 PetscReal oshift; /* old shift applied, prevent to recompute the IJacobian */ 68 } AppCtx; 69 70 /* 71 User-defined routines 72 */ 73 extern PetscErrorCode InitialConditions(Vec, AppCtx *); 74 extern PetscErrorCode RHSMatrixHeat(TS, PetscReal, Vec, Mat, Mat, void *); 75 extern PetscErrorCode IFunctionHeat(TS, PetscReal, Vec, Vec, Vec, void *); 76 extern PetscErrorCode IJacobianHeat(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *); 77 extern PetscErrorCode Monitor(TS, PetscInt, PetscReal, Vec, void *); 78 extern PetscErrorCode ExactSolution(PetscReal, Vec, AppCtx *); 79 80 int main(int argc, char **argv) { 81 AppCtx appctx; /* user-defined application context */ 82 TS ts; /* timestepping context */ 83 Mat A; /* matrix data structure */ 84 Vec u; /* approximate solution vector */ 85 PetscReal time_total_max = 100.0; /* default max total time */ 86 PetscInt time_steps_max = 100; /* default max timesteps */ 87 PetscDraw draw; /* drawing context */ 88 PetscInt steps, m; 89 PetscMPIInt size; 90 PetscReal dt; 91 PetscBool flg, flg_string; 92 93 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 94 Initialize program and set problem parameters 95 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 96 97 PetscFunctionBeginUser; 98 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 99 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 100 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!"); 101 102 m = 60; 103 PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL)); 104 PetscCall(PetscOptionsHasName(NULL, NULL, "-debug", &appctx.debug)); 105 flg_string = PETSC_FALSE; 106 PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_string_viewer", &flg_string, NULL)); 107 108 appctx.m = m; 109 appctx.h = 1.0 / (m - 1.0); 110 appctx.norm_2 = 0.0; 111 appctx.norm_max = 0.0; 112 113 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Solving a linear TS problem on 1 processor\n")); 114 115 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 116 Create vector data structures 117 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 118 119 /* 120 Create vector data structures for approximate and exact solutions 121 */ 122 PetscCall(VecCreateSeq(PETSC_COMM_SELF, m, &u)); 123 PetscCall(VecDuplicate(u, &appctx.solution)); 124 125 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 126 Set up displays to show graphs of the solution and error 127 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 128 129 PetscCall(PetscViewerDrawOpen(PETSC_COMM_SELF, 0, "", 80, 380, 400, 160, &appctx.viewer1)); 130 PetscCall(PetscViewerDrawGetDraw(appctx.viewer1, 0, &draw)); 131 PetscCall(PetscDrawSetDoubleBuffer(draw)); 132 PetscCall(PetscViewerDrawOpen(PETSC_COMM_SELF, 0, "", 80, 0, 400, 160, &appctx.viewer2)); 133 PetscCall(PetscViewerDrawGetDraw(appctx.viewer2, 0, &draw)); 134 PetscCall(PetscDrawSetDoubleBuffer(draw)); 135 136 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 137 Create timestepping solver context 138 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 139 140 PetscCall(TSCreate(PETSC_COMM_SELF, &ts)); 141 PetscCall(TSSetProblemType(ts, TS_LINEAR)); 142 143 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 144 Set optional user-defined monitoring routine 145 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 146 147 if (!flg_string) PetscCall(TSMonitorSet(ts, Monitor, &appctx, NULL)); 148 149 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 150 151 Create matrix data structure; set matrix evaluation routine. 152 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 153 154 PetscCall(MatCreate(PETSC_COMM_SELF, &A)); 155 PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m, m)); 156 PetscCall(MatSetFromOptions(A)); 157 PetscCall(MatSetUp(A)); 158 159 flg = PETSC_FALSE; 160 PetscCall(PetscOptionsGetBool(NULL, NULL, "-use_ifunc", &flg, NULL)); 161 if (!flg) { 162 appctx.A = NULL; 163 PetscCall(PetscOptionsGetBool(NULL, NULL, "-time_dependent_rhs", &flg, NULL)); 164 if (flg) { 165 /* 166 For linear problems with a time-dependent f(u,t) in the equation 167 u_t = f(u,t), the user provides the discretized right-hand-side 168 as a time-dependent matrix. 169 */ 170 PetscCall(TSSetRHSFunction(ts, NULL, TSComputeRHSFunctionLinear, &appctx)); 171 PetscCall(TSSetRHSJacobian(ts, A, A, RHSMatrixHeat, &appctx)); 172 } else { 173 /* 174 For linear problems with a time-independent f(u) in the equation 175 u_t = f(u), the user provides the discretized right-hand-side 176 as a matrix only once, and then sets the special Jacobian evaluation 177 routine TSComputeRHSJacobianConstant() which will NOT recompute the Jacobian. 178 */ 179 PetscCall(RHSMatrixHeat(ts, 0.0, u, A, A, &appctx)); 180 PetscCall(TSSetRHSFunction(ts, NULL, TSComputeRHSFunctionLinear, &appctx)); 181 PetscCall(TSSetRHSJacobian(ts, A, A, TSComputeRHSJacobianConstant, &appctx)); 182 } 183 } else { 184 Mat J; 185 186 PetscCall(RHSMatrixHeat(ts, 0.0, u, A, A, &appctx)); 187 PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &J)); 188 PetscCall(TSSetIFunction(ts, NULL, IFunctionHeat, &appctx)); 189 PetscCall(TSSetIJacobian(ts, J, J, IJacobianHeat, &appctx)); 190 PetscCall(MatDestroy(&J)); 191 192 PetscCall(PetscObjectReference((PetscObject)A)); 193 appctx.A = A; 194 appctx.oshift = PETSC_MIN_REAL; 195 } 196 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 197 Set solution vector and initial timestep 198 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 199 200 dt = appctx.h * appctx.h / 2.0; 201 PetscCall(TSSetTimeStep(ts, dt)); 202 203 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 204 Customize timestepping solver: 205 - Set the solution method to be the Backward Euler method. 206 - Set timestepping duration info 207 Then set runtime options, which can override these defaults. 208 For example, 209 -ts_max_steps <maxsteps> -ts_max_time <maxtime> 210 to override the defaults set by TSSetMaxSteps()/TSSetMaxTime(). 211 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 212 213 PetscCall(TSSetMaxSteps(ts, time_steps_max)); 214 PetscCall(TSSetMaxTime(ts, time_total_max)); 215 PetscCall(TSSetExactFinalTime(ts, TS_EXACTFINALTIME_STEPOVER)); 216 PetscCall(TSSetFromOptions(ts)); 217 218 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 219 Solve the problem 220 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 221 222 /* 223 Evaluate initial conditions 224 */ 225 PetscCall(InitialConditions(u, &appctx)); 226 227 /* 228 Run the timestepping solver 229 */ 230 PetscCall(TSSolve(ts, u)); 231 PetscCall(TSGetStepNumber(ts, &steps)); 232 233 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 234 View timestepping solver info 235 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 236 237 PetscCall(PetscPrintf(PETSC_COMM_SELF, "avg. error (2 norm) = %g, avg. error (max norm) = %g\n", (double)(appctx.norm_2 / steps), (double)(appctx.norm_max / steps))); 238 if (!flg_string) { 239 PetscCall(TSView(ts, PETSC_VIEWER_STDOUT_SELF)); 240 } else { 241 PetscViewer stringviewer; 242 char string[512]; 243 const char *outstring; 244 245 PetscCall(PetscViewerStringOpen(PETSC_COMM_WORLD, string, sizeof(string), &stringviewer)); 246 PetscCall(TSView(ts, stringviewer)); 247 PetscCall(PetscViewerStringGetStringRead(stringviewer, &outstring, NULL)); 248 PetscCheck((char *)outstring == (char *)string, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "String returned from viewer does not equal original string"); 249 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Output from string viewer:%s\n", outstring)); 250 PetscCall(PetscViewerDestroy(&stringviewer)); 251 } 252 253 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 254 Free work space. All PETSc objects should be destroyed when they 255 are no longer needed. 256 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 257 258 PetscCall(TSDestroy(&ts)); 259 PetscCall(MatDestroy(&A)); 260 PetscCall(VecDestroy(&u)); 261 PetscCall(PetscViewerDestroy(&appctx.viewer1)); 262 PetscCall(PetscViewerDestroy(&appctx.viewer2)); 263 PetscCall(VecDestroy(&appctx.solution)); 264 PetscCall(MatDestroy(&appctx.A)); 265 266 /* 267 Always call PetscFinalize() before exiting a program. This routine 268 - finalizes the PETSc libraries as well as MPI 269 - provides summary and diagnostic information if certain runtime 270 options are chosen (e.g., -log_view). 271 */ 272 PetscCall(PetscFinalize()); 273 return 0; 274 } 275 /* --------------------------------------------------------------------- */ 276 /* 277 InitialConditions - Computes the solution at the initial time. 278 279 Input Parameter: 280 u - uninitialized solution vector (global) 281 appctx - user-defined application context 282 283 Output Parameter: 284 u - vector with solution at initial time (global) 285 */ 286 PetscErrorCode InitialConditions(Vec u, AppCtx *appctx) { 287 PetscScalar *u_localptr, h = appctx->h; 288 PetscInt i; 289 290 /* 291 Get a pointer to vector data. 292 - For default PETSc vectors, VecGetArray() returns a pointer to 293 the data array. Otherwise, the routine is implementation dependent. 294 - You MUST call VecRestoreArray() when you no longer need access to 295 the array. 296 - Note that the Fortran interface to VecGetArray() differs from the 297 C version. See the users manual for details. 298 */ 299 PetscCall(VecGetArrayWrite(u, &u_localptr)); 300 301 /* 302 We initialize the solution array by simply writing the solution 303 directly into the array locations. Alternatively, we could use 304 VecSetValues() or VecSetValuesLocal(). 305 */ 306 for (i = 0; i < appctx->m; i++) u_localptr[i] = PetscSinScalar(PETSC_PI * i * 6. * h) + 3. * PetscSinScalar(PETSC_PI * i * 2. * h); 307 308 /* 309 Restore vector 310 */ 311 PetscCall(VecRestoreArrayWrite(u, &u_localptr)); 312 313 /* 314 Print debugging information if desired 315 */ 316 if (appctx->debug) { 317 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Initial guess vector\n")); 318 PetscCall(VecView(u, PETSC_VIEWER_STDOUT_SELF)); 319 } 320 321 return 0; 322 } 323 /* --------------------------------------------------------------------- */ 324 /* 325 ExactSolution - Computes the exact solution at a given time. 326 327 Input Parameters: 328 t - current time 329 solution - vector in which exact solution will be computed 330 appctx - user-defined application context 331 332 Output Parameter: 333 solution - vector with the newly computed exact solution 334 */ 335 PetscErrorCode ExactSolution(PetscReal t, Vec solution, AppCtx *appctx) { 336 PetscScalar *s_localptr, h = appctx->h, ex1, ex2, sc1, sc2, tc = t; 337 PetscInt i; 338 339 /* 340 Get a pointer to vector data. 341 */ 342 PetscCall(VecGetArrayWrite(solution, &s_localptr)); 343 344 /* 345 Simply write the solution directly into the array locations. 346 Alternatively, we culd use VecSetValues() or VecSetValuesLocal(). 347 */ 348 ex1 = PetscExpScalar(-36. * PETSC_PI * PETSC_PI * tc); 349 ex2 = PetscExpScalar(-4. * PETSC_PI * PETSC_PI * tc); 350 sc1 = PETSC_PI * 6. * h; 351 sc2 = PETSC_PI * 2. * h; 352 for (i = 0; i < appctx->m; i++) s_localptr[i] = PetscSinScalar(sc1 * (PetscReal)i) * ex1 + 3. * PetscSinScalar(sc2 * (PetscReal)i) * ex2; 353 354 /* 355 Restore vector 356 */ 357 PetscCall(VecRestoreArrayWrite(solution, &s_localptr)); 358 return 0; 359 } 360 /* --------------------------------------------------------------------- */ 361 /* 362 Monitor - User-provided routine to monitor the solution computed at 363 each timestep. This example plots the solution and computes the 364 error in two different norms. 365 366 This example also demonstrates changing the timestep via TSSetTimeStep(). 367 368 Input Parameters: 369 ts - the timestep context 370 step - the count of the current step (with 0 meaning the 371 initial condition) 372 time - the current time 373 u - the solution at this timestep 374 ctx - the user-provided context for this monitoring routine. 375 In this case we use the application context which contains 376 information about the problem size, workspace and the exact 377 solution. 378 */ 379 PetscErrorCode Monitor(TS ts, PetscInt step, PetscReal time, Vec u, void *ctx) { 380 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 381 PetscReal norm_2, norm_max, dt, dttol; 382 383 /* 384 View a graph of the current iterate 385 */ 386 PetscCall(VecView(u, appctx->viewer2)); 387 388 /* 389 Compute the exact solution 390 */ 391 PetscCall(ExactSolution(time, appctx->solution, appctx)); 392 393 /* 394 Print debugging information if desired 395 */ 396 if (appctx->debug) { 397 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Computed solution vector\n")); 398 PetscCall(VecView(u, PETSC_VIEWER_STDOUT_SELF)); 399 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Exact solution vector\n")); 400 PetscCall(VecView(appctx->solution, PETSC_VIEWER_STDOUT_SELF)); 401 } 402 403 /* 404 Compute the 2-norm and max-norm of the error 405 */ 406 PetscCall(VecAXPY(appctx->solution, -1.0, u)); 407 PetscCall(VecNorm(appctx->solution, NORM_2, &norm_2)); 408 norm_2 = PetscSqrtReal(appctx->h) * norm_2; 409 PetscCall(VecNorm(appctx->solution, NORM_MAX, &norm_max)); 410 411 PetscCall(TSGetTimeStep(ts, &dt)); 412 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep %3" PetscInt_FMT ": step size = %g, time = %g, 2-norm error = %g, max norm error = %g\n", step, (double)dt, (double)time, (double)norm_2, (double)norm_max)); 413 414 appctx->norm_2 += norm_2; 415 appctx->norm_max += norm_max; 416 417 dttol = .0001; 418 PetscCall(PetscOptionsGetReal(NULL, NULL, "-dttol", &dttol, NULL)); 419 if (dt < dttol) { 420 dt *= .999; 421 PetscCall(TSSetTimeStep(ts, dt)); 422 } 423 424 /* 425 View a graph of the error 426 */ 427 PetscCall(VecView(appctx->solution, appctx->viewer1)); 428 429 /* 430 Print debugging information if desired 431 */ 432 if (appctx->debug) { 433 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Error vector\n")); 434 PetscCall(VecView(appctx->solution, PETSC_VIEWER_STDOUT_SELF)); 435 } 436 437 return 0; 438 } 439 /* --------------------------------------------------------------------- */ 440 /* 441 RHSMatrixHeat - User-provided routine to compute the right-hand-side 442 matrix for the heat equation. 443 444 Input Parameters: 445 ts - the TS context 446 t - current time 447 global_in - global input vector 448 dummy - optional user-defined context, as set by TSetRHSJacobian() 449 450 Output Parameters: 451 AA - Jacobian matrix 452 BB - optionally different preconditioning matrix 453 str - flag indicating matrix structure 454 455 Notes: 456 Recall that MatSetValues() uses 0-based row and column numbers 457 in Fortran as well as in C. 458 */ 459 PetscErrorCode RHSMatrixHeat(TS ts, PetscReal t, Vec X, Mat AA, Mat BB, void *ctx) { 460 Mat A = AA; /* Jacobian matrix */ 461 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 462 PetscInt mstart = 0; 463 PetscInt mend = appctx->m; 464 PetscInt i, idx[3]; 465 PetscScalar v[3], stwo = -2. / (appctx->h * appctx->h), sone = -.5 * stwo; 466 467 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 468 Compute entries for the locally owned part of the matrix 469 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 470 /* 471 Set matrix rows corresponding to boundary data 472 */ 473 474 mstart = 0; 475 v[0] = 1.0; 476 PetscCall(MatSetValues(A, 1, &mstart, 1, &mstart, v, INSERT_VALUES)); 477 mstart++; 478 479 mend--; 480 v[0] = 1.0; 481 PetscCall(MatSetValues(A, 1, &mend, 1, &mend, v, INSERT_VALUES)); 482 483 /* 484 Set matrix rows corresponding to interior data. We construct the 485 matrix one row at a time. 486 */ 487 v[0] = sone; 488 v[1] = stwo; 489 v[2] = sone; 490 for (i = mstart; i < mend; i++) { 491 idx[0] = i - 1; 492 idx[1] = i; 493 idx[2] = i + 1; 494 PetscCall(MatSetValues(A, 1, &i, 3, idx, v, INSERT_VALUES)); 495 } 496 497 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 498 Complete the matrix assembly process and set some options 499 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 500 /* 501 Assemble matrix, using the 2-step process: 502 MatAssemblyBegin(), MatAssemblyEnd() 503 Computations can be done while messages are in transition 504 by placing code between these two statements. 505 */ 506 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 507 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 508 509 /* 510 Set and option to indicate that we will never add a new nonzero location 511 to the matrix. If we do, it will generate an error. 512 */ 513 PetscCall(MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_TRUE)); 514 515 return 0; 516 } 517 518 PetscErrorCode IFunctionHeat(TS ts, PetscReal t, Vec X, Vec Xdot, Vec r, void *ctx) { 519 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 520 521 PetscCall(MatMult(appctx->A, X, r)); 522 PetscCall(VecAYPX(r, -1.0, Xdot)); 523 return 0; 524 } 525 526 PetscErrorCode IJacobianHeat(TS ts, PetscReal t, Vec X, Vec Xdot, PetscReal s, Mat A, Mat B, void *ctx) { 527 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 528 529 if (appctx->oshift == s) return 0; 530 PetscCall(MatCopy(appctx->A, A, SAME_NONZERO_PATTERN)); 531 PetscCall(MatScale(A, -1)); 532 PetscCall(MatShift(A, s)); 533 PetscCall(MatCopy(A, B, SAME_NONZERO_PATTERN)); 534 appctx->oshift = s; 535 return 0; 536 } 537 538 /*TEST 539 540 test: 541 args: -nox -ts_type ssp -ts_dt 0.0005 542 543 test: 544 suffix: 2 545 args: -nox -ts_type ssp -ts_dt 0.0005 -time_dependent_rhs 1 546 547 test: 548 suffix: 3 549 args: -nox -ts_type rosw -ts_max_steps 3 -ksp_converged_reason 550 filter: sed "s/ATOL/RTOL/g" 551 requires: !single 552 553 test: 554 suffix: 4 555 args: -nox -ts_type beuler -ts_max_steps 3 -ksp_converged_reason 556 filter: sed "s/ATOL/RTOL/g" 557 558 test: 559 suffix: 5 560 args: -nox -ts_type beuler -ts_max_steps 3 -ksp_converged_reason -time_dependent_rhs 561 filter: sed "s/ATOL/RTOL/g" 562 563 test: 564 requires: !single 565 suffix: pod_guess 566 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type pod -pc_type none -ksp_converged_reason 567 568 test: 569 requires: !single 570 suffix: pod_guess_Ainner 571 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type pod -ksp_guess_pod_Ainner -pc_type none -ksp_converged_reason 572 573 test: 574 requires: !single 575 suffix: fischer_guess 576 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type fischer -pc_type none -ksp_converged_reason 577 578 test: 579 requires: !single 580 suffix: fischer_guess_2 581 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type fischer -ksp_guess_fischer_model 2,10 -pc_type none -ksp_converged_reason 582 583 test: 584 requires: !single 585 suffix: fischer_guess_3 586 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type fischer -ksp_guess_fischer_model 3,10 -pc_type none -ksp_converged_reason 587 588 test: 589 requires: !single 590 suffix: stringview 591 args: -nox -ts_type rosw -test_string_viewer 592 593 test: 594 requires: !single 595 suffix: stringview_euler 596 args: -nox -ts_type euler -test_string_viewer 597 598 TEST*/ 599