DISPATCH
download_mod.f90
1 !===============================================================================
2 !> download_link: takes care of downloads to link%task
3 !> same: called for patches on the same level
4 !> differ: called for patches on different levels
5 !>
6 !> Since we keep parent patches, there can be any number of levels overlapping
7 !> in the guard zones of a patch. Therefore, to avoid wasting cycles on
8 !> computing guard cell values for all levels, we sort the nbor list in order
9 !> of decreasing levels, and check if a cell is already filled before computing
10 !> a value from an nbor patch.
11 !>
12 !> The "filled" array stores the level of the value stored in the cell, and since
13 !> staggering can vary for different variables, there needs to be an array for
14 !> each variable.
15 !===============================================================================
17  USE iso_fortran_env, only: int8
18  USE io_mod
19  USE io_unit_mod
20  USE omp_mod
21  USE omp_timer_mod
22  USE mpi_mod
23  USE bits_mod
24  USE kinds_mod
25  USE mesh_mod
26  USE trace_mod
27  USE task_mod
28  USE patch_mod
29  USE link_mod
30  USE shared_mod
31  USE interpolate_mod
32  USE lagrange_mod
33  USE validate_mod
34  USE omp_lock_mod
35  USE timer_mod
36  USE guard_zones_mod
37  USE remesh_mod
38  implicit none
39  private
40  type, public:: download_t
41  integer(8):: n_same=0, n_differ=0
42  integer:: order_interpolator
43  logical:: check_filled=.false.
44  logical:: link_locks=.false.
45  procedure(interpolator_interface), pointer, nopass:: interpolate => null()
46  procedure(same_linear) , pointer :: same => null()
47  contains
48  procedure:: init
49  procedure:: download_link
50  procedure, nopass:: dnload_range
51  procedure:: prolong
52  procedure:: different
53  procedure:: different_prolong
54  procedure:: test
55  procedure, nopass:: set_verbose
56  end type
57  integer(kind=int8), dimension(:,:,:,:), allocatable:: filled
58  integer(kind=4), dimension(:,:,:), allocatable:: src
59  !$omp threadprivate (filled, src)
60  integer:: verbose=0
61  integer:: id_debug=0
62  integer:: use_different=0
63  logical:: sorted=.false.
64  logical:: use_locks=.true.
65  logical:: detailed_timer=.false.
66  character(len=16):: use_restrict='', use_prolong='', use_same=''
67  type(download_t), public:: download
68 CONTAINS
69 
70 !===============================================================================
71 !> Initialise downloading and interpolating options.
72 !===============================================================================
73 SUBROUTINE init (self)
74  class(download_t):: self
75  integer:: iostat
76  integer, save:: order_interpolator=-1, order_time=1
77  logical, save:: check_filled=.false.
78  logical, save:: first_time=.true.
79  namelist /download_params/ order_interpolator, order_time, check_filled, &
80  use_prolong, use_restrict, use_same, sorted, use_locks, id_debug, &
81  detailed_timer, verbose
82  character(len=120):: id = &
83  '$Id: bbf1b560f7d4227e33c65e71682b111f384e61b8 $ interpolation/download_mod.f90'
84  !-----------------------------------------------------------------------------
85  call trace%begin ('download_t%init')
86  call trace%print_id (id)
87  !$omp critical (download_cr)
88  if (first_time) then
89  first_time = .false.
90  rewind(io%input)
91  read (io%input, download_params, iostat=iostat)
92  if (io%master) write (*, download_params)
93  call guard_zones%init
94  call remesh%init
95  call lagrange%init
96  end if
97  if (order_time==1) then
98  self%same => same_linear
99  else
100  self%same => same_lagrange
101  end if
102  self%check_filled = check_filled
103  self%order_interpolator = order_interpolator
104  self%interpolate => selectinterpolator(order_interpolator)
105  !$omp end critical (download_cr)
106  call trace%end()
107 END SUBROUTINE init
108 
109 !===============================================================================
110 !> Select the best method for each case. So far, this only supports static
111 !> patches at the same level, or static patches one level up or down
112 !===============================================================================
113 SUBROUTINE download_link (self, link, all_cells, only, all_times)
114  class(download_t):: self
115  class(link_t), pointer:: link
116  logical, optional:: all_cells, all_times
117  integer, optional:: only
118  real(kind=KindScalarVar), dimension(:,:,:), pointer:: d
119  !.............................................................................
120  class(link_t), pointer:: nbor, sorted_head, nbors
121  class(patch_t), pointer:: target, source
122  integer:: edges, l(3), u(3), ix, iy, iz, jt(2)
123  integer, save:: itimer=0
124  real:: pt(2)
125  character(len=64):: info
126  logical:: all_cells_l, lock_it
127  integer:: n_not
128  !-----------------------------------------------------------------------------
129  call trace%begin('download_t%download_link', itimer=itimer)
130  if (.not.associated(link)) then
131  call io%abort('download_t%download_link: link not associated')
132  else if (.not.associated(link%task)) then
133  call io%abort('download_t%download_link: link%task not associated')
134  end if
135  !-----------------------------------------------------------------------------
136  ! Guard against download twice, which could happen if an extras (e.g. RT) task
137  ! requests dnload before the parent task
138  !-----------------------------------------------------------------------------
139  target => task2patch(link%task)
140  call target%log ('download', 2)
141  if (target%dnload_time >= target%time .or. &
142  target%is_set(bits%no_download)) then
143  call trace%end (itimer)
144  return
145  end if
146  target%dnload_time = target%time
147  if (target%is_clear(bits%no_density)) then
148  d => target%mem(:,:,:,target%idx%d,target%it,1)
149  call validate%check (target%link, d, 'before download')
150  if (target%id == id_debug .or. verbose > 0) then
151  write(io_unit%log,*) 'download_t%download_link: id, time =', &
152  target%id, target%time
153  end if
154  end if
155  if (.not.associated(link%nbor)) then
156  sorted = .false.
157  end if
158  !-----------------------------------------------------------------------------
159  ! Allocate and initialise array to track filling of cells.
160  ! A value of `-1` defines unfilled; all other values denote the *level* of the
161  ! patch that filled that cell.
162  !-----------------------------------------------------------------------------
163  if (present(all_cells)) then
164  all_cells_l = all_cells
165  else if (target%quality == 1024) then
166  all_cells_l = .true.
167  else
168  all_cells_l = .false.
169  end if
170  if (self%check_filled) then
171  allocate (filled(target%gn(1),target%gn(2),target%gn(3),target%nv))
172  if (verbose > 0) &
173  allocate (src(target%gn(1),target%gn(2),target%gn(3)))
174  l = target%li
175  u = target%ui
176  filled = -1
177  if (.not.all_cells_l) then
178  filled(l(1):u(1),l(2):u(2),l(3):u(3),:) = target%level
179  if (verbose > 0) &
180  src(l(1):u(1),l(2):u(2),l(3):u(3)) = target%id
181  end if
182  end if
183  if (target%is_clear(bits%no_density) .and. .not.all_cells_l) &
184  call target%check_density ("before download")
185  if (target%is_periodic()) call target%periodic_grid
186  !-----------------------------------------------------------------------------
187  ! Use an nbor list sorted by decreasing level, recording in each cell the
188  ! level that supplied values, and skipping nbor contributions from lower levels.
189  ! To avoid dead-lock, make a copy of it, and release the lock
190  !-----------------------------------------------------------------------------
191  nbor => link%nbors_by_level
192  !-----------------------------------------------------------------------------
193  ! Insert a value only if "filled" indicates that one is missing. With the
194  ! nbor list sorted in decreasing level order, values are only computed where
195  ! actually needed
196  !-----------------------------------------------------------------------------
197  do while (associated(nbor))
198  if (.not.associated(nbor%task)) then
199  write (io_unit%log,*) 'WARNING: nbor task not associated'
200  flush (io_unit%log)
201  nbor => nbor%next
202  cycle
203  end if
204  source => task2patch(nbor%task)
205  call source%log ('download', 3)
206  if (source%quality == 1024) then
207  all_cells_l = .true.
208  end if
209  if (verbose > 1) then
210  write(io_unit%log,'(f12.6,2x,a,2(i6,i3),2x,2l1,i4)') &
211  wallclock(), 'target, level, source, level, ND, n_needed =', &
212  target%id, target%level, source%id, source%level, &
213  nbor%needed, nbor%download, source%n_needed
214  flush (io_unit%log)
215  end if
216  !-------------------------------------------------------------------------
217  ! If the source task is being actively updated, and the slot could become
218  ! the one after source%new, lock the task mem, so either it cannot change
219  ! during the guard zone loading or, if the lock is already set, access to
220  ! the source will have to wait until after the source update.
221  !-------------------------------------------------------------------------
222  call source%time_interval (target%time, jt, pt)
223  call source%mem_lock(jt(1))%set ('download_link1')
224  if (jt(2) /= jt(1)) &
225  call source%mem_lock(jt(2))%set ('download_link2')
226  if (verbose > 0) &
227  write (io_unit%log,'(a,2i3,1x,2i3,2x,2f6.3,f12.6,2x,10f12.6)') &
228  'download_link: it, new, jt, pt, target%time, source%t =', &
229  source%it, source%new, jt, pt, target%time, source%t
230  if (verbose > 1) then
231  write(io_unit%log,'(f12.6,i6,i3,3x,a,i6,i3,2x,2l1)') wallclock(), target%id, target%level, &
232  'download_link: source, level =', &
233  source%id, source%level, nbor%download, lock_it
234  flush (io_unit%log)
235  end if
236  !---------------------------------------------------------------------------
237  ! Skip nbor not marked for download
238  !---------------------------------------------------------------------------
239  if (nbor%download) then
240  if (source%nt == 1) &
241  self%same => same_linear
242  if (trim(target%kind) /= 'rt_solver') then
243  if (source%is_clear(bits%no_download)) then
244  !----------------------------------------------------------------------
245  ! Default same interpolation when target and source levels agree
246  !----------------------------------------------------------------------
247  if (source%level==target%level) then
248  select case (trim(use_same))
249  case ('different')
250  call different (self, target, source, jt, pt, all_cells_l, only=only)
251  case ('guard_zones')
252  call guard_zones%method (target, source, jt, pt)
253  case ('guard_same')
254  call guard_zones%method (target, source, jt, pt, same=.true.)
255  case default
256  call self%same (target, source, jt, pt, all_cells_l, only=only)
257  end select
258  !----------------------------------------------------------------------
259  ! Prolong with one level difference and all_cells
260  !----------------------------------------------------------------------
261  else if (source%level==target%level-1 .and. present(all_cells)) then
262  call remesh%prolong (target, source, jt, pt, all_cells)
263  !----------------------------------------------------------------------
264  ! Prolong cases
265  !----------------------------------------------------------------------
266  else if ((source%level-target%level) <= -1) then
267  select case (trim(use_prolong))
268  case ('different')
269  call different (self, target, source, jt, pt, all_cells_l, only=only)
270  case ('prolong')
271  call prolong (self, target, source, jt, pt, all_cells_l, only=only)
272  case ('guard_zones')
273  call guard_zones%method (target, source, jt, pt)
274  case ('remesh')
275  call remesh%prolong (target, source, jt, pt, all_cells_l)
276  case default
277  call different_prolong (self, target, source, jt, pt, all_cells_l, only=only)
278  end select
279  !----------------------------------------------------------------------
280  ! Prolong cases with all_cell
281  !----------------------------------------------------------------------
282  else if (source%level < target%level .and. present(all_cells)) then
283  call different_prolong (self, target, source, jt, pt, all_cells_l, only=only)
284  !----------------------------------------------------------------------
285  ! Restrict cases
286  !----------------------------------------------------------------------
287  else if ((source%level-target%level) == 1) then
288  select case (trim(use_restrict))
289  case ('different')
290  call different (self, target, source, jt, pt, all_cells_l, only=only)
291  case ('diff_pro')
292  call different_prolong (self, target, source, jt, pt, all_cells_l, only=only)
293  case ('guard_zones')
294  call guard_zones%method (target, source, jt, pt)
295  case default
296  call different_restrict (self, target, source, jt, pt, only=only)
297  end select
298  end if
299  end if
300  if (target%id==id_debug) then
301  write(io_unit%log,'(a,i6,1p,4e12.5)') 'source: id,min,max,fmin,fmaxval:', &
302  source%id, &
303  minval(source%mem(:,:,:,source%idx%d,source%it,1)), &
304  maxval(source%mem(:,:,:,source%idx%d,source%it,1)), &
305  source%fminval(source%idx%d), &
306  source%fmaxval(source%idx%d)
307  write(io_unit%log,'(a,i6,1p,4e12.5)') 'target: id,min,max,fmin,fmaxval:', &
308  target%id, &
309  minval(target%mem(:,:,:,target%idx%d,target%it,1)), &
310  maxval(target%mem(:,:,:,target%idx%d,target%it,1)), &
311  target%fminval(target%idx%d), &
312  target%fmaxval(target%idx%d)
313  end if
314  else
315  if (trim(source%kind)=='rt_solver') then
316  if (abs(source%level-target%level) < 2) then
317  if (source%level==target%level) then
318  call same_rt (self, target, source, jt, pt)
319  else
320  call different_rt (self, target, source, jt, pt)
321  end if
322  end if
323  end if
324  end if
325  end if
326  if (jt(2) /= jt(1)) &
327  call source%mem_lock(jt(2))%unset ('download_link2')
328  call source%mem_lock(jt(1))%unset ('download_link1')
329  if (verbose > 1) then
330  write(io_unit%log,'(f12.6,i6,i3,3x,a,i6,i3,2x,2l1)') wallclock(), target%id, target%level, &
331  'download_link: source, level =', source%id, source%level
332  flush (io_unit%log)
333  end if
334  nbor => nbor%next
335  end do
336  if (present(only)) then
337  if (only==target%idx%d) then
338  call target%check_density (" after download")
339  end if
340  else
341  call target%check_density (" after download")
342  end if
343  if (self%check_filled) then
344  if (all_cells_l) then
345  l = target%mesh%li
346  u = target%mesh%ui
347  else
348  l = target%mesh%lb
349  u = target%mesh%ub
350  end if
351  n_not = 0
352  if (any(filled(l(1):u(1),l(2):u(2),l(3):u(3),:) == -1)) then
353  !-------------------------------------------------------------------------
354  ! Since only the parent patch is used, this condition is expected, and OK
355  !-------------------------------------------------------------------------
356  if (verbose > 1) then
357  if (all_cells_l) then
358  write (io_unit%log,'(a,i6)') &
359  'INFO: some internal cells were not filled, iD =', target%id
360  else
361  write (io_unit%log,'(a,i6)') &
362  'INFO: some guard zone cells were not filled, iD =', target%id
363  end if
364  do iz=l(3),u(3)
365  do iy=l(2),u(2)
366  do ix=l(1),u(1)
367  if (filled(ix,iy,iz,1) == -1) then
368  n_not = n_not+1
369  write(io_unit%log,'(a,3i3)') 'not filled:', ix, iy, iz
370  end if
371  end do
372  end do
373  end do
374  if (all_cells_l) then
375  write (io%output,'(i6,2x,a,i6,i4)') &
376  n_not, 'internal cells were not filled: id, level =', &
377  target%id, target%level
378  else
379  write (io%output,'(i6,2x,a,i6,i4)') &
380  n_not, 'guard zone cells were not filled: id, level =', &
381  target%id, target%level
382  end if
383  end if
384  end if
385  associate(task => link%task)
386  select type (task)
387  class is (patch_t)
388  task%not_filled = n_not
389  end select
390  end associate
391  deallocate (filled)
392  if (verbose > 0) &
393  deallocate (src)
394  end if
395  !-----------------------------------------------------------------------------
396  call guard_zones%info()
397  if (target%is_clear(bits%no_density)) &
398  call validate%check (target%link, d, ' after download')
399  call trace%end (itimer)
400 END SUBROUTINE download_link
401 
402 !===============================================================================
403 !> Find the target index ranges that correspond to the overlap region, for
404 !> download from source (source has higher, same, or lower level than target).
405 !>
406 !> The first index, l(i), returned is the index of the first target point that
407 !> is inside the authoritative region; i.e., it must lie to the right of the
408 !> source point r(li).
409 !> FIXME: To be validated/improved for moving meshes
410 !===============================================================================
411 SUBROUTINE dnload_range (target, source, iv, l, u, offset)
412  class(patch_t), pointer :: target, source
413  integer :: iv, l(3), u(3), offset(3), li, ui
414  optional :: offset
415  !.............................................................................
416  real(8), dimension(3) :: h, dist, la, ua
417  integer :: i
418  class(mesh_t), pointer :: sm, tm
419  logical :: debug
420  !-----------------------------------------------------------------------------
421  call trace%begin ('download_t%dnload_range', 2)
422  !-----------------------------------------------------------------------------
423  ! First take care of mapping the distance between the source and target
424  ! centers, periodically, so this can be ignored in the next steps.
425  !-----------------------------------------------------------------------------
426  dist = source%position - target%position
427  if (target%id == id_debug) write (io_unit%log,*) 'dist(1)', dist
428  where (source%periodic)
429  dist = modulo(dist+0.5_8*source%mesh%b,source%mesh%b) - 0.5_8*source%mesh%b
430  end where
431  if (target%id == id_debug) write (io_unit%log,*) 'dist(2)', dist
432  !-----------------------------------------------------------------------------
433  ! Add the distances from the source center to the start and end of the
434  ! authoritative domain of each variable.
435  !-----------------------------------------------------------------------------
436  debug = source%level /= target%level
437  la = dist - source%size/2d0
438  ua = dist + source%size/2d0
439  if (target%id == id_debug) write (io_unit%log,*) 'la,ua(1)', la, ua
440  if (source%level /= target%level) then
441  do i=1,3
442  sm => source%mesh(i)
443  tm => target%mesh(i)
444  la(i) = la(i) + sm%h(iv)*sm%d - tm%h(iv)*tm%d
445  ua(i) = ua(i) + sm%h(iv)*sm%d - tm%h(iv)*tm%d
446  end do
447  end if
448  if (target%id == id_debug) write (io_unit%log,*) 'la,ua(2)', la, ua
449  l = ceiling(la/target%mesh%d + target%mesh%o)
450  u = floor(ua/target%mesh%d + target%mesh%o)
451  if (target%id == id_debug) write (io_unit%log,*) 'l,u(1)', l, u
452  !-----------------------------------------------------------------------------
453  ! Make sure the range is legal. The u index can become smaller than lb, which
454  ! means there is no overlap. Likewise, the l index can become larger than ub
455  !-----------------------------------------------------------------------------
456  l = max(l,target%mesh%lb)
457  u = min(u,target%mesh%ub)
458  if (target%id == id_debug) write (io_unit%log,*) 'l,u(2)', l, u
459  !-----------------------------------------------------------------------------
460  ! Account for collapsed dimensions
461  !-----------------------------------------------------------------------------
462  where (target%n == 1)
463  l = min(l,target%mesh%ub)
464  u = max(u,target%mesh%lb)
465  end where
466  if (target%id == id_debug) write (io_unit%log,*) 'l,u(3)', l, u
467  !-----------------------------------------------------------------------------
468  ! For same level patch, the integer offset is what should be added to the
469  ! target index to get the source index. For aligned meshes, offset=0.
470  !-----------------------------------------------------------------------------
471  if (present(offset)) then
472  if (source%level==target%level) then
473  offset = nint(-dist/target%ds + source%offset - target%offset)
474  else
475  write(io_unit%log,*) &
476  'ERROR: fixed offset not valid for different levels', &
477  target%id, target%level, source%id, source%level
478  call mpi%abort('dnload_range')
479  end if
480  if (verbose > 1 .or. io_unit%do_validate .or. target%id==id_debug) &
481  write(io_unit%log,1) &
482  target%id, target%level, iv, source%id, source%level, l, u, offset
483  1 format("download_t%dnload_range: id,lv,iv:",i6,2i3, &
484  " source,lv,l,u,o:",i6,i3,3(2x,3i4))
485  if (target%id == id_debug) write (io_unit%log,*) 'l,u,o', l, u, offset
486  else
487  if (verbose > 1 .or. io_unit%do_validate .or. target%id==id_debug) &
488  write(io_unit%log,2) &
489  target%id, target%level, iv, source%id, source%level, l, u
490  2 format("download_t%dnload_range: id,lv,iv:",i6,2i3, &
491  " source,lv,l,u:",i6,i3,3(2x,3i4))
492  end if
493  call trace%end()
494 contains
495 !===============================================================================
496 !> Local function that returns the index in target index space, given a position
497 !> p, in relative target coordinates.
498 !===============================================================================
499  function target_index (p, roundup, nearest) result (ii)
500  real(8) :: p(3)
501  logical, optional :: roundup, nearest
502  integer :: ii(3)
503  !---------------------------------------------------------------------------
504  p = p/target%mesh%d + target%offset ! float index
505  if (present(roundup)) then
506  ii = ceiling(p) ! round up
507  else if (present(nearest)) then
508  ii = nint(p) ! nearest int
509  else
510  ii = p + 0.0001 ! integer part
511  end if
512  end function target_index
513 END SUBROUTINE dnload_range
514 
515 !===============================================================================
516 !> Copy guard zone info from source patch at the same level
517 !>
518 !> same
519 !> dnload_range
520 !> time_interval
521 !===============================================================================
522 SUBROUTINE same_linear (self, target, source, jt, pt, all_cells, only, all_times)
523  class(download_t):: self
524  class(patch_t), pointer :: target, source
525  logical, optional :: all_cells, all_times
526  integer, optional :: only
527  !.............................................................................
528  integer :: l(3), u(3), o(3), ii(3), jt(2), i
529  real(8) :: position(3), xx, dist(3), start
530  real :: pt(2), fmin
531  integer :: ix, iy, iz, iv, iv1, iv2, li(3), ui(3), cells
532  integer, save :: itimer=0
533  logical :: all_cells_l, no
534  !-----------------------------------------------------------------------------
535  call trace%begin('download_t%same_linear', itimer=itimer)
536  if (detailed_timer) then
537  cells = 0
538  start = wallclock()
539  end if
540  if (target%id == id_debug) &
541  write(io_unit%log,'(a,1p,e14.5,i7,e14.5)') &
542  'dbg download_t%same_linear: time, source =', &
543  target%time, source%id, source%time
544  if (present(only)) then
545  iv1 = only
546  iv2 = only
547  else
548  iv1 = 1
549  iv2 = target%nv
550  end if
551  if (present(all_cells)) then
552  all_cells_l = all_cells
553  else
554  all_cells_l = .false.
555  end if
556  if (verbose > 1) then
557  write(io_unit%log,'(a,2i7,2i4)') &
558  ' same_linear: target, source, iv1, iv2 =', target%id, source%id, iv1, iv2
559  end if
560  if (iv2 > source%nv) then
561  write(stderr,*) 'iv loop limit:', iv1, iv2
562  write(stderr,*) 'source: id, nv =', source%id, source%nv
563  write(stderr,*) 'source kind ',source%kind, source%position
564  write(stderr,*) 'target kind ',target%kind, target%position
565  write(stderr,*) 'target: id, nv =', target%id, target%nv
566  call io%abort("download::same_linear::iv2 > source%nv")
567  end if
568  do iv=iv1, iv2
569  if (verbose > 1 .and. source%unsigned(iv)) then
570  fmin = min(minval(source%mem(:,:,:,iv,jt(1),1)), &
571  minval(source%mem(:,:,:,iv,jt(2),1)))
572  if (fmin <= 0.0) then
573  write(io_unit%log,*) 'iv, jt, pt / minval, minloc =', iv, jt, pt
574  write(io_unit%log,*) minval(source%mem(:,:,:,iv,jt(1),1)), &
575  minloc(source%mem(:,:,:,iv,jt(1),1))
576  write(io_unit%log,*) minval(source%mem(:,:,:,iv,jt(2),1)), &
577  minloc(source%mem(:,:,:,iv,jt(2),1))
578  flush(io_unit%log)
579  end if
580  end if
581  if (iv == target%idx%phi) then
582  if (verbose > 1) write(io_unit%log,*) &
583  target%id, ' using same level potential from', source%id
584  end if
585  if (source%quality > 1000.) then
586  if (iv==target%idx%d .or. (iv >= target%idx%px .and. iv <= target%idx%pz)) then
587  all_cells_l = .true.
588  write (io_unit%log,*) 'dnload all cells for target, source =', &
589  target%id, source%id, source%quality, source%level
590  end if
591  end if
592  ! -- `dphi` is part of `mem` when `time_order` = 0; do not download it.
593  if (iv == target%idx%dphi) cycle
594  li = target%mesh%li
595  ui = target%mesh%ui
596  call dnload_range (target, source, iv, l, u, o)
597  if (target%id==id_debug .or. verbose > 1) &
598  write(io_unit%log,'(2i7," a l:",3i3,3x,"u:",3i3,3x,"o:",3i3)') &
599  target%id,source%id,l,u,o
600  l = max(l,source%mesh%lb-o)
601  u = min(u,source%mesh%ub-o)
602  if (target%id==id_debug .or. verbose > 1) then
603  write(io_unit%log,'(2i7," b l:",3i3,3x,"u:",3i3,3x,"o:",3i3)') &
604  target%id,source%id,l,u,o
605  flush(io_unit%log)
606  end if
607  if (product(max(0,u-l+1))>0) then
608  !-------------------------------------------------------------------------
609  ! Extend the pick-up region to fill the guard zones also in the BC regions
610  !-------------------------------------------------------------------------
611  do i=1,3
612  !if (source%mesh(i)%lower_boundary) l(i) = target%mesh(i)%lb
613  !if (source%mesh(i)%upper_boundary) u(i) = target%mesh(i)%ub
614  !if (i==2) then
615  !if (source%boundaries%is_set (bits%yl)) l(i) = target%mesh(i)%lb
616  !if (source%boundaries%is_set (bits%yu)) u(i) = target%mesh(i)%ub
617  !end if
618  ! if (target%mesh(i)%h(iv) /= 0.0) li(i) = li(i)+1
619  ! In no-no-mans-land, the `li`-th cell of variables that are staggered in
620  ! direction `i` lie just outside the patch edge and should therefore be
621  ! interpolated. Currently only implemented for ZEUS solvers; TBD if it
622  ! should be applied to STAGGER solvers.
623  if (target%kind(1:4) == 'zeus' .and. &
624  target%mesh(i)%h(iv) == -0.5 .and. &
625  ui(i) > 1) &
626  li(i) = li(i)+1
627  end do
628  ! for staggered quantities, include the ub+1 zone for the appropriate
629  ! components note: p* and b* are face-centred in the same direction as the
630  ! component, while emfs are face-centred in both directions *orthogonal*
631  ! to the component direction.
632  if (target%kind(1:4) == 'zeus' .and. &
633  iv >= target%idx%px .and. iv <= target%idx%pz) then
634  do i=1,3
635  if (iv == target%idx%px+i-1 .and. target%n(i) > 1) u(i) = u(i) + 1
636  end do
637  end if
638  if (target%kind(1:4) == 'zeus' .and. &
639  iv >= target%idx%bx .and. &
640  iv <= target%idx%bz) then
641  do i=1,3
642  if (iv == target%idx%bx+i-1 .and. target%n(i) > 1) u(i) = u(i) + 1
643  end do
644  end if
645  if (target%kind(1:4) == 'zeus' .and. &
646  iv >= target%idx%ex .and. &
647  iv <= target%idx%ez) then
648  do i=1,3
649  if (iv /= target%idx%ex+i-1 .and. target%n(i) > 1) u(i) = u(i) + 1
650  end do
651  end if
652  if (target%id==id_debug) &
653  write(io_unit%log,1) source%id, l, u, o, jt, pt
654  1 format("download_t%same:", i6," l,u,o:",3(3i4,2x)," jt,pt:", &
655  2i4,2x,2f6.2)
656  !-------------------------------------------------------------------------
657  ! Special case, with nt==1
658  !-------------------------------------------------------------------------
659  cells = product(u-l+1)
660  if (source%nt == 1) then
661  if (all_cells_l) then
662  do iz=l(3),u(3)
663  do iy=l(2),u(2)
664  do ix=l(1),u(1)
665  target%mem(ix,iy,iz,iv,target%it,1) = &
666  source%mem(ix+o(1),iy+o(2),iz+o(3),iv,1,1)
667  end do
668  if (self%check_filled) then
669  do ix=l(1),u(1)
670  filled(ix,iy,iz,iv) = source%level
671  end do
672  if (verbose > 0) then
673  do ix=l(1),u(1)
674  src(ix,iy,iz) = source%id
675  end do
676  end if
677  end if
678  end do
679  end do
680  else
681  do iz=l(3),u(3)
682  do iy=l(2),u(2)
683  do ix=l(1),u(1)
684  ii = [ix,iy,iz]
685  if (ix+o(1) > size(source%mem,1)) then
686  print '(5(3i4,2x))', l, u, o, shape(source%mem)
687  end if
688  target%mem(ix,iy,iz,iv,target%it,1) = &
689  merge(target%mem(ix,iy,iz,iv,target%it,1), &
690  source%mem(ix+o(1),iy+o(2),iz+o(3),iv,1,1), &
691  all((ii>=li).and.(ii<=ui)))
692  end do
693  if (self%check_filled) then
694  do ix=l(1),u(1)
695  ii = [ix,iy,iz]
696  no = all((ii>=li).and.(ii<=ui))
697  filled(ix,iy,iz,iv) = merge(filled(ix,iy,iz,iv), &
698  int(source%level,kind=int8), no)
699  end do
700  if (verbose > 0) then
701  do ix=l(1),u(1)
702  ii = [ix,iy,iz]
703  no = all((ii>=li).and.(ii<=ui))
704  src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
705  end do
706  end if
707  end if
708  end do
709  end do
710  end if
711  !-------------------------------------------------------------------------
712  ! Unsigned values; interpolated in the log
713  !-------------------------------------------------------------------------
714  else if (source%unsigned(iv)) then
715  !-----------------------------------------------------------------------
716  ! Fill all cells, even interior ones, when all_cells==.true., but only
717  ! for cells that are not already filled. The filling is expensive in
718  ! this case, so use an if test.
719  !-----------------------------------------------------------------------
720  if (all_cells_l) then
721  do iz=l(3),u(3)
722  do iy=l(2),u(2)
723  if (self%check_filled) then
724  do ix=l(1),u(1)
725  if (filled(ix,iy,iz,iv) < source%level) then
726  target%mem(ix,iy,iz,iv,target%it,1) = exp( &
727  pt(1)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1)) + &
728  pt(2)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1)))
729  filled(ix,iy,iz,iv) = source%level
730  if (verbose > 0) &
731  src(ix,iy,iz) = source%id
732  else
733  write (io_unit%log,*)'WARNING: prevented overwrite', target%id, &
734  source%id, ix,iy,iz,o,iv, filled(ix,iy,iz,iv), src(ix,iy,iz)
735  end if
736  end do
737  else
738  do ix=l(1),u(1)
739  target%mem(ix,iy,iz,iv,target%it,1) = exp( &
740  pt(1)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1)) + &
741  pt(2)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1)))
742  end do
743  end if
744  end do
745  end do
746  !-----------------------------------------------------------------------
747  ! Make sure to change only guard cells, when all_cells==.false.. Here
748  ! the calculation is cheap, so we use a merge to allow vectorization.
749  !-----------------------------------------------------------------------
750  else
751  do iz=l(3),u(3)
752  do iy=l(2),u(2)
753  do ix=l(1),u(1)
754  ii = [ix,iy,iz]
755 ! if (source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1) <= 0.0) then
756 ! write (io_unit%log,9) wallclock(), &
757 ! l,u,ix,iy,iz,o,iv,source%it,source%new,1,jt, &
758 ! source%id,source%time,target%id,target%time,source%iit
759 ! 9 format(f12.6,6(3i4,2x),2(i6,f12.6),7i3)
760 ! end if
761 ! if (source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1) <= 0.0) then
762 ! write (io_unit%log,9) wallclock(), &
763 ! l,u,ix,iy,iz,o,iv,source%it,source%new,2,jt, &
764 ! source%id,source%time,target%id,target%time,source%iit
765 ! end if
766  target%mem(ix,iy,iz,iv,target%it,1) = &
767  merge(target%mem(ix,iy,iz,iv,target%it,1), &
768  exp(pt(1)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1)) + &
769  pt(2)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1))), &
770  all((ii>=li).and.(ii<=ui)))
771  end do
772  if (self%check_filled) then
773  do ix=l(1),u(1)
774  ii = [ix,iy,iz]
775  no = all((ii>=li).and.(ii<=ui))
776  filled(ix,iy,iz,iv) = merge(filled(ix,iy,iz,iv), &
777  int(source%level,kind=int8), no)
778  end do
779  if (verbose > 0) then
780  do ix=l(1),u(1)
781  ii = [ix,iy,iz]
782  no = all((ii>=li).and.(ii<=ui))
783  src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
784  end do
785  end if
786  end if
787  end do
788  end do
789  end if
790  !-------------------------------------------------------------------------
791  ! Signed values
792  !-------------------------------------------------------------------------
793  else
794  !-----------------------------------------------------------------------
795  ! All cells
796  !-----------------------------------------------------------------------
797  if (all_cells_l) then
798  do iz=l(3),u(3)
799  do iy=l(2),u(2)
800  !-----------------------------------------------------------------
801  ! Prevent overwriting better values
802  !-----------------------------------------------------------------
803  if (self%check_filled) then
804  do ix=l(1),u(1)
805  no = filled(ix,iy,iz,iv) > source%level
806  target%mem(ix,iy,iz,iv,target%it,1) = &
807  merge(target%mem(ix,iy,iz,iv,target%it,1), &
808  pt(1)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1) + &
809  pt(2)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1), &
810  no)
811  filled(ix,iy,iz,iv) = merge(filled(ix,iy,iz,iv), &
812  int(source%level, kind=int8), no)
813  end do
814  if (verbose > 0) then
815  do ix=l(1),u(1)
816  no = filled(ix,iy,iz,iv) > source%level
817  src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
818  end do
819  end if
820  !-----------------------------------------------------------------
821  ! No questions asked
822  !-----------------------------------------------------------------
823  else
824  do ix=l(1),u(1)
825  target%mem(ix,iy,iz,iv,target%it,1) = &
826  pt(1)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1) + &
827  pt(2)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1)
828  end do
829  end if
830  end do
831  end do
832  !-----------------------------------------------------------------------
833  ! Only guard cells. Prefer same level, so don't check for levels
834  !-----------------------------------------------------------------------
835  else
836  do iz=l(3),u(3)
837  do iy=l(2),u(2)
838  do ix=l(1),u(1)
839  ii = [ix,iy,iz]
840  target%mem(ix,iy,iz,iv,target%it,1) = &
841  merge(target%mem(ix,iy,iz,iv,target%it,1), &
842  pt(1)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1) + &
843  pt(2)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1), &
844  all((ii>=li).and.(ii<=ui)))
845  end do
846  if (self%check_filled) then
847  do ix=l(1),u(1)
848  ii = [ix,iy,iz]
849  no = all((ii>=li).and.(ii<=ui))
850  filled(ix,iy,iz,iv) = &
851  merge(filled(ix,iy,iz,iv), int(source%level,kind=int8), no)
852  end do
853  if (verbose > 0) then
854  do ix=l(1),u(1)
855  ii = [ix,iy,iz]
856  no = all((ii>=li).and.(ii<=ui))
857  src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
858  end do
859  end if
860  end if
861  end do
862  end do
863  end if
864  end if
865  if (verbose > 0 .or. target%id==id_debug) then
866  if (iv == 1 .or. verbose > 2) &
867  write (io_unit%log,2) target%id, target%level, source%id, source%level, &
868  'same_line picked up', &
869  product(u-l+1), ' values, l, u, mind, pt =', l, u, &
870  minval(target%mem(l(1):u(1),l(2):u(2),l(3):u(3),1,target%it,1)), pt
871  2 format(2(i6,i3),2x,a,i6,a,2(2x,3i4),1p,e14.4,2x,0p2f6.3)
872  end if
873  end if
874  end do
875  if (detailed_timer) then
876  if (cells > 0) then
877  write (stdout,'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus/cell")') &
878  "download_t%dsame_linear: target, level, source, level,", &
879  target%id, target%level, source%id, source%level, cells, &
880  (wallclock()-start)/max(1,cells)
881  else
882  write (stdout,'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus")') &
883  "download_t%same_linear: target, level, source, level,", &
884  target%id, target%level, source%id, source%level, cells, wallclock()-start
885  end if
886  end if
887  call trace%end (itimer)
888 END SUBROUTINE same_linear
889 
890 !===============================================================================
891 !> Copy guard zone info from source patch at the same level
892 !>
893 !> same
894 !> dnload_range
895 !> time_interval
896 !===============================================================================
897 SUBROUTINE same_lagrange (self, target, source, jj, pt, all_cells, only, all_times)
898  class(download_t):: self
899  class(patch_t), pointer :: target, source
900  logical, optional :: all_cells, all_times
901  integer, optional :: only
902  integer :: jj(2)
903  real :: pt(2)
904  !.............................................................................
905  integer :: l(3), u(3), o(3), ii(3), jt, i
906  real(8) :: position(3), xx, dist(3)
907  integer :: ix, iy, iz, iv, iv1, iv2, li(3), ui(3)
908  integer, save :: itimer=0, nprint=20
909  logical :: all_cells_l, filter
910  integer, parameter :: order=2
911  real(8) :: w(order+1)
912  integer :: j0, j1, iit(source%nt-1)
913  real(8) :: times(source%nt-1)
914  !-----------------------------------------------------------------------------
915  call trace%begin('download_t%same_lagrange', itimer=itimer)
916  if (target%id == id_debug) &
917  write(io_unit%log,*) 'dbg download_t%same_lagrange: time =', target%time
918  !-----------------------------------------------------------------------------
919  ! Compute weights for Lagrange interpolation
920  !-----------------------------------------------------------------------------
921  call source%timeslots (iit, times)
922  call lagrange%sequence_weights (target%time, times, j0, j1, w, order)
923  if (verbose > 2) &
924  write(io_unit%log,'(a,2i4,1p,10e14.6)') &
925  'same_lagrange: j0, j1, time, times, w:', &
926  j0, j1, target%time, times(j0:j1), w
927  if (present(only)) then
928  iv1 = only
929  iv2 = only
930  else
931  iv1 = 1
932  iv2 = target%nv
933  end if
934  if (present(all_cells)) then
935  all_cells_l = all_cells
936  else
937  all_cells_l = .false.
938  end if
939  do iv=iv1, iv2
940  if (iv == target%idx%phi) then
941  if (verbose>1) write(io_unit%log,*) &
942  target%id, ' using same level potential from', source%id
943  end if
944  if (source%quality > 1000.) then
945  if (iv==target%idx%d .or. &
946  (iv >= target%idx%px .and. iv <= target%idx%pz)) then
947  all_cells_l = .true.
948  write (io_unit%log,*) '(2)dnload all cells for target, source =', &
949  target%id, source%id, source%quality, source%level
950  end if
951  end if
952  ! -- `dphi` is part of `mem` when `time_order` = 0; do not download it.
953  if (iv == target%idx%dphi) cycle
954  li = target%mesh%li
955  ui = target%mesh%ui
956  call dnload_range (target, source, iv, l, u, o)
957  if (target%id==id_debug) &
958  write(io_unit%log,'(2i7," a l:",3i3,3x,"u:",3i3,3x,"o:",3i3)') &
959  target%id,source%id,l,u,o
960  l = max(l,source%mesh%lb-o)
961  u = min(u,source%mesh%ub-o)
962  if (target%id==id_debug) &
963  write(io_unit%log,'(2i7," b l:",3i3,3x,"u:",3i3,3x,"o:",3i3)') &
964  target%id,source%id,l,u,o
965  if (product(max(0,u-l+1))>0) then
966  !-------------------------------------------------------------------------
967  ! Extend the pick-up region to fill the guard zones also in the BC regions
968  !-------------------------------------------------------------------------
969  do i=1,3
970  if (source%mesh(i)%lower_boundary) l(i) = target%mesh(i)%lb
971  if (source%mesh(i)%upper_boundary) u(i) = target%mesh(i)%ub
972  end do
973  !-------------------------------------------------------------------------
974  ! Check if we need to filter away internal points
975  !-------------------------------------------------------------------------
976  filter = .false.
977  if (.not. all_cells_l) then
978  do iz=l(3),u(3)
979  do iy=l(2),u(2)
980  do ix=l(1),u(1)
981  ii = [ix,iy,iz]
982  if (all((ii>=li).and.(ii<=ui))) then
983  filter = .true.
984  exit
985  end if
986  end do
987  if (filter) exit
988  end do
989  if (filter) exit
990  end do
991  end if
992  if (filter) then
993  if (nprint > 0) then
994  write(io_unit%log,*) &
995 'WARNING: attempt to set internal cells in download_t%same_lagrange, using filter!'
996  nprint = nprint - 1
997  if (nprint==0) then
998  write(io_unit%log,*) &
999  'target id, position =', target%id, target%position
1000  write(io_unit%log,*) &
1001  'source id, position =', source%id, source%position
1002  write(io_unit%log,*) 'l , u =', l, u
1003  write(io_unit%log,*) 'gn, o =', target%gn, o
1004  end if
1005  end if
1006  end if
1007  !---------------------------------------------------------------------------
1008  ! Unsigned => use interpolation in the log
1009  !---------------------------------------------------------------------------
1010  if (source%unsigned(iv)) then
1011  do iz=l(3),u(3)
1012  do iy=l(2),u(2)
1013  if (filter) then
1014  do ix=l(1),u(1)
1015  ii = [ix,iy,iz]
1016  if (all_cells_l .or. any((ii<li).or.(ii>ui))) then
1017  do jt=j0,j1
1018  if (jt==j0) then
1019  target%mem(ix,iy,iz,iv,target%it,1) = &
1020  w(1+jt-j0)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1021  else
1022  target%mem(ix,iy,iz,iv,target%it,1) = &
1023  target%mem(ix,iy,iz,iv,target%it,1) + &
1024  w(1+jt-j0)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1025  end if
1026  end do
1027  target%mem(ix,iy,iz,iv,target%it,1) = &
1028  exp(target%mem(ix,iy,iz,iv,target%it,1))
1029  if (self%check_filled) then
1030  filled(ix,iy,iz,iv) = source%level
1031  if (verbose > 0) &
1032  src(ix,iy,iz) = source%id
1033  end if
1034  end if
1035  end do
1036  else
1037  do jt=j0,j1
1038  if (jt==j0) then
1039  do ix=l(1),u(1)
1040  target%mem(ix,iy,iz,iv,target%it,1) = &
1041  w(1+jt-j0)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1042  end do
1043  else
1044  do ix=l(1),u(1)
1045  target%mem(ix,iy,iz,iv,target%it,1) = &
1046  target%mem(ix,iy,iz,iv,target%it,1) + &
1047  w(1+jt-j0)*log(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1048  end do
1049  end if
1050  if (self%check_filled) then
1051  do ix=l(1),u(1)
1052  target%mem(ix,iy,iz,iv,target%it,1) = &
1053  exp(target%mem(ix,iy,iz,iv,target%it,1))
1054  filled(ix,iy,iz,iv) = source%level
1055  end do
1056  if (verbose > 0) then
1057  do ix=l(1),u(1)
1058  src(ix,iy,iz) = source%id
1059  end do
1060  end if
1061  else
1062  do ix=l(1),u(1)
1063  target%mem(ix,iy,iz,iv,target%it,1) = &
1064  exp(target%mem(ix,iy,iz,iv,target%it,1))
1065  end do
1066  end if
1067  end do
1068  end if
1069  end do
1070  end do
1071  !---------------------------------------------------------------------------
1072  ! Use linear interpolation
1073  !---------------------------------------------------------------------------
1074  else
1075  do iz=l(3),u(3)
1076  do iy=l(2),u(2)
1077  if (filter) then
1078  do ix=l(1),u(1)
1079  ii = [ix,iy,iz]
1080  if (all_cells_l .or. any((ii<li).or.(ii>ui))) then
1081  do jt=j0,j1
1082  if (jt==j0) then
1083  target%mem(ix,iy,iz,iv,target%it,1) = &
1084  w(1+jt-j0)*(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1085  else
1086  target%mem(ix,iy,iz,iv,target%it,1) = &
1087  target%mem(ix,iy,iz,iv,target%it,1) + &
1088  w(1+jt-j0)*(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1089  end if
1090  end do
1091  if (self%check_filled) then
1092  filled(ix,iy,iz,iv) = source%level
1093  if (verbose > 0) &
1094  src(ix,iy,iz) = source%id
1095  end if
1096  end if
1097  end do
1098  else
1099  do jt=j0,j1
1100  if (jt==j0) then
1101  do ix=l(1),u(1)
1102  target%mem(ix,iy,iz,iv,target%it,1) = &
1103  w(1+jt-j0)*(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1104  end do
1105  else
1106  do ix=l(1),u(1)
1107  target%mem(ix,iy,iz,iv,target%it,1) = &
1108  target%mem(ix,iy,iz,iv,target%it,1) + &
1109  w(1+jt-j0)*(source%mem(ix+o(1),iy+o(2),iz+o(3),iv,iit(jt),1))
1110  end do
1111  end if
1112  if (self%check_filled) then
1113  do ix=l(1),u(1)
1114  filled(ix,iy,iz,iv) = source%level
1115  end do
1116  if (verbose > 0) then
1117  do ix=l(1),u(1)
1118  src(ix,iy,iz) = source%id
1119  end do
1120  end if
1121  end if
1122  end do
1123  end if
1124  end do
1125  end do
1126  end if
1127  if ((verbose>2.or.target%id==id_debug).and.product(max(0,u-l+1))>0) &
1128  write(io_unit%log,'(i6,4(2x,a,i6),1p,2e12.4)') target%id, &
1129  'same picked up', product(max(0,u-l+1)), &
1130  'values of', product(target%n), &
1131  'from', source%id, &
1132  ' iv, min, max =', iv, &
1133  minval(target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,target%it,1)), &
1134  maxval(target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,target%it,1))
1135  end if
1136  end do
1137  call trace%end (itimer)
1138 END SUBROUTINE same_lagrange
1139 
1140 !===============================================================================
1141 !> Same as above, but for RT_solver type patches
1142 !===============================================================================
1143 SUBROUTINE same_rt (self, target, source, jt, pt)
1144  class(download_t):: self
1145  class(patch_t), pointer :: target, source
1146  !integer(int8), dimension(:,:,:,:), allocatable:: filled
1147  !.............................................................................
1148  integer :: l(3), u(3), o(3), ii(3), jt(2), li(3), ui(3)
1149  real(8) :: position(3)
1150  real :: pt(2)
1151  integer :: ix, iy, iz, iv
1152  integer, save :: itimer=0
1153  !-----------------------------------------------------------------------------
1154  call trace%begin('download_t%same_rt', itimer=itimer)
1155  call dnload_range (target, source, 1, l, u, o)
1156  li = target%mesh%li
1157  ui = target%mesh%ui
1158  do iv=1,target%nv
1159  do iz=l(3),u(3)
1160  do iy=l(2),u(2)
1161  do ix=l(1),u(1)
1162  ii = [ix,iy,iz]
1163  if (self%check_filled) then
1164  if (filled(ix,iy,iz,iv) > source%level) cycle
1165  end if
1166  if (all((ii>=li).and.(ii<=ui))) cycle
1167  target%mem(ix,iy,iz,iv,target%it,1) = &
1168  pt(1)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(1),1) + &
1169  pt(2)*source%mem(ix+o(1),iy+o(2),iz+o(3),iv,jt(2),1)
1170  if (self%check_filled) then
1171  filled(ix,iy,iz,iv) = source%level
1172  if (verbose > 0) &
1173  src(ix,iy,iz) = source%id
1174  end if
1175  end do
1176  end do
1177  end do
1178  end do
1179  call trace%end (itimer)
1180 END SUBROUTINE same_rt
1181 
1182 !===============================================================================
1183 !> Interpolate guard zone info from source patch at a different level
1184 !===============================================================================
1185 SUBROUTINE different (self, target, source, jt, pt, all_cells, only)
1186  class(download_t):: self
1187  class(patch_t), pointer :: target, source
1188  logical, optional:: all_cells
1189  integer, optional:: only
1190  !.............................................................................
1191  integer :: ix, iy, iz, iv, jt(2), iv1, iv2, i, cells
1192  integer, dimension(3) :: l, u, ii, li, ui, jj
1193  real(8), dimension(3) :: pos
1194  logical :: valid
1195  real :: pt(2), pp(3)
1196  real(kind=KindScalarVar):: qint
1197  logical :: all_cells_l
1198  real(8) :: start
1199  integer, save :: itimer=0
1200  !-----------------------------------------------------------------------------
1201  call trace%begin('download_t%different', itimer=itimer)
1202  if (detailed_timer) then
1203  cells = 0
1204  start = wallclock()
1205  end if
1206  if (present(only)) then
1207  iv1 = only
1208  iv2 = only
1209  else
1210  iv1 = 1
1211  iv2 = target%nv
1212  end if
1213  if (present(all_cells)) then
1214  all_cells_l = all_cells
1215  else
1216  all_cells_l = .false.
1217  end if
1218  do iv=iv1,iv2
1219  call dnload_range (target, source, iv, l, u)
1220  li = target%mesh%li
1221  ui = target%mesh%ui
1222  all_cells_l = all_cells_l .or. all((l >= li) .and. (u<=ui))
1223  !---------------------------------------------------------------------------
1224  ! For no-no-mans-land meshes, use guard zone values for the li index, to
1225  ! make fluxes up/down symmetric
1226  !---------------------------------------------------------------------------
1227 ! do i=1,3
1228 ! if (target%mesh(i)%h(iv) /= 0.0) li(i) = li(i)+1
1229 ! end do
1230  if (verbose > 1 .or. target%id==id_debug) &
1231  write (io_unit%log,1) &
1232  source%id, (source%position-target%position)/target%size, l, u, all_cells_l
1233  1 format("DBG download_t%different:", i6," displ:",3f11.6," l,u:",2(3i4,2x),l3)
1234  associate(q1 => source%mem(:,:,:,:,jt(1),1), q2 => source%mem(:,:,:,:,jt(2),1))
1235  do iz=l(3),u(3)
1236  do iy=l(2),u(2)
1237  do ix=l(1),u(1)
1238  ii=[ix,iy,iz]
1239  if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0)) then
1240  if (iv==1) cells = cells+1
1241  if (verbose > 1 .and. target%id==id_debug) &
1242  write (io_unit%log,*) &
1243  target%id, 'DBG interpolate_t%different: using', &
1244  iv,ii, present(all_cells), (ii-li) < 0
1245  pos = target%myposition (ii, iv)
1246  call source%index_space (pos, iv, jj, pp)
1247  if (source%pervolume(iv)) then
1248  call source%interpolate_specific (target, ii, iv, jt, pt)
1249  else
1250  call source%interpolate (target, ii, iv, jt, pt)
1251  end if
1252  if (.not.all_cells_l) then
1253  if (all((ii>=li .and. ii<=ui))) then
1254  write (io_unit%log,*)'INTERNAL BREACH!', ii
1255  end if
1256  end if
1257  else
1258  if (verbose > 1 .or. target%id==id_debug) &
1259  write (io_unit%log,'(i6,i3,2x,a,i6,i2,i3,2x,3i4,L2,3L1)') &
1260  target%id, target%level, 'DBG interpolate_t%different: skipping', &
1261  source%id, source%level, iv,ii, present(all_cells), (ii-li) < 0
1262  cycle
1263  end if
1264  end do
1265  end do
1266  end do
1267  end associate
1268  if (verbose > 0 .or. target%id==id_debug) then
1269  if (iv == 1 .or. verbose > 2) &
1270  write (io_unit%log,2) target%id, target%level, source%id, source%level, &
1271  'different picked up', &
1272  product(u-l+1), ' values, l, u, mind, pt =', l, u, &
1273  minval(target%mem(l(1):u(1),l(2):u(2),l(3):u(3),1,target%it,1)), pt
1274  2 format(2(i6,i3),2x,a,i6,a,2(2x,3i4),1p,e14.4,2x,0p2f6.3)
1275  end if
1276  end do
1277  if (detailed_timer) then
1278  if (cells > 0) then
1279  write (stdout,'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus/cell")') &
1280  "download_t%different: target, level, source, level,", &
1281  target%id, target%level, source%id, source%level, cells, &
1282  (wallclock()-start)/max(1,cells)
1283  else
1284  write (stdout,'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus")') &
1285  "download_t%different: target, level, source, level,", &
1286  target%id, target%level, source%id, source%level, cells, wallclock()-start
1287  end if
1288  end if
1289  call trace%end (itimer)
1290 END SUBROUTINE different
1291 
1292 !===============================================================================
1293 !> Interpolate guard zone info from source patch at a lower level
1294 !===============================================================================
1295 SUBROUTINE different_prolong (self, target, source, jt, pt, all_cells, only)
1296  class(download_t):: self
1297  class(patch_t), pointer :: target, source
1298  logical, optional:: all_cells
1299  integer, optional:: only
1300  !.............................................................................
1301  integer :: ntot, nskip, cells, nint, nperv, nexp
1302  integer :: ix, iy, iz, iv, jt(2), iv1, iv2, i, n, srcid, m(3)
1303  integer, dimension(3) :: l, u, ii, li, ui, jj
1304  real(8), dimension(3) :: pos
1305  real(4), dimension(3) :: ot, os, ratio
1306  real(8) :: start
1307  logical :: valid
1308  real :: pt(2), pp(3)
1309  real(kind=KindScalarVar):: qint, dint
1310  logical :: all_cells_l
1311  real(kind=KindScalarVar), pointer:: q1(:,:,:), q2(:,:,:), d1(:,:,:), d2(:,:,:)
1312  real(kind=KindScalarVar), pointer:: q(:,:,:,:)
1313  real(kind=KindScalarVar), allocatable:: dtmp(:,:,:)
1314  integer, save :: itimer=0
1315  !-----------------------------------------------------------------------------
1316  call trace%begin('download_t%different_prolong', itimer=itimer)
1317  start = wallclock()
1318  if (present(only)) then
1319  iv1 = only
1320  iv2 = only
1321  else
1322  iv1 = 1
1323  iv2 = target%nv
1324  end if
1325  if (present(all_cells)) then
1326  all_cells_l = all_cells
1327  else
1328  all_cells_l = .false.
1329  end if
1330  if (verbose > 1 .or. target%id == id_debug) &
1331  write(io_unit%log,'(2(a,i7,i4,1p,3e11.3,4x),l3)') &
1332  'DBG different_prolong; target:', target%id, target%level, target%ds, &
1333  'source:', source%id, source%level, source%ds, all_cells_l
1334  li = target%mesh%li
1335  ui = target%mesh%ui
1336  ot = target%mesh%o
1337  os = target%position - source%position
1338  os = modulo(os + 0.5_8*target%mesh%b, target%mesh%b) - 0.5_8*target%mesh%b
1339  os = source%mesh%o + os/source%mesh%d
1340  ratio = 2.0
1341  ntot = 0
1342  nint = 0
1343  nexp = 0
1344  nskip = 0
1345  nperv = 0
1346  cells = 0
1347  allocate (dtmp(target%gn(1),target%gn(2),target%gn(3)))
1348  d1 => source%mem(:,:,:,source%idx%d,jt(1),1)
1349  d2 => source%mem(:,:,:,source%idx%d,jt(2),1)
1350  do iv=iv1,iv2
1351  ! -- `dphi` is part of `mem` when `time_order` = 0; do not download it.
1352  if (iv == target%idx%dphi) cycle
1353  call dnload_range (target, source, iv, l, u)
1354  !---------------------------------------------------------------------------
1355  ! For no-no-mans-land meshes, use guard zone values for the li index, to
1356  ! make fluxes up/down symmetric
1357  !---------------------------------------------------------------------------
1358 ! do i=1,3
1359 ! if (target%mesh(i)%h(iv) /= 0.0) li(i) = li(i)+1
1360 ! end do
1361  q1 => source%mem(:,:,:,iv,jt(1),1)
1362  q2 => source%mem(:,:,:,iv,jt(2),1)
1363  q => source%mem(:,:,:,iv,:,1)
1364  all_cells_l = all_cells_l .or. all((l >= li) .and. (u <= ui))
1365  if (verbose > 0 .or. target%id==id_debug) &
1366  write (io_unit%log,1) &
1367  source%id, source%distance(target)/(0.5*(target%size+source%size)), &
1368  l, u, all_cells_l
1369  1 format("download_t%different_prolong: source =", i6," displ:",3f11.6, &
1370  " l,u:",2(3i4,2x)," all_cells:",l2)
1371  n = 0
1372  do iz=l(3),u(3)
1373  do iy=l(2),u(2)
1374  do ix=l(1),u(1)
1375  ntot = ntot+1
1376  ii=[ix,iy,iz]
1377  if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0)) then
1378  if (self%check_filled) then
1379  !-----------------------------------------------------------------
1380  ! If either the cell is already filled with higher level data,
1381  ! or filled with target level data (which is preferred), cycle.
1382  !-----------------------------------------------------------------
1383  if (filled(ix,iy,iz,iv) > source%level .or. &
1384  filled(ix,iy,iz,iv) == target%level) then
1385  if (allocated(src)) then
1386  srcid = src(ix,iy,iz)
1387  else
1388  srcid = -1
1389  end if
1390  if (verbose > 2 .or. verbose > 1 .and. target%id == id_debug) &
1391  write (io_unit%log,'(a,2(i6,i3),i4,2x,3i4,2x,i4,i6)') &
1392  'different_prolong: skipping', target%id, target%level, &
1393  source%id, source%level, iv,ii, filled(ix,iy,iz,iv), srcid
1394  nskip = nskip+1
1395  cycle
1396  end if
1397  end if
1398  if (iv==1) cells = cells+1
1399  if (verbose > 1 .and. target%id == id_debug) &
1400  write (io_unit%log,*) target%id, 'different_prolong: using', &
1401  iv,ii, present(all_cells), (ii-li) < 0
1402  pos = target%myposition (ii, iv)
1403  call source%index_space (pos, iv, jj, pp)
1404  !call source%index_space_of (target, ii, iv, jj, pp)
1405  if (source%unsigned(iv)) then
1406  nexp = nexp+1
1407  qint = exp(interpolator%four_d_log (q, jj(1), jj(2), jj(3), jt, [pp,pt(2)]))
1408  if (iv==target%idx%d) &
1409  dtmp(ix,iy,iz) = qint
1410  else if (source%pervolume(iv)) then
1411  nperv = nperv+1
1412  qint = interpolator%four_d_pv (q, d1, d2, jj(1), jj(2), jj(3), jt, [pp,pt(2)])
1413  qint = qint * dtmp(ix,iy,iz)
1414  else
1415  qint = interpolator%four_d (q, jj(1), jj(2), jj(3), jt, [pp,pt(2)])
1416  end if
1417  n = n + 1
1418  !if (all((ii >= li .and. ii <= ui)) .and. .not. all_cells_l) then
1419  ! write (io_unit%log,*) 'INTERNAL BREACH!', ii
1420  !end if
1421  target%mem(ix,iy,iz,iv,target%it,1) = qint
1422  if (self%check_filled) filled(ix,iy,iz,iv) = source%level
1423  else
1424  nint = nint+1
1425  if (verbose > 1 .or. target%id==id_debug) &
1426  write (io_unit%log,'(i6,i3,2x,a,i6,i2,i3,2x,3i4,L2,2(1x,3L1))') &
1427  target%id, target%level, 'different_prolong: skipping internal', &
1428  source%id, source%level, iv, ii, all_cells_l, (ii-li) < 0, (ii-ui) > 0
1429  cycle
1430  end if
1431  end do
1432  end do
1433  end do
1434  if (verbose > 0 .or. target%id==id_debug) then
1435  if (iv == 1 .or. verbose > 2) &
1436  write (io_unit%log,2) target%id, target%level, source%id, source%level, &
1437  'diff_prol picked up', &
1438  cells, nskip, nint, ' values, l, u, mind, pt =', l, u, &
1439  minval(target%mem(l(1):u(1),l(2):u(2),l(3):u(3),1,target%it,1)), pt
1440  2 format(2(i6,i3),2x,a,3i6,a,2(2x,3i4),1p,e14.4,2x,0p2f6.3)
1441  end if
1442  end do
1443  deallocate (dtmp)
1444  if (detailed_timer) then
1445  if (cells > 0) then
1446  write (stdout,'(a,6p,2(2x,i6,i3),",",i6," cells,",f9.3," mus/cell")') &
1447  "download_t%different_prolong: target, level, source, level =", &
1448  target%id, target%level, source%id, source%level, cells, &
1449  (wallclock()-start)/cells
1450  else
1451  write (stdout,'(a,2(2x,i6,i3),2x,3f8.3,",",i6," cells,",6p,f9.3," mus")') &
1452  "download_t%different_prolong: target,lvl,source,lvl,rel=", &
1453  target%id, target%level, source%id, source%level, &
1454  source%distance(target)/(0.5d0*(target%size+source%size)), &
1455  cells, wallclock()-start
1456  end if
1457  end if
1458  call trace%end (itimer)
1459 END SUBROUTINE different_prolong
1460 
1461 !> =============================================================================
1462 !> Interpolate guard zone info from source patch at a lower level
1463 !> =============================================================================
1464 SUBROUTINE prolong (self, target, source, jt, pt, all_cells, only)
1465  class(download_t):: self
1466  class(patch_t), pointer :: target, source
1467  logical, optional:: all_cells
1468  integer, optional:: only
1469  !.............................................................................
1470  integer :: ntot, nskip, cells, nint, nperv, nexp
1471  integer :: ix, iy, iz, iv, jt(2), iv1, iv2, i, m, n, srcid
1472  integer, dimension(3) :: l, u, ii, li, ui, jj, jjmin, jjmax, iimin, iimax
1473  real(8), dimension(3) :: pos
1474  real(4), dimension(3) :: ot, os, ratio
1475  real(8) :: start
1476  logical :: valid
1477  integer, save :: itimer=0
1478  real :: pt(2), pp(3)
1479  real(kind=KindScalarVar):: qint, dint
1480  logical :: all_cells_l
1481  real(kind=KindScalarVar), allocatable:: sbuf(:,:,:,:,:), tbuf(:), pbuf(:,:)
1482  real(kind=KindScalarVar), allocatable:: sdbuf(:,:,:,:,:), tdbuf(:)
1483  integer :: mv(target%nv)
1484  !-----------------------------------------------------------------------------
1485  call trace%begin('download_t%prolong', itimer=itimer)
1486  start = wallclock()
1487  if (present(only)) then
1488  iv1 = only
1489  iv2 = only
1490  else
1491  iv1 = 1
1492  iv2 = target%nv
1493  end if
1494  if (present(all_cells)) then
1495  all_cells_l = all_cells
1496  else
1497  all_cells_l = .false.
1498  end if
1499  if (verbose > 1 .or. target%id == id_debug) &
1500  write(io_unit%log,'(2(a,i7,i4,1p,3e11.3,4x),l3)') &
1501  'DBG prolong; target:', target%id, target%level, target%ds, &
1502  'source:', source%id, source%level, source%ds, all_cells_l
1503  li = target%mesh%li
1504  ui = target%mesh%ui
1505  call dnload_range (target, source, target%idx%d, l, u)
1506  n = product(u-l+2)
1507  allocate ( sbuf(n,2,2,2,2))
1508  allocate (sdbuf(n,2,2,2,2))
1509  allocate ( pbuf(n,4 ))
1510  allocate ( tbuf(n ))
1511  allocate (tdbuf(n ))
1512  ot = target%mesh%o
1513  os = source%mesh%o
1514  ratio = 2.0
1515  ntot = 0
1516  nint = 0
1517  nexp = 0
1518  nskip = 0
1519  nperv = 0
1520  cells = 0
1521  do iv=iv1,iv2
1522  ! `dphi` is part of `mem` when `time_order` = 0; do not download it.
1523  if (iv == target%idx%dphi) cycle
1524  call dnload_range (target, source, iv, l, u)
1525  !---------------------------------------------------------------------------
1526  ! For no-no-mans-land meshes, use guard zone values for the li index, to
1527  ! make fluxes up/down symmetric
1528  !---------------------------------------------------------------------------
1529 ! do i=1,3
1530 ! if (target%mesh(i)%h(iv) /= 0.0) li(i) = li(i)+1
1531 ! end do
1532  if (verbose > 1 .or. target%id==id_debug) &
1533  write (io_unit%log,1) &
1534  source%id, (source%position-target%position)/target%size, l, u
1535  1 format("download_t%different_prolong:", i6," displ:",3f11.6," l,u:", &
1536  2(3i4,2x))
1537  all_cells_l = all_cells_l .or. all((l >= li) .and. (u <= ui))
1538  n = 0
1539  do iz=l(3),u(3)
1540  do iy=l(2),u(2)
1541  do ix=l(1),u(1)
1542  ntot = ntot+1
1543  ii=[ix,iy,iz]
1544  if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0)) then
1545  !if (self%check_filled) then
1546  ! !-----------------------------------------------------------------
1547  ! ! If either the cell is already filled with higher level data,
1548  ! ! or filled with target level data (which is preferred), cycle.
1549  ! !-----------------------------------------------------------------
1550  ! if (filled(ix,iy,iz,iv) > source%level .or. &
1551  ! filled(ix,iy,iz,iv) == target%level) then
1552  ! if (allocated(src)) then
1553  ! srcid = src(ix,iy,iz)
1554  ! else
1555  ! srcid = -1
1556  ! end if
1557  ! nskip = nskip+1
1558  ! cycle
1559  ! end if
1560  !end if
1561  if (iv==1) cells = cells+1
1562  pos = target%myposition (ii, iv)
1563  call source%index_space (pos, iv, jj, pp)
1564  !pp = (ii-ot)*ratio + os
1565  !jj = pp
1566  !pp = pp - jj
1567  n = n + 1
1568  !if (n==1) then
1569  ! iimin = ii
1570  ! iimax = ii
1571  ! jjmin = jj
1572  ! jjmax = jj
1573  !else
1574  ! iimin = min(ii,iimin)
1575  ! iimax = max(ii,iimax)
1576  ! jjmin = min(jj,jjmin)
1577  ! jjmax = max(jj,jjmax)
1578  !end if
1579  pbuf(n,1:3) = pp
1580  pbuf(n, 4) = pt(1)
1581  sbuf(n,:,:,:,:) = &
1582  source%mem(jj(1):jj(1)+1,jj(2):jj(2)+1,jj(3):jj(3)+1,iv,jt(:),1)
1583  !if (self%check_filled) filled(ix,iy,iz,iv) = source%level
1584  !if (iv==1) &
1585  !write (io_unit%output,'(a,i4,2(3i4,2x),2x,4f7.3)') &
1586  !'prol', n, ii, jj, pbuf(n,:)
1587  else
1588  nint = nint+1
1589  !if (verbose > 1 .and. target%id==id_debug) &
1590  !write (io_unit%output,'(i6,i3,2x,a,i6,i2,i3,2x,3i4,L2,2(1x,3L1))') &
1591  ! target%id, target%level, 'prolong: skipping internal', &
1592  ! source%id, source%level, iv, ii, all_cells_l, (ii-li) < 0,
1593  ! (ii-ui) > 0
1594  cycle
1595  end if
1596  end do
1597  end do
1598  end do
1599  m = n
1600  mv(iv) = m
1601  call mem_interp
1602  n = 0
1603  do iz=l(3),u(3)
1604  do iy=l(2),u(2)
1605  do ix=l(1),u(1)
1606  ntot = ntot+1
1607  ii=[ix,iy,iz]
1608  if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0)) then
1609  n = n + 1
1610  target%mem(ix,iy,iz,iv,target%it,1) = tbuf(n)
1611  end if
1612  end do
1613  end do
1614  end do
1615  end do
1616  !print '(a,2(2x,3i4))', 'prolong: iimin, iimax =', iimin, iimax
1617  !print '(a,2(2x,3i4))', 'prolong: jjmin, jjmax =', jjmin, jjmax
1618  deallocate (sbuf, pbuf, tbuf, sdbuf, tdbuf)
1619  if (detailed_timer) then
1620  if (cells > 0) then
1621  write (stdout,'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus/cell")') &
1622  "download_t%prolong: target, level, source, level,", &
1623  target%id, target%level, source%id, source%level, cells, &
1624  (wallclock()-start)/max(1,cells)
1625  else
1626  write (stdout,'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus")') &
1627  "download_t%prolong: target, level, source, level,", &
1628  target%id, target%level, source%id, source%level, cells, wallclock()-start
1629  end if
1630  end if
1631  call trace%end (itimer)
1632 contains
1633  subroutine mem_interp
1634  real:: q1, q2, q3, q4
1635  integer:: j
1636  !print *, iv, m, target%pervolume(iv), target%unsigned(iv)
1637  !---------------------------------------------------------------------------
1638  ! Save density corner values
1639  !---------------------------------------------------------------------------
1640  if (iv == 1) then
1641  sdbuf = sbuf
1642  end if
1643  !---------------------------------------------------------------------------
1644  ! Use density corner values for per-volume variables
1645  !---------------------------------------------------------------------------
1646  if (target%pervolume(iv)) then
1647  sbuf(1:m,:,:,:,:) = sbuf(1:m,:,:,:,:)/sdbuf(1:m,:,:,:,:)
1648  end if
1649  !---------------------------------------------------------------------------
1650  ! Take log of unsigned variables
1651  !---------------------------------------------------------------------------
1652  if (target%unsigned(iv)) then
1653  sbuf(1:m,:,:,:,:) = log(sbuf(1:m,:,:,:,:))
1654  end if
1655  !---------------------------------------------------------------------------
1656  ! Interpolate, vectorized
1657  !---------------------------------------------------------------------------
1658  do j=1,m
1659  associate(p1=>pbuf(j,1), p2=>pbuf(j,2), p3=>pbuf(j,3), p4=>pbuf(j,4))
1660  q1 = 1.0-p1
1661  q2 = 1.0-p2
1662  q3 = 1.0-p3
1663  q4 = 1.0-p4
1664  tbuf(j) = q4*(q3*(q2*(q1*sbuf(j,1,1,1,1)+p1*sbuf(j,2,1,1,1)) + &
1665  p2*(q1*sbuf(j,1,2,1,1)+p1*sbuf(j,2,2,1,1))) + &
1666  p3*(q2*(q1*sbuf(j,1,1,2,1)+p1*sbuf(j,2,1,2,1)) + &
1667  p2*(q1*sbuf(j,1,2,2,1)+p1*sbuf(j,2,2,2,1)))) + &
1668  p4*(q3*(q2*(q1*sbuf(j,1,1,1,2)+p1*sbuf(j,2,1,1,2)) + &
1669  p2*(q1*sbuf(j,1,2,1,2)+p1*sbuf(j,2,2,1,2))) + &
1670  p3*(q2*(q1*sbuf(j,1,1,2,2)+p1*sbuf(j,2,1,2,2)) + &
1671  p2*(q1*sbuf(j,1,2,2,2)+p1*sbuf(j,2,2,2,2))))
1672  end associate
1673  end do
1674  !---------------------------------------------------------------------------
1675  ! Take exp of unsigned variables
1676  !---------------------------------------------------------------------------
1677  if (target%unsigned(iv)) then
1678  tbuf(1:m) = exp(tbuf(1:m))
1679  end if
1680  !---------------------------------------------------------------------------
1681  ! Save density values
1682  !---------------------------------------------------------------------------
1683  if (iv == 1) then
1684  tdbuf = tbuf
1685  end if
1686  !---------------------------------------------------------------------------
1687  ! Use density values for per-volume variables
1688  !---------------------------------------------------------------------------
1689  if (target%pervolume(iv)) then
1690  tbuf(1:m) = tbuf(1:m)*tdbuf(1:m)
1691  end if
1692  end subroutine mem_interp
1693 END SUBROUTINE prolong
1694 
1695 !===============================================================================
1696 !> Same as above, but for RT_solver type patches
1697 !===============================================================================
1698 SUBROUTINE different_rt (self, target, source, jt, pt)
1699  class(download_t):: self
1700  class(patch_t), pointer :: target, source
1701  !.............................................................................
1702  integer :: ix, iy, iz, iv, jt(2)
1703  integer, dimension(3) :: l, u, ii, lb, ub, li, ui
1704  real(8), dimension(3) :: pos
1705  logical :: valid
1706  integer, save :: itimer=0
1707  real :: pt(2)
1708  !-----------------------------------------------------------------------------
1709  call trace%begin('download_t%different_rt', itimer=itimer)
1710  lb = target%mesh%lb
1711  ub = target%mesh%ub
1712  call dnload_range (target, source, 1, l, u)
1713  li = target%mesh%li
1714  ui = target%mesh%ui
1715  do iv=1,target%nv
1716  do iz=l(3),u(3)
1717  do iy=l(2),u(2)
1718  do ix=l(1),u(1)
1719  ii=[ix,iy,iz]
1720  if (self%check_filled) then
1721  if (filled(ix,iy,iz,iv) > source%level) cycle
1722  end if
1723  if (all((ii>=li).and.(ii<=ui))) cycle
1724  call source%interpolate (target, ii, iv, jt, pt)
1725  if (self%check_filled) filled(ix,iy,iz,iv) = source%level
1726  end do
1727  end do
1728  end do
1729  end do
1730  call trace%end (itimer)
1731 END SUBROUTINE different_rt
1732 
1733 !===============================================================================
1734 !> Restrict guard zone info from source patch at a higher level
1735 !===============================================================================
1736 SUBROUTINE different_restrict (self, target, source, jt, pt, only)
1737  class(download_t):: self
1738  class(patch_t), pointer :: target, source
1739  !integer(int8), dimension(:,:,:,:), allocatable:: filled
1740  integer, optional:: only
1741  integer :: ix, iy, iz, iv, jt(2), iv1, iv2, i, jx, jy, jz
1742  integer, dimension(3) :: l, u, ii, li, ui, jj
1743  real(8), dimension(3) :: target_pos, source_pos, target_upper_edges, &
1744  target_lower_edges, source_upper_edges, &
1745  source_lower_edges, le, ue
1746  logical :: valid
1747  integer, save :: itimer=0
1748  real :: pt(2), pp(3)
1749  real(kind=KindScalarVar):: qint
1750  logical :: qmask(3)
1751  real(8) :: start
1752  integer :: cells
1753  real(8):: rratio(3), dvsource, dvtarget, dwsource, dwtarget, filling_factor, &
1754  accumulate, filling_total, dwtargeti
1755  integer:: ilook(3)
1756  !-----------------------------------------------------------------------------
1757  call trace%begin('download_t%different_restrict', itimer=itimer)
1758  if (detailed_timer) then
1759  cells = 0
1760  start = wallclock()
1761  end if
1762  rratio = target%ds / source%ds
1763  ilook = merge(floor(rratio)/2 + 1, 0, source%n>1)
1764  dvsource = product(source%ds) ! only valid for Cartesian
1765  dvtarget = product(target%ds) ! onlt valid for Cartesian
1766  if (present(only)) then
1767  iv1 = only
1768  iv2 = only
1769  else
1770  iv1 = 1
1771  iv2 = target%nv
1772  end if
1773  li = target%mesh%li
1774  ui = target%mesh%ui
1775  do iv=iv1,iv2
1776  ! -- `dphi` is part of `mem` when `time_order` = 0; do not download it.
1777  if (iv == target%idx%dphi) cycle
1778  ! -- logical mask to track if variable is averaged over a volume or an area
1779  qmask(:) = .true.
1780  if (iv == target%idx%bx) then
1781  dwsource = source%ds(2) * source%ds(3) ! only valid for Cartesian
1782  dwtarget = target%ds(2) * target%ds(3) ! only valid for Cartesian
1783  qmask(1) = .false.
1784  else if (iv == target%idx%by) then
1785  dwsource = source%ds(3) * source%ds(1) ! only valid for Cartesian
1786  dwtarget = target%ds(3) * target%ds(1) ! only valid for Cartesian
1787  qmask(2) = .false.
1788  else if (iv == target%idx%bz) then
1789  dwsource = source%ds(1) * source%ds(2) ! only valid for Cartesian
1790  dwtarget = target%ds(1) * target%ds(2) ! only valid for Cartesian
1791  qmask(3) = .false.
1792  else
1793  dwsource = dvsource
1794  dwtarget = dvtarget
1795  end if
1796  dwtargeti = 1.0 / dwtarget
1797  call dnload_range (target, source, iv, l, u)
1798  associate(q1 => source%mem(:,:,:,iv,jt(1),1), &
1799  q2 => source%mem(:,:,:,iv,jt(2),1))
1800  do iz=l(3),u(3)
1801  do iy=l(2),u(2)
1802  do ix=l(1),u(1)
1803  ii = [ix,iy,iz]
1804  if (iv==1) cells = cells+1
1805  ! prefer same level for guard zone filling
1806  if (self%check_filled) then
1807  if (filled(ix,iy,iz,iv) == target%level .and. &
1808  (any(ii < target%li) .or. any(ii > target%ui))) cycle
1809  end if
1810  target_pos = target%myposition (ii, iv)
1811  target_lower_edges = target_pos - target%ds * 0.499999_8
1812  target_upper_edges = target_pos + target%ds * 0.499999_8
1813  ! skip cells that are not *entirely* covered by the finer patch
1814  if (any(target_lower_edges < source%llc_nat)) cycle
1815  if (any(target_upper_edges > (source%llc_nat+source%size))) cycle
1816  call source%index_space (target_pos, iv, jj, pp)
1817  accumulate = 0.0
1818  filling_total = 0.0
1819  do jz=jj(3)-ilook(3),jj(3)+ilook(3)
1820  do jy=jj(2)-ilook(2),jj(2)+ilook(2)
1821  do jx=jj(1)-ilook(1),jj(1)+ilook(1)
1822  ! skip ghost cells of source
1823  if (any([jx,jy,jz] < source%li .or. [jx,jy,jz] > source%ui)) cycle
1824  ! find intersection between source and target cells
1825  source_pos = source%myposition([jx,jy,jz],iv)
1826  source_lower_edges = source_pos - source%ds * 0.499999_8
1827  source_upper_edges = source_pos + source%ds * 0.499999_8
1828  ue = min(source_upper_edges, target_upper_edges)
1829  le = max(source_lower_edges, target_lower_edges)
1830  ! filling factor; if zero, then no overlap.
1831  ! the logical mask is used to produce an area filling factor for
1832  ! mag. fields (but a volume filling factor otherwise)
1833  filling_factor = product(max(ue - le, 0.0),mask=qmask)
1834  filling_total = filling_total + filling_factor
1835  ! `accumulate` should be multipled by `dWsource`, but
1836  ! `filling_factor` should be divided by `dWsource`, so they cancel.
1837  if (filling_factor > 0.0) then
1838  accumulate = accumulate &
1839  + (q1(jx,jy,jz) * pt(1) + q2(jx,jy,jz) * pt(2)) * &
1840  filling_factor
1841  end if
1842  end do
1843  end do
1844  end do
1845  ! replace a portion of the current cell with values from the finer patch
1846  if (filling_total > 0.0) then
1847  target%mem(ix,iy,iz,iv,target%it,1) = accumulate * dwtargeti &
1848  + target%mem(ix,iy,iz,iv,target%it,1) &
1849  * (1.0 - filling_total * dwtargeti)
1850  if (self%check_filled) then
1851  filled(ix,iy,iz,iv) = source%level
1852  end if
1853  end if
1854  end do
1855  end do
1856  end do
1857  end associate
1858  if (verbose > 0 .or. target%id==id_debug) then
1859  if (iv == 1 .or. verbose > 2) &
1860  write (io_unit%log,2) target%id, target%level, source%id, source%level, &
1861  'diff_rest picked up', &
1862  product(u-l+1), ' values, l, u, mind, pt =', l, u, &
1863  minval(target%mem(l(1):u(1),l(2):u(2),l(3):u(3),1,target%it,1)), pt
1864  2 format(2(i6,i3),2x,a,i6,a,2(2x,3i4),1p,e14.4,2x,0p2f6.3)
1865  end if
1866  end do
1867  if (detailed_timer) then
1868  write (stdout,'(a,6p,f6.3," mus/cell")') &
1869  "download_t%different_restrict:", (wallclock()-start)/cells
1870  end if
1871  call trace%end (itimer)
1872 END SUBROUTINE different_restrict
1873 
1874 !===============================================================================
1875 !> Restrict guard zone info from source patch at a higher level
1876 !===============================================================================
1877 SUBROUTINE restrict_centered (self, target, source, jt, pt, only)
1878  class(download_t):: self
1879  class(patch_t), pointer :: target, source
1880  !integer(int8), dimension(:,:,:,:), allocatable:: filled
1881  integer, optional:: only
1882  integer :: ix, iy, iz, iv, jt(2), iv1, iv2, i, jx, jy, jz
1883  integer, dimension(3) :: l, u, ii, li, ui, jj
1884  real(8), dimension(3) :: target_pos, source_pos, target_upper_edges, &
1885  target_lower_edges, source_upper_edges, &
1886  source_lower_edges, le, ue
1887  logical :: valid
1888  integer, save :: itimer=0
1889  real :: pt(2), pp(3)
1890  real(kind=KindScalarVar):: qint
1891  logical :: qmask(3)
1892  real(8):: rratio(3), dvsource, dvtarget, dwsource, dwtarget, filling_factor, &
1893  accumulate(target%nv), filling_total, dwtargeti
1894  integer:: ilook(3)
1895  !-----------------------------------------------------------------------------
1896  call trace%begin('download_t%different_restrict', itimer=itimer)
1897  rratio = target%ds / source%ds
1898  ilook = merge(floor(rratio)/2 + 1, 0, source%n>1)
1899  dvsource = product(source%ds) ! only valid for Cartesian
1900  dvtarget = product(target%ds) ! onlt valid for Cartesian
1901  if (present(only)) then
1902  iv1 = only
1903  iv2 = only
1904  else
1905  iv1 = 1
1906  iv2 = target%nv
1907  end if
1908  li = target%mesh%li
1909  ui = target%mesh%ui
1910  dwsource = dvsource
1911  dwtarget = dvtarget
1912  dwtargeti = 1.0 / dwtarget
1913  call dnload_range (target, source, 1, l, u)
1914  associate(q1 => source%mem(:,:,:,:,jt(1),1), q2 => source%mem(:,:,:,:,jt(2),1))
1915  do iz=l(3),u(3)
1916  do iy=l(2),u(2)
1917  do ix=l(1),u(1)
1918  ii = [ix,iy,iz]
1919  ! prefer same level for guard zone filling
1920  if (self%check_filled) then
1921  if (filled(ix,iy,iz,iv) == target%level .and. &
1922  (any(ii < target%li) .or. any(ii > target%ui))) cycle
1923  end if
1924  target_pos = target%myposition (ii, iv)
1925  target_lower_edges = target_pos - target%ds * 0.499999_8
1926  target_upper_edges = target_pos + target%ds * 0.499999_8
1927  ! skip cells that are not *entirely* covered by the finer patch
1928  if (any(target_lower_edges < source%llc_nat)) cycle
1929  if (any(target_upper_edges > (source%llc_nat+source%size))) cycle
1930  call source%index_space (target_pos, iv, jj, pp)
1931  accumulate = 0.0
1932  filling_total = 0.0
1933  do jz=jj(3)-ilook(3),jj(3)+ilook(3)
1934  do jy=jj(2)-ilook(2),jj(2)+ilook(2)
1935  do jx=jj(1)-ilook(1),jj(1)+ilook(1)
1936  ! skip ghost cells of source
1937  if (any([jx,jy,jz] < source%li .or. [jx,jy,jz] > source%ui)) cycle
1938  ! find intersection between source and target cells
1939  source_pos = source%myposition([jx,jy,jz],iv)
1940  source_lower_edges = source_pos - source%ds * 0.499999_8
1941  source_upper_edges = source_pos + source%ds * 0.499999_8
1942  ue = min(source_upper_edges, target_upper_edges)
1943  le = max(source_lower_edges, target_lower_edges)
1944  ! filling factor; if zero, then no overlap.
1945  ! the logical mask is used to produce an area filling factor for
1946  ! mag. fields (but a volume filling factor otherwise)
1947  filling_factor = product(max(ue - le, 0.0),mask=qmask)
1948  filling_total = filling_total + filling_factor
1949  ! `accumulate` should be multipled by `dWsource`. but
1950  ! `filling_factor` should be divided by `dWsource`, so they cancel.
1951  if (filling_factor > 0.0) then
1952  do iv=iv1,iv2
1953  accumulate(iv) = accumulate(iv) &
1954  + (q1(jx,jy,jz,iv) * pt(1) + q2(jx,jy,jz,iv) &
1955  * pt(2)) * filling_factor
1956  end do
1957  end if
1958  end do
1959  end do
1960  end do
1961  ! replace a portion of the current cell with values from the finer patch
1962  if (filling_total > 0.0) then
1963  do iv=iv1,iv2
1964  target%mem(ix,iy,iz,iv,target%it,1) = accumulate(iv) * dwtargeti &
1965  + target%mem(ix,iy,iz,iv,target%it,1) &
1966  * (1.0 - filling_total * dwtargeti)
1967  end do
1968  if (self%check_filled) then
1969  filled(ix,iy,iz,iv) = source%level
1970  end if
1971  end if
1972  end do
1973  end do
1974  end do
1975  end associate
1976  call trace%end (itimer)
1977 END SUBROUTINE restrict_centered
1978 
1979 !===============================================================================
1980 SUBROUTINE set_verbose (new)
1981  integer:: new
1982  verbose = new
1983 END SUBROUTINE set_verbose
1984 
1985 !===============================================================================
1986 !> Initialise downloading and interpolating options.
1987 !===============================================================================
1988 SUBROUTINE test (self)
1989  class(download_t):: self
1990  !.............................................................................
1991  class(patch_t), pointer :: source(:), target, src
1992  integer, dimension(3) :: l, u, o, ii
1993  integer :: i, j, k, iv, dir, id_start, id_end, jt(2)
1994  real :: pt(2)
1995  real(8), dimension(3) :: pos
1996  real(4), dimension(3) :: pp
1997  class(mesh_t), pointer :: mesh
1998  logical :: ok, allok
1999  !-----------------------------------------------------------------------------
2000  return
2001  call trace%begin ('download_t%test')
2002  id_start = mpi%id%update(0)
2003  self%interpolate => selectinterpolator(-1)
2004  !-----------------------------------------------------------------------------
2005  ! Generate meshes with three resolutions, differing with a factor of three,
2006  ! with and w/o staggering
2007  !-----------------------------------------------------------------------------
2008  allocate (target)
2009  call target%setup (size=0.1_8*[1,1,1], n=[16,16,16], nv=2, nt=2, nw=1)
2010  filled = -1
2011  !-----------------------------------------------------------------------------
2012  ! Loop over the three resolutions for source patches
2013  !-----------------------------------------------------------------------------
2014  allocate (source(3))
2015  do i=1,3
2016  call source(i)%setup (size=0.1_8*[1,1,1]*(shared%f_zoom)**(i-2), &
2017  n=[16,16,16], nv=2, nt=2, nw=1)
2018  end do
2019  !---------------------------------------------------------------------------
2020  ! Loop over source size
2021  !---------------------------------------------------------------------------
2022  target%level = 2
2023  do i=1,3
2024  source(i)%level = i
2025  !---------------------------------------------------------------------------
2026  ! Loop over no-stagger / stagger
2027  !---------------------------------------------------------------------------
2028  do iv=1,2
2029  do dir=1,3
2030  mesh => source(i)%mesh(dir)
2031  mesh%h(iv) = -0.5*(iv-1)
2032  mesh => target%mesh(dir)
2033  mesh%h(iv) = -0.5*(iv-1)
2034  end do
2035  if (io%verbose>0 .or. io_unit%do_validate) then
2036  mesh => source(i)%mesh(1)
2037  write (io_unit%log,*) 'source size, stagger:', source(i)%size(1), mesh%h(iv)
2038  end if
2039  !-------------------------------------------------------------------------
2040  ! Loop over placement of the source mesh: ol, il, c, ir, or
2041  !-------------------------------------------------------------------------
2042  src => source(i)
2043  do j=1,5
2044  call set_position (src, target, j)
2045  if (i==2.and.(j==2.or.j==4)) cycle
2046  if (src%level==target%level) then
2047  call dnload_range (target, src, iv, l, u, o)
2048  else
2049  call dnload_range (target, src, iv, l, u)
2050  end if
2051  do k=l(1),u(1)
2052  do dir=1,3
2053  mesh => target%mesh(dir)
2054  pos(dir) = target%position(dir) + mesh%r(k) + mesh%h(iv)*target%ds(dir)
2055  end do
2056  call source(i)%index_space (pos, iv, ii, pp)
2057  if (io%verbose>0 .or. io_unit%do_validate) then
2058  if (src%level==target%level) then
2059  print 1, i, iv, j, "k:", k, "o:", o(1), "ii:", ii(1), "pp:", pp(1)
2060  1 format("download_t%test",3i3,3(2x,a,i4),2x,a,f5.2)
2061  else
2062  print 2, i, iv, j, "k:", k, "o: N/A", "ii:", ii(1), "pp:", pp(1)
2063  2 format("download_t%test",3i3,2x,a,i4,2x,a,2x,a,i4,2x,a,f5.2)
2064  end if
2065  end if
2066  end do
2067  if (io%verbose>0 .or. io_unit%do_validate) &
2068  write (io_unit%log,*)' '
2069  end do
2070  end do
2071  end do
2072  !-----------------------------------------------------------------------------
2073  ! Test actual download values
2074  !-----------------------------------------------------------------------------
2075  allok = .true.
2076  do i=1,3
2077  src => source(i)
2078  do j=1,5
2079  call set_position (src, target, j)
2080  if (i==2.and.(j==2.or.j==4)) cycle
2081  if (io%verbose>1) &
2082  print "('source size:',i2,3x,'position:',i2,2x,'levels:',2i2)", i, j, &
2083  src%level, target%level
2084  call set_values (src)
2085  call set_values (target)
2086  if (src%level==target%level) then
2087  call src%time_interval (target%time, jt, pt)
2088  call self%same (target, src, jt, pt, all_cells=.false.)
2089  if (self%check_filled) then
2090  if (any(filled < 0)) then
2091  write (io_unit%log,*) target%id, &
2092  'WARNING: some cells not filled by download_t%same'
2093  end if
2094  end if
2095  else
2096  !call different_prolong (self, target, src, all_cells=.false.)
2097  call src%time_interval (target%time, jt, pt)
2098  call prolong (self, target, src, jt, pt, all_cells=.false.)
2099  end if
2100  call check_values (target, ok)
2101  allok = allok.and.ok
2102  end do
2103  end do
2104  if (io%master) write (io_unit%log,*)io%hl
2105  if (allok) then
2106  if (io%master) then
2107  write (io_unit%log,*)'download_t%test passed'
2108  write (io_unit%log,*)io%hl
2109  end if
2110  else
2111  write (io_unit%log,*)mpi%rank,'download_t%test failed'
2112  if (io%master) write (io_unit%log,*)io%hl
2113  call mpi%abort('download_t%test::check_values')
2114  end if
2115  !-----------------------------------------------------------------------------
2116  ! Clean up
2117  !-----------------------------------------------------------------------------
2118  deallocate (filled)
2119  call target%dealloc
2120  do i=1,3
2121  call source(i)%dealloc
2122  end do
2123  !-----------------------------------------------------------------------------
2124  ! Set the MPI-global id counter back to what it was
2125  !-----------------------------------------------------------------------------
2126  id_end = mpi%id%update(0)
2127  id_start = mpi%id%update (id_start-id_end)
2128  deallocate (target, source)
2129  self%interpolate => selectinterpolator(self%order_interpolator)
2130  call trace%end()
2131 contains
2132 !===============================================================================
2133 !> Set position systematically
2134 !===============================================================================
2135  subroutine set_position (source, target, j)
2136  type(patch_t) :: source, target
2137  integer :: j
2138  !---------------------------------------------------------------------------
2139  select case (j)
2140  case (1)
2141  source%position = target%position - 0.5_8*(source%size+target%size) + &
2142  2.*abs(source%ds-target%ds)
2143  case (2)
2144  source%position = target%position - 0.5_8*(source%size+target%size) + &
2145  2.*abs(source%ds-target%ds)
2146  case (3)
2147  source%position = target%position
2148  case (4)
2149  source%position = target%position + 0.5_8*(source%size+target%size) - &
2150  2.*abs(source%ds-target%ds)
2151  case (5)
2152  source%position = target%position + 0.5_8*(source%size+target%size) - &
2153  2.*abs(source%ds-target%ds)
2154  end select
2155  source%mesh%p = source%position
2156  source%t = 0d0; source%time=0d0; source%iit=1
2157  target%t = 0d0; target%time=0d0; target%iit=1
2158 end subroutine set_position
2159 !===============================================================================
2160 !> Set some values to recover
2161 !===============================================================================
2162 subroutine set_values (patch)
2163  type(patch_t) :: patch
2164  class(mesh_t), pointer :: m1, m2, m3
2165  integer :: ix, iy, iz, iv
2166  real :: x, y, z, v
2167  !-----------------------------------------------------------------------------
2168  m1 => patch%mesh(1)
2169  m2 => patch%mesh(2)
2170  m3 => patch%mesh(3)
2171  do iv=1,2
2172  patch%mem(:,:,:,iv,:,:) = 0.0
2173  do iz=m3%li,m3%ui
2174  z = m3%p + m3%r(iz) + m3%h(iv)*m3%d
2175  do iy=m2%li,m2%ui
2176  y = m2%p + m2%r(iy) + m2%h(iv)*m2%d
2177  do ix=m1%li,m1%ui
2178  x = m1%p + m1%r(ix) + m1%h(iv)*m1%d
2179  v = x + y*2.0 + z*3.0
2180  patch%mem(ix,iy,iz,iv,:,:) = v
2181  end do
2182  end do
2183  end do
2184  end do
2185 end subroutine set_values
2186 !===============================================================================
2187 !> Check values, including in guard zones
2188 !===============================================================================
2189 subroutine check_values (patch, ok)
2190  type(patch_t) :: patch
2191  class(mesh_t), pointer :: m1, m2, m3
2192  logical :: ok
2193  integer :: ix, iy, iz, iv, n(3)
2194  real :: x, y, z, v, w, eps
2195  !-----------------------------------------------------------------------------
2196  m1 => patch%mesh(1)
2197  m2 => patch%mesh(2)
2198  m3 => patch%mesh(3)
2199  eps=0.1*max(m1%d,m2%d,m3%d)
2200  n = 0
2201  do iv=1,2
2202  do iz=m3%lb,m3%ub
2203  z = m3%p + m3%r(iz) + m3%h(iv)*m3%d
2204  do iy=m2%lb,m2%ub
2205  y = m2%p + m2%r(iy) + m2%h(iv)*m2%d
2206  do ix=m1%lb,m1%ub
2207  x = m1%p + m1%r(ix) + m1%h(iv)*m1%d
2208  v = x + y*2.0 + z*3.0
2209  w = patch%mem(ix,iy,iz,iv,1,1)
2210  if (w == 0.0) then
2211  n(1) = n(1)+1
2212  else
2213  if (abs(w-v) > eps) then
2214  print "('error: ix, iy, iz, iv, value:',4i4,1p,2e15.5)", &
2215  ix, iy, iz, iv, w, v
2216  n(3) = n(3)+1
2217  else
2218  n(2) = n(2)+1
2219  end if
2220  end if
2221  end do
2222  end do
2223  end do
2224  end do
2225  if (io%verbose>0) &
2226  write (io_unit%log,*)'check_values: n_zero, n_ok, n_nok =', n
2227  if (n(3) == 0) then
2228  ok = .true.
2229  else
2230  ok = .false.
2231  write (io_unit%log,*)'check_values: n_zero, n_ok, n_nok =', n
2232  end if
2233 end subroutine check_values
2234 !===============================================================================
2235 END SUBROUTINE test
2236 
2237 END MODULE
Each thread uses a private timer data type, with arrays for start time and total time for each regist...
Definition: timer_mod.f90:11
download_link: takes care of downloads to linktask same: called for patches on the same level differ:...
Support tic/toc timing, as in MATLAB, and accurate wallclock() function. The timing is generally much...
Generic validation module. The general idea is to be able to compare two runs at critical points in t...
This module handles remeshing, eg from lower resolution to higher resolution, with special attention ...
Definition: remesh_mod.f90:6
Template module for patches, which adds pointers to memory and mesh, and number of dimensions and var...
Definition: patch_mod.f90:6
This module is a placeholder for shared information, solving the problem of making a piece of informa...
Definition: shared_mod.f90:7
Module for Lagrange interpolation.
Definition: lagrange_mod.f90:4
Template module for mesh.
Definition: mesh_mod.f90:4
Definition: io_mod.f90:4
Optimized restrict operations, intended for guard zones (not conservative)
The lock module uses nested locks, to allow versatile use of locks, where a procedure may want to mak...
Template module for tasks.
Definition: task_mod.f90:4