2D AMRCLAW
check_hdf.f
Go to the documentation of this file.
1 c
2 c ---------------------------------------------------------
3 c
4  subroutine check(nsteps,time,nvar,naux)
5 c
6 c :::::::::::::::::::::: CHECK ::::::::::::::::::::::::::::::::;
8 c :::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::;
9 
10  use amr_module
11  implicit double precision (a-h,o-z)
12 
13  character*16 chkname
14 c
15 c # HDF: Declare variables that describe datasets and HDF files.
16 c
17  integer sd_id
18 c
19 c # HDF: Arrays which will hold output arrays.
20  dimension iqout(15), qout(4)
21 c
22 c # HDF: Declare external HDF functions
23 c
24  integer sfstart, sfend
25  external sfstart, sfend
26 c
27 c # HDF: Set up HDF constants
28 c
29  integer DFACC_CREATE
30  parameter(dfacc_create = 4)
31 
32  integer SUCCEED, FAIL
33  parameter(succeed = 0, fail = -1)
34 c
35 c ### make the file name showing the time step
36 c
37  chkname = 'fort.chk'
38  & // char(ichar('0') + mod(nsteps/1000,10))
39  & // char(ichar('0') + mod(nsteps/100,10))
40  & // char(ichar('0') + mod(nsteps/10,10))
41  & // char(ichar('0') + mod(nsteps,10))
42  & // '.hdf'
43 c
44 c # HDF: create hdf file.
45 c
46  sd_id = sfstart(chkname, dfacc_create)
47  if (sd_id.eq.fail) THEN
48  WRITE(*,*) 'Failed to create HDF file (call to sfstart)'
49  stop
50  end if
51 c
52 c # HDF: create array of integers for storage.
53 c
54 c # HDF: DFNT_INT32
55  iqout(1) = lenmax
56  iqout(2) = lendim
57  iqout(3) = memsize
58  iqout(4) = lenf
59  iqout(5) = ibuff
60  iqout(6) = mstart
61  iqout(7) = ndfree
62  iqout(8) = lfine
63  iqout(9) = iorder
64  iqout(10) = mxnest
65  iqout(11) = kcheck
66  iqout(12) = nsteps
67  iqout(13) = matlabu
68  iqout(14) = lentot
69  iqout(15) = ndfree_bnd
70  call dump_integer_vector(sd_id,15,'iqout ',iqout)
71  call dump_integer_vector(sd_id,maxlv,'icheck',icheck)
72  call dump_integer_vector(sd_id,maxlv,'lstart',lstart)
73  call dump_integer_vector(sd_id,maxlv,'newstl',newstl)
74  call dump_integer_vector(sd_id,maxlv,'listsp',listsp)
75  call dump_integer_vector(sd_id,maxlv,'intratx',intratx)
76  call dump_integer_vector(sd_id,maxlv,'intraty',intraty)
77  call dump_integer_vector(sd_id,maxlv,'kratio',kratio)
78  call dump_integer_vector(sd_id,maxlv,'iregsz',iregsz)
79  call dump_integer_vector(sd_id,maxlv,'jregsz',jregsz)
80 
81  call dump_integer_array(sd_id,lfdim,2,'lfree ',lfree)
82  call dump_integer_array(sd_id,rsize,maxgr,'node ',node)
83 
84 c # HDF: DFNT_FLOAT64
85 
86  qout(1) = tol
87  qout(2) = time
88  qout(3) = evol
89  qout(4) = rvol
90  call dump_double_vector(sd_id,4,'qout ',qout)
91  call dump_double_vector(sd_id,maxlv,'hxposs',hxposs)
92  call dump_double_vector(sd_id,maxlv,'hyposs',hyposs)
93  call dump_double_vector(sd_id,maxlv,'possk ',possk)
94  call dump_double_vector(sd_id,10,'rvoll ',rvoll)
95  call dump_double_vector(sd_id,lendim,'alloc ',alloc)
96 
97  call dump_double_array(sd_id,rsize,maxgr,'rnode ',rnode)
98 c
99 c # HDF: Close HDF file.
100 c
101  istat = sfend(sd_id)
102  if (istat.eq.fail) then
103  WRITE(*,*) 'Failed to close SDS (call to sfend)'
104  stop
105  end if
106 
107  return
108  end
109 
110 c=================================================================
111 c # HDF: NOTE THAT ALL DOUBLE PRECISION HDF CHECKPOINTING ROUTINES
112 c # FOR BOTH WRITING AND READING ARE INCLUDED HERE, ALTHOUGH
113 c # THE DOUBLE PRECISION READING ROUTINES ARE CALLED IN
114 c # restrt_hdf.f. THIS IS DONE TO PREVENT COMPILER ERRORS
115 c # SINCE THE SAME HDF ROUTINES ARE USED TO READ BOTH INTEGER
116 c # AND DOUBLE PRECISION ARRAYS.
117 c=================================================================
118  subroutine dump_double_vector(sd_id,idims,qname,out)
119  implicit double precision (a-h,o-z)
120  character*6 qname
121 c
122 c # HDF: Declare variables that describe datasets and HDF files.
123 c
124  integer sd_id, sds_id
125  dimension out(idims), istart(1), istride(1), iedges(1), idim(1)
126 c
127 c # HDF: Declare external HDF functions
128 c
129  integer sfcreate, sfwdata, sfscompress, sfendacc
130  external sfcreate, sfwdata, sfscompress, sfendacc
131 c
132 c # HDF: Set up HDF constants
133 c
134  integer DFNT_FLOAT64, DFNT_INT32
135  parameter(dfnt_float64 = 6, dfnt_int32 = 24)
136 
137  integer SUCCEED, FAIL
138  parameter(succeed = 0, fail = -1)
139 c
140 c # HDF: Set up compression constants for HDF file.
141 c
142  integer COMP_CODE_DEFLATE, DEFLATE_LEVEL
143  parameter(comp_code_deflate = 4, deflate_level = 6)
144 c
145 c # HDF: Create double vector in HDF file.
146 c # NOTE THAT THE RANK MUST BE ONE.
147 c
148  irank = 1
149  idim(1) = idims
150  sds_id = sfcreate(sd_id,qname,dfnt_float64,irank,idim)
151  if (sds_id.eq.fail) THEN
152  WRITE(*,*) 'Failed to create variable ', qname,
153  & ' in restart HDF file'
154  WRITE(*,*) '(call to sfcreate in check_hdf.f)'
155  stop
156  end if
157 c
158 c # HDF: Set up dimensions for double vector.
159 c
160  istart(1) = 0
161  iedges(1) = idims
162  istride(1) = 1
163 c
164 c # HDF: set compression mode and write data to hdf file.
165 c
166  istat=sfscompress(sds_id,comp_code_deflate,deflate_level)
167  istat = sfwdata(sds_id,istart,istride,iedges,out)
168  if (istat.eq.fail) THEN
169  WRITE(*,*) 'Failed to write variable ', qname,
170  & ' in restart HDF file'
171  WRITE(*,*) '(call to sfwdata in check_hdf.f)'
172  stop
173  end if
174 c
175 c # HDF: End access to double vector in HDF file.
176 c
177  istat = sfendacc(sds_id)
178  if (istat.eq.fail) THEN
179  WRITE(*,*) 'Failed to end access to variable ', qname,
180  & ' in restart HDF file'
181  WRITE(*,*) '(call to sfendacc in check_hdf.f)'
182  stop
183  end if
184 
185  return
186  end
187 
188 c==========================================================
189  subroutine dump_double_array(sd_id,idim1,idim2,qname,out)
190  implicit double precision (a-h,o-z)
191  character*6 qname
192 c
193 c # HDF: Declare variables that describe datasets and HDF files.
194 c
195  integer sd_id, sds_id
196  dimension idims(2), istart(2), istride(2), iedges(2)
197  dimension out(idim1,idim2)
198 c
199 c # HDF: Declare external HDF functions
200 c
201  integer sfcreate, sfwdata, sfscompress, sfendacc
202  external sfcreate, sfwdata, sfscompress, sfendacc
203 c
204 c # HDF: Set up HDF constants
205 c
206  integer DFNT_FLOAT64, DFNT_INT32
207  parameter(dfnt_float64 = 6, dfnt_int32 = 24)
208 
209  integer SUCCEED, FAIL
210  parameter(succeed = 0, fail = -1)
211 c
212 c # HDF: Set up compression constants for HDF file.
213 c
214  integer COMP_CODE_DEFLATE, DEFLATE_LEVEL
215  parameter(comp_code_deflate = 4, deflate_level = 6)
216 c
217 c # HDF: Create double array in HDF file.
218 c # NOTE THAT THE RANK MUST BE TWO.
219 c
220  irank = 2
221  idims(1) = idim1
222  idims(2) = idim2
223  sds_id = sfcreate(sd_id,qname,dfnt_float64,irank,idims)
224  if (sds_id.eq.fail) THEN
225  WRITE(*,*) 'Failed to create variable ', qname,
226  & ' in restart HDF file'
227  WRITE(*,*) '(call to sfcreate in check_hdf.f)'
228  stop
229  end if
230 c
231 c # HDF: Set up dimensions for double array.
232 c
233  istart(1) = 0
234  istart(2) = 0
235  iedges(1) = idims(1)
236  iedges(2) = idims(2)
237  istride(1) = 1
238  istride(2) = 1
239 c
240 c # HDF: set compression mode and write data to hdf file.
241 c
242  istat=sfscompress(sds_id,comp_code_deflate,deflate_level)
243  istat = sfwdata(sds_id,istart,istride,iedges,out)
244  if (istat.eq.fail) THEN
245  WRITE(*,*) 'Failed to write variable ', qname,
246  & ' in restart HDF file'
247  WRITE(*,*) '(call to sfwdata in check_hdf.f)'
248  stop
249  end if
250 c
251 c # HDF: End access to double array in HDF file.
252 c
253  istat = sfendacc(sds_id)
254  if (istat.eq.fail) THEN
255  WRITE(*,*) 'Failed to end access to variable ', qname,
256  & ' in restart HDF file'
257  WRITE(*,*) '(call to sfendacc in check_hdf.f)'
258  stop
259  end if
260 
261  return
262  end
263 
264 
265 c==========================================================
266  subroutine read_double_vector(sd_id,idims,index,qname,out)
267  implicit double precision (a-h,o-z)
268  character*6 qname
269 c
270 c # HDF: Declare variables that describe datasets and HDF files.
271 c
272  integer sd_id, sds_id
273  dimension out(idims), istart(1), istride(1), iedges(1)
274 c
275 c # HDF: Declare external HDF functions
276 c
277  integer sfcreate, sfrdata, sfselect, sfendacc
278  external sfcreate, sfrdata, sfselect, sfendacc
279 c
280 c # HDF: Set up HDF constants
281 c
282  integer SUCCEED, FAIL
283  parameter(succeed = 0, fail = -1)
284 c
285 c # HDF: Select dataset in HDF file.
286 c
287  sds_id = sfselect(sd_id,index)
288  if (sds_id.eq.fail) THEN
289  WRITE(*,*) 'Failed to select data set for variable ', qname,
290  & ' in restart HDF file'
291  WRITE(*,*) '(call to sfselect in restrt_hdf.f)'
292  stop
293  end if
294 c
295 c # HDF: Set up dimensions for double vector.
296 c
297  istart(1) = 0
298  iedges(1) = idims
299  istride(1) = 1
300 c
301 c # HDF: read double vector from hdf file.
302 c
303  istat = sfrdata(sds_id,istart,istride,iedges,out)
304  if (istat.eq.fail) THEN
305  WRITE(*,*) 'Failed to read variable ', qname,
306  & ' from restart HDF file'
307  WRITE(*,*) '(call to sfrdata in restrt_hdf.f)'
308  stop
309  end if
310 c
311 c # HDF: End access to double vector in HDF file.
312 c
313  istat = sfendacc(sds_id)
314  if (istat.eq.fail) THEN
315  WRITE(*,*) 'Failed to end access to variable ', qname,
316  & ' in restart HDF file'
317  WRITE(*,*) '(call to sfendacc in restrt_hdf.f)'
318  stop
319  end if
320 
321  return
322  end
323 
324 c==========================================================
325  subroutine read_double_array(sd_id,idim1,idim2,index,qname,out)
326  implicit double precision (a-h,o-z)
327  character*6 qname
328 c
329 c # HDF: Declare variables that describe datasets and HDF files.
330 c
331  integer sd_id, sds_id
332  dimension istart(2), istride(2), iedges(2)
333  dimension out(idim1,idim2)
334 c
335 c # HDF: Declare external HDF functions
336 c
337  integer sfcreate, sfrdata, sfselect, sfendacc
338  external sfcreate, sfrdata, sfselect, sfendacc
339 c
340 c # HDF: Set up HDF constants
341 c
342  integer SUCCEED, FAIL
343  parameter(succeed = 0, fail = -1)
344 c
345 c # HDF: Select dataset in HDF file.
346 c
347  sds_id = sfselect(sd_id,index)
348  if (sds_id.eq.fail) THEN
349  WRITE(*,*) 'Failed to select data set for variable ', qname,
350  & ' in restart HDF file'
351  WRITE(*,*) '(call to sfselect in restrt_hdf.f)'
352  stop
353  end if
354 c
355 c # HDF: Set up dimensions for double array.
356 c
357  istart(1) = 0
358  istart(2) = 0
359  iedges(1) = idim1
360  iedges(2) = idim2
361  istride(1) = 1
362  istride(2) = 1
363 c
364 c # HDF: read double array from hdf file.
365 c
366  istat = sfrdata(sds_id,istart,istride,iedges,out)
367  if (istat.eq.fail) THEN
368  WRITE(*,*) 'Failed to read variable ', qname,
369  & ' from restart HDF file'
370  WRITE(*,*) '(call to sfrdata in restrt_hdf.f)'
371  stop
372  end if
373 c
374 c # HDF: End access to double array in HDF file.
375 c
376  istat = sfendacc(sds_id)
377  if (istat.eq.fail) THEN
378  WRITE(*,*) 'Failed to end access to variable ', qname,
379  & ' in restart HDF file'
380  WRITE(*,*) '(call to sfendacc in restrt_hdf.f)'
381  stop
382  end if
383 
384  return
385  end
integer, dimension(maxlv) kratio
Definition: amr_module.f90:198
integer, parameter rsize
Definition: amr_module.f90:30
integer lendim
Definition: amr_module.f90:247
subroutine read_double_array(sd_id, idim1, idim2, index, qname, out)
Definition: check_hdf.f:326
integer, parameter maxlv
Definition: amr_module.f90:174
integer ndfree
Definition: amr_module.f90:198
subroutine dump_double_array(sd_id, idim1, idim2, qname, out)
Definition: check_hdf.f:190
integer, parameter maxgr
Definition: amr_module.f90:173
real(kind=8) evol
Definition: amr_module.f90:237
integer lenmax
Definition: amr_module.f90:247
integer, dimension(maxlv) iregsz
Definition: amr_module.f90:198
real(kind=8), dimension(maxlv) hyposs
Definition: amr_module.f90:193
real(kind=8), dimension(maxlv) hxposs
Definition: amr_module.f90:193
integer, dimension(maxlv) newstl
Definition: amr_module.f90:198
subroutine dump_double_vector(sd_id, idims, qname, out)
Definition: check_hdf.f:119
integer, dimension(maxlv) jregsz
Definition: amr_module.f90:198
integer, dimension(nsize, maxgr) node
Definition: amr_module.f90:198
integer, parameter lfdim
Definition: amr_module.f90:224
integer memsize
Definition: amr_module.f90:219
real(kind=8), dimension(rsize, maxgr) rnode
Definition: amr_module.f90:193
integer ndfree_bnd
Definition: amr_module.f90:198
integer lentot
Definition: amr_module.f90:247
real(kind=8) tol
Definition: amr_module.f90:197
subroutine check(nsteps, time, nvar, naux)
Definition: check.f:5
subroutine dump_integer_array(sd_id, idim1, idim2, qname, iout)
Definition: restrt_hdf.f:433
integer, dimension(maxlv) icheck
Definition: amr_module.f90:198
integer matlabu
Definition: amr_module.f90:283
integer iorder
Definition: amr_module.f90:198
integer, dimension(maxlv) lstart
Definition: amr_module.f90:198
real(kind=8) rvol
Definition: amr_module.f90:237
integer ibuff
Definition: amr_module.f90:198
subroutine read_double_vector(sd_id, idims, index, qname, out)
Definition: check_hdf.f:267
integer, dimension(maxlv) intraty
Definition: amr_module.f90:198
real(kind=8), dimension(maxlv) rvoll
Definition: amr_module.f90:237
real(kind=8), dimension(maxlv) possk
Definition: amr_module.f90:193
integer, dimension(lfdim, 2) lfree
Definition: amr_module.f90:225
integer, dimension(maxlv) listsp
Definition: amr_module.f90:198
integer lenf
Definition: amr_module.f90:225
integer, dimension(maxlv) intratx
Definition: amr_module.f90:198
integer mxnest
Definition: amr_module.f90:198
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
subroutine dump_integer_vector(sd_id, idims, qname, iout)
Definition: restrt_hdf.f:362
integer mstart
Definition: amr_module.f90:198
integer kcheck
Definition: amr_module.f90:198
integer lfine
Definition: amr_module.f90:198
real(kind=8), dimension(:), allocatable alloc
Definition: amr_module.f90:218