19 class(
patch_t),
pointer::
target => null()
20 class(
patch_t),
pointer:: source => null()
24 integer,
dimension(:,:),
pointer :: t_l => null(), t_u => null()
25 integer,
dimension(:,:),
pointer :: s_l => null(), s_u => null()
26 logical :: first_time = .true.
27 logical :: rt_target = .false.
28 logical :: all_cells = .false.
29 procedure(dnload_from_same),
pointer :: dnload => null()
30 integer :: int_order = -1
35 procedure:: check_in_nbors
36 procedure:: check_in_target
37 procedure:: dnload_from_higher
38 procedure:: dnload_from_lower
39 procedure:: dnload_from_same
41 procedure:: init_higher
42 procedure:: init_lower
44 procedure:: interpolate
45 procedure:: interpolate_unsigned
46 procedure:: interpolate_one
47 procedure:: interpolate_unsigned_one
48 procedure:: remove_hub
49 procedure:: select_dnloader
57 SUBROUTINE init (self, link, all_cells)
59 class(
link_t),
pointer :: link
60 logical,
optional :: all_cells
62 class(
patch_t),
pointer :: patch
64 call trace%begin (
'data_hub_t%init')
65 if (self%first_time)
then 66 patch => task2patch(link%task)
68 call self%select_dnloader(link, all_cells)
69 self%first_time = .false.
77 SUBROUTINE select_dnloader (self, link, all_cells)
79 class(
link_t),
pointer :: link
80 logical,
optional :: all_cells
82 class(
link_t),
pointer :: nbor
83 class(
patch_t),
pointer :: source, patch
88 do while (
associated(nbor))
89 source => task2patch(nbor%task)
90 if (nbor%download)
then 91 if (trim(patch%kind)/=
'rt_solver')
then 92 if ((source%is_clear(bits%no_download)))
then 93 if (abs(source%level-patch%level) < 2)
then 94 call self%add_hub (patch, source, hub, all_cells)
95 if (source%level == patch%level)
then 96 hub%dnload => dnload_from_same
98 else if (source%level > patch%level)
then 99 hub%dnload => dnload_from_higher
102 hub%dnload => dnload_from_lower
107 else if (trim(source%kind)==
'rt_solver')
then 108 if (abs(source%level-patch%level) < 2)
then 109 call self%add_hub (patch, source, hub, all_cells)
110 hub%rt_target = .true.
111 if (source%level == patch%level)
then 112 hub%dnload => dnload_from_same
114 else if (source%level > patch%level)
then 115 hub%dnload => dnload_from_higher
118 hub%dnload => dnload_from_lower
126 END SUBROUTINE select_dnloader
134 SUBROUTINE add_hub (self, target, source, hub, all_cells)
136 class(
patch_t),
pointer ::
target, source
138 logical,
optional:: all_cells
146 hub%level = source%level
148 if (
present(all_cells)) hub%all_cells = all_cells
149 if (.not.
associated(self%head))
then 155 do while (
associated(next).and.(.not.found))
156 if (.not.
associated(next%source,hub%source))
then 157 if (next%source%level > hub%source%level)
then 160 if (
associated(prev))
then 171 if (
associated(prev))
then 174 call io%abort(
"data_hub_t::add_hub:: location not found?")
178 END SUBROUTINE add_hub
184 SUBROUTINE remove_hub (self, id, clear)
186 integer,
optional :: id
187 logical,
optional :: clear
192 if (
present(id))
then 194 if (
associated(self%head))
then 196 do while(
associated(this).and.(.not.cleared))
197 if (this%id == id)
then 198 if (
associated(this%prev))
then 199 this%prev%next => this%next
200 this%next%prev => this%prev
201 deallocate(this%s_l, this%s_u)
202 deallocate(this%t_l, this%t_u)
206 if (
associated(this%next))
then 207 self%head => this%next
209 self%first_time = .true.
211 deallocate(this%s_l, this%s_u)
212 deallocate(this%t_l, this%t_u)
221 if (
present(clear))
then 223 do while(
associated(this))
225 deallocate(this%s_l, this%s_u)
226 deallocate(this%t_l, this%t_u)
231 END SUBROUTINE remove_hub
239 SUBROUTINE init_same (self)
242 class(
mesh_t),
pointer ::m
246 dist = self%source%distance(self%target)
247 allocate (self%t_l(3,1))
248 allocate (self%t_u(3,1))
249 allocate (self%s_l(3,1))
250 allocate (self%s_u(3,1))
252 m => self%target%mesh(i)
253 if (dist(i) > self%source%ds(i))
then 256 self%s_l(i,1) = self%t_l(i,1) - m%n
257 self%s_u(i,1) = self%t_u(i,1) - m%n
258 else if (dist(i) < -self%source%ds(i))
then 261 self%s_l(i,1) = self%t_l(i,1) + m%n
262 self%s_u(i,1) = self%t_u(i,1) + m%n
266 self%s_l(i,1) = self%t_l(i,1)
267 self%s_u(i,1) = self%t_u(i,1)
270 END SUBROUTINE init_same
277 SUBROUTINE dnload_from_same (self, only)
279 integer,
optional :: only
281 integer :: jt(2), iv, sl(3), su(3), tl(3), tu(3), iv1, iv2
283 integer,
save:: itimer=0
285 call trace%begin (
'data_hub_t%dnload_same',itimer=itimer)
286 call self%source%time_interval (self%target%time, jt, pt)
287 sl = self%s_l(:,1); su = self%s_u(:,1)
288 tl = self%t_l(:,1); tu = self%t_u(:,1)
289 if (self%rt_target)
then 290 self%target%mem(tl(1):tu(1),tl(2):tu(2),tl(3):tu(3),:,self%target%it,1)= &
291 pt(1) *self%source%mem(sl(1):su(1),sl(2):su(2),sl(3):su(3),:,jt(1),1) + &
292 pt(2) *self%source%mem(sl(1):su(1),sl(2):su(2),sl(3):su(3),:,jt(2),1)
294 if (
present(only))
then 302 if (self%source%unsigned(iv))
then 303 self%target%mem(tl(1):tu(1),tl(2):tu(2),tl(3):tu(3),iv,self%target%it,1) = exp &
304 (pt(1) *log(self%source%mem(sl(1):su(1),sl(2):su(2),sl(3):su(3),iv,jt(1),1)) + &
305 pt(2) *log(self%source%mem(sl(1):su(1),sl(2):su(2),sl(3):su(3),iv,jt(2),1)))
307 self%target%mem(tl(1):tu(1),tl(2):tu(2),tl(3):tu(3),iv,self%target%it,1)= &
308 pt(1) *self%source%mem(sl(1):su(1),sl(2):su(2),sl(3):su(3),iv,jt(1),1) + &
309 pt(2) *self%source%mem(sl(1):su(1),sl(2):su(2),sl(3):su(3),iv,jt(2),1)
313 call trace%end(itimer)
314 END SUBROUTINE dnload_from_same
320 SUBROUTINE init_higher (self)
323 class(
mesh_t),
pointer:: m, m2
324 class(
patch_t),
pointer:: source, target
326 integer :: iv, i, j(3), nv
327 real(8) :: eps, pl(3), pu(3), dp
329 source => self%source
330 target => self%target
331 if (self%rt_target)
then 332 allocate (self%t_l(3,1))
333 allocate (self%t_u(3,1))
336 allocate (self%t_l(3,
target%nv))
337 allocate (self%t_u(3,
target%nv))
347 pl(i) = m%r(m%li) - eps
348 pu(i) = m%r(m%ui) + eps
353 self%t_l(:,iv) =
target%index_only2 (source, pl, iv, roundup=.true.)
354 self%t_u(:,iv) =
target%index_only2 (source, pu, iv)
358 if (.not.self%all_cells)
then 362 self%t_l(i,iv) = max(m%lb,min(m%uo,self%t_l(i,iv)))
363 self%t_u(i,iv) = max(m%lo,min(m%ub,self%t_u(i,iv)))
364 dp = round(m%p-m2%p,12)+m%h(iv)*(m%d-m2%d)
365 if (self%t_u(i,iv)>m%lo)
then 366 if (m%r(m%li)+dp > m2%r(m2%ui)) &
367 self%t_u(i,iv) = min(self%t_u(i,iv),m%lo)
369 if (self%t_l(i,iv)<m%uo)
then 370 if (m%r(m%ui)+dp < m2%r(m2%li)) &
371 self%t_l(i,iv) = max(self%t_l(i,iv),m%uo)
373 if (self%t_u(i,iv)>m%ui)
then 374 if (m%r(m%uo)+dp > m2%uf)
then 375 self%t_u(i,iv) = min(self%t_u(i,iv),m%ui)
377 self%t_u(i,iv) = max(self%t_l(i,iv),m%ub)
380 if (self%t_l(i,iv)<m%li)
then 381 if (m%r(m%lo)+dp < m2%lf)
then 382 self%t_l(i,iv) = max(self%t_l(i,iv),m%li)
384 self%t_l(i,iv) = min(self%t_l(i,iv),m%lb)
390 END SUBROUTINE init_higher
399 SUBROUTINE dnload_from_higher (self, only)
401 integer,
optional :: only
403 integer :: jt(2), iv, tl(3), tu(3), i, j, k, iv1, iv2
404 integer :: d, l(3), u(3), ii
405 real(kind=KindScalarVar),
dimension(:,:,:,:,:),
pointer :: mem
406 real(kind=KindScalarVar),
dimension(:,:,:,:),
pointer :: cube
407 real(kind=KindScalarVar) :: tcell
408 integer,
save:: itimer=0
410 class(
patch_t),
pointer :: source, target
411 class(
mesh_t),
pointer :: m
414 call trace%begin (
'data_hub_t%dload_higher',itimer=itimer)
415 source => self%source
416 target => self%target
417 allocate(cube(2,2,2,2))
418 call source%time_interval (
target%time, jt, pt)
419 if (self%rt_target)
then 420 tl = self%t_l(:,1); tu = self%t_u(:,1)
430 if ((self%all_cells) .or. &
431 any((l-
target%mesh%li) < 0 .or. (l-
target%mesh%ui) > 0))
then 432 call source%index_space2(
target,pos, 1, l, p)
434 mem => source%mem(l(1):u(1),l(2):u(2),l(3):u(3),:,:,1)
436 cube(:,:,:,1) = source%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1)
437 cube(:,:,:,2) = source%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1)
438 call self%interpolate(cube, pt, p,tcell)
439 target%mem(i,j,k,iv,
target%it,1) = tcell
446 if (
present(only))
then 453 d = self%source%idx%d
455 tl = self%t_l(:,iv); tu = self%t_u(:,iv)
464 if ((self%all_cells) .or. &
465 any((l-
target%mesh%li) < 0 .or. (l-
target%mesh%ui) > 0))
then 466 call source%index_space2(
target,pos, iv, l, p)
468 mem => source%mem(l(1):u(1),l(2):u(2),l(3):u(3),:,:,1)
469 if (source%pervolume(iv))
then 470 cube(:,:,:,1) = source%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1) / &
471 source%mem(l(1):u(1),l(2):u(2),l(3):u(3), d,jt(1),1)
472 cube(:,:,:,2) = source%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1) / &
473 source%mem(l(1):u(1),l(2):u(2),l(3):u(3), d,jt(2),1)
475 cube(:,:,:,1) = source%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1)
476 cube(:,:,:,2) = source%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1)
478 if (source%unsigned(iv))
then 479 call self%interpolate_unsigned(cube, pt, p,tcell)
480 target%mem(i,j,k,iv,
target%it,1) = tcell
481 else if (source%pervolume(iv))
then 482 call self%interpolate(cube, pt, p,tcell)
483 target%mem(i,j,k,iv,
target%it,1) = tcell *
target%mem(i,j,k,d,
target%it,1)
485 call self%interpolate(cube, pt, p,tcell)
486 target%mem(i,j,k,iv,
target%it,1) = tcell
495 call trace%end(itimer)
496 END SUBROUTINE dnload_from_higher
502 SUBROUTINE init_lower (self)
505 class(
mesh_t),
pointer:: m, m2
506 class(
patch_t),
pointer:: source, target
507 integer :: iv, i, j(3), nv
508 real(8) :: eps, pl(3), pu(3), dp
510 source => self%source
511 target => self%target
513 if (self%rt_target)
then 514 allocate (self%t_l(3,1))
515 allocate (self%t_u(3,1))
516 allocate (self%s_l(3,1))
517 allocate (self%s_u(3,1))
520 allocate (self%t_l(3,
target%nv))
521 allocate (self%t_u(3,
target%nv))
522 allocate (self%s_l(3,
target%nv))
523 allocate (self%s_u(3,
target%nv))
538 self%t_l(:,iv) =
target%index_only2 (source, pl, iv, roundup = .true.)
539 self%t_u(:,iv) =
target%index_only2 (source, pu, iv)
543 if (.not.self%all_cells)
then 547 self%t_l(i,iv) = max(m%lb,min(m%uo,self%t_l(i,iv)))
548 self%t_u(i,iv) = max(m%lo,min(m%ub,self%t_u(i,iv)))
549 dp = round(m%p-m2%p,12)+m%h(iv)*(m%d-m2%d)
550 if (self%t_u(i,iv)>m%lo)
then 551 if (m%r(m%li)+dp > m2%r(m2%ui)) &
552 self%t_u(i,iv) = min(self%t_u(i,iv),m%lo)
554 if (self%t_l(i,iv)<m%uo)
then 555 if (m%r(m%ui)+dp < m2%r(m2%li)) &
556 self%t_l(i,iv) = max(self%t_l(i,iv),m%uo)
558 if (self%t_u(i,iv)>m%ui)
then 559 if (m%r(m%uo)+dp > m2%uf)
then 560 self%t_u(i,iv) = min(self%t_u(i,iv),m%ui)
562 self%t_u(i,iv) = max(self%t_l(i,iv),m%ub)
565 if (self%t_l(i,iv)<m%li)
then 566 if (m%r(m%lo)+dp < m2%lf)
then 567 self%t_l(i,iv) = max(self%t_l(i,iv),m%li)
569 self%t_l(i,iv) = min(self%t_l(i,iv),m%lb)
579 pl(i) = m%r(self%t_l(i,iv))
580 pu(i) = m%r(self%t_u(i,iv))
582 self%s_l(:,iv) = source%index_only2 (
target,pl, iv)
583 self%s_u(:,iv) = source%index_only2 (
target,pu, iv, roundup=.true.)
584 if (.not.self%all_cells)
then 587 self%s_l(i,iv) = max(m%lo,min(m%ui,self%s_l(i,iv)))
588 self%s_u(i,iv) = min(m%uo,max(m%li,self%s_u(i,iv)))
592 END SUBROUTINE init_lower
604 SUBROUTINE dnload_from_lower (self, only)
606 integer,
optional :: only
608 integer :: jt(2),jt2(2), iv, sl(3), su(3), tl(3), tu(3), i, j, k
609 integer :: d, cl(3), ii, l(3), u(3), id, n(3), iv1, iv2
610 class(
mesh_t),
pointer :: m
611 class(
patch_t),
pointer:: source, target
612 real(kind=KindScalarVar),
dimension(:,:,:,:,:),
pointer :: mem
613 real(kind=KindScalarVar),
dimension(:,:,:),
pointer :: scr
614 real(kind=KindScalarVar),
dimension(:,:,:,:),
pointer :: rt_scr, rt_cube
615 real(kind=KindScalarVar),
dimension(:,:,:),
pointer :: cube
616 real(kind=KindScalarVar),
pointer :: scr3
617 real(kind=KindScalarVar),
dimension(:),
pointer :: rt_scr3
618 real(kind=KindScalarVar) :: tcell
619 real :: pt(2), pt2(2), p(3)
620 real(8) :: pos(3), time, pos2(3)
621 logical :: outside, found
622 integer,
save:: itimer=0
624 call trace%begin (
'data_hub_t%dnload_lower',itimer=itimer)
625 source => self%source
626 target => self%target
628 mem => self%source%mem(:,:,:,:,:,1)
629 if (
present(only))
then 642 call source%time_interval (
target%time, jt2, pt2)
645 if (self%rt_target)
then 646 sl = self%s_l(:,1); su = self%s_u(:,1)
647 tl = self%t_l(:,1); tu = self%t_u(:,1)
649 allocate (rt_scr(n(1), n(2), n(3),
target%nv))
653 cl = [sl(1)+i-1,sl(2)+j-1,sl(3)+k-1]
656 if ((cl(ii)<source%li(ii)).or.(cl(ii)>source%ui(ii))) outside = .true.
662 pos(ii) = m%r(cl(ii))
663 pos2(ii) = m%r(cl(ii))+m%p+m%h(1)*m%d
665 rt_scr3 => rt_scr(i,j,k,:)
666 call self%check_in_target(pos, 1, time, jt, pt, id, found, rt_scr=rt_scr3)
668 call self%check_in_nbors (pos, 1, time, jt, pt, id, found, rt_scr=rt_scr3)
672 call io%abort(
"data_hub::dnload_lower:: cell not found")
675 rt_scr(i,j,k,:) = pt2(1) * mem(cl(1),cl(2),cl(3),:, jt2(1)) + &
676 pt2(2) * mem(cl(1),cl(2),cl(3),:, jt2(2))
689 if ((self%all_cells) .or. &
690 any((cl-
target%mesh%li) < 0 .or. (cl-
target%mesh%ui) > 0))
then 693 pos(ii) = m%r(cl(ii))
695 call source%index_space2(
target,pos, 1, cl, p)
699 cube => rt_scr(l(1):u(1),l(2):u(2),l(3):u(3), iv)
700 call self%interpolate_one(cube, p,tcell)
701 target%mem(i,j,k,iv,
target%it,1) = tcell
710 sl = self%s_l(:,iv); su = self%s_u(:,iv)
711 tl = self%t_l(:,iv); tu = self%t_u(:,iv)
713 allocate (scr(n(1), n(2), n(3)))
717 cl = [sl(1)+i-1,sl(2)+j-1,sl(3)+k-1]
720 if ((cl(ii)<source%li(ii)).or.(cl(ii)>source%ui(ii))) outside = .true.
726 pos(ii) = m%r(cl(ii))
727 pos2(ii) = m%r(cl(ii))+m%p+m%h(iv)*m%d
730 call self%check_in_target(pos, iv, time, jt, pt, id, found, scr=scr3)
732 call self%check_in_nbors (pos, iv, time, jt, pt, id, found, scr=scr3)
736 call io%abort(
"data_hub::dnload_lower:: cell not found")
739 if (source%unsigned(iv))
then 740 scr(i,j,k) = exp(pt2(1) * log(mem(cl(1),cl(2),cl(3),iv, jt2(1))) + &
741 pt2(2) * log(mem(cl(1),cl(2),cl(3),iv, jt2(2))))
742 else if (source%pervolume(iv))
then 743 scr(i,j,k) = pt2(1) * (mem(cl(1),cl(2),cl(3),iv, jt2(1)) / &
744 mem(cl(1),cl(2),cl(3), d, jt2(1))) + &
745 pt2(2) * (mem(cl(1),cl(2),cl(3),iv, jt2(2)) / &
746 mem(cl(1),cl(2),cl(3), d, jt2(2)))
748 scr(i,j,k) = pt2(1) * mem(cl(1),cl(2),cl(3),iv, jt2(1)) + &
749 pt2(2) * mem(cl(1),cl(2),cl(3),iv, jt2(2))
763 if ((self%all_cells) .or. &
764 any((cl-
target%mesh%li) < 0 .or. (cl-
target%mesh%ui) > 0))
then 767 pos(ii) = m%r(cl(ii))
769 call source%index_space2(
target,pos, iv, cl, p)
772 cube => scr(l(1):u(1),l(2):u(2),l(3):u(3))
773 if (source%unsigned(iv))
then 774 call self%interpolate_unsigned_one(cube, p,tcell)
775 target%mem(i,j,k,iv,
target%it,1) = tcell
776 else if (self%source%pervolume(iv))
then 777 call self%interpolate_one(cube, p,tcell)
778 target%mem(i,j,k,iv,
target%it,1) = tcell *
target%mem(i,j,k,d,
target%it,1)
780 call self%interpolate_one(cube, p,tcell)
781 target%mem(i,j,k,iv,
target%it,1) = tcell
790 call trace%end(itimer)
791 END SUBROUTINE dnload_from_lower
798 SUBROUTINE check_in_target (self, pos, iv, time, jt, pt, id, found, scr, rt_scr)
800 real(8) :: pos(3), time
801 real(kind=KindScalarVar),
pointer,
optional :: scr
802 real(kind=KindScalarVar),
dimension(:),
pointer,
optional :: rt_scr
803 integer :: iv, jt(2), id
807 integer :: sl(3), su(3), tl(3), tu(3), i, j, k, aa, bb, cc
808 integer :: d, cl(3), ii, l(3), u(3), v
810 real(8) :: eps = 0.001_8
811 class(
mesh_t),
pointer :: m
812 class(
patch_t),
pointer:: source, target
813 real(kind=KindScalarVar),
dimension(:,:,:,:),
pointer :: cube
816 source => self%source
817 target => self%target
820 call target%index_space2(source,pos,iv, cl, p)
821 inside = all(cl>=
target%mesh%li .and. cl<=
target%mesh%ui)
826 if (id /=
target%id)
then 827 call target%time_interval (time, jt, pt)
830 allocate(cube(2,2,2,2))
833 if (self%rt_target)
then 835 cube(:,:,:,1) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),v,jt(1),1)
836 cube(:,:,:,2) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),v,jt(2),1)
837 call self%interpolate(cube, pt, p,rt_scr(v))
840 if (source%unsigned(iv))
then 841 cube(:,:,:,1) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1)
842 cube(:,:,:,2) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1)
843 call self%interpolate_unsigned(cube, pt, p,scr)
844 else if (source%pervolume(iv))
then 845 cube(:,:,:,1) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1) / &
846 target%mem(l(1):u(1),l(2):u(2),l(3):u(3), d,jt(1),1)
847 cube(:,:,:,2) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1) / &
848 target%mem(l(1):u(1),l(2):u(2),l(3):u(3), d,jt(2),1)
849 call self%interpolate(cube, pt, p,scr)
851 cube(:,:,:,1) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1)
852 cube(:,:,:,2) =
target%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1)
853 call self%interpolate(cube, pt, p,scr)
859 END SUBROUTINE check_in_target
869 SUBROUTINE check_in_nbors (self, pos, iv, time, jt, pt, id, found, scr, rt_scr)
871 real(8) :: pos(3), time
872 real(kind=KindScalarVar),
pointer,
optional :: scr
873 real(kind=KindScalarVar),
dimension(:),
pointer,
optional :: rt_scr
878 integer :: iv, sl(3), su(3), tl(3), tu(3), i, j, k
879 integer :: d, cl(3), ii, l(3), u(3), v
881 real(8) :: eps = 0.001_8
882 class(
mesh_t),
pointer :: m
883 class(
link_t),
pointer :: nbor
884 class(
patch_t),
pointer:: source,
target, nbpatch
885 real(kind=KindScalarVar),
dimension(:,:,:,:),
pointer :: cube
888 source => self%source
889 target => self%target
891 nbor => source%link%nbor
892 do while((
associated(nbor)).and..not.found)
893 if ((nbor%task%id /= source%id).and.(nbor%task%id /=
target%id))
then 894 nbpatch => task2patch(nbor%task)
895 if (abs(nbpatch%level-source%level)<2)
then 896 if (nbpatch%level == source%level)
then 897 cl = nbpatch%index_only2 (source,pos, iv, closest=.true.)
898 inside = all(cl>=nbpatch%mesh%li .and. cl<=nbpatch%mesh%ui)
900 if (id /=nbpatch%id)
then 901 call nbpatch%time_interval (time, jt, pt)
904 if(self%rt_target)
then 906 pt(1) * nbpatch%mem(cl(1),cl(2),cl(3),:,jt(1),1) + &
907 pt(2) * nbpatch%mem(cl(1),cl(2),cl(3),:,jt(2),1)
910 if (source%unsigned(iv))
then 912 pt(1) * log(nbpatch%mem(cl(1),cl(2),cl(3),iv,jt(1),1)) + &
913 pt(2) * log(nbpatch%mem(cl(1),cl(2),cl(3),iv,jt(2),1)))
914 else if (source%pervolume(iv))
then 916 pt(1) * (nbpatch%mem(cl(1),cl(2),cl(3),iv,jt(1),1) / &
917 nbpatch%mem(cl(1),cl(2),cl(3), d,jt(1),1))+ &
918 pt(2) * (nbpatch%mem(cl(1),cl(2),cl(3),iv,jt(2),1) / &
919 nbpatch%mem(cl(1),cl(2),cl(3), d,jt(2),1))
922 pt(1) * nbpatch%mem(cl(1),cl(2),cl(3),iv,jt(1),1) + &
923 pt(2) * nbpatch%mem(cl(1),cl(2),cl(3),iv,jt(2),1)
929 call nbpatch%index_space2(source,pos,iv, cl, p)
930 inside = all(cl>=nbpatch%mesh%li .and. cl<=nbpatch%mesh%ui)
932 if (id /=nbpatch%id)
then 933 call nbpatch%time_interval (time, jt, pt)
936 allocate(cube(2,2,2,2))
939 if (self%rt_target)
then 941 cube(:,:,:,1) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),v,jt(1),1)
942 cube(:,:,:,2) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),v,jt(2),1)
943 call self%interpolate(cube, pt, p,rt_scr(v))
946 if (source%unsigned(iv))
then 947 cube(:,:,:,1) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1)
948 cube(:,:,:,2) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1)
949 call self%interpolate_unsigned(cube, pt, p,scr)
950 else if (source%pervolume(iv))
then 951 cube(:,:,:,1) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1) / &
952 nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3), d,jt(1),1)
953 cube(:,:,:,2) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1) / &
954 nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3), d,jt(2),1)
955 call self%interpolate(cube, pt, p,scr)
957 cube(:,:,:,1) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(1),1)
958 cube(:,:,:,2) = nbpatch%mem(l(1):u(1),l(2):u(2),l(3):u(3),iv,jt(2),1)
959 call self%interpolate(cube, pt, p,scr)
970 END SUBROUTINE check_in_nbors
975 SUBROUTINE update (self, link, all_cells, only)
977 class(
link_t),
pointer :: link
978 logical,
optional :: all_cells
979 integer,
optional :: only
982 integer,
save:: itimer=0
985 if (self%first_time)
call self%init(link, all_cells)
988 do while(
associated(hub))
989 call hub%dnload(only=only)
993 END SUBROUTINE update
995 SUBROUTINE interpolate (self, src, pt, p, res)
997 real(kind=KindScalarVar),
dimension(:,:,:,:):: src
999 real(kind=KindScalarVar) :: res
1001 if (any(p < 0.0).or.any(p>1.0))
then 1003 call io%abort(
"p negative")
1006 pt(1)*((1.-p(3))*((1.-p(2))*((1.-p(1))*src(1,1,1,1) + &
1007 ( p(1))*src(2,1,1,1)) + &
1008 ( p(2))*((1.-p(1))*src(1,2,1,1) + &
1009 ( p(1))*src(2,2,1,1))) + &
1010 ( p(3))*((1.-p(2))*((1.-p(1))*src(1,1,2,1) + &
1011 ( p(1))*src(2,1,2,1)) + &
1012 ( p(2))*((1.-p(1))*src(1,2,2,1) + &
1013 ( p(1))*src(2,2,2,1)))) + &
1014 pt(2)*((1.-p(3))*((1.-p(2))*((1.-p(1))*src(1,1,1,2) + &
1015 ( p(1))*src(2,1,1,2)) + &
1016 ( p(2))*((1.-p(1))*src(1,2,1,2) + &
1017 ( p(1))*src(2,2,1,2))) + &
1018 ( p(3))*((1.-p(2))*((1.-p(1))*src(1,1,2,2) + &
1019 ( p(1))*src(2,1,2,2)) + &
1020 ( p(2))*((1.-p(1))*src(1,2,2,2) + &
1021 ( p(1))*src(2,2,2,2))))
1022 END SUBROUTINE interpolate
1024 SUBROUTINE interpolate_unsigned (self, src, pt, p, res)
1026 real(kind=KindScalarVar),
dimension(:,:,:,:):: src
1028 real(kind=KindScalarVar) :: res
1030 if (any(p < 0.0).or.any(p>1.0))
then 1032 call io%abort(
"p negative")
1035 pt(1)*((1.-p(3))*((1.-p(2))*((1.-p(1))*log(src(1,1,1,1)) + &
1036 ( p(1))*log(src(2,1,1,1))) + &
1037 ( p(2))*((1.-p(1))*log(src(1,2,1,1)) + &
1038 ( p(1))*log(src(2,2,1,1)))) + &
1039 ( p(3))*((1.-p(2))*((1.-p(1))*log(src(1,1,2,1)) + &
1040 ( p(1))*log(src(2,1,2,1))) + &
1041 ( p(2))*((1.-p(1))*log(src(1,2,2,1)) + &
1042 ( p(1))*log(src(2,2,2,1))))) + &
1043 pt(2)*((1.-p(3))*((1.-p(2))*((1.-p(1))*log(src(1,1,1,2)) + &
1044 ( p(1))*log(src(2,1,1,2))) + &
1045 ( p(2))*((1.-p(1))*log(src(1,2,1,2)) + &
1046 ( p(1))*log(src(2,2,1,2)))) + &
1047 ( p(3))*((1.-p(2))*((1.-p(1))*log(src(1,1,2,2)) + &
1048 ( p(1))*log(src(2,1,2,2))) + &
1049 ( p(2))*((1.-p(1))*log(src(1,2,2,2)) + &
1050 ( p(1))*log(src(2,2,2,2))))))
1051 END SUBROUTINE interpolate_unsigned
1053 SUBROUTINE interpolate_one (self, src, p, res)
1055 real(kind=KindScalarVar),
dimension(:,:,:):: src
1057 real(kind=KindScalarVar) :: res
1059 if (any(p < 0.0).or.any(p>1.0))
then 1061 call io%abort(
"p negative")
1064 ((1.-p(3))*((1.-p(2))*((1.-p(1))*src(1,1,1) + &
1065 ( p(1))*src(2,1,1)) + &
1066 ( p(2))*((1.-p(1))*src(1,2,1) + &
1067 ( p(1))*src(2,2,1))) + &
1068 ( p(3))*((1.-p(2))*((1.-p(1))*src(1,1,2) + &
1069 ( p(1))*src(2,1,2)) + &
1070 ( p(2))*((1.-p(1))*src(1,2,2) + &
1071 ( p(1))*src(2,2,2))))
1072 END SUBROUTINE interpolate_one
1074 SUBROUTINE interpolate_unsigned_one (self, src, p, res)
1076 real(kind=KindScalarVar),
dimension(:,:,:):: src
1078 real(kind=KindScalarVar) :: res
1080 if (any(p < 0.0).or.any(p>1.0))
then 1082 call io%abort(
"p negative")
1085 ((1.-p(3))*((1.-p(2))*((1.-p(1))*log(src(1,1,1)) + &
1086 ( p(1))*log(src(2,1,1))) + &
1087 ( p(2))*((1.-p(1))*log(src(1,2,1)) + &
1088 ( p(1))*log(src(2,2,1)))) + &
1089 ( p(3))*((1.-p(2))*((1.-p(1))*log(src(1,1,2)) + &
1090 ( p(1))*log(src(2,1,2))) + &
1091 ( p(2))*((1.-p(1))*log(src(1,2,2)) + &
1092 ( p(1))*log(src(2,2,2)))))
1093 END SUBROUTINE interpolate_unsigned_one
A communications hub between patches. Three methods are used: download_same, download_higher and down...
Template module for patches, which adds pointers to memory and mesh, and number of dimensions and var...
Template module for mesh.
Module with list handling for generic class task_t objects.
Template module for tasks.