1 static char help[] = "Solves a simple time-dependent linear PDE (the heat equation).\n\ 2 Input parameters include:\n\ 3 -m <points>, where <points> = number of grid points\n\ 4 -time_dependent_rhs : Treat the problem as having a time-dependent right-hand side\n\ 5 -use_ifunc : Use IFunction/IJacobian interface\n\ 6 -debug : Activate debugging printouts\n\ 7 -nox : Deactivate x-window graphics\n\n"; 8 9 /* ------------------------------------------------------------------------ 10 11 This program solves the one-dimensional heat equation (also called the 12 diffusion equation), 13 u_t = u_xx, 14 on the domain 0 <= x <= 1, with the boundary conditions 15 u(t,0) = 0, u(t,1) = 0, 16 and the initial condition 17 u(0,x) = sin(6*pi*x) + 3*sin(2*pi*x). 18 This is a linear, second-order, parabolic equation. 19 20 We discretize the right-hand side using finite differences with 21 uniform grid spacing h: 22 u_xx = (u_{i+1} - 2u_{i} + u_{i-1})/(h^2) 23 We then demonstrate time evolution using the various TS methods by 24 running the program via 25 ex3 -ts_type <timestepping solver> 26 27 We compare the approximate solution with the exact solution, given by 28 u_exact(x,t) = exp(-36*pi*pi*t) * sin(6*pi*x) + 29 3*exp(-4*pi*pi*t) * sin(2*pi*x) 30 31 Notes: 32 This code demonstrates the TS solver interface to two variants of 33 linear problems, u_t = f(u,t), namely 34 - time-dependent f: f(u,t) is a function of t 35 - time-independent f: f(u,t) is simply f(u) 36 37 The parallel version of this code is ts/tutorials/ex4.c 38 39 ------------------------------------------------------------------------- */ 40 41 /* 42 Include "petscts.h" so that we can use TS solvers. Note that this file 43 automatically includes: 44 petscsys.h - base PETSc routines petscvec.h - vectors 45 petscmat.h - matrices 46 petscis.h - index sets petscksp.h - Krylov subspace methods 47 petscviewer.h - viewers petscpc.h - preconditioners 48 petscksp.h - linear solvers petscsnes.h - nonlinear solvers 49 */ 50 51 #include <petscts.h> 52 #include <petscdraw.h> 53 54 /* 55 User-defined application context - contains data needed by the 56 application-provided call-back routines. 57 */ 58 typedef struct { 59 Vec solution; /* global exact solution vector */ 60 PetscInt m; /* total number of grid points */ 61 PetscReal h; /* mesh width h = 1/(m-1) */ 62 PetscBool debug; /* flag (1 indicates activation of debugging printouts) */ 63 PetscViewer viewer1, viewer2; /* viewers for the solution and error */ 64 PetscReal norm_2, norm_max; /* error norms */ 65 Mat A; /* RHS mat, used with IFunction interface */ 66 PetscReal oshift; /* old shift applied, prevent to recompute the IJacobian */ 67 } AppCtx; 68 69 /* 70 User-defined routines 71 */ 72 extern PetscErrorCode InitialConditions(Vec, AppCtx *); 73 extern PetscErrorCode RHSMatrixHeat(TS, PetscReal, Vec, Mat, Mat, void *); 74 extern PetscErrorCode IFunctionHeat(TS, PetscReal, Vec, Vec, Vec, void *); 75 extern PetscErrorCode IJacobianHeat(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *); 76 extern PetscErrorCode Monitor(TS, PetscInt, PetscReal, Vec, void *); 77 extern PetscErrorCode ExactSolution(PetscReal, Vec, AppCtx *); 78 79 int main(int argc, char **argv) 80 { 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 { 288 PetscScalar *u_localptr, h = appctx->h; 289 PetscInt i; 290 291 PetscFunctionBeginUser; 292 /* 293 Get a pointer to vector data. 294 - For default PETSc vectors, VecGetArray() returns a pointer to 295 the data array. Otherwise, the routine is implementation dependent. 296 - You MUST call VecRestoreArray() when you no longer need access to 297 the array. 298 - Note that the Fortran interface to VecGetArray() differs from the 299 C version. See the users manual for details. 300 */ 301 PetscCall(VecGetArrayWrite(u, &u_localptr)); 302 303 /* 304 We initialize the solution array by simply writing the solution 305 directly into the array locations. Alternatively, we could use 306 VecSetValues() or VecSetValuesLocal(). 307 */ 308 for (i = 0; i < appctx->m; i++) u_localptr[i] = PetscSinScalar(PETSC_PI * i * 6. * h) + 3. * PetscSinScalar(PETSC_PI * i * 2. * h); 309 310 /* 311 Restore vector 312 */ 313 PetscCall(VecRestoreArrayWrite(u, &u_localptr)); 314 315 /* 316 Print debugging information if desired 317 */ 318 if (appctx->debug) { 319 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Initial guess vector\n")); 320 PetscCall(VecView(u, PETSC_VIEWER_STDOUT_SELF)); 321 } 322 323 PetscFunctionReturn(PETSC_SUCCESS); 324 } 325 /* --------------------------------------------------------------------- */ 326 /* 327 ExactSolution - Computes the exact solution at a given time. 328 329 Input Parameters: 330 t - current time 331 solution - vector in which exact solution will be computed 332 appctx - user-defined application context 333 334 Output Parameter: 335 solution - vector with the newly computed exact solution 336 */ 337 PetscErrorCode ExactSolution(PetscReal t, Vec solution, AppCtx *appctx) 338 { 339 PetscScalar *s_localptr, h = appctx->h, ex1, ex2, sc1, sc2, tc = t; 340 PetscInt i; 341 342 PetscFunctionBeginUser; 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; 355 sc2 = PETSC_PI * 2. * h; 356 for (i = 0; i < appctx->m; i++) s_localptr[i] = PetscSinScalar(sc1 * (PetscReal)i) * ex1 + 3. * PetscSinScalar(sc2 * (PetscReal)i) * ex2; 357 358 /* 359 Restore vector 360 */ 361 PetscCall(VecRestoreArrayWrite(solution, &s_localptr)); 362 PetscFunctionReturn(PETSC_SUCCESS); 363 } 364 /* --------------------------------------------------------------------- */ 365 /* 366 Monitor - User-provided routine to monitor the solution computed at 367 each timestep. This example plots the solution and computes the 368 error in two different norms. 369 370 This example also demonstrates changing the timestep via TSSetTimeStep(). 371 372 Input Parameters: 373 ts - the timestep context 374 step - the count of the current step (with 0 meaning the 375 initial condition) 376 time - the current time 377 u - the solution at this timestep 378 ctx - the user-provided context for this monitoring routine. 379 In this case we use the application context which contains 380 information about the problem size, workspace and the exact 381 solution. 382 */ 383 PetscErrorCode Monitor(TS ts, PetscInt step, PetscReal time, Vec u, void *ctx) 384 { 385 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 386 PetscReal norm_2, norm_max, dt, dttol; 387 388 PetscFunctionBeginUser; 389 /* 390 View a graph of the current iterate 391 */ 392 PetscCall(VecView(u, appctx->viewer2)); 393 394 /* 395 Compute the exact solution 396 */ 397 PetscCall(ExactSolution(time, appctx->solution, appctx)); 398 399 /* 400 Print debugging information if desired 401 */ 402 if (appctx->debug) { 403 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Computed solution vector\n")); 404 PetscCall(VecView(u, PETSC_VIEWER_STDOUT_SELF)); 405 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Exact solution vector\n")); 406 PetscCall(VecView(appctx->solution, PETSC_VIEWER_STDOUT_SELF)); 407 } 408 409 /* 410 Compute the 2-norm and max-norm of the error 411 */ 412 PetscCall(VecAXPY(appctx->solution, -1.0, u)); 413 PetscCall(VecNorm(appctx->solution, NORM_2, &norm_2)); 414 norm_2 = PetscSqrtReal(appctx->h) * norm_2; 415 PetscCall(VecNorm(appctx->solution, NORM_MAX, &norm_max)); 416 417 PetscCall(TSGetTimeStep(ts, &dt)); 418 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)); 419 420 appctx->norm_2 += norm_2; 421 appctx->norm_max += norm_max; 422 423 dttol = .0001; 424 PetscCall(PetscOptionsGetReal(NULL, NULL, "-dttol", &dttol, NULL)); 425 if (dt < dttol) { 426 dt *= .999; 427 PetscCall(TSSetTimeStep(ts, dt)); 428 } 429 430 /* 431 View a graph of the error 432 */ 433 PetscCall(VecView(appctx->solution, appctx->viewer1)); 434 435 /* 436 Print debugging information if desired 437 */ 438 if (appctx->debug) { 439 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Error vector\n")); 440 PetscCall(VecView(appctx->solution, PETSC_VIEWER_STDOUT_SELF)); 441 } 442 443 PetscFunctionReturn(PETSC_SUCCESS); 444 } 445 /* --------------------------------------------------------------------- */ 446 /* 447 RHSMatrixHeat - User-provided routine to compute the right-hand-side 448 matrix for the heat equation. 449 450 Input Parameters: 451 ts - the TS context 452 t - current time 453 global_in - global input vector 454 dummy - optional user-defined context, as set by TSetRHSJacobian() 455 456 Output Parameters: 457 AA - Jacobian matrix 458 BB - optionally different preconditioning matrix 459 str - flag indicating matrix structure 460 461 Notes: 462 Recall that MatSetValues() uses 0-based row and column numbers 463 in Fortran as well as in C. 464 */ 465 PetscErrorCode RHSMatrixHeat(TS ts, PetscReal t, Vec X, Mat AA, Mat BB, void *ctx) 466 { 467 Mat A = AA; /* Jacobian matrix */ 468 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 469 PetscInt mstart = 0; 470 PetscInt mend = appctx->m; 471 PetscInt i, idx[3]; 472 PetscScalar v[3], stwo = -2. / (appctx->h * appctx->h), sone = -.5 * stwo; 473 474 PetscFunctionBeginUser; 475 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 476 Compute entries for the locally owned part of the matrix 477 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 478 /* 479 Set matrix rows corresponding to boundary data 480 */ 481 482 mstart = 0; 483 v[0] = 1.0; 484 PetscCall(MatSetValues(A, 1, &mstart, 1, &mstart, v, INSERT_VALUES)); 485 mstart++; 486 487 mend--; 488 v[0] = 1.0; 489 PetscCall(MatSetValues(A, 1, &mend, 1, &mend, v, INSERT_VALUES)); 490 491 /* 492 Set matrix rows corresponding to interior data. We construct the 493 matrix one row at a time. 494 */ 495 v[0] = sone; 496 v[1] = stwo; 497 v[2] = sone; 498 for (i = mstart; i < mend; i++) { 499 idx[0] = i - 1; 500 idx[1] = i; 501 idx[2] = i + 1; 502 PetscCall(MatSetValues(A, 1, &i, 3, idx, v, INSERT_VALUES)); 503 } 504 505 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 506 Complete the matrix assembly process and set some options 507 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 508 /* 509 Assemble matrix, using the 2-step process: 510 MatAssemblyBegin(), MatAssemblyEnd() 511 Computations can be done while messages are in transition 512 by placing code between these two statements. 513 */ 514 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 515 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 516 517 /* 518 Set and option to indicate that we will never add a new nonzero location 519 to the matrix. If we do, it will generate an error. 520 */ 521 PetscCall(MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_TRUE)); 522 523 PetscFunctionReturn(PETSC_SUCCESS); 524 } 525 526 PetscErrorCode IFunctionHeat(TS ts, PetscReal t, Vec X, Vec Xdot, Vec r, void *ctx) 527 { 528 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 529 530 PetscFunctionBeginUser; 531 PetscCall(MatMult(appctx->A, X, r)); 532 PetscCall(VecAYPX(r, -1.0, Xdot)); 533 PetscFunctionReturn(PETSC_SUCCESS); 534 } 535 536 PetscErrorCode IJacobianHeat(TS ts, PetscReal t, Vec X, Vec Xdot, PetscReal s, Mat A, Mat B, void *ctx) 537 { 538 AppCtx *appctx = (AppCtx *)ctx; /* user-defined application context */ 539 540 PetscFunctionBeginUser; 541 if (appctx->oshift == s) PetscFunctionReturn(PETSC_SUCCESS); 542 PetscCall(MatCopy(appctx->A, A, SAME_NONZERO_PATTERN)); 543 PetscCall(MatScale(A, -1)); 544 PetscCall(MatShift(A, s)); 545 PetscCall(MatCopy(A, B, SAME_NONZERO_PATTERN)); 546 appctx->oshift = s; 547 PetscFunctionReturn(PETSC_SUCCESS); 548 } 549 550 /*TEST 551 552 test: 553 args: -nox -ts_type ssp -ts_dt 0.0005 554 555 test: 556 suffix: 2 557 args: -nox -ts_type ssp -ts_dt 0.0005 -time_dependent_rhs 1 558 559 test: 560 suffix: 3 561 args: -nox -ts_type rosw -ts_max_steps 3 -ksp_converged_reason 562 filter: sed "s/ATOL/RTOL/g" 563 requires: !single 564 565 test: 566 suffix: 4 567 args: -nox -ts_type beuler -ts_max_steps 3 -ksp_converged_reason 568 filter: sed "s/ATOL/RTOL/g" 569 570 test: 571 suffix: 5 572 args: -nox -ts_type beuler -ts_max_steps 3 -ksp_converged_reason -time_dependent_rhs 573 filter: sed "s/ATOL/RTOL/g" 574 575 test: 576 requires: !single 577 suffix: pod_guess 578 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type pod -pc_type none -ksp_converged_reason 579 580 test: 581 requires: !single 582 suffix: pod_guess_Ainner 583 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 584 585 test: 586 requires: !single 587 suffix: fischer_guess 588 args: -nox -ts_type beuler -use_ifunc -ts_dt 0.0005 -ksp_guess_type fischer -pc_type none -ksp_converged_reason 589 590 test: 591 requires: !single 592 suffix: fischer_guess_2 593 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 594 595 test: 596 requires: !single 597 suffix: fischer_guess_3 598 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 599 600 test: 601 requires: !single 602 suffix: stringview 603 args: -nox -ts_type rosw -test_string_viewer 604 605 test: 606 requires: !single 607 suffix: stringview_euler 608 args: -nox -ts_type euler -test_string_viewer 609 610 TEST*/ 611