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 { 82 AppCtx appctx; /* user-defined application context */ 83 TS ts; /* timestepping context */ 84 Mat A; /* matrix data structure */ 85 Vec u; /* approximate solution vector */ 86 PetscReal time_total_max = 100.0; /* default max total time */ 87 PetscInt time_steps_max = 100; /* default max timesteps */ 88 PetscDraw draw; /* drawing context */ 89 PetscInt steps,m; 90 PetscMPIInt size; 91 PetscReal dt; 92 PetscBool flg,flg_string; 93 94 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 95 Initialize program and set problem parameters 96 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 97 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) { 148 PetscCall(TSMonitorSet(ts,Monitor,&appctx,NULL)); 149 } 150 151 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 152 153 Create matrix data structure; set matrix evaluation routine. 154 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 155 156 PetscCall(MatCreate(PETSC_COMM_SELF,&A)); 157 PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,m)); 158 PetscCall(MatSetFromOptions(A)); 159 PetscCall(MatSetUp(A)); 160 161 flg = PETSC_FALSE; 162 PetscCall(PetscOptionsGetBool(NULL,NULL,"-use_ifunc",&flg,NULL)); 163 if (!flg) { 164 appctx.A = NULL; 165 PetscCall(PetscOptionsGetBool(NULL,NULL,"-time_dependent_rhs",&flg,NULL)); 166 if (flg) { 167 /* 168 For linear problems with a time-dependent f(u,t) in the equation 169 u_t = f(u,t), the user provides the discretized right-hand-side 170 as a time-dependent matrix. 171 */ 172 PetscCall(TSSetRHSFunction(ts,NULL,TSComputeRHSFunctionLinear,&appctx)); 173 PetscCall(TSSetRHSJacobian(ts,A,A,RHSMatrixHeat,&appctx)); 174 } else { 175 /* 176 For linear problems with a time-independent f(u) in the equation 177 u_t = f(u), the user provides the discretized right-hand-side 178 as a matrix only once, and then sets the special Jacobian evaluation 179 routine TSComputeRHSJacobianConstant() which will NOT recompute the Jacobian. 180 */ 181 PetscCall(RHSMatrixHeat(ts,0.0,u,A,A,&appctx)); 182 PetscCall(TSSetRHSFunction(ts,NULL,TSComputeRHSFunctionLinear,&appctx)); 183 PetscCall(TSSetRHSJacobian(ts,A,A,TSComputeRHSJacobianConstant,&appctx)); 184 } 185 } else { 186 Mat J; 187 188 PetscCall(RHSMatrixHeat(ts,0.0,u,A,A,&appctx)); 189 PetscCall(MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,&J)); 190 PetscCall(TSSetIFunction(ts,NULL,IFunctionHeat,&appctx)); 191 PetscCall(TSSetIJacobian(ts,J,J,IJacobianHeat,&appctx)); 192 PetscCall(MatDestroy(&J)); 193 194 PetscCall(PetscObjectReference((PetscObject)A)); 195 appctx.A = A; 196 appctx.oshift = PETSC_MIN_REAL; 197 } 198 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 199 Set solution vector and initial timestep 200 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 201 202 dt = appctx.h*appctx.h/2.0; 203 PetscCall(TSSetTimeStep(ts,dt)); 204 205 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 206 Customize timestepping solver: 207 - Set the solution method to be the Backward Euler method. 208 - Set timestepping duration info 209 Then set runtime options, which can override these defaults. 210 For example, 211 -ts_max_steps <maxsteps> -ts_max_time <maxtime> 212 to override the defaults set by TSSetMaxSteps()/TSSetMaxTime(). 213 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 214 215 PetscCall(TSSetMaxSteps(ts,time_steps_max)); 216 PetscCall(TSSetMaxTime(ts,time_total_max)); 217 PetscCall(TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER)); 218 PetscCall(TSSetFromOptions(ts)); 219 220 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 221 Solve the problem 222 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 223 224 /* 225 Evaluate initial conditions 226 */ 227 PetscCall(InitialConditions(u,&appctx)); 228 229 /* 230 Run the timestepping solver 231 */ 232 PetscCall(TSSolve(ts,u)); 233 PetscCall(TSGetStepNumber(ts,&steps)); 234 235 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 236 View timestepping solver info 237 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 238 239 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))); 240 if (!flg_string) { 241 PetscCall(TSView(ts,PETSC_VIEWER_STDOUT_SELF)); 242 } else { 243 PetscViewer stringviewer; 244 char string[512]; 245 const char *outstring; 246 247 PetscCall(PetscViewerStringOpen(PETSC_COMM_WORLD,string,sizeof(string),&stringviewer)); 248 PetscCall(TSView(ts,stringviewer)); 249 PetscCall(PetscViewerStringGetStringRead(stringviewer,&outstring,NULL)); 250 PetscCheck((char*)outstring == (char*)string,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"String returned from viewer does not equal original string"); 251 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Output from string viewer:%s\n",outstring)); 252 PetscCall(PetscViewerDestroy(&stringviewer)); 253 } 254 255 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 256 Free work space. All PETSc objects should be destroyed when they 257 are no longer needed. 258 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 259 260 PetscCall(TSDestroy(&ts)); 261 PetscCall(MatDestroy(&A)); 262 PetscCall(VecDestroy(&u)); 263 PetscCall(PetscViewerDestroy(&appctx.viewer1)); 264 PetscCall(PetscViewerDestroy(&appctx.viewer2)); 265 PetscCall(VecDestroy(&appctx.solution)); 266 PetscCall(MatDestroy(&appctx.A)); 267 268 /* 269 Always call PetscFinalize() before exiting a program. This routine 270 - finalizes the PETSc libraries as well as MPI 271 - provides summary and diagnostic information if certain runtime 272 options are chosen (e.g., -log_view). 273 */ 274 PetscCall(PetscFinalize()); 275 return 0; 276 } 277 /* --------------------------------------------------------------------- */ 278 /* 279 InitialConditions - Computes the solution at the initial time. 280 281 Input Parameter: 282 u - uninitialized solution vector (global) 283 appctx - user-defined application context 284 285 Output Parameter: 286 u - vector with solution at initial time (global) 287 */ 288 PetscErrorCode InitialConditions(Vec u,AppCtx *appctx) 289 { 290 PetscScalar *u_localptr,h = appctx->h; 291 PetscInt i; 292 293 /* 294 Get a pointer to vector data. 295 - For default PETSc vectors, VecGetArray() returns a pointer to 296 the data array. Otherwise, the routine is implementation dependent. 297 - You MUST call VecRestoreArray() when you no longer need access to 298 the array. 299 - Note that the Fortran interface to VecGetArray() differs from the 300 C version. See the users manual for details. 301 */ 302 PetscCall(VecGetArrayWrite(u,&u_localptr)); 303 304 /* 305 We initialize the solution array by simply writing the solution 306 directly into the array locations. Alternatively, we could use 307 VecSetValues() or VecSetValuesLocal(). 308 */ 309 for (i=0; i<appctx->m; i++) u_localptr[i] = PetscSinScalar(PETSC_PI*i*6.*h) + 3.*PetscSinScalar(PETSC_PI*i*2.*h); 310 311 /* 312 Restore vector 313 */ 314 PetscCall(VecRestoreArrayWrite(u,&u_localptr)); 315 316 /* 317 Print debugging information if desired 318 */ 319 if (appctx->debug) { 320 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Initial guess vector\n")); 321 PetscCall(VecView(u,PETSC_VIEWER_STDOUT_SELF)); 322 } 323 324 return 0; 325 } 326 /* --------------------------------------------------------------------- */ 327 /* 328 ExactSolution - Computes the exact solution at a given time. 329 330 Input Parameters: 331 t - current time 332 solution - vector in which exact solution will be computed 333 appctx - user-defined application context 334 335 Output Parameter: 336 solution - vector with the newly computed exact solution 337 */ 338 PetscErrorCode ExactSolution(PetscReal t,Vec solution,AppCtx *appctx) 339 { 340 PetscScalar *s_localptr,h = appctx->h,ex1,ex2,sc1,sc2,tc = t; 341 PetscInt i; 342 343 /* 344 Get a pointer to vector data. 345 */ 346 PetscCall(VecGetArrayWrite(solution,&s_localptr)); 347 348 /* 349 Simply write the solution directly into the array locations. 350 Alternatively, we culd use VecSetValues() or VecSetValuesLocal(). 351 */ 352 ex1 = PetscExpScalar(-36.*PETSC_PI*PETSC_PI*tc); 353 ex2 = PetscExpScalar(-4.*PETSC_PI*PETSC_PI*tc); 354 sc1 = PETSC_PI*6.*h; sc2 = PETSC_PI*2.*h; 355 for (i=0; i<appctx->m; i++) s_localptr[i] = PetscSinScalar(sc1*(PetscReal)i)*ex1 + 3.*PetscSinScalar(sc2*(PetscReal)i)*ex2; 356 357 /* 358 Restore vector 359 */ 360 PetscCall(VecRestoreArrayWrite(solution,&s_localptr)); 361 return 0; 362 } 363 /* --------------------------------------------------------------------- */ 364 /* 365 Monitor - User-provided routine to monitor the solution computed at 366 each timestep. This example plots the solution and computes the 367 error in two different norms. 368 369 This example also demonstrates changing the timestep via TSSetTimeStep(). 370 371 Input Parameters: 372 ts - the timestep context 373 step - the count of the current step (with 0 meaning the 374 initial condition) 375 time - the current time 376 u - the solution at this timestep 377 ctx - the user-provided context for this monitoring routine. 378 In this case we use the application context which contains 379 information about the problem size, workspace and the exact 380 solution. 381 */ 382 PetscErrorCode Monitor(TS ts,PetscInt step,PetscReal time,Vec u,void *ctx) 383 { 384 AppCtx *appctx = (AppCtx*) ctx; /* user-defined application context */ 385 PetscReal norm_2,norm_max,dt,dttol; 386 387 /* 388 View a graph of the current iterate 389 */ 390 PetscCall(VecView(u,appctx->viewer2)); 391 392 /* 393 Compute the exact solution 394 */ 395 PetscCall(ExactSolution(time,appctx->solution,appctx)); 396 397 /* 398 Print debugging information if desired 399 */ 400 if (appctx->debug) { 401 PetscCall(PetscPrintf(PETSC_COMM_SELF,"Computed solution vector\n")); 402 PetscCall(VecView(u,PETSC_VIEWER_STDOUT_SELF)); 403 PetscCall(PetscPrintf(PETSC_COMM_SELF,"Exact solution vector\n")); 404 PetscCall(VecView(appctx->solution,PETSC_VIEWER_STDOUT_SELF)); 405 } 406 407 /* 408 Compute the 2-norm and max-norm of the error 409 */ 410 PetscCall(VecAXPY(appctx->solution,-1.0,u)); 411 PetscCall(VecNorm(appctx->solution,NORM_2,&norm_2)); 412 norm_2 = PetscSqrtReal(appctx->h)*norm_2; 413 PetscCall(VecNorm(appctx->solution,NORM_MAX,&norm_max)); 414 415 PetscCall(TSGetTimeStep(ts,&dt)); 416 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)); 417 418 appctx->norm_2 += norm_2; 419 appctx->norm_max += norm_max; 420 421 dttol = .0001; 422 PetscCall(PetscOptionsGetReal(NULL,NULL,"-dttol",&dttol,NULL)); 423 if (dt < dttol) { 424 dt *= .999; 425 PetscCall(TSSetTimeStep(ts,dt)); 426 } 427 428 /* 429 View a graph of the error 430 */ 431 PetscCall(VecView(appctx->solution,appctx->viewer1)); 432 433 /* 434 Print debugging information if desired 435 */ 436 if (appctx->debug) { 437 PetscCall(PetscPrintf(PETSC_COMM_SELF,"Error vector\n")); 438 PetscCall(VecView(appctx->solution,PETSC_VIEWER_STDOUT_SELF)); 439 } 440 441 return 0; 442 } 443 /* --------------------------------------------------------------------- */ 444 /* 445 RHSMatrixHeat - User-provided routine to compute the right-hand-side 446 matrix for the heat equation. 447 448 Input Parameters: 449 ts - the TS context 450 t - current time 451 global_in - global input vector 452 dummy - optional user-defined context, as set by TSetRHSJacobian() 453 454 Output Parameters: 455 AA - Jacobian matrix 456 BB - optionally different preconditioning matrix 457 str - flag indicating matrix structure 458 459 Notes: 460 Recall that MatSetValues() uses 0-based row and column numbers 461 in Fortran as well as in C. 462 */ 463 PetscErrorCode RHSMatrixHeat(TS ts,PetscReal t,Vec X,Mat AA,Mat BB,void *ctx) 464 { 465 Mat A = AA; /* Jacobian matrix */ 466 AppCtx *appctx = (AppCtx*)ctx; /* user-defined application context */ 467 PetscInt mstart = 0; 468 PetscInt mend = appctx->m; 469 PetscInt i,idx[3]; 470 PetscScalar v[3],stwo = -2./(appctx->h*appctx->h),sone = -.5*stwo; 471 472 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 473 Compute entries for the locally owned part of the matrix 474 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 475 /* 476 Set matrix rows corresponding to boundary data 477 */ 478 479 mstart = 0; 480 v[0] = 1.0; 481 PetscCall(MatSetValues(A,1,&mstart,1,&mstart,v,INSERT_VALUES)); 482 mstart++; 483 484 mend--; 485 v[0] = 1.0; 486 PetscCall(MatSetValues(A,1,&mend,1,&mend,v,INSERT_VALUES)); 487 488 /* 489 Set matrix rows corresponding to interior data. We construct the 490 matrix one row at a time. 491 */ 492 v[0] = sone; v[1] = stwo; v[2] = sone; 493 for (i=mstart; i<mend; i++) { 494 idx[0] = i-1; idx[1] = i; idx[2] = i+1; 495 PetscCall(MatSetValues(A,1,&i,3,idx,v,INSERT_VALUES)); 496 } 497 498 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 499 Complete the matrix assembly process and set some options 500 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 501 /* 502 Assemble matrix, using the 2-step process: 503 MatAssemblyBegin(), MatAssemblyEnd() 504 Computations can be done while messages are in transition 505 by placing code between these two statements. 506 */ 507 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 508 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 509 510 /* 511 Set and option to indicate that we will never add a new nonzero location 512 to the matrix. If we do, it will generate an error. 513 */ 514 PetscCall(MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_TRUE)); 515 516 return 0; 517 } 518 519 PetscErrorCode IFunctionHeat(TS ts,PetscReal t,Vec X,Vec Xdot,Vec r,void *ctx) 520 { 521 AppCtx *appctx = (AppCtx*)ctx; /* user-defined application context */ 522 523 PetscCall(MatMult(appctx->A,X,r)); 524 PetscCall(VecAYPX(r,-1.0,Xdot)); 525 return 0; 526 } 527 528 PetscErrorCode IJacobianHeat(TS ts,PetscReal t,Vec X,Vec Xdot,PetscReal s,Mat A,Mat B,void *ctx) 529 { 530 AppCtx *appctx = (AppCtx*)ctx; /* user-defined application context */ 531 532 if (appctx->oshift == s) return 0; 533 PetscCall(MatCopy(appctx->A,A,SAME_NONZERO_PATTERN)); 534 PetscCall(MatScale(A,-1)); 535 PetscCall(MatShift(A,s)); 536 PetscCall(MatCopy(A,B,SAME_NONZERO_PATTERN)); 537 appctx->oshift = s; 538 return 0; 539 } 540 541 /*TEST 542 543 test: 544 args: -nox -ts_type ssp -ts_dt 0.0005 545 546 test: 547 suffix: 2 548 args: -nox -ts_type ssp -ts_dt 0.0005 -time_dependent_rhs 1 549 550 test: 551 suffix: 3 552 args: -nox -ts_type rosw -ts_max_steps 3 -ksp_converged_reason 553 filter: sed "s/ATOL/RTOL/g" 554 requires: !single 555 556 test: 557 suffix: 4 558 args: -nox -ts_type beuler -ts_max_steps 3 -ksp_converged_reason 559 filter: sed "s/ATOL/RTOL/g" 560 561 test: 562 suffix: 5 563 args: -nox -ts_type beuler -ts_max_steps 3 -ksp_converged_reason -time_dependent_rhs 564 filter: sed "s/ATOL/RTOL/g" 565 566 test: 567 requires: !single 568 suffix: pod_guess 569 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type pod -pc_type none -ksp_converged_reason 570 571 test: 572 requires: !single 573 suffix: pod_guess_Ainner 574 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 575 576 test: 577 requires: !single 578 suffix: fischer_guess 579 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type fischer -pc_type none -ksp_converged_reason 580 581 test: 582 requires: !single 583 suffix: fischer_guess_2 584 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 585 586 test: 587 requires: !single 588 suffix: fischer_guess_3 589 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 590 591 test: 592 requires: !single 593 suffix: stringview 594 args: -nox -ts_type rosw -test_string_viewer 595 596 test: 597 requires: !single 598 suffix: stringview_euler 599 args: -nox -ts_type euler -test_string_viewer 600 601 TEST*/ 602