17 USE iso_fortran_env
, only: int8
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()
49 procedure:: download_link
50 procedure,
nopass:: dnload_range
53 procedure:: different_prolong
55 procedure,
nopass:: set_verbose
57 integer(kind=int8),
dimension(:,:,:,:),
allocatable:: filled
58 integer(kind=4),
dimension(:,:,:),
allocatable:: src
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=
'' 73 SUBROUTINE init (self)
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' 85 call trace%begin (
'download_t%init')
86 call trace%print_id (id)
91 read (io%input, download_params, iostat=iostat)
92 if (io%master)
write (*, download_params)
97 if (order_time==1)
then 98 self%same => same_linear
100 self%same => same_lagrange
102 self%check_filled = check_filled
103 self%order_interpolator = order_interpolator
104 self%interpolate => selectinterpolator(order_interpolator)
113 SUBROUTINE download_link (self, link, all_cells, only, all_times)
115 class(
link_t),
pointer:: link
116 logical,
optional:: all_cells, all_times
117 integer,
optional:: only
118 real(kind=KindScalarVar),
dimension(:,:,:),
pointer:: d
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
125 character(len=64):: info
126 logical:: all_cells_l, lock_it
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')
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)
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
155 if (.not.
associated(link%nbor))
then 163 if (
present(all_cells))
then 164 all_cells_l = all_cells
165 else if (
target%quality == 1024)
then 168 all_cells_l = .false.
170 if (self%check_filled)
then 171 allocate (filled(
target%gn(1),
target%gn(2),
target%gn(3),
target%nv))
173 allocate (src(
target%gn(1),
target%gn(2),
target%gn(3)))
177 if (.not.all_cells_l)
then 178 filled(l(1):u(1),l(2):u(2),l(3):u(3),:) =
target%level
180 src(l(1):u(1),l(2):u(2),l(3):u(3)) =
target%id
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
191 nbor => link%nbors_by_level
197 do while (
associated(nbor))
198 if (.not.
associated(nbor%task))
then 199 write (io_unit%log,*)
'WARNING: nbor task not associated' 204 source => task2patch(nbor%task)
205 call source%log (
'download', 3)
206 if (source%quality == 1024)
then 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
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')
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
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 247 if (source%level==
target%level)
then 248 select case (trim(use_same))
250 call different (self,
target, source, jt, pt, all_cells_l, only=only)
252 call guard_zones%method (
target, source, jt, pt)
254 call guard_zones%method (
target, source, jt, pt, same=.true.)
256 call self%same (
target, source, jt, pt, all_cells_l, only=only)
261 else if (source%level==
target%level-1 .and.
present(all_cells))
then 262 call remesh%prolong (
target, source, jt, pt, all_cells)
266 else if ((source%level-
target%level) <= -1)
then 267 select case (trim(use_prolong))
269 call different (self,
target, source, jt, pt, all_cells_l, only=only)
271 call prolong (self,
target, source, jt, pt, all_cells_l, only=only)
273 call guard_zones%method (
target, source, jt, pt)
275 call remesh%prolong (
target, source, jt, pt, all_cells_l)
277 call different_prolong (self,
target, source, jt, pt, all_cells_l, only=only)
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)
287 else if ((source%level-
target%level) == 1)
then 288 select case (trim(use_restrict))
290 call different (self,
target, source, jt, pt, all_cells_l, only=only)
292 call different_prolong (self,
target, source, jt, pt, all_cells_l, only=only)
294 call guard_zones%method (
target, source, jt, pt)
296 call different_restrict (self,
target, source, jt, pt, only=only)
300 if (
target%id==id_debug)
then 301 write(io_unit%log,
'(a,i6,1p,4e12.5)')
'source: id,min,max,fmin,fmaxval:', &
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:', &
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)
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)
320 call different_rt (self,
target, source, jt, pt)
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
336 if (
present(only))
then 337 if (only==
target%idx%d)
then 338 call target%check_density (
" after download")
341 call target%check_density (
" after download")
343 if (self%check_filled)
then 344 if (all_cells_l)
then 352 if (any(filled(l(1):u(1),l(2):u(2),l(3):u(3),:) == -1))
then 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
361 write (io_unit%log,
'(a,i6)') &
362 'INFO: some guard zone cells were not filled, iD =',
target%id
367 if (filled(ix,iy,iz,1) == -1)
then 369 write(io_unit%log,
'(a,3i3)')
'not filled:', ix, iy, iz
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
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
385 associate(task => link%task)
388 task%not_filled = n_not
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
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
416 real(8),
dimension(3) :: h, dist, la, ua
418 class(mesh_t),
pointer :: sm, tm
421 call trace%begin (
'download_t%dnload_range', 2)
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
431 if (
target%id == id_debug)
write (io_unit%log,*)
'dist(2)', dist
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 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
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
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
462 where (
target%n == 1)
463 l = min(l,
target%mesh%ub)
464 u = max(u,
target%mesh%lb)
466 if (
target%id == id_debug)
write (io_unit%log,*)
'l,u(3)', l, u
471 if (
present(offset))
then 472 if (source%level==
target%level)
then 473 offset = nint(-dist/
target%ds + source%offset -
target%offset)
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')
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
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))
499 function target_index (p, roundup, nearest)
result (ii)
501 logical,
optional :: roundup, nearest
504 p = p/
target%mesh%d +
target%offset
505 if (
present(roundup))
then 507 else if (
present(nearest))
then 512 end function target_index
513 END SUBROUTINE dnload_range
522 SUBROUTINE same_linear (self, target, source, jt, pt, all_cells, only, all_times)
524 class(patch_t),
pointer ::
target, source
525 logical,
optional :: all_cells, all_times
526 integer,
optional :: only
528 integer :: l(3), u(3), o(3), ii(3), jt(2), i
529 real(8) :: position(3), xx, dist(3), start
531 integer :: ix, iy, iz, iv, iv1, iv2, li(3), ui(3), cells
532 integer,
save :: itimer=0
533 logical :: all_cells_l, no
535 call trace%begin(
'download_t%same_linear', itimer=itimer)
536 if (detailed_timer)
then 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 551 if (
present(all_cells))
then 552 all_cells_l = all_cells
554 all_cells_l = .false.
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
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")
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))
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
585 if (source%quality > 1000.)
then 586 if (iv==
target%idx%d .or. (iv >=
target%idx%px .and. iv <=
target%idx%pz))
then 588 write (io_unit%log,*)
'dnload all cells for target, source =', &
589 target%id, source%id, source%quality, source%level
593 if (iv ==
target%idx%dphi) cycle
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
607 if (product(max(0,u-l+1))>0)
then 623 if (
target%kind(1:4) ==
'zeus' .and. &
624 target%mesh(i)%h(iv) == -0.5 .and. &
632 if (
target%kind(1:4) ==
'zeus' .and. &
633 iv >=
target%idx%px .and. iv <=
target%idx%pz)
then 635 if (iv ==
target%idx%px+i-1 .and.
target%n(i) > 1) u(i) = u(i) + 1
638 if (
target%kind(1:4) ==
'zeus' .and. &
639 iv >=
target%idx%bx .and. &
640 iv <=
target%idx%bz)
then 642 if (iv ==
target%idx%bx+i-1 .and.
target%n(i) > 1) u(i) = u(i) + 1
645 if (
target%kind(1:4) ==
'zeus' .and. &
646 iv >=
target%idx%ex .and. &
647 iv <=
target%idx%ez)
then 649 if (iv /=
target%idx%ex+i-1 .and.
target%n(i) > 1) u(i) = u(i) + 1
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:", &
659 cells = product(u-l+1)
660 if (source%nt == 1)
then 661 if (all_cells_l)
then 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)
668 if (self%check_filled)
then 670 filled(ix,iy,iz,iv) = source%level
672 if (verbose > 0)
then 674 src(ix,iy,iz) = source%id
685 if (ix+o(1) >
size(source%mem,1))
then 686 print
'(5(3i4,2x))', l, u, o, shape(source%mem)
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)))
693 if (self%check_filled)
then 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)
700 if (verbose > 0)
then 703 no = all((ii>=li).and.(ii<=ui))
704 src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
714 else if (source%unsigned(iv))
then 720 if (all_cells_l)
then 723 if (self%check_filled)
then 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
731 src(ix,iy,iz) = source%id
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)
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)))
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)))
772 if (self%check_filled)
then 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)
779 if (verbose > 0)
then 782 no = all((ii>=li).and.(ii<=ui))
783 src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
797 if (all_cells_l)
then 803 if (self%check_filled)
then 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), &
811 filled(ix,iy,iz,iv) = merge(filled(ix,iy,iz,iv), &
812 int(source%level, kind=int8), no)
814 if (verbose > 0)
then 816 no = filled(ix,iy,iz,iv) > source%level
817 src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
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)
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)))
846 if (self%check_filled)
then 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)
853 if (verbose > 0)
then 856 no = all((ii>=li).and.(ii<=ui))
857 src(ix,iy,iz) = merge(src(ix,iy,iz), source%id, no)
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)
875 if (detailed_timer)
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)
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
887 call trace%end (itimer)
888 END SUBROUTINE same_linear
897 SUBROUTINE same_lagrange (self, target, source, jj, pt, all_cells, only, all_times)
899 class(patch_t),
pointer ::
target, source
900 logical,
optional :: all_cells, all_times
901 integer,
optional :: only
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)
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
921 call source%timeslots (iit, times)
922 call lagrange%sequence_weights (
target%time, times, j0, j1, w, order)
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 934 if (
present(all_cells))
then 935 all_cells_l = all_cells
937 all_cells_l = .false.
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
944 if (source%quality > 1000.)
then 945 if (iv==
target%idx%d .or. &
946 (iv >=
target%idx%px .and. iv <=
target%idx%pz))
then 948 write (io_unit%log,*)
'(2)dnload all cells for target, source =', &
949 target%id, source%id, source%quality, source%level
953 if (iv ==
target%idx%dphi) cycle
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 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
977 if (.not. all_cells_l)
then 982 if (all((ii>=li).and.(ii<=ui)))
then 994 write(io_unit%log,*) &
995 'WARNING: attempt to set internal cells in download_t%same_lagrange, using filter!' 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
1010 if (source%unsigned(iv))
then 1016 if (all_cells_l .or. any((ii<li).or.(ii>ui)))
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))
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))
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
1032 src(ix,iy,iz) = source%id
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))
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))
1050 if (self%check_filled)
then 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
1056 if (verbose > 0)
then 1058 src(ix,iy,iz) = source%id
1063 target%mem(ix,iy,iz,iv,
target%it,1) = &
1064 exp(
target%mem(ix,iy,iz,iv,
target%it,1))
1080 if (all_cells_l .or. any((ii<li).or.(ii>ui)))
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))
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))
1091 if (self%check_filled)
then 1092 filled(ix,iy,iz,iv) = source%level
1094 src(ix,iy,iz) = source%id
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))
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))
1112 if (self%check_filled)
then 1114 filled(ix,iy,iz,iv) = source%level
1116 if (verbose > 0)
then 1118 src(ix,iy,iz) = source%id
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))
1137 call trace%end (itimer)
1138 END SUBROUTINE same_lagrange
1143 SUBROUTINE same_rt (self, target, source, jt, pt)
1145 class(patch_t),
pointer ::
target, source
1148 integer :: l(3), u(3), o(3), ii(3), jt(2), li(3), ui(3)
1149 real(8) :: position(3)
1151 integer :: ix, iy, iz, iv
1152 integer,
save :: itimer=0
1154 call trace%begin(
'download_t%same_rt', itimer=itimer)
1155 call dnload_range (
target, source, 1, l, u, o)
1163 if (self%check_filled)
then 1164 if (filled(ix,iy,iz,iv) > source%level) cycle
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
1173 src(ix,iy,iz) = source%id
1179 call trace%end (itimer)
1180 END SUBROUTINE same_rt
1185 SUBROUTINE different (self, target, source, jt, pt, all_cells, only)
1187 class(patch_t),
pointer ::
target, source
1188 logical,
optional:: all_cells
1189 integer,
optional:: only
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
1195 real :: pt(2), pp(3)
1196 real(kind=KindScalarVar):: qint
1197 logical :: all_cells_l
1199 integer,
save :: itimer=0
1201 call trace%begin(
'download_t%different', itimer=itimer)
1202 if (detailed_timer)
then 1206 if (
present(only))
then 1213 if (
present(all_cells))
then 1214 all_cells_l = all_cells
1216 all_cells_l = .false.
1219 call dnload_range (
target, source, iv, l, u)
1222 all_cells_l = all_cells_l .or. all((l >= li) .and. (u<=ui))
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))
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)
1250 call source%interpolate (
target, ii, iv, jt, pt)
1252 if (.not.all_cells_l)
then 1253 if (all((ii>=li .and. ii<=ui)))
then 1254 write (io_unit%log,*)
'INTERNAL BREACH!', ii
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
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)
1277 if (detailed_timer)
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)
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
1289 call trace%end (itimer)
1290 END SUBROUTINE different
1295 SUBROUTINE different_prolong (self, target, source, jt, pt, all_cells, only)
1297 class(patch_t),
pointer ::
target, source
1298 logical,
optional:: all_cells
1299 integer,
optional:: only
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
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
1316 call trace%begin(
'download_t%different_prolong', itimer=itimer)
1318 if (
present(only))
then 1325 if (
present(all_cells))
then 1326 all_cells_l = all_cells
1328 all_cells_l = .false.
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
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
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)
1352 if (iv ==
target%idx%dphi) cycle
1353 call dnload_range (
target, source, iv, l, u)
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)), &
1369 1
format(
"download_t%different_prolong: source =", i6,
" displ:",3f11.6, &
1370 " l,u:",2(3i4,2x),
" all_cells:",l2)
1377 if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0))
then 1378 if (self%check_filled)
then 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)
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
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)
1405 if (source%unsigned(iv))
then 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 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)
1415 qint = interpolator%four_d (q, jj(1), jj(2), jj(3), jt, [pp,pt(2)])
1421 target%mem(ix,iy,iz,iv,
target%it,1) = qint
1422 if (self%check_filled) filled(ix,iy,iz,iv) = source%level
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
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)
1444 if (detailed_timer)
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
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
1458 call trace%end (itimer)
1459 END SUBROUTINE different_prolong
1464 SUBROUTINE prolong (self, target, source, jt, pt, all_cells, only)
1466 class(patch_t),
pointer ::
target, source
1467 logical,
optional:: all_cells
1468 integer,
optional:: only
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
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)
1485 call trace%begin(
'download_t%prolong', itimer=itimer)
1487 if (
present(only))
then 1494 if (
present(all_cells))
then 1495 all_cells_l = all_cells
1497 all_cells_l = .false.
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
1505 call dnload_range (
target, source,
target%idx%d, l, u)
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 ))
1523 if (iv ==
target%idx%dphi) cycle
1524 call dnload_range (
target, source, iv, l, u)
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:", &
1537 all_cells_l = all_cells_l .or. all((l >= li) .and. (u <= ui))
1544 if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0))
then 1561 if (iv==1) cells = cells+1
1562 pos =
target%myposition (ii, iv)
1563 call source%index_space (pos, iv, jj, pp)
1582 source%mem(jj(1):jj(1)+1,jj(2):jj(2)+1,jj(3):jj(3)+1,iv,jt(:),1)
1608 if (all_cells_l .or. any((ii-li) < 0 .or. (ii-ui) > 0))
then 1610 target%mem(ix,iy,iz,iv,
target%it,1) = tbuf(n)
1618 deallocate (sbuf, pbuf, tbuf, sdbuf, tdbuf)
1619 if (detailed_timer)
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)
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
1631 call trace%end (itimer)
1633 subroutine mem_interp
1634 real:: q1, q2, q3, q4
1646 if (
target%pervolume(iv))
then 1647 sbuf(1:m,:,:,:,:) = sbuf(1:m,:,:,:,:)/sdbuf(1:m,:,:,:,:)
1652 if (
target%unsigned(iv))
then 1653 sbuf(1:m,:,:,:,:) = log(sbuf(1:m,:,:,:,:))
1659 associate(p1=>pbuf(j,1), p2=>pbuf(j,2), p3=>pbuf(j,3), p4=>pbuf(j,4))
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))))
1677 if (
target%unsigned(iv))
then 1678 tbuf(1:m) = exp(tbuf(1:m))
1689 if (
target%pervolume(iv))
then 1690 tbuf(1:m) = tbuf(1:m)*tdbuf(1:m)
1692 end subroutine mem_interp
1693 END SUBROUTINE prolong
1698 SUBROUTINE different_rt (self, target, source, jt, pt)
1700 class(patch_t),
pointer ::
target, source
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
1706 integer,
save :: itimer=0
1709 call trace%begin(
'download_t%different_rt', itimer=itimer)
1712 call dnload_range (
target, source, 1, l, u)
1720 if (self%check_filled)
then 1721 if (filled(ix,iy,iz,iv) > source%level) cycle
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
1730 call trace%end (itimer)
1731 END SUBROUTINE different_rt
1736 SUBROUTINE different_restrict (self, target, source, jt, pt, only)
1738 class(patch_t),
pointer ::
target, source
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
1747 integer,
save :: itimer=0
1748 real :: pt(2), pp(3)
1749 real(kind=KindScalarVar):: qint
1753 real(8):: rratio(3), dvsource, dvtarget, dwsource, dwtarget, filling_factor, &
1754 accumulate, filling_total, dwtargeti
1757 call trace%begin(
'download_t%different_restrict', itimer=itimer)
1758 if (detailed_timer)
then 1762 rratio =
target%ds / source%ds
1763 ilook = merge(floor(rratio)/2 + 1, 0, source%n>1)
1764 dvsource = product(source%ds)
1765 dvtarget = product(
target%ds)
1766 if (
present(only))
then 1777 if (iv ==
target%idx%dphi) cycle
1780 if (iv ==
target%idx%bx)
then 1781 dwsource = source%ds(2) * source%ds(3)
1782 dwtarget =
target%ds(2) *
target%ds(3)
1784 else if (iv ==
target%idx%by)
then 1785 dwsource = source%ds(3) * source%ds(1)
1786 dwtarget =
target%ds(3) *
target%ds(1)
1788 else if (iv ==
target%idx%bz)
then 1789 dwsource = source%ds(1) * source%ds(2)
1790 dwtarget =
target%ds(1) *
target%ds(2)
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))
1804 if (iv==1) cells = cells+1
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
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
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)
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)
1823 if (any([jx,jy,jz] < source%li .or. [jx,jy,jz] > source%ui)) cycle
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)
1833 filling_factor = product(max(ue - le, 0.0),mask=qmask)
1834 filling_total = filling_total + filling_factor
1837 if (filling_factor > 0.0)
then 1838 accumulate = accumulate &
1839 + (q1(jx,jy,jz) * pt(1) + q2(jx,jy,jz) * pt(2)) * &
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
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)
1867 if (detailed_timer)
then 1868 write (stdout,
'(a,6p,f6.3," mus/cell")') &
1869 "download_t%different_restrict:", (wallclock()-start)/cells
1871 call trace%end (itimer)
1872 END SUBROUTINE different_restrict
1877 SUBROUTINE restrict_centered (self, target, source, jt, pt, only)
1879 class(patch_t),
pointer ::
target, source
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
1888 integer,
save :: itimer=0
1889 real :: pt(2), pp(3)
1890 real(kind=KindScalarVar):: qint
1892 real(8):: rratio(3), dvsource, dvtarget, dwsource, dwtarget, filling_factor, &
1893 accumulate(target%nv), filling_total, dwtargeti
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)
1900 dvtarget = product(
target%ds)
1901 if (
present(only))
then 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))
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
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
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)
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)
1937 if (any([jx,jy,jz] < source%li .or. [jx,jy,jz] > source%ui)) cycle
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)
1947 filling_factor = product(max(ue - le, 0.0),mask=qmask)
1948 filling_total = filling_total + filling_factor
1951 if (filling_factor > 0.0)
then 1953 accumulate(iv) = accumulate(iv) &
1954 + (q1(jx,jy,jz,iv) * pt(1) + q2(jx,jy,jz,iv) &
1955 * pt(2)) * filling_factor
1962 if (filling_total > 0.0)
then 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)
1968 if (self%check_filled)
then 1969 filled(ix,iy,iz,iv) = source%level
1976 call trace%end (itimer)
1977 END SUBROUTINE restrict_centered
1980 SUBROUTINE set_verbose (new)
1983 END SUBROUTINE set_verbose
1988 SUBROUTINE test (self)
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)
1995 real(8),
dimension(3) :: pos
1996 real(4),
dimension(3) :: pp
1997 class(mesh_t),
pointer :: mesh
1998 logical :: ok, allok
2001 call trace%begin (
'download_t%test')
2002 id_start = mpi%id%update(0)
2003 self%interpolate => selectinterpolator(-1)
2009 call target%setup (size=0.1_8*[1,1,1], n=[16,16,16], nv=2, nt=2, nw=1)
2014 allocate (source(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)
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)
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)
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)
2049 call dnload_range (
target, src, iv, l, u)
2053 mesh =>
target%mesh(dir)
2054 pos(dir) =
target%position(dir) + mesh%r(k) + mesh%h(iv)*
target%ds(dir)
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)
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)
2067 if (io%verbose>0 .or. io_unit%do_validate) &
2068 write (io_unit%log,*)
' ' 2079 call set_position (src,
target, j)
2080 if (i==2.and.(j==2.or.j==4)) cycle
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' 2097 call src%time_interval (
target%time, jt, pt)
2098 call prolong (self,
target, src, jt, pt, all_cells=.false.)
2100 call check_values (
target, ok)
2101 allok = allok.and.ok
2104 if (io%master)
write (io_unit%log,*)io%hl
2107 write (io_unit%log,*)
'download_t%test passed' 2108 write (io_unit%log,*)io%hl
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')
2121 call source(i)%dealloc
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)
2135 subroutine set_position (source, target, j)
2136 type(patch_t) :: source, target
2141 source%position =
target%position - 0.5_8*(source%size+
target%size) + &
2142 2.*abs(source%ds-
target%ds)
2144 source%position =
target%position - 0.5_8*(source%size+
target%size) + &
2145 2.*abs(source%ds-
target%ds)
2147 source%position =
target%position
2149 source%position =
target%position + 0.5_8*(source%size+
target%size) - &
2150 2.*abs(source%ds-
target%ds)
2152 source%position =
target%position + 0.5_8*(source%size+
target%size) - &
2153 2.*abs(source%ds-
target%ds)
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
2162 subroutine set_values (patch)
2163 type(patch_t) :: patch
2164 class(mesh_t),
pointer :: m1, m2, m3
2165 integer :: ix, iy, iz, iv
2172 patch%mem(:,:,:,iv,:,:) = 0.0
2174 z = m3%p + m3%r(iz) + m3%h(iv)*m3%d
2176 y = m2%p + m2%r(iy) + m2%h(iv)*m2%d
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
2185 end subroutine set_values
2189 subroutine check_values (patch, ok)
2190 type(patch_t) :: patch
2191 class(mesh_t),
pointer :: m1, m2, m3
2193 integer :: ix, iy, iz, iv, n(3)
2194 real :: x, y, z, v, w, eps
2199 eps=0.1*max(m1%d,m2%d,m3%d)
2203 z = m3%p + m3%r(iz) + m3%h(iv)*m3%d
2205 y = m2%p + m2%r(iy) + m2%h(iv)*m2%d
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)
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
2226 write (io_unit%log,*)
'check_values: n_zero, n_ok, n_nok =', n
2231 write (io_unit%log,*)
'check_values: n_zero, n_ok, n_nok =', n
2233 end subroutine check_values
Each thread uses a private timer data type, with arrays for start time and total time for each regist...
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 ...
Template module for patches, which adds pointers to memory and mesh, and number of dimensions and var...
This module is a placeholder for shared information, solving the problem of making a piece of informa...
Module for Lagrange interpolation.
Template module for mesh.
Module with list handling for generic class task_t objects.
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.