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