1 /* 2 The Problem: 3 Solve the convection-diffusion equation: 4 5 u_t+a*(u_x+u_y)=epsilon*(u_xx+u_yy) 6 u=0 at x=0, y=0 7 u_x=0 at x=1 8 u_y=0 at y=1 9 u = exp(-20.0*(pow(x-0.5,2.0)+pow(y-0.5,2.0))) at t=0 10 11 This program tests the routine of computing the Jacobian by the 12 finite difference method as well as PETSc. 13 14 */ 15 16 static char help[] = "Solve the convection-diffusion equation. \n\n"; 17 18 #include <petscts.h> 19 20 typedef struct 21 { 22 PetscInt m; /* the number of mesh points in x-direction */ 23 PetscInt n; /* the number of mesh points in y-direction */ 24 PetscReal dx; /* the grid space in x-direction */ 25 PetscReal dy; /* the grid space in y-direction */ 26 PetscReal a; /* the convection coefficient */ 27 PetscReal epsilon; /* the diffusion coefficient */ 28 PetscReal tfinal; 29 } Data; 30 31 extern PetscErrorCode Monitor(TS,PetscInt,PetscReal,Vec,void*); 32 extern PetscErrorCode Initial(Vec,void*); 33 extern PetscErrorCode RHSFunction(TS,PetscReal,Vec,Vec,void*); 34 extern PetscErrorCode RHSJacobian(TS,PetscReal,Vec,Mat,Mat,void*); 35 extern PetscErrorCode PostStep(TS); 36 37 int main(int argc,char **argv) 38 { 39 PetscErrorCode ierr; 40 PetscInt time_steps=100,iout,NOUT=1; 41 Vec global; 42 PetscReal dt,ftime,ftime_original; 43 TS ts; 44 PetscViewer viewfile; 45 Mat J = 0; 46 Vec x; 47 Data data; 48 PetscInt mn; 49 PetscBool flg; 50 MatColoring mc; 51 ISColoring iscoloring; 52 MatFDColoring matfdcoloring = 0; 53 PetscBool fd_jacobian_coloring = PETSC_FALSE; 54 SNES snes; 55 KSP ksp; 56 PC pc; 57 58 ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 59 60 /* set data */ 61 data.m = 9; 62 data.n = 9; 63 data.a = 1.0; 64 data.epsilon = 0.1; 65 data.dx = 1.0/(data.m+1.0); 66 data.dy = 1.0/(data.n+1.0); 67 mn = (data.m)*(data.n); 68 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-time",&time_steps,NULL)); 69 70 /* set initial conditions */ 71 CHKERRQ(VecCreate(PETSC_COMM_WORLD,&global)); 72 CHKERRQ(VecSetSizes(global,PETSC_DECIDE,mn)); 73 CHKERRQ(VecSetFromOptions(global)); 74 CHKERRQ(Initial(global,&data)); 75 CHKERRQ(VecDuplicate(global,&x)); 76 77 /* create timestep context */ 78 CHKERRQ(TSCreate(PETSC_COMM_WORLD,&ts)); 79 CHKERRQ(TSMonitorSet(ts,Monitor,&data,NULL)); 80 CHKERRQ(TSSetType(ts,TSEULER)); 81 dt = 0.1; 82 ftime_original = data.tfinal = 1.0; 83 84 CHKERRQ(TSSetTimeStep(ts,dt)); 85 CHKERRQ(TSSetMaxSteps(ts,time_steps)); 86 CHKERRQ(TSSetMaxTime(ts,ftime_original)); 87 CHKERRQ(TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER)); 88 CHKERRQ(TSSetSolution(ts,global)); 89 90 /* set user provided RHSFunction and RHSJacobian */ 91 CHKERRQ(TSSetRHSFunction(ts,NULL,RHSFunction,&data)); 92 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&J)); 93 CHKERRQ(MatSetSizes(J,PETSC_DECIDE,PETSC_DECIDE,mn,mn)); 94 CHKERRQ(MatSetFromOptions(J)); 95 CHKERRQ(MatSeqAIJSetPreallocation(J,5,NULL)); 96 CHKERRQ(MatMPIAIJSetPreallocation(J,5,NULL,5,NULL)); 97 98 CHKERRQ(PetscOptionsHasName(NULL,NULL,"-ts_fd",&flg)); 99 if (!flg) { 100 CHKERRQ(TSSetRHSJacobian(ts,J,J,RHSJacobian,&data)); 101 } else { 102 CHKERRQ(TSGetSNES(ts,&snes)); 103 CHKERRQ(PetscOptionsHasName(NULL,NULL,"-fd_color",&fd_jacobian_coloring)); 104 if (fd_jacobian_coloring) { /* Use finite differences with coloring */ 105 /* Get data structure of J */ 106 PetscBool pc_diagonal; 107 CHKERRQ(PetscOptionsHasName(NULL,NULL,"-pc_diagonal",&pc_diagonal)); 108 if (pc_diagonal) { /* the preconditioner of J is a diagonal matrix */ 109 PetscInt rstart,rend,i; 110 PetscScalar zero=0.0; 111 CHKERRQ(MatGetOwnershipRange(J,&rstart,&rend)); 112 for (i=rstart; i<rend; i++) { 113 CHKERRQ(MatSetValues(J,1,&i,1,&i,&zero,INSERT_VALUES)); 114 } 115 CHKERRQ(MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY)); 116 CHKERRQ(MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY)); 117 } else { 118 /* Fill the structure using the expensive SNESComputeJacobianDefault. Temporarily set up the TS so we can call this function */ 119 CHKERRQ(TSSetType(ts,TSBEULER)); 120 CHKERRQ(TSSetUp(ts)); 121 CHKERRQ(SNESComputeJacobianDefault(snes,x,J,J,ts)); 122 } 123 124 /* create coloring context */ 125 CHKERRQ(MatColoringCreate(J,&mc)); 126 CHKERRQ(MatColoringSetType(mc,MATCOLORINGSL)); 127 CHKERRQ(MatColoringSetFromOptions(mc)); 128 CHKERRQ(MatColoringApply(mc,&iscoloring)); 129 CHKERRQ(MatColoringDestroy(&mc)); 130 CHKERRQ(MatFDColoringCreate(J,iscoloring,&matfdcoloring)); 131 CHKERRQ(MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode (*)(void))SNESTSFormFunction,ts)); 132 CHKERRQ(MatFDColoringSetFromOptions(matfdcoloring)); 133 CHKERRQ(MatFDColoringSetUp(J,iscoloring,matfdcoloring)); 134 CHKERRQ(SNESSetJacobian(snes,J,J,SNESComputeJacobianDefaultColor,matfdcoloring)); 135 CHKERRQ(ISColoringDestroy(&iscoloring)); 136 } else { /* Use finite differences (slow) */ 137 CHKERRQ(SNESSetJacobian(snes,J,J,SNESComputeJacobianDefault,NULL)); 138 } 139 } 140 141 /* Pick up a Petsc preconditioner */ 142 /* one can always set method or preconditioner during the run time */ 143 CHKERRQ(TSGetSNES(ts,&snes)); 144 CHKERRQ(SNESGetKSP(snes,&ksp)); 145 CHKERRQ(KSPGetPC(ksp,&pc)); 146 CHKERRQ(PCSetType(pc,PCJACOBI)); 147 CHKERRQ(TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER)); 148 149 CHKERRQ(TSSetFromOptions(ts)); 150 CHKERRQ(TSSetUp(ts)); 151 152 /* Test TSSetPostStep() */ 153 CHKERRQ(PetscOptionsHasName(NULL,NULL,"-test_PostStep",&flg)); 154 if (flg) { 155 CHKERRQ(TSSetPostStep(ts,PostStep)); 156 } 157 158 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-NOUT",&NOUT,NULL)); 159 for (iout=1; iout<=NOUT; iout++) { 160 CHKERRQ(TSSetMaxSteps(ts,time_steps)); 161 CHKERRQ(TSSetMaxTime(ts,iout*ftime_original/NOUT)); 162 CHKERRQ(TSSolve(ts,global)); 163 CHKERRQ(TSGetSolveTime(ts,&ftime)); 164 CHKERRQ(TSSetTime(ts,ftime)); 165 CHKERRQ(TSSetTimeStep(ts,dt)); 166 } 167 /* Interpolate solution at tfinal */ 168 CHKERRQ(TSGetSolution(ts,&global)); 169 CHKERRQ(TSInterpolate(ts,ftime_original,global)); 170 171 CHKERRQ(PetscOptionsHasName(NULL,NULL,"-matlab_view",&flg)); 172 if (flg) { /* print solution into a MATLAB file */ 173 CHKERRQ(PetscViewerASCIIOpen(PETSC_COMM_WORLD,"out.m",&viewfile)); 174 CHKERRQ(PetscViewerPushFormat(viewfile,PETSC_VIEWER_ASCII_MATLAB)); 175 CHKERRQ(VecView(global,viewfile)); 176 CHKERRQ(PetscViewerPopFormat(viewfile)); 177 CHKERRQ(PetscViewerDestroy(&viewfile)); 178 } 179 180 /* free the memories */ 181 CHKERRQ(TSDestroy(&ts)); 182 CHKERRQ(VecDestroy(&global)); 183 CHKERRQ(VecDestroy(&x)); 184 CHKERRQ(MatDestroy(&J)); 185 if (fd_jacobian_coloring) CHKERRQ(MatFDColoringDestroy(&matfdcoloring)); 186 ierr = PetscFinalize(); 187 return ierr; 188 } 189 190 /* -------------------------------------------------------------------*/ 191 /* the initial function */ 192 PetscReal f_ini(PetscReal x,PetscReal y) 193 { 194 PetscReal f; 195 196 f=PetscExpReal(-20.0*(PetscPowRealInt(x-0.5,2)+PetscPowRealInt(y-0.5,2))); 197 return f; 198 } 199 200 PetscErrorCode Initial(Vec global,void *ctx) 201 { 202 Data *data = (Data*)ctx; 203 PetscInt m,row,col; 204 PetscReal x,y,dx,dy; 205 PetscScalar *localptr; 206 PetscInt i,mybase,myend,locsize; 207 208 PetscFunctionBeginUser; 209 /* make the local copies of parameters */ 210 m = data->m; 211 dx = data->dx; 212 dy = data->dy; 213 214 /* determine starting point of each processor */ 215 CHKERRQ(VecGetOwnershipRange(global,&mybase,&myend)); 216 CHKERRQ(VecGetLocalSize(global,&locsize)); 217 218 /* Initialize the array */ 219 CHKERRQ(VecGetArrayWrite(global,&localptr)); 220 221 for (i=0; i<locsize; i++) { 222 row = 1+(mybase+i)-((mybase+i)/m)*m; 223 col = (mybase+i)/m+1; 224 x = dx*row; 225 y = dy*col; 226 localptr[i] = f_ini(x,y); 227 } 228 229 CHKERRQ(VecRestoreArrayWrite(global,&localptr)); 230 PetscFunctionReturn(0); 231 } 232 233 PetscErrorCode Monitor(TS ts,PetscInt step,PetscReal time,Vec global,void *ctx) 234 { 235 VecScatter scatter; 236 IS from,to; 237 PetscInt i,n,*idx,nsteps,maxsteps; 238 Vec tmp_vec; 239 const PetscScalar *tmp; 240 241 PetscFunctionBeginUser; 242 CHKERRQ(TSGetStepNumber(ts,&nsteps)); 243 /* display output at selected time steps */ 244 CHKERRQ(TSGetMaxSteps(ts, &maxsteps)); 245 if (nsteps % 10 != 0) PetscFunctionReturn(0); 246 247 /* Get the size of the vector */ 248 CHKERRQ(VecGetSize(global,&n)); 249 250 /* Set the index sets */ 251 CHKERRQ(PetscMalloc1(n,&idx)); 252 for (i=0; i<n; i++) idx[i]=i; 253 254 /* Create local sequential vectors */ 255 CHKERRQ(VecCreateSeq(PETSC_COMM_SELF,n,&tmp_vec)); 256 257 /* Create scatter context */ 258 CHKERRQ(ISCreateGeneral(PETSC_COMM_SELF,n,idx,PETSC_COPY_VALUES,&from)); 259 CHKERRQ(ISCreateGeneral(PETSC_COMM_SELF,n,idx,PETSC_COPY_VALUES,&to)); 260 CHKERRQ(VecScatterCreate(global,from,tmp_vec,to,&scatter)); 261 CHKERRQ(VecScatterBegin(scatter,global,tmp_vec,INSERT_VALUES,SCATTER_FORWARD)); 262 CHKERRQ(VecScatterEnd(scatter,global,tmp_vec,INSERT_VALUES,SCATTER_FORWARD)); 263 264 CHKERRQ(VecGetArrayRead(tmp_vec,&tmp)); 265 CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"At t[%D] =%14.2e u= %14.2e at the center \n",nsteps,(double)time,(double)PetscRealPart(tmp[n/2]))); 266 CHKERRQ(VecRestoreArrayRead(tmp_vec,&tmp)); 267 268 CHKERRQ(PetscFree(idx)); 269 CHKERRQ(ISDestroy(&from)); 270 CHKERRQ(ISDestroy(&to)); 271 CHKERRQ(VecScatterDestroy(&scatter)); 272 CHKERRQ(VecDestroy(&tmp_vec)); 273 PetscFunctionReturn(0); 274 } 275 276 PetscErrorCode RHSJacobian(TS ts,PetscReal t,Vec x,Mat A,Mat BB,void *ptr) 277 { 278 Data *data = (Data*)ptr; 279 PetscScalar v[5]; 280 PetscInt idx[5],i,j,row; 281 PetscInt m,n,mn; 282 PetscReal dx,dy,a,epsilon,xc,xl,xr,yl,yr; 283 284 PetscFunctionBeginUser; 285 m = data->m; 286 n = data->n; 287 mn = m*n; 288 dx = data->dx; 289 dy = data->dy; 290 a = data->a; 291 epsilon = data->epsilon; 292 293 xc = -2.0*epsilon*(1.0/(dx*dx)+1.0/(dy*dy)); 294 xl = 0.5*a/dx+epsilon/(dx*dx); 295 xr = -0.5*a/dx+epsilon/(dx*dx); 296 yl = 0.5*a/dy+epsilon/(dy*dy); 297 yr = -0.5*a/dy+epsilon/(dy*dy); 298 299 row = 0; 300 v[0] = xc; v[1] = xr; v[2] = yr; 301 idx[0] = 0; idx[1] = 2; idx[2] = m; 302 CHKERRQ(MatSetValues(A,1,&row,3,idx,v,INSERT_VALUES)); 303 304 row = m-1; 305 v[0] = 2.0*xl; v[1] = xc; v[2] = yr; 306 idx[0] = m-2; idx[1] = m-1; idx[2] = m-1+m; 307 CHKERRQ(MatSetValues(A,1,&row,3,idx,v,INSERT_VALUES)); 308 309 for (i=1; i<m-1; i++) { 310 row = i; 311 v[0] = xl; v[1] = xc; v[2] = xr; v[3] = yr; 312 idx[0] = i-1; idx[1] = i; idx[2] = i+1; idx[3] = i+m; 313 CHKERRQ(MatSetValues(A,1,&row,4,idx,v,INSERT_VALUES)); 314 } 315 316 for (j=1; j<n-1; j++) { 317 row = j*m; 318 v[0] = xc; v[1] = xr; v[2] = yl; v[3] = yr; 319 idx[0] = j*m; idx[1] = j*m; idx[2] = j*m-m; idx[3] = j*m+m; 320 CHKERRQ(MatSetValues(A,1,&row,4,idx,v,INSERT_VALUES)); 321 322 row = j*m+m-1; 323 v[0] = xc; v[1] = 2.0*xl; v[2] = yl; v[3] = yr; 324 idx[0] = j*m+m-1; idx[1] = j*m+m-1-1; idx[2] = j*m+m-1-m; idx[3] = j*m+m-1+m; 325 CHKERRQ(MatSetValues(A,1,&row,4,idx,v,INSERT_VALUES)); 326 327 for (i=1; i<m-1; i++) { 328 row = j*m+i; 329 v[0] = xc; v[1] = xl; v[2] = xr; v[3] = yl; v[4]=yr; 330 idx[0] = j*m+i; idx[1] = j*m+i-1; idx[2] = j*m+i+1; idx[3] = j*m+i-m; 331 idx[4] = j*m+i+m; 332 CHKERRQ(MatSetValues(A,1,&row,5,idx,v,INSERT_VALUES)); 333 } 334 } 335 336 row = mn-m; 337 v[0] = xc; v[1] = xr; v[2] = 2.0*yl; 338 idx[0] = mn-m; idx[1] = mn-m+1; idx[2] = mn-m-m; 339 CHKERRQ(MatSetValues(A,1,&row,3,idx,v,INSERT_VALUES)); 340 341 row = mn-1; 342 v[0] = xc; v[1] = 2.0*xl; v[2] = 2.0*yl; 343 idx[0] = mn-1; idx[1] = mn-2; idx[2] = mn-1-m; 344 CHKERRQ(MatSetValues(A,1,&i,3,idx,v,INSERT_VALUES)); 345 346 for (i=1; i<m-1; i++) { 347 row = mn-m+i; 348 v[0] = xl; v[1] = xc; v[2] = xr; v[3] = 2.0*yl; 349 idx[0] = mn-m+i-1; idx[1] = mn-m+i; idx[2] = mn-m+i+1; idx[3] = mn-m+i-m; 350 CHKERRQ(MatSetValues(A,1,&row,4,idx,v,INSERT_VALUES)); 351 } 352 353 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 354 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 355 356 PetscFunctionReturn(0); 357 } 358 359 /* globalout = -a*(u_x+u_y) + epsilon*(u_xx+u_yy) */ 360 PetscErrorCode RHSFunction(TS ts,PetscReal t,Vec globalin,Vec globalout,void *ctx) 361 { 362 Data *data = (Data*)ctx; 363 PetscInt m,n,mn; 364 PetscReal dx,dy; 365 PetscReal xc,xl,xr,yl,yr; 366 PetscReal a,epsilon; 367 PetscScalar *outptr; 368 const PetscScalar *inptr; 369 PetscInt i,j,len; 370 IS from,to; 371 PetscInt *idx; 372 VecScatter scatter; 373 Vec tmp_in,tmp_out; 374 375 PetscFunctionBeginUser; 376 m = data->m; 377 n = data->n; 378 mn = m*n; 379 dx = data->dx; 380 dy = data->dy; 381 a = data->a; 382 epsilon = data->epsilon; 383 384 xc = -2.0*epsilon*(1.0/(dx*dx)+1.0/(dy*dy)); 385 xl = 0.5*a/dx+epsilon/(dx*dx); 386 xr = -0.5*a/dx+epsilon/(dx*dx); 387 yl = 0.5*a/dy+epsilon/(dy*dy); 388 yr = -0.5*a/dy+epsilon/(dy*dy); 389 390 /* Get the length of parallel vector */ 391 CHKERRQ(VecGetSize(globalin,&len)); 392 393 /* Set the index sets */ 394 CHKERRQ(PetscMalloc1(len,&idx)); 395 for (i=0; i<len; i++) idx[i]=i; 396 397 /* Create local sequential vectors */ 398 CHKERRQ(VecCreateSeq(PETSC_COMM_SELF,len,&tmp_in)); 399 CHKERRQ(VecDuplicate(tmp_in,&tmp_out)); 400 401 /* Create scatter context */ 402 CHKERRQ(ISCreateGeneral(PETSC_COMM_SELF,len,idx,PETSC_COPY_VALUES,&from)); 403 CHKERRQ(ISCreateGeneral(PETSC_COMM_SELF,len,idx,PETSC_COPY_VALUES,&to)); 404 CHKERRQ(VecScatterCreate(globalin,from,tmp_in,to,&scatter)); 405 CHKERRQ(VecScatterBegin(scatter,globalin,tmp_in,INSERT_VALUES,SCATTER_FORWARD)); 406 CHKERRQ(VecScatterEnd(scatter,globalin,tmp_in,INSERT_VALUES,SCATTER_FORWARD)); 407 CHKERRQ(VecScatterDestroy(&scatter)); 408 409 /*Extract income array - include ghost points */ 410 CHKERRQ(VecGetArrayRead(tmp_in,&inptr)); 411 412 /* Extract outcome array*/ 413 CHKERRQ(VecGetArrayWrite(tmp_out,&outptr)); 414 415 outptr[0] = xc*inptr[0]+xr*inptr[1]+yr*inptr[m]; 416 outptr[m-1] = 2.0*xl*inptr[m-2]+xc*inptr[m-1]+yr*inptr[m-1+m]; 417 for (i=1; i<m-1; i++) { 418 outptr[i] = xc*inptr[i]+xl*inptr[i-1]+xr*inptr[i+1]+yr*inptr[i+m]; 419 } 420 421 for (j=1; j<n-1; j++) { 422 outptr[j*m] = xc*inptr[j*m]+xr*inptr[j*m+1]+ yl*inptr[j*m-m]+yr*inptr[j*m+m]; 423 outptr[j*m+m-1] = xc*inptr[j*m+m-1]+2.0*xl*inptr[j*m+m-1-1]+ yl*inptr[j*m+m-1-m]+yr*inptr[j*m+m-1+m]; 424 for (i=1; i<m-1; i++) { 425 outptr[j*m+i] = xc*inptr[j*m+i]+xl*inptr[j*m+i-1]+xr*inptr[j*m+i+1]+yl*inptr[j*m+i-m]+yr*inptr[j*m+i+m]; 426 } 427 } 428 429 outptr[mn-m] = xc*inptr[mn-m]+xr*inptr[mn-m+1]+2.0*yl*inptr[mn-m-m]; 430 outptr[mn-1] = 2.0*xl*inptr[mn-2]+xc*inptr[mn-1]+2.0*yl*inptr[mn-1-m]; 431 for (i=1; i<m-1; i++) { 432 outptr[mn-m+i] = xc*inptr[mn-m+i]+xl*inptr[mn-m+i-1]+xr*inptr[mn-m+i+1]+2*yl*inptr[mn-m+i-m]; 433 } 434 435 CHKERRQ(VecRestoreArrayRead(tmp_in,&inptr)); 436 CHKERRQ(VecRestoreArrayWrite(tmp_out,&outptr)); 437 438 CHKERRQ(VecScatterCreate(tmp_out,from,globalout,to,&scatter)); 439 CHKERRQ(VecScatterBegin(scatter,tmp_out,globalout,INSERT_VALUES,SCATTER_FORWARD)); 440 CHKERRQ(VecScatterEnd(scatter,tmp_out,globalout,INSERT_VALUES,SCATTER_FORWARD)); 441 442 /* Destroy idx aand scatter */ 443 CHKERRQ(VecDestroy(&tmp_in)); 444 CHKERRQ(VecDestroy(&tmp_out)); 445 CHKERRQ(ISDestroy(&from)); 446 CHKERRQ(ISDestroy(&to)); 447 CHKERRQ(VecScatterDestroy(&scatter)); 448 449 CHKERRQ(PetscFree(idx)); 450 PetscFunctionReturn(0); 451 } 452 453 PetscErrorCode PostStep(TS ts) 454 { 455 PetscReal t; 456 457 PetscFunctionBeginUser; 458 CHKERRQ(TSGetTime(ts,&t)); 459 CHKERRQ(PetscPrintf(PETSC_COMM_SELF," PostStep, t: %g\n",(double)t)); 460 PetscFunctionReturn(0); 461 } 462 463 /*TEST 464 465 test: 466 args: -ts_fd -ts_type beuler 467 output_file: output/ex4.out 468 469 test: 470 suffix: 2 471 args: -ts_fd -ts_type beuler 472 nsize: 2 473 output_file: output/ex4.out 474 475 test: 476 suffix: 3 477 args: -ts_fd -ts_type cn 478 479 test: 480 suffix: 4 481 args: -ts_fd -ts_type cn 482 output_file: output/ex4_3.out 483 nsize: 2 484 485 test: 486 suffix: 5 487 args: -ts_type beuler -ts_fd -fd_color -mat_coloring_type sl 488 output_file: output/ex4.out 489 490 test: 491 suffix: 6 492 args: -ts_type beuler -ts_fd -fd_color -mat_coloring_type sl 493 output_file: output/ex4.out 494 nsize: 2 495 496 test: 497 suffix: 7 498 requires: !single 499 args: -ts_fd -ts_type beuler -test_PostStep -ts_dt .1 500 501 test: 502 suffix: 8 503 requires: !single 504 args: -ts_type rk -ts_rk_type 5dp -ts_dt .01 -ts_adapt_type none -ts_view 505 506 TEST*/ 507