2D AMRCLAW
Data Types | Functions/Subroutines | Variables
gauges_module Module Reference

Data Types

type  gauge_type
 

Functions/Subroutines

subroutine set_gauges (restart, num_eqn, num_aux, fname)
 
subroutine setbestsrc ()
 
subroutine update_gauges (q, aux, xlow, ylow, num_eqn, mitot, mjtot, num_aux, mptr)
 
subroutine print_gauges_and_reset_nextloc (gauge_num)
 

Variables

logical, private module_setup = .false.
 
integer, parameter outgaugeunit = 89
 
integer num_gauges
 
integer, parameter max_buffer = 1000
 
type(gauge_type), dimension(:), allocatable gauges
 
integer, dimension(:), allocatable mbestsrc
 
integer, dimension(:), allocatable mbestorder
 
integer, dimension(:), allocatable igauge
 
integer, dimension(:), allocatable mbestg1
 
integer, dimension(:), allocatable mbestg2
 

Function/Subroutine Documentation

◆ print_gauges_and_reset_nextloc()

subroutine gauges_module::print_gauges_and_reset_nextloc ( integer, intent(in)  gauge_num)

Definition at line 499 of file gauges_module.f90.

References gauges, and outgaugeunit.

Referenced by check(), tick(), and update_gauges().

499  ! Write out gauge data for the gauge specified
500 
501  implicit none
502 
503  ! Input
504  integer, intent(in) :: gauge_num
505 
506  ! Locals
507  integer :: j, k, myunit
508  integer :: omp_get_thread_num, mythread
509  character(len=32) :: out_format
510 
511  ! Open unit dependent on thread number
512  mythread = 0
513 !$ mythread = omp_get_thread_num()
514  myunit = outgaugeunit + mythread
515 
516  ! ASCII output
517  if (gauges(gauge_num)%file_format == 1) then
518  ! Construct output format based on number of output variables and
519  ! request format
520  write(out_format, "(A7, i2, A6, A1)") "(i5.2,", &
521  gauges(gauge_num)%num_out_vars + 1, gauges(gauge_num)%display_format, ")"
522 
523  open(unit=myunit, file=gauges(gauge_num)%file_name, status='old', &
524  position='append', form='formatted')
525 
526  ! Loop through gauge's buffer writing out all available data. Also
527  ! reset buffer_index back to beginning of buffer since we are emptying
528  ! the buffer here
529  do j = 1, gauges(gauge_num)%buffer_index - 1
530  write(myunit, out_format) gauges(gauge_num)%level(j), &
531  (gauges(gauge_num)%data(k, j), k=1, gauges(gauge_num)%num_out_vars + 1)
532  end do
533  gauges(gauge_num)%buffer_index = 1
534 
535  ! close file
536  close(myunit)
537  else
538  print *, "Unhandled file format ", gauges(gauge_num)%file_format
539  stop
540  end if
541 
Here is the caller graph for this function:

◆ set_gauges()

subroutine gauges_module::set_gauges ( logical, intent(in)  restart,
integer  num_eqn,
integer  num_aux,
character(len=*), intent(in), optional  fname 
)

Definition at line 90 of file gauges_module.f90.

References gauges, max_buffer, amr_module::maxgr, mbestg1, mbestg2, mbestorder, mbestsrc, module_setup, num_gauges, opendatafile(), and outgaugeunit.

Referenced by amr2().

90 
91  use amr_module, only: maxgr
92  use utility_module, only: get_value_count
93 
94  implicit none
95 
96  ! Input
97  logical, intent(in) :: restart
98  integer :: num_eqn, num_aux
99  character(len=*), intent(in), optional :: fname
100 
101  ! Locals
102  integer :: i, n, index
103  integer :: num, pos, digit
104  integer, parameter :: unit = 7
105  character(len=128) :: header_1
106  character(len=40) :: q_column, aux_column
107 
108  if (.not. module_setup) then
109 
110  ! Open file
111  if (present(fname)) then
112  call opendatafile(unit, fname)
113  else
114  call opendatafile(unit, 'gauges.data')
115  endif
116 
117  read(unit, *) num_gauges
118  allocate(gauges(num_gauges))
119 
120  ! Initialize gauge source data
121  allocate(mbestsrc(num_gauges), mbestorder(num_gauges))
122  allocate(mbestg1(maxgr), mbestg2(maxgr))
123  mbestsrc = 0
124 
125  ! Original gauge information
126  do i=1,num_gauges
127  read(unit, *) gauges(i)%gauge_num, gauges(i)%x, gauges(i)%y, &
128  gauges(i)%t_start, gauges(i)%t_end
129  gauges(i)%buffer_index = 1
130  gauges(i)%last_time = gauges(i)%t_start
131  enddo
132 
133  ! Read in output formats
134  read(unit, *)
135  read(unit, *)
136  read(unit, *) (gauges(i)%file_format, i=1, num_gauges)
137  read(unit, *)
138  read(unit, *)
139  read(unit, *) (gauges(i)%display_format, i=1, num_gauges)
140  read(unit, *)
141  read(unit, *)
142  read(unit, *) (gauges(i)%min_time_increment, i=1, num_gauges)
143 
144  ! Read in q fields
145  read(unit, *)
146  read(unit, *)
147  do i = 1, num_gauges
148  allocate(gauges(i)%q_out_vars(num_eqn))
149  read(unit, *) gauges(i)%q_out_vars
150 
151  ! Count number of vars to be output
152  gauges(i)%num_out_vars = 0
153  do n = 1, size(gauges(i)%q_out_vars, 1)
154  if (gauges(i)%q_out_vars(n)) then
155  gauges(i)%num_out_vars = gauges(i)%num_out_vars + 1
156  end if
157  end do
158  end do
159 
160  ! Read in aux fields
161  if (num_aux > 0) then
162  read(unit, *)
163  read(unit, *)
164  do i = 1, num_gauges
165  allocate(gauges(i)%aux_out_vars(num_aux))
166  read(unit, *) gauges(i)%aux_out_vars
167 
168  ! Count number of vars to be output
169  do n = 1, size(gauges(i)%aux_out_vars, 1)
170  if (gauges(i)%aux_out_vars(n)) then
171  gauges(i)%num_out_vars = gauges(i)%num_out_vars + 1
172  end if
173  end do
174  end do
175  end if
176 
177  close(unit)
178  ! Done reading =====================================================
179 
180  ! Allocate data buffer
181  do i = 1, num_gauges
182  allocate(gauges(i)%data(gauges(i)%num_out_vars + 1, max_buffer))
183  end do
184 
185  ! Create gauge output files
186  do i = 1, num_gauges
187  gauges(i)%file_name = 'gaugexxxxx.txt'
188  num = gauges(i)%gauge_num
189  do pos = 10, 6, -1
190  digit = mod(num,10)
191  gauges(i)%file_name(pos:pos) = char(ichar('0') + digit)
192  num = num / 10
193  end do
194 
195  ! Handle restart
196  if (restart) then
197  open(unit=outgaugeunit, file=gauges(i)%file_name, &
198  status='old', position='append', form='formatted')
199  else
200  open(unit=outgaugeunit, file=gauges(i)%file_name, &
201  status='unknown', position='append', form='formatted')
202  rewind outgaugeunit
203 
204  ! Write header
205  header_1 = "('# gauge_id= ',i5,' " // &
206  "location=( ',1e17.10,' ',1e17.10,' ) " // &
207  "num_var= ',i2)"
208  write(outgaugeunit, header_1) gauges(i)%gauge_num, &
209  gauges(i)%x, &
210  gauges(i)%y, &
211  gauges(i)%num_out_vars
212 
213  ! Construct column labels
214  index = 0
215  q_column = "["
216  do n=1, size(gauges(i)%q_out_vars, 1)
217  if (gauges(i)%q_out_vars(n)) then
218  write(q_column(3 * index + 2:4 + 3 * index), "(i3)") n
219  index = index + 1
220  end if
221  end do
222  q_column(3 * index + 2:4 + 3 * index) = "]"
223 
224  aux_column = "["
225  index = 0
226  if (allocated(gauges(i)%aux_out_vars)) then
227  do n=1, size(gauges(i)%aux_out_vars, 1)
228  if (gauges(i)%aux_out_vars(n)) then
229  write(aux_column(3 * index + 2:4 + 3 * index), "(i3)") n
230  index = index + 1
231  end if
232  end do
233  end if
234  aux_column(3 * index + 2:4 + 3 * index) = "]"
235 
236  write(outgaugeunit, "(a,a,a,a)") "# level, time, q", &
237  trim(q_column), ", aux", &
238  trim(aux_column)
239  endif
240 
241  close(outgaugeunit)
242 
243  end do
244 
245  module_setup = .true.
246  end if
247 
integer, parameter maxgr
Definition: amr_module.f90:173
subroutine opendatafile(iunit, fname)
Definition: opendatafile.f:3
logical, private module_setup
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
Here is the call graph for this function:
Here is the caller graph for this function:

◆ setbestsrc()

subroutine gauges_module::setbestsrc ( )

Definition at line 255 of file gauges_module.f90.

References amr_module::cornxhi, amr_module::cornxlo, amr_module::cornyhi, amr_module::cornylo, gauges, amr_module::levelptr, amr_module::lfine, amr_module::lstart, mbestg1, mbestg2, mbestorder, mbestsrc, amr_module::node, num_gauges, qsorti(), and amr_module::rnode.

Referenced by tick().

255 !
256 ! Called every time grids change, to set the best source grid patch
257 ! for each gauge, i.e. the finest level patch that includes the gauge.
258 !
259 ! lbase is grid level that didn't change, but since fine
260 ! grid may have disappeared, we still have to look starting
261 ! at coarsest level 1.
262 !
263  use amr_module
264  implicit none
265 
266  integer :: lev, mptr, i, k1, ki
267 
268 !
269 ! ## set source grid for each loc from coarsest level to finest.
270 ! ## that way finest src grid left and old ones overwritten
271 ! ## this code uses fact that grids do not overlap
272 
273 ! # for debugging, initialize sources to 0 then check that all set
274  mbestsrc = 0
275 
276  do lev = 1, lfine
277  mptr = lstart(lev)
278  do
279  do i = 1, num_gauges
280  if ((gauges(i)%x >= rnode(cornxlo,mptr)) .and. &
281  (gauges(i)%x <= rnode(cornxhi,mptr)) .and. &
282  (gauges(i)%y >= rnode(cornylo,mptr)) .and. &
283  (gauges(i)%y <= rnode(cornyhi,mptr)) ) then
284  mbestsrc(i) = mptr
285  end if
286  end do
287  mptr = node(levelptr, mptr)
288  if (mptr == 0) exit
289  end do
290  end do
291 
292 
293  do i = 1, num_gauges
294  if (mbestsrc(i) .eq. 0) &
295  print *, "ERROR in setting grid src for gauge data", i
296  end do
297 
298  ! Sort the source arrays for easy testing during integration
299  call qsorti(mbestorder,num_gauges,mbestsrc)
300 
301 ! After sorting,
302 ! mbestsrc(mbestorder(i)) = grid index to be used for gauge i
303 ! and mbestsrc(mbestorder(i)) is non-decreasing as i=1,2,..., num_gauges
304 
305 ! write(6,*) '+++ mbestorder: ',mbestorder
306 ! write(6,*) '+++ mbestsrc: ',mbestsrc
307 
308 ! Figure out the set of gauges that should be handled on each grid:
309 ! after loop below, grid k should handle gauges numbered
310 ! mbestorder(i) for i = mbestg1(k), mbestg1(k)+1, ..., mbestg2(k)
311 ! This will be used for looping in print_gauges subroutine.
312 
313  ! initialize arrays to default indicating grids that contain no gauges:
314  mbestg1 = 0
315  mbestg2 = 0
316 
317  k1 = 0
318  do i=1,num_gauges
319  ki = mbestsrc(mbestorder(i))
320  if (ki > k1) then
321  ! new grid number seen for first time in list
322  if (k1 > 0) then
323  ! mark end of gauges seen by previous grid
324  mbestg2(k1) = i-1
325 ! write(6,*) '+++ k1, mbestg2(k1): ',k1,mbestg2(k1)
326  endif
327  mbestg1(ki) = i
328 ! write(6,*) '+++ ki, mbestg1(ki): ',ki,mbestg1(ki)
329  endif
330  k1 = ki
331  enddo
332  if (num_gauges > 0) then
333  ! finalize
334  mbestg2(ki) = num_gauges
335 ! write(6,*) '+++ ki, mbestg2(ki): ',ki,mbestg2(ki)
336  endif
integer, parameter cornxlo
x-coordinate of the left border of this grid
Definition: amr_module.f90:143
subroutine qsorti(ORD, N, A)
Definition: quick_sort1.f:23
integer, dimension(nsize, maxgr) node
Definition: amr_module.f90:198
real(kind=8), dimension(rsize, maxgr) rnode
Definition: amr_module.f90:193
integer, dimension(maxlv) lstart
Definition: amr_module.f90:198
integer, parameter cornylo
y-coordinate of the lower border of this grid
Definition: amr_module.f90:145
integer, parameter levelptr
node number (index) of next grid on the same level
Definition: amr_module.f90:35
integer, parameter cornxhi
x-coordinate of the right border of this grid
Definition: amr_module.f90:147
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
integer lfine
Definition: amr_module.f90:198
integer, parameter cornyhi
y-coordinate of the upper border of this grid
Definition: amr_module.f90:149
Here is the call graph for this function:
Here is the caller graph for this function:

◆ update_gauges()

subroutine gauges_module::update_gauges ( real(kind=8), dimension(num_eqn, mitot, mjtot), intent(in)  q,
real(kind=8), dimension(num_aux, mitot, mjtot), intent(in)  aux,
real(kind=8), intent(in)  xlow,
real(kind=8), intent(in)  ylow,
integer, intent(in)  num_eqn,
integer, intent(in)  mitot,
integer, intent(in)  mjtot,
integer, intent(in)  num_aux,
integer, intent(in)  mptr 
)

Definition at line 344 of file gauges_module.f90.

References gauges, amr_module::hxposs, amr_module::hyposs, max_buffer, amr_module::maxaux, amr_module::maxvar, mbestg1, mbestg2, mbestorder, mbestsrc, amr_module::nestlevel, amr_module::nghost, amr_module::node, num_gauges, print_gauges_and_reset_nextloc(), amr_module::rnode, and amr_module::timemult.

Referenced by par_advanc().

344  !
345  ! This routine is called each time step for each grid patch, to output
346  ! gauge values for all gauges for which this patch is the best one to
347  ! use (i.e. at the finest refinement level).
348 
349  ! It is called after ghost cells have been filled from adjacent grids
350  ! at the same level, so bilinear interpolation can be used to
351  ! to compute values at any gauge location that is covered by this grid.
352 
353  ! The grid patch is designated by mptr.
354  ! We only want to set gauges i for which mbestsrc(i) == mptr.
355  ! The array mbestsrc is reset after each regridding to indicate which
356  ! grid patch is best to use for each gauge.
357 
358  ! This is a refactoring of dumpgauge.f from Clawpack 5.2
359  ! Loops over only the gauges to be handled by this grid, as specified
360  ! by indices from mbestg1(mptr) to mbestg2(mptr)
361 
363  use amr_module, only: maxaux, hxposs, hyposs
364 
365  implicit none
366 
367  ! Input
368  integer, intent(in) :: num_eqn, mitot, mjtot, num_aux, mptr
369  real(kind=8), intent(in) :: q(num_eqn, mitot, mjtot)
370  real(kind=8), intent(in) :: aux(num_aux, mitot, mjtot)
371  real(kind=8), intent(in) :: xlow, ylow
372 
373  ! Locals
374  real(kind=8) :: var(maxvar + maxaux)
375  real(kind=8) :: xcent, ycent, xoff, yoff, tgrid, hx, hy
376  integer :: i, j, i1, i2, iindex, jindex, n, ii, index, level, var_index
377 
378  ! No gauges to record, exit
379  if (num_gauges == 0) then
380  return
381  endif
382 
383  i1 = mbestg1(mptr)
384  i2 = mbestg2(mptr)
385 
386  if (i1 == 0) then
387  ! no gauges to be handled by this grid
388  return
389  endif
390 
391  ! Grid info
392  tgrid = rnode(timemult, mptr)
393  level = node(nestlevel, mptr)
394  hx = hxposs(level)
395  hy = hyposs(level)
396 
397  ! Main Gauge Loop ======================================================
398  do i = i1, i2
399  ii = mbestorder(i)
400  if (mptr /= mbestsrc(ii)) then
401  print *, '*** should not happen... i, ii, mbestsrc(ii), mptr:'
402  print *, i, ii, mbestsrc(ii), mptr
403  stop
404  endif
405  if (tgrid < gauges(ii)%t_start .or. tgrid > gauges(ii)%t_end) then
406  cycle
407  endif
408  ! Minimum increment
409  ! TODO Maybe always allow last time output recording?
410  if (tgrid - gauges(ii)%last_time < gauges(ii)%min_time_increment) then
411  cycle
412  end if
413 
414  ! Compute indexing and bilinear interpolant weights
415  ! Note: changed 0.5 to 0.5d0 etc.
416  iindex = int(.5d0 + (gauges(ii)%x - xlow) / hx)
417  jindex = int(.5d0 + (gauges(ii)%y - ylow) / hy)
418  if ((iindex < nghost .or. iindex > mitot-nghost) .or. &
419  (jindex < nghost .or. jindex > mjtot-nghost)) then
420  print *, "ERROR in output of Gauge Data "
421  end if
422  xcent = xlow + (iindex - 0.5d0) * hx
423  ycent = ylow + (jindex - 0.5d0) * hy
424  xoff = (gauges(ii)%x - xcent) / hx
425  yoff = (gauges(ii)%y - ycent) / hy
426 
427  ! Gauge interpolation seems to work, so error test is commented out.
428  ! For debugging, use the code below...
429  ! Note: we expect 0 <= xoff, yoff <= 1 but if gauge is exactly
430  ! at center of cell these might be off by rounding error
431 
432  !if (xoff .lt. -1.d-4 .or. xoff .gt. 1.0001d0 .or. &
433  ! yoff .lt. -1.d-4 .or. yoff .gt. 1.0001d0) then
434  ! write(6,*) "*** print_gauges: Interpolation problem at gauge ",&
435  ! igauge(ii)
436  ! write(6,*) " xoff,yoff: ", xoff,yoff
437  !endif
438 
439  ! Bilinear interpolation
440  var_index = 0
441  do n = 1, size(gauges(ii)%q_out_vars, 1)
442  if (gauges(ii)%q_out_vars(n)) then
443  var_index = var_index + 1
444  var(var_index) = (1.d0 - xoff) * (1.d0 - yoff) * q(n, iindex, jindex) &
445  + xoff * (1.d0 - yoff) * q(n, iindex + 1, jindex) &
446  + (1.d0 - xoff) * yoff * q(n, iindex, jindex + 1) &
447  + xoff * yoff * q(n, iindex + 1, jindex + 1)
448  endif
449  end do
450 
451  if (allocated(gauges(ii)%aux_out_vars)) then
452  do n = 1, size(gauges(ii)%aux_out_vars, 1)
453  if (gauges(ii)%aux_out_vars(n)) then
454  var_index = var_index + 1
455  var(var_index) = (1.d0 - xoff) * (1.d0 - yoff) * aux(n, iindex, jindex) &
456  + xoff * (1.d0 - yoff) * aux(n, iindex + 1, jindex) &
457  + (1.d0 - xoff) * yoff * aux(n, iindex, jindex + 1) &
458  + xoff * yoff * aux(n, iindex + 1, jindex + 1)
459  endif
460  end do
461  end if
462 
463  ! Check to make sure we grabbed all the values
464  if (gauges(ii)%num_out_vars /= var_index) then
465  print *, gauges(ii)%num_out_vars, var_index
466  print *, gauges(ii)%q_out_vars
467  print *, gauges(ii)%aux_out_vars
468  stop "Somehow we did not grab all the values we wanted..."
469  end if
470 
471  ! Zero out tiny values to prevent underflow problems
472  do j = 1, gauges(ii)%num_out_vars
473  if (abs(var(j)) < 1d-90) var(j) = 0.d0
474  end do
475 
476  ! save info for this time
477  index = gauges(ii)%buffer_index
478 
479  gauges(ii)%level(index) = level
480  gauges(ii)%data(1,index) = tgrid
481  do j = 1, gauges(ii)%num_out_vars
482  gauges(ii)%data(1 + j, index) = var(j)
483  end do
484 
485  gauges(ii)%buffer_index = index + 1
486  if (gauges(ii)%buffer_index > max_buffer) then
488  endif
489 
490  gauges(ii)%last_time = tgrid
491 
492  end do ! End of gauge loop =============================================
493 
integer, parameter timemult
current simulation time on this grid
Definition: amr_module.f90:151
integer, parameter maxvar
Definition: amr_module.f90:183
subroutine print_gauges_and_reset_nextloc(gauge_num)
real(kind=8), dimension(maxlv) hyposs
Definition: amr_module.f90:193
real(kind=8), dimension(maxlv) hxposs
Definition: amr_module.f90:193
integer, dimension(nsize, maxgr) node
Definition: amr_module.f90:198
integer, parameter nestlevel
AMR level of the grid.
Definition: amr_module.f90:44
real(kind=8), dimension(rsize, maxgr) rnode
Definition: amr_module.f90:193
integer, parameter maxaux
Definition: amr_module.f90:184
integer nghost
Definition: amr_module.f90:232
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ gauges

type(gauge_type), dimension(:), allocatable gauges_module::gauges

Definition at line 81 of file gauges_module.f90.

Referenced by print_gauges_and_reset_nextloc(), set_gauges(), setbestsrc(), and update_gauges().

81  type(gauge_type), allocatable :: gauges(:)

◆ igauge

integer, dimension(:), allocatable gauges_module::igauge

Definition at line 84 of file gauges_module.f90.

Referenced by setgauges().

◆ max_buffer

integer, parameter gauges_module::max_buffer = 1000

Definition at line 49 of file gauges_module.f90.

Referenced by set_gauges(), and update_gauges().

49  integer, parameter :: max_buffer = 1000

◆ mbestg1

integer, dimension(:), allocatable gauges_module::mbestg1

Definition at line 84 of file gauges_module.f90.

Referenced by set_gauges(), setbestsrc(), and update_gauges().

◆ mbestg2

integer, dimension(:), allocatable gauges_module::mbestg2

Definition at line 84 of file gauges_module.f90.

Referenced by set_gauges(), setbestsrc(), and update_gauges().

◆ mbestorder

integer, dimension(:), allocatable gauges_module::mbestorder

Definition at line 84 of file gauges_module.f90.

Referenced by set_gauges(), setbestsrc(), and update_gauges().

◆ mbestsrc

integer, dimension(:), allocatable gauges_module::mbestsrc

Definition at line 84 of file gauges_module.f90.

Referenced by set_gauges(), setbestsrc(), setgauges(), and update_gauges().

84  integer, allocatable, dimension(:) :: mbestsrc, mbestorder, &
85  igauge, mbestg1, mbestg2

◆ module_setup

logical, private gauges_module::module_setup = .false.
private

Definition at line 44 of file gauges_module.f90.

Referenced by set_gauges().

44  logical, private :: module_setup = .false.
logical, private module_setup

◆ num_gauges

integer gauges_module::num_gauges

Definition at line 47 of file gauges_module.f90.

Referenced by amr2(), check(), par_advanc(), set_gauges(), setbestsrc(), tick(), and update_gauges().

47  integer :: num_gauges

◆ outgaugeunit

integer, parameter gauges_module::outgaugeunit = 89

Definition at line 46 of file gauges_module.f90.

Referenced by print_gauges_and_reset_nextloc(), set_gauges(), and setgauges().

46  integer, parameter :: outgaugeunit = 89