Actual source code: ex22f_mf.F90
1: ! Time-dependent advection-reaction PDE in 1d. Demonstrates IMEX methods
2: !
3: ! u_t + a1*u_x = -k1*u + k2*v + s1
4: ! v_t + a2*v_x = k1*u - k2*v + s2
5: ! 0 < x < 1;
6: ! a1 = 1, k1 = 10^6, s1 = 0,
7: ! a2 = 0, k2 = 2*k1, s2 = 1
8: !
9: ! Initial conditions:
10: ! u(x,0) = 1 + s2*x
11: ! v(x,0) = k0/k1*u(x,0) + s1/k1
12: !
13: ! Upstream boundary conditions:
14: ! u(0,t) = 1-sin(12*t)^4
15: !
17: module PETScShiftMod
18: #include <petsc/finclude/petscts.h>
19: use petscts
20: PetscScalar::PETSC_SHIFT
21: TS::tscontext
22: Mat::Jmat
23: PetscReal::MFuser(6)
24: end module PETScShiftMod
26: program main
27: use PETScShiftMod
28: use petscdmda
29: implicit none
31: !
32: ! Create an application context to contain data needed by the
33: ! application-provided call-back routines, FormJacobian() and
34: ! FormFunction(). We use a double precision array with six
35: ! entries, two for each problem parameter a, k, s.
36: !
37: PetscReal user(6)
38: integer user_a,user_k,user_s
39: parameter (user_a = 0,user_k = 2,user_s = 4)
41: external FormRHSFunction,FormIFunction
42: external FormInitialSolution
43: external FormIJacobian
44: external MyMult,FormIJacobianMF
46: TS ts
47: Vec X
48: Mat J
49: PetscInt mx
50: PetscBool OptionSaveToDisk
51: PetscErrorCode ierr
52: DM da
53: PetscReal ftime,dt
54: PetscReal one,pone
55: PetscInt im11,i2
56: PetscBool flg
58: PetscInt xs,xe,gxs,gxe,dof,gdof
59: PetscScalar shell_shift
60: Mat A
62: im11 = 11
63: i2 = 2
64: one = 1.0
65: pone = one / 10
67: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
68: if (ierr .ne. 0) then
69: print*,'PetscInitialize failed'
70: stop
71: endif
73: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: ! Create distributed array (DMDA) to manage parallel grid and vectors
75: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
76: call DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,im11,i2,i2,PETSC_NULL_INTEGER,da,ierr);CHKERRA(ierr)
77: call DMSetFromOptions(da,ierr);CHKERRA(ierr)
78: call DMSetUp(da,ierr);CHKERRA(ierr)
80: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
81: ! Extract global vectors from DMDA;
82: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
83: call DMCreateGlobalVector(da,X,ierr);CHKERRA(ierr)
85: ! Initialize user application context
86: ! Use zero-based indexing for command line parameters to match ex22.c
87: user(user_a+1) = 1.0
88: call PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-a0',user(user_a+1),flg,ierr);CHKERRA(ierr)
89: user(user_a+2) = 0.0
90: call PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-a1',user(user_a+2),flg,ierr);CHKERRA(ierr)
91: user(user_k+1) = 1000000.0
92: call PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-k0', user(user_k+1),flg,ierr);CHKERRA(ierr)
93: user(user_k+2) = 2*user(user_k+1)
94: call PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-k1', user(user_k+2),flg,ierr);CHKERRA(ierr)
95: user(user_s+1) = 0.0
96: call PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s0',user(user_s+1),flg,ierr);CHKERRA(ierr)
97: user(user_s+2) = 1.0
98: call PetscOptionsGetReal(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s1',user(user_s+2),flg,ierr);CHKERRA(ierr)
100: OptionSaveToDisk=.FALSE.
101: call PetscOptionsGetBool(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-sdisk',OptionSaveToDisk,flg,ierr);CHKERRA(ierr)
102: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
103: ! Create timestepping solver context
104: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
105: call TSCreate(PETSC_COMM_WORLD,ts,ierr);CHKERRA(ierr)
106: tscontext=ts
107: call TSSetDM(ts,da,ierr);CHKERRA(ierr)
108: call TSSetType(ts,TSARKIMEX,ierr);CHKERRA(ierr)
109: call TSSetRHSFunction(ts,PETSC_NULL_VEC,FormRHSFunction,user,ierr);CHKERRA(ierr)
111: ! - - - - - - - - -- - - - -
112: ! Matrix free setup
113: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr);CHKERRA(ierr)
114: dof=i2*(xe-xs+1)
115: gdof=i2*(gxe-gxs+1)
116: call MatCreateShell(PETSC_COMM_WORLD,dof,dof,gdof,gdof,shell_shift,A,ierr);CHKERRA(ierr)
118: call MatShellSetOperation(A,MATOP_MULT,MyMult,ierr);CHKERRA(ierr)
119: ! - - - - - - - - - - - -
121: call TSSetIFunction(ts,PETSC_NULL_VEC,FormIFunction,user,ierr);CHKERRA(ierr)
122: call DMSetMatType(da,MATAIJ,ierr);CHKERRA(ierr)
123: call DMCreateMatrix(da,J,ierr);CHKERRA(ierr)
125: Jmat=J
127: call TSSetIJacobian(ts,J,J,FormIJacobian,user,ierr);CHKERRA(ierr)
128: call TSSetIJacobian(ts,A,A,FormIJacobianMF,user,ierr);CHKERRA(ierr)
130: ftime = 1.0
131: call TSSetMaxTime(ts,ftime,ierr);CHKERRA(ierr)
132: call TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER,ierr);CHKERRA(ierr)
134: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
135: ! Set initial conditions
136: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
137: call FormInitialSolution(ts,X,user,ierr);CHKERRA(ierr)
138: call TSSetSolution(ts,X,ierr);CHKERRA(ierr)
139: call VecGetSize(X,mx,ierr);CHKERRA(ierr)
140: ! Advective CFL, I don't know why it needs so much safety factor.
141: dt = pone * max(user(user_a+1),user(user_a+2)) / mx;
142: call TSSetTimeStep(ts,dt,ierr);CHKERRA(ierr)
144: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
145: ! Set runtime options
146: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
147: call TSSetFromOptions(ts,ierr);CHKERRA(ierr)
149: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
150: ! Solve nonlinear system
151: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
152: call TSSolve(ts,X,ierr);CHKERRA(ierr)
154: if (OptionSaveToDisk) then
155: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr);CHKERRA(ierr)
156: dof=i2*(xe-xs+1)
157: gdof=i2*(gxe-gxs+1)
158: call SaveSolutionToDisk(da,X,gdof,xs,xe);CHKERRA(ierr)
159: end if
161: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
162: ! Free work space.
163: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
164: call MatDestroy(A,ierr);CHKERRA(ierr)
165: call MatDestroy(J,ierr);CHKERRA(ierr)
166: call VecDestroy(X,ierr);CHKERRA(ierr)
167: call TSDestroy(ts,ierr);CHKERRA(ierr)
168: call DMDestroy(da,ierr);CHKERRA(ierr)
169: call PetscFinalize(ierr)
170: end program main
172: ! Small helper to extract the layout, result uses 1-based indexing.
173: subroutine GetLayout(da,mx,xs,xe,gxs,gxe,ierr)
174: use petscdmda
175: implicit none
177: DM da
178: PetscInt mx,xs,xe,gxs,gxe
179: PetscErrorCode ierr
180: PetscInt xm,gxm
181: call DMDAGetInfo(da,PETSC_NULL_INTEGER,mx,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER, &
182: PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr);CHKERRQ(ierr)
183: call DMDAGetCorners(da,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,xm,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr);CHKERRQ(ierr)
184: call DMDAGetGhostCorners(da,gxs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,gxm,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr);CHKERRQ(ierr)
185: xs = xs + 1
186: gxs = gxs + 1
187: xe = xs + xm - 1
188: gxe = gxs + gxm - 1
189: end subroutine GetLayout
191: subroutine FormIFunctionLocal(mx,xs,xe,gxs,gxe,x,xdot,f,a,k,s,ierr)
192: implicit none
193: PetscInt mx,xs,xe,gxs,gxe
194: PetscScalar x(2,xs:xe)
195: PetscScalar xdot(2,xs:xe)
196: PetscScalar f(2,xs:xe)
197: PetscReal a(2),k(2),s(2)
198: PetscErrorCode ierr
199: PetscInt i
200: do i = xs,xe
201: f(1,i) = xdot(1,i) + k(1)*x(1,i) - k(2)*x(2,i) - s(1)
202: f(2,i) = xdot(2,i) - k(1)*x(1,i) + k(2)*x(2,i) - s(2)
203: end do
204: end subroutine FormIFunctionLocal
206: subroutine FormIFunction(ts,t,X,Xdot,F,user,ierr)
207: use petscdmda
208: use petscts
209: implicit none
211: TS ts
212: PetscReal t
213: Vec X,Xdot,F
214: PetscReal user(6)
215: PetscErrorCode ierr
216: integer user_a,user_k,user_s
217: parameter (user_a = 1,user_k = 3,user_s = 5)
219: DM da
220: PetscInt mx,xs,xe,gxs,gxe
221: PetscOffset ixx,ixxdot,iff
222: PetscScalar xx(0:1),xxdot(0:1),ff(0:1)
224: call TSGetDM(ts,da,ierr);CHKERRQ(ierr)
225: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr);CHKERRQ(ierr)
227: ! Get access to vector data
228: call VecGetArrayRead(X,xx,ixx,ierr);CHKERRQ(ierr)
229: call VecGetArrayRead(Xdot,xxdot,ixxdot,ierr);CHKERRQ(ierr)
230: call VecGetArray(F,ff,iff,ierr);CHKERRQ(ierr)
232: call FormIFunctionLocal(mx,xs,xe,gxs,gxe,xx(ixx),xxdot(ixxdot),ff(iff),user(user_a),user(user_k),user(user_s),ierr);CHKERRQ(ierr)
234: call VecRestoreArrayRead(X,xx,ixx,ierr);CHKERRQ(ierr)
235: call VecRestoreArrayRead(Xdot,xxdot,ixxdot,ierr);CHKERRQ(ierr)
236: call VecRestoreArray(F,ff,iff,ierr);CHKERRQ(ierr)
237: end subroutine FormIFunction
239: subroutine FormRHSFunctionLocal(mx,xs,xe,gxs,gxe,t,x,f,a,k,s,ierr)
240: implicit none
241: PetscInt mx,xs,xe,gxs,gxe
242: PetscReal t
243: PetscScalar x(2,gxs:gxe),f(2,xs:xe)
244: PetscReal a(2),k(2),s(2)
245: PetscErrorCode ierr
246: PetscInt i,j
247: PetscReal hx,u0t(2)
248: PetscReal one,two,three,four,six,twelve
249: PetscReal half,third,twothird,sixth
250: PetscReal twelfth
252: one = 1.0
253: two = 2.0
254: three = 3.0
255: four = 4.0
256: six = 6.0
257: twelve = 12.0
258: hx = one / mx
259: u0t(1) = one - sin(twelve*t)**four
260: u0t(2) = 0.0
261: half = one/two
262: third = one / three
263: twothird = two / three
264: sixth = one / six
265: twelfth = one / twelve
266: do i = xs,xe
267: do j = 1,2
268: if (i .eq. 1) then
269: f(j,i) = a(j)/hx*(third*u0t(j) + half*x(j,i) - x(j,i+1) + sixth*x(j,i+2))
270: else if (i .eq. 2) then
271: f(j,i) = a(j)/hx*(-twelfth*u0t(j) + twothird*x(j,i-1) - twothird*x(j,i+1) + twelfth*x(j,i+2))
272: else if (i .eq. mx-1) then
273: f(j,i) = a(j)/hx*(-sixth*x(j,i-2) + x(j,i-1) - half*x(j,i) -third*x(j,i+1))
274: else if (i .eq. mx) then
275: f(j,i) = a(j)/hx*(-x(j,i) + x(j,i-1))
276: else
277: f(j,i) = a(j)/hx*(-twelfth*x(j,i-2) + twothird*x(j,i-1) - twothird*x(j,i+1) + twelfth*x(j,i+2))
278: end if
279: end do
280: end do
282: #ifdef EXPLICIT_INTEGRATOR22
283: do i = xs,xe
284: f(1,i) = f(1,i) -( k(1)*x(1,i) - k(2)*x(2,i) - s(1))
285: f(2,i) = f(2,i) -(- k(1)*x(1,i) + k(2)*x(2,i) - s(2))
286: end do
287: #endif
289: end subroutine FormRHSFunctionLocal
291: subroutine FormRHSFunction(ts,t,X,F,user,ierr)
292: use petscts
293: use petscdmda
294: implicit none
296: TS ts
297: PetscReal t
298: Vec X,F
299: PetscReal user(6)
300: PetscErrorCode ierr
301: integer user_a,user_k,user_s
302: parameter (user_a = 1,user_k = 3,user_s = 5)
303: DM da
304: Vec Xloc
305: PetscInt mx,xs,xe,gxs,gxe
306: PetscOffset ixx,iff
307: PetscScalar xx(0:1),ff(0:1)
309: call TSGetDM(ts,da,ierr);CHKERRQ(ierr)
310: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr);CHKERRQ(ierr)
312: ! Scatter ghost points to local vector,using the 2-step process
313: ! DMGlobalToLocalBegin(),DMGlobalToLocalEnd().
314: ! By placing code between these two statements, computations can be
315: ! done while messages are in transition.
316: call DMGetLocalVector(da,Xloc,ierr);CHKERRQ(ierr)
317: call DMGlobalToLocalBegin(da,X,INSERT_VALUES,Xloc,ierr);CHKERRQ(ierr)
318: call DMGlobalToLocalEnd(da,X,INSERT_VALUES,Xloc,ierr);CHKERRQ(ierr)
320: ! Get access to vector data
321: call VecGetArrayRead(Xloc,xx,ixx,ierr);CHKERRQ(ierr)
322: call VecGetArray(F,ff,iff,ierr);CHKERRQ(ierr)
324: call FormRHSFunctionLocal(mx,xs,xe,gxs,gxe,t,xx(ixx),ff(iff),user(user_a),user(user_k),user(user_s),ierr);CHKERRQ(ierr)
326: call VecRestoreArrayRead(Xloc,xx,ixx,ierr);CHKERRQ(ierr)
327: call VecRestoreArray(F,ff,iff,ierr);CHKERRQ(ierr)
328: call DMRestoreLocalVector(da,Xloc,ierr);CHKERRQ(ierr)
329: end subroutine FormRHSFunction
331: ! ---------------------------------------------------------------------
332: !
333: ! IJacobian - Compute IJacobian = dF/dU + shift*dF/dUdot
334: !
335: subroutine FormIJacobian(ts,t,X,Xdot,shift,J,Jpre,user,ierr)
336: use petscts
337: use petscdmda
338: implicit none
340: TS ts
341: PetscReal t,shift
342: Vec X,Xdot
343: Mat J,Jpre
344: PetscReal user(6)
345: PetscErrorCode ierr
346: integer user_a,user_k,user_s
347: parameter (user_a = 0,user_k = 2,user_s = 4)
349: DM da
350: PetscInt mx,xs,xe,gxs,gxe
351: PetscInt i,i1,row,col
352: PetscReal k1,k2;
353: PetscScalar val(4)
355: call TSGetDM(ts,da,ierr);CHKERRQ(ierr)
356: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr);CHKERRQ(ierr)
358: i1 = 1
359: k1 = user(user_k+1)
360: k2 = user(user_k+2)
361: do i = xs,xe
362: row = i-gxs
363: col = i-gxs
364: val(1) = shift + k1
365: val(2) = -k2
366: val(3) = -k1
367: val(4) = shift + k2
368: call MatSetValuesBlockedLocal(Jpre,i1,row,i1,col,val,INSERT_VALUES,ierr);CHKERRQ(ierr)
369: end do
370: call MatAssemblyBegin(Jpre,MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
371: call MatAssemblyEnd(Jpre,MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
372: if (J /= Jpre) then
373: call MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
374: call MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
375: end if
376: end subroutine FormIJacobian
378: subroutine FormInitialSolutionLocal(mx,xs,xe,gxs,gxe,x,a,k,s,ierr)
379: implicit none
380: PetscInt mx,xs,xe,gxs,gxe
381: PetscScalar x(2,xs:xe)
382: PetscReal a(2),k(2),s(2)
383: PetscErrorCode ierr
385: PetscInt i
386: PetscReal one,hx,r,ik
387: one = 1.0
388: hx = one / mx
389: do i=xs,xe
390: r = i*hx
391: if (k(2) .ne. 0.0) then
392: ik = one/k(2)
393: else
394: ik = one
395: end if
396: x(1,i) = one + s(2)*r
397: x(2,i) = k(1)*ik*x(1,i) + s(2)*ik
398: end do
399: end subroutine FormInitialSolutionLocal
401: subroutine FormInitialSolution(ts,X,user,ierr)
402: use petscts
403: use petscdmda
404: implicit none
406: TS ts
407: Vec X
408: PetscReal user(6)
409: PetscErrorCode ierr
410: integer user_a,user_k,user_s
411: parameter (user_a = 1,user_k = 3,user_s = 5)
413: DM da
414: PetscInt mx,xs,xe,gxs,gxe
415: PetscOffset ixx
416: PetscScalar xx(0:1)
418: call TSGetDM(ts,da,ierr);CHKERRQ(ierr)
419: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr);CHKERRQ(ierr)
421: ! Get access to vector data
422: call VecGetArray(X,xx,ixx,ierr);CHKERRQ(ierr)
424: call FormInitialSolutionLocal(mx,xs,xe,gxs,gxe,xx(ixx),user(user_a),user(user_k),user(user_s),ierr);CHKERRQ(ierr)
426: call VecRestoreArray(X,xx,ixx,ierr);CHKERRQ(ierr)
427: end subroutine FormInitialSolution
429: ! ---------------------------------------------------------------------
430: !
431: ! IJacobian - Compute IJacobian = dF/dU + shift*dF/dUdot
432: !
433: subroutine FormIJacobianMF(ts,t,X,Xdot,shift,J,Jpre,user,ierr)
434: use PETScShiftMod
435: implicit none
436: TS ts
437: PetscReal t,shift
438: Vec X,Xdot
439: Mat J,Jpre
440: PetscReal user(6)
441: PetscErrorCode ierr
443: ! call MatShellSetContext(J,shift,ierr)
444: PETSC_SHIFT=shift
445: MFuser=user
447: end subroutine FormIJacobianMF
449: ! -------------------------------------------------------------------
450: !
451: ! MyMult - user provided matrix multiply
452: !
453: ! Input Parameters:
454: !. X - input vector
455: !
456: ! Output Parameter:
457: !. F - function vector
458: !
459: subroutine MyMult(A,X,F,ierr)
460: use PETScShiftMod
461: implicit none
463: Mat A
464: Vec X,F
466: PetscErrorCode ierr
467: PetscScalar shift
469: ! Mat J,Jpre
471: PetscReal user(6)
473: integer user_a,user_k,user_s
474: parameter (user_a = 0,user_k = 2,user_s = 4)
476: DM da
477: PetscInt mx,xs,xe,gxs,gxe
478: PetscInt i,i1,row,col
479: PetscReal k1,k2;
480: PetscScalar val(4)
482: !call MatShellGetContext(A,shift,ierr)
483: shift=PETSC_SHIFT
484: user=MFuser
486: call TSGetDM(tscontext,da,ierr)
487: call GetLayout(da,mx,xs,xe,gxs,gxe,ierr)
489: i1 = 1
490: k1 = user(user_k+1)
491: k2 = user(user_k+2)
493: do i = xs,xe
494: row = i-gxs
495: col = i-gxs
496: val(1) = shift + k1
497: val(2) = -k2
498: val(3) = -k1
499: val(4) = shift + k2
500: call MatSetValuesBlockedLocal(Jmat,i1,row,i1,col,val,INSERT_VALUES,ierr)
501: end do
503: ! call MatAssemblyBegin(Jpre,MAT_FINAL_ASSEMBLY,ierr)
504: ! call MatAssemblyEnd(Jpre,MAT_FINAL_ASSEMBLY,ierr)
505: ! if (J /= Jpre) then
506: call MatAssemblyBegin(Jmat,MAT_FINAL_ASSEMBLY,ierr)
507: call MatAssemblyEnd(Jmat,MAT_FINAL_ASSEMBLY,ierr)
508: ! end if
510: call MatMult(Jmat,X,F,ierr)
512: return
513: end subroutine MyMult
515: !
516: subroutine SaveSolutionToDisk(da,X,gdof,xs,xe)
517: use petscdmda
518: implicit none
520: Vec X
521: DM da
522: PetscInt xs,xe,two
523: PetscInt gdof,i
524: PetscErrorCode ierr
525: PetscOffset ixx
526: PetscScalar data2(2,xs:xe),data(gdof)
527: PetscScalar xx(0:1)
529: call VecGetArrayRead(X,xx,ixx,ierr)
531: two = 2
532: data2=reshape(xx(ixx:ixx+gdof),(/two,xe-xs+1/))
533: data=reshape(data2,(/gdof/))
534: open(1020,file='solution_out_ex22f_mf.txt')
535: do i=1,gdof
536: write(1020,'(e24.16,1x)') data(i)
537: end do
538: close(1020)
540: call VecRestoreArrayRead(X,xx,ixx,ierr)
541: end subroutine SaveSolutionToDisk
543: !/*TEST
544: !
545: ! test:
546: ! args: -da_grid_x 200 -ts_arkimex_type 4
547: ! requires: !single
548: ! output_file: output/ex22f_mf_1.out
549: !
550: !TEST*/