28 procedure(restrict_simple),
pointer:: method => null()
31 procedure:: restrict_simple
32 procedure:: different_opt
33 procedure:: restrict_vect
34 procedure,
nopass:: info
36 real(8):: update_t(3)=0d0
38 logical:: detailed_timer=.false.
45 SUBROUTINE init (self)
47 character(len=6):: method=
'opt' 48 namelist /guard_zone_params/ method, detailed_timer, verbose
50 logical,
save:: first_time=.true.
52 call trace%begin (
'guard_zones_t%init')
57 read (io%input, guard_zone_params, iostat=iostat)
58 write (io%output, guard_zone_params)
63 select case (trim(method))
65 self%method => restrict_simple
67 self%method => different_opt
69 self%method => restrict_vect
71 call io%abort (
'guard zone method unknown')
80 SUBROUTINE restrict_simple (self, target, source, jt, pt, same)
82 class(
patch_t),
pointer::
target, source
85 logical,
optional:: same
87 real,
dimension(:,:,:),
pointer:: tmem, smem, td, sd
88 class(
mesh_t),
pointer:: mesh
89 integer:: i1, i2, i3, iv, jv, tl(3), tu(3), sl(3), su(3)
90 real(8):: tlc(3), tuc(3), slc(3), suc(3), dist(3)
91 real(4):: pl(3), pu(3)
93 integer,
save:: itimer=0
97 call trace%begin (
'guard_zones_t%restrict_simple', itimer=itimer)
100 call source%time_interval (
target%time, jt, pt)
106 tlc =
target%myposition (
target%mesh%lb, jv)
107 tuc =
target%myposition (
target%mesh%ub, jv)
111 sl = source%index_only (tlc, jv, nowrap=.true.)
112 su = source%index_only (tuc, jv, nowrap=.true.)
113 sl = min(max(sl,source%mesh%li),source%mesh%ui-1)
114 su = min(max(su,source%mesh%li),source%mesh%ui-1)
118 slc = source%myposition(sl, jv)
119 suc = source%myposition(su, jv)
120 tl =
target%index_only (slc, jv, nowrap=.true., roundup=.true.)
121 tu =
target%index_only (suc, jv, nowrap=.true., roundup=.true.)
122 tl = min(max(tl,
target%mesh%lb),
target%mesh%ub)
123 tu = min(max(tu,
target%mesh%lb),
target%mesh%ub)
125 write(io_unit%log,
'(a,i3,2(2x,3i4),2(2x,3f7.3))') &
126 'guard_zones_t%restrict_simple1: iv, tl, tu =', iv, tl, tu
131 do i3=tl(3),tu(3);
do i2=tl(2),tu(2);
do i1=tl(1),tu(1)
132 call interpolate (source,
target, [i1,i2,i3], iv, jt, pt)
133 end do; end do; end do
135 if (detailed_timer)
then 136 used = wallclock()-used
137 write (stdout,
'(a,6p,f6.3," mus/cell")') &
138 "guard_zone_t%restrict_simple:", used/product(tu-tl+1)
140 call trace%end (itimer)
141 END SUBROUTINE restrict_simple
147 SUBROUTINE different_opt (self, target, source, jt, pt, same)
149 class(patch_t),
pointer::
target, source
152 logical,
optional:: same
154 real,
parameter:: eps=1e-4
155 real,
dimension(:,:,:),
pointer:: tmem, smem, td, sd
156 class(mesh_t),
pointer:: mesh
157 integer:: i, i1, i2, i3, iv, jv, tl(3), tu(3), sl(3), su(3), cells
159 real(8):: tlc, tuc, slc, suc, dist, box, tc
160 class(mesh_t),
pointer:: tm, sm
162 integer,
save:: itimer=0
166 call trace%begin (
'guard_zones_t%different_opt', itimer=itimer)
167 if (detailed_timer)
then 170 call source%time_interval (
target%time, jt, pt)
181 if (
target%periodic(i)) &
182 dist = modulo(dist+0.5d0*box,box) - 0.5d0*box
191 slc = sm%r(sm%li) + dist
192 suc = sm%r(sm%ui) + dist
202 tl(i) = ceiling(tlc*tc + tm%o - eps)
203 tu(i) = ceiling(tuc*tc + tm%o - eps)
205 write(io_unit%log,
'(a,i3,2x,5f8.3,2x,2i4)') &
206 'guard_zones_t%different_opt: i,dist,tlc,tuc,slc,suc,tl,tu =', &
207 i, dist, tlc, tuc, slc, suc, tl(i), tu(i)
213 if (
present(same))
then 214 call interpolate (source,
target, [i1,i2,i3], iv, jt, pt, pi=pi, set=.true.)
215 do i3=tl(3),tu(3);
do i2=tl(2),tu(2);
do i1=tl(1),tu(1)
216 call interpolate (source,
target, [i1,i2,i3], iv, jt, pt, pi=pi)
217 end do; end do; end do
219 do i3=tl(3),tu(3);
do i2=tl(2),tu(2);
do i1=tl(1),tu(1)
220 call interpolate (source,
target, [i1,i2,i3], iv, jt, pt)
221 end do; end do; end do
224 if (detailed_timer)
then 225 cells = product(tu-tl+1)
226 write (stdout,
'(a,6p,2(2x,i6,i3),i6," cells,",f9.3," mus/cell")') &
227 "guard_zones_t%different_opt: target, level, source, level,", &
228 target%id,
target%level, source%id, source%level, cells, &
229 (wallclock()-start)/max(1,cells)
231 call trace%end (itimer)
232 END SUBROUTINE different_opt
237 SUBROUTINE restrict_vect (self, target, source, jt, pt, same)
239 class(patch_t),
pointer::
target, source
242 logical,
optional:: same
244 real,
dimension(:,:,:),
pointer:: tmem, smem, td, sd
245 class(mesh_t),
pointer:: mesh
246 integer:: i, i1, i2, i3, j, k, iv, tl(3), tu(3), sl(3), su(3), j0(3)
247 real(8):: tlc(3), tuc(3), slc(3), suc(3)
248 real(4):: p(3), p1, p2, p3, w(8)
250 integer,
save:: itimer=0
254 call trace%begin (
'guard_zones_t%restrict_vect', itimer=itimer)
255 if (detailed_timer) &
260 call source%time_interval (
target%time, jt, pt)
262 write(io_unit%log,
'(a,2(2x,i6,i3,2x,3f7.3),2x,3f7.3)') &
263 'target, level, lc, source, level, lc, relative =', &
264 target%id,
target%level,
target%llc_cart, &
265 source%id, source%level, source%llc_cart, &
266 source%distance(
target)/(0.5d0*(
target%size+source%size))
272 tlc =
target%myposition (
target%mesh%lb, iv)
273 tuc =
target%myposition (
target%mesh%ub, iv)
279 sl = source%index_only (tlc, 1, nowrap=.true.)
280 su = source%index_only (tuc, 1, nowrap=.true.)
281 sl = min(max(sl,source%mesh%li),source%mesh%ui-1)
282 su = min(max(su,source%mesh%li),source%mesh%ui-1)
287 slc = source%myposition(sl, iv)
288 suc = source%myposition(su, iv)
298 tl =
target%index_only (slc, 1, roundup=.true., nowrap=.true.)
299 tu =
target%index_only (suc, 1, roundup=.true., nowrap=.true.)
304 tl = min(max(tl,
target%mesh%lb),
target%mesh%ub)
305 tu = min(max(tu,
target%mesh%lb),
target%mesh%ub)
306 if (verbose > 0)
then 307 write(io_unit%log,1) &
308 'guard_zones_t%restrict_vect: iv, tlc, sl, slc, tl =', &
310 1
format(a,i3,2(2(2x,3f8.4,2x,3i4)))
311 write(io_unit%log,1) &
312 'guard_zones_t%restrict_vect: iv, tuc, su, suc, tu =', &
319 call source%index_space (
target%myposition(tl), 1, j0, p)
322 p3 = merge(p(3), 1.0-p(3), i3==1)
324 p2 = merge(p(2), 1.0-p(2), i2==1)
326 p1 = merge(p(1), 1.0-p(1), i1==1)
334 write(io_unit%log,
'(a,i3,4(2x,3i4))') &
335 'guard_zones_t%restrict_vect: iv, tl, tu, j0 =', &
340 tmem =>
target%mem(:,:,:,iv,
target%it,1)
341 td =>
target%mem(:,:,:,
target%idx%d,
target%it,1)
343 smem => source%mem(:,:,:,iv,jt(k),1)
344 sd => source%mem(:,:,:,source%idx%d,jt(k),1)
345 call update (
target, tmem, td, tl, tu, iv, source, smem, sd, j0, k, &
349 if (detailed_timer)
then 350 used = wallclock()-used
351 write (stdout,
'(a,6p,f6.3," mus/cell")') &
352 "guard_zone_t%restrict_vect:", used/product(tu-tl+1)
354 call trace%end (itimer)
355 END SUBROUTINE restrict_vect
370 SUBROUTINE update (target, tmem, td, l, u, iv, source, smem, sd, j0, k, &
372 class(patch_t),
pointer::
target, source
373 real,
dimension(:,:,:),
contiguous,
intent(out):: tmem
374 real,
dimension(:,:,:),
contiguous,
intent(in):: smem
375 real,
dimension(:,:,:),
contiguous,
intent(inout):: td, sd
376 integer:: l(3), u(3), iv, j0(3), k, verbose
380 real,
allocatable:: f(:,:), tot(:)
381 integer:: i, j, n, i1 , i2, i3, j1, j2, j3, jj(3)
383 call trace%begin (
'guard_zones_t%update')
384 if (detailed_timer)
then 388 allocate (f(n,8), tot(n))
398 if (
target%pervolume(iv))
then 404 if (any(jj < source%li .or. jj > source%ui))
then 405 write(stderr ,1) iv, jj, j0, l, u
406 write(io_unit%log,1) iv, jj, j0, l, u
407 1
format(
"guard_zones_t%update: ERROR iv,jj,j0,l,u =",i3,4(2x,3i4))
410 f(i,1) = smem(j1 ,j2 ,j3 )/sd(j1 ,j2 ,j3 )
411 f(i,2) = smem(j1+1,j2 ,j3 )/sd(j1+1,j2 ,j3 )
412 f(i,3) = smem(j1 ,j2+1,j3 )/sd(j1 ,j2+1,j3 )
413 f(i,4) = smem(j1+1,j2+1,j3 )/sd(j1+1,j2+1,j3 )
414 f(i,5) = smem(j1 ,j2 ,j3+1)/sd(j1 ,j2 ,j3+1)
415 f(i,6) = smem(j1+1,j2 ,j3+1)/sd(j1+1,j2 ,j3+1)
416 f(i,7) = smem(j1 ,j2+1,j3+1)/sd(j1 ,j2+1,j3+1)
417 f(i,8) = smem(j1+1,j2+1,j3+1)/sd(j1+1,j2+1,j3+1)
425 if (any(jj < source%li .or. jj > source%ui))
then 426 write(stderr ,1) iv, jj, j0, l, u
427 write(io_unit%log,1) iv, jj, j0, l, u
430 f(i,1) = smem(j1 ,j2 ,j3 )
431 f(i,2) = smem(j1+1,j2 ,j3 )
432 f(i,3) = smem(j1 ,j2+1,j3 )
433 f(i,4) = smem(j1+1,j2+1,j3 )
434 f(i,5) = smem(j1 ,j2 ,j3+1)
435 f(i,6) = smem(j1+1,j2 ,j3+1)
436 f(i,7) = smem(j1 ,j2+1,j3+1)
437 f(i,8) = smem(j1+1,j2+1,j3+1)
445 if (detailed_timer)
then 452 if (
target%unsigned(iv))
then 455 tot(i) = tot(i) + pt*w(j)*log(f(i,j))
462 tot(i) = tot(i) + pt*w(j)*f(i,j)
466 if (detailed_timer)
then 478 tmem(i1,i2,i3) = tot(i)
483 tmem(i1,i2,i3) = tmem(i1,i2,i3) + tot(i)
491 if (
target%pervolume(iv) .and. k==2)
then 495 tmem(i1,i2,i3) = tmem(i1,i2,i3)*td(i1,i2,i3)
501 if (detailed_timer)
then 505 update_t(i) = update_t(i) + (wc(i+1)-wc(i))
509 END SUBROUTINE update
516 SUBROUTINE interpolate (source, target, ii, iv, jt, pt, pi, set)
517 class(patch_t) :: source
518 class(patch_t),
pointer :: target
519 integer :: ii(3), iv, jt(2)
521 real,
optional :: pi(:)
522 logical,
optional :: set
527 integer,
save :: itimer=0
529 if (
present(set))
then 530 pi(01) = pt(1)*(1.-p(3))*(1.-p(2))*(1.0-p(1))
531 pi(02) = pt(1)*(1.-p(3))*(1.-p(2))*( p(1))
532 pi(03) = pt(1)*(1.-p(3))*( p(2))*(1.0-p(1))
533 pi(04) = pt(1)*(1.-p(3))*( p(2))*( p(1))
534 pi(05) = pt(1)*( p(3))*(1.-p(2))*(1.0-p(1))
535 pi(06) = pt(1)*( p(3))*(1.-p(2))*( p(1))
536 pi(07) = pt(1)*( p(3))*( p(2))*(1.0-p(1))
537 pi(08) = pt(1)*( p(3))*( p(2))*( p(1))
538 pi(09) = pt(2)*(1.-p(3))*(1.-p(2))*(1.0-p(1))
539 pi(10) = pt(2)*(1.-p(3))*(1.-p(2))*( p(1))
540 pi(11) = pt(2)*(1.-p(3))*( p(2))*(1.0-p(1))
541 pi(12) = pt(2)*(1.-p(3))*( p(2))*( p(1))
542 pi(13) = pt(2)*( p(3))*(1.-p(2))*(1.0-p(1))
543 pi(14) = pt(2)*( p(3))*(1.-p(2))*( p(1))
544 pi(15) = pt(2)*( p(3))*( p(2))*(1.0-p(1))
545 pi(16) = pt(2)*( p(3))*( p(2))*( p(1))
546 else if (source%pervolume(iv))
then 547 call interpolate_specific (source,
target, ii, iv, jt, pt)
549 pos =
target%myposition (ii, iv)
550 call source%index_space (pos, iv, j, p)
555 if (source%unsigned(iv))
then 556 if (
present(pi))
then 557 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) = exp( &
558 pi(01)*log(source%mem(j(1) ,j(2) ,j(3) ,iv,jt(1),1)) + &
559 pi(02)*log(source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(1),1)) + &
560 pi(03)*log(source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(1),1)) + &
561 pi(04)*log(source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(1),1)) + &
562 pi(05)*log(source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(1),1)) + &
563 pi(06)*log(source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(1),1)) + &
564 pi(07)*log(source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(1),1)) + &
565 pi(08)*log(source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(1),1)) + &
566 pi(09)*log(source%mem(j(1) ,j(2) ,j(3) ,iv,jt(2),1)) + &
567 pi(10)*log(source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(2),1)) + &
568 pi(11)*log(source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(2),1)) + &
569 pi(12)*log(source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(2),1)) + &
570 pi(13)*log(source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(2),1)) + &
571 pi(14)*log(source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(2),1)) + &
572 pi(15)*log(source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(2),1)) + &
573 pi(16)*log(source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(2),1)))
575 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) = exp( &
576 pt(1)*((1.-p(3))*((1.-p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2) ,j(3) ,iv,jt(1),1)) + &
577 ( p(1))*log(source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(1),1))) + &
578 ( p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(1),1)) + &
579 ( p(1))*log(source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(1),1)))) + &
580 ( p(3))*((1.-p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(1),1)) + &
581 ( p(1))*log(source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(1),1))) + &
582 ( p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(1),1)) + &
583 ( p(1))*log(source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(1),1))))) + &
584 pt(2)*((1.-p(3))*((1.-p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2) ,j(3) ,iv,jt(2),1)) + &
585 ( p(1))*log(source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(2),1))) + &
586 ( p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(2),1)) + &
587 ( p(1))*log(source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(2),1)))) + &
588 ( p(3))*((1.-p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(2),1)) + &
589 ( p(1))*log(source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(2),1))) + &
590 ( p(2))*((1.0-p(1))*log(source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(2),1)) + &
591 ( p(1))*log(source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(2),1))))))
594 if (
present(pi))
then 595 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) = &
596 pi(01)*source%mem(j(1) ,j(2) ,j(3) ,iv,jt(1),1) + &
597 pi(02)*source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(1),1) + &
598 pi(03)*source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(1),1) + &
599 pi(04)*source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(1),1) + &
600 pi(05)*source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(1),1) + &
601 pi(06)*source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(1),1) + &
602 pi(07)*source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(1),1) + &
603 pi(08)*source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(1),1) + &
604 pi(09)*source%mem(j(1) ,j(2) ,j(3) ,iv,jt(2),1) + &
605 pi(10)*source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(2),1) + &
606 pi(11)*source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(2),1) + &
607 pi(12)*source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(2),1) + &
608 pi(13)*source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(2),1) + &
609 pi(14)*source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(2),1) + &
610 pi(15)*source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(2),1) + &
611 pi(16)*source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(2),1)
613 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) = &
614 pt(1)*((1.-p(3))*((1.-p(2))*((1.0-p(1))*source%mem(j(1) ,j(2) ,j(3) ,iv,jt(1),1) + &
615 ( p(1))*source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(1),1)) + &
616 ( p(2))*((1.0-p(1))*source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(1),1) + &
617 ( p(1))*source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(1),1))) + &
618 ( p(3))*((1.-p(2))*((1.0-p(1))*source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(1),1) + &
619 ( p(1))*source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(1),1)) + &
620 ( p(2))*((1.0-p(1))*source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(1),1) + &
621 ( p(1))*source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(1),1)))) + &
622 pt(2)*((1.-p(3))*((1.-p(2))*((1.0-p(1))*source%mem(j(1) ,j(2) ,j(3) ,iv,jt(2),1) + &
623 ( p(1))*source%mem(j(1)+1,j(2) ,j(3) ,iv,jt(2),1)) + &
624 ( p(2))*((1.0-p(1))*source%mem(j(1) ,j(2)+1,j(3) ,iv,jt(2),1) + &
625 ( p(1))*source%mem(j(1)+1,j(2)+1,j(3) ,iv,jt(2),1))) + &
626 ( p(3))*((1.-p(2))*((1.0-p(1))*source%mem(j(1) ,j(2) ,j(3)+1,iv,jt(2),1) + &
627 ( p(1))*source%mem(j(1)+1,j(2) ,j(3)+1,iv,jt(2),1)) + &
628 ( p(2))*((1.0-p(1))*source%mem(j(1) ,j(2)+1,j(3)+1,iv,jt(2),1) + &
629 ( p(1))*source%mem(j(1)+1,j(2)+1,j(3)+1,iv,jt(2),1))))
633 if (
target%id==io%id_debug .and. io%verbose>2) &
634 write(io%output,1)
target%id, ii, iv, source%id, j, jt, p, pt, &
635 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1)
636 1
format(i6,1x,3i3,i4,2x,
"DBG patch_t%interpolate src:",i6, &
637 " j,p;",3i4,1x,2i4,3f6.2,1x,2f6.2,
" out:",1p,e15.6)
639 END SUBROUTINE interpolate
646 SUBROUTINE interpolate_specific (source, target, ii, iv, jt, pt)
647 class(patch_t) :: source
648 class(patch_t),
pointer :: target
649 integer :: ii(3), iv, jt(2)
656 integer,
save :: itimer=0
658 pos =
target%myposition (ii, iv)
659 call source%index_space (pos, iv, j, p)
664 associate(id=>source%idx%d)
665 tmp(:,:,:,:) = source%mem(j(1):j(1)+1,j(2):j(2)+1,j(3):j(3)+1,iv,jt(1:2),1) &
666 / source%mem(j(1):j(1)+1,j(2):j(2)+1,j(3):j(3)+1,id,jt(1:2),1)
667 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) = &
668 pt(1)*((1.-p(3))*((1.-p(2))*((1.0-p(1))*tmp(1,1,1,1) + &
669 ( p(1))*tmp(2,1,1,1)) + &
670 ( p(2))*((1.0-p(1))*tmp(1,2,1,1) + &
671 ( p(1))*tmp(2,2,1,1))) + &
672 ( p(3))*((1.-p(2))*((1.0-p(1))*tmp(1,1,2,1) + &
673 ( p(1))*tmp(2,1,2,1)) + &
674 ( p(2))*((1.0-p(1))*tmp(1,2,2,1) + &
675 ( p(1))*tmp(2,2,2,1)))) + &
676 pt(2)*((1.-p(3))*((1.-p(2))*((1.0-p(1))*tmp(1,1,1,2) + &
677 ( p(1))*tmp(2,1,1,2)) + &
678 ( p(2))*((1.0-p(1))*tmp(1,2,1,2) + &
679 ( p(1))*tmp(2,2,1,2))) + &
680 ( p(3))*((1.-p(2))*((1.0-p(1))*tmp(1,1,2,2) + &
681 ( p(1))*tmp(2,1,2,2)) + &
682 ( p(2))*((1.0-p(1))*tmp(1,2,2,2) + &
683 ( p(1))*tmp(2,2,2,2))))
687 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) = &
688 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1) * &
689 target%mem(ii(1),ii(2),ii(3),id,
target%it,1)
691 if (
target%id==io%id_debug .and. io%verbose>2) &
692 write(io%output,1)
target%id, ii, iv, source%id, j, jt, p, pt, &
693 target%mem(ii(1),ii(2),ii(3),iv,
target%it,1)
694 1
format(i6,1x,3i3,i4,2x,
"DBG patch_t%interpolate_specific src:", &
695 i6,
" j,p;",3i4,1x,2i4,3f6.2,1x,2f6.2,
" out:",1p,e15.6)
697 END SUBROUTINE interpolate_specific
700 if (verbose > 0)
then 701 write (stdout,*)
'guard_zones_t%info: update-t =', update_t
Support tic/toc timing, as in MATLAB, and accurate wallclock() function. The timing is generally much...
Template module for patches, which adds pointers to memory and mesh, and number of dimensions and var...
Template module for mesh.
Optimized restrict operations, intended for guard zones (not conservative)