module funnel ! ! highly specific routines for serial i/o on parallel computer: ! to "funnel" arrays from multiple processors to a single processor ! ! (c) Copyright 1991 to 1998 by Michael A. Beer, William D. Dorland, ! P. B. Snyder, Q. P. Liu, and Gregory W. Hammett. ALL RIGHTS RESERVED. ! use mp, only: barrier, send, receive, proc0 use gryffin_layouts use itg_data, only: nspecies use convert implicit none private ! I would call the subroutine just "funnel", but apparently it has to be ! different than the module name?! public :: funnel20 ! "funnel data 2 processor 0" interface funnel20 module procedure funnel4 module procedure funnel4cr module procedure rfunnel4 module procedure funnel3 module procedure funnel3cr module procedure rfunnel3 module procedure funnel2 module procedure funnel2cr module procedure rfunnel2 module procedure funnel1 end interface public :: mfunnel3 contains subroutine funnel4(a, b) ! funnel a 4-D complex distributed array "a" to a ! 4-D array "b" on a single processor (processor 0) complex, dimension(:,m_low:,n_low:,:) :: a complex, dimension(:,:,:,:) :: b integer :: i, m, n do i=1,nspecies do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(:,m,n,i) = a(:,m,n,i) else call receive(b(:,m,n,i), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(:,m,n,i), 0) endif call barrier enddo enddo enddo end subroutine funnel4 subroutine funnel4cr(a, ri) ! funnel a 4-D complex distributed array "a" to a ! 4-D (plus real-imaginary) array "ri" on a single processor complex, dimension(:,m_low:,n_low:,:) :: a real, dimension(:,:,:,:,:) :: ri complex, dimension(size(ri,2),size(ri,3),size(ri,4),size(ri,5)) :: b integer :: i, m, n do i=1,nspecies do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(:,m,n,i) = a(:,m,n,i) else call receive(b(:,m,n,i), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(:,m,n,i), 0) endif call barrier enddo enddo enddo if(proc0) call c2r(b,ri) end subroutine funnel4cr subroutine rfunnel4(a, b, lmin, lmax) real, dimension(:,m_low:,n_low:,:) :: a real, dimension(:,:,:,:) :: b integer :: lmin, lmax integer :: i, m, n do i=1,nspecies do n = m_lo%nl_world, m_lo%nu_world do m = m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(lmin:lmax,m,n,i) = a(lmin:lmax,m,n,i) else call receive(b(lmin:lmax,m,n,i), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(lmin:lmax,m,n,i), 0) endif call barrier enddo enddo enddo end subroutine rfunnel4 subroutine funnel3(a, b, lmin, lmax) complex, dimension(:,m_low:,n_low:) :: a complex, dimension(:,:,:) :: b integer :: lmin, lmax integer :: m, n do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(lmin:lmax,m,n) = a(lmin:lmax,m,n) else call receive(b(lmin:lmax,m,n), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(lmin:lmax,m,n), 0) endif call barrier enddo enddo end subroutine funnel3 subroutine funnel3cr(a, ri, lmin, lmax) complex, dimension(:,m_low:,n_low:) :: a real, dimension(:,:,:,:) :: ri complex, dimension(size(ri,2),size(ri,3),size(ri,4)) :: b integer :: lmin, lmax integer :: m, n do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(lmin:lmax,m,n) = a(lmin:lmax,m,n) else call receive(b(lmin:lmax,m,n), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(lmin:lmax,m,n), 0) endif call barrier enddo enddo if(proc0) call c2r(b,ri) end subroutine funnel3cr subroutine rfunnel3(a, b, lmin, lmax) real, dimension(:,m_low:,n_low:) :: a real, dimension(:,:,:) :: b integer :: lmin, lmax integer :: m, n do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(lmin:lmax,m,n) = a(lmin:lmax,m,n) else call receive(b(lmin:lmax,m,n), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(lmin:lmax,m,n), 0) endif call barrier enddo enddo end subroutine rfunnel3 subroutine mfunnel3(a, b, n1, n2) real, dimension(m_low:,n_low:,:) :: a real, dimension(:,:,:) :: b integer :: n1, n2 integer :: i, m, n do i=n1,n2 do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(m,n,i) = a(m,n,i) else call receive(b(m,n,i), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(m,n,i), 0) endif call barrier enddo enddo enddo end subroutine mfunnel3 subroutine funnel2(a, b) complex, dimension(m_low:,n_low:) :: a complex, dimension(:,:) :: b integer :: m, n do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(m,n) = a(m,n) else call receive(b(m,n), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(m,n), 0) endif call barrier enddo enddo end subroutine funnel2 subroutine funnel2cr(a, ri) complex, dimension(m_low:,n_low:) :: a real, dimension(:,:,:) :: ri complex, dimension(size(ri,2),size(ri,3)) :: b integer :: m, n do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(m,n) = a(m,n) else call receive(b(m,n), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(m,n), 0) endif call barrier enddo enddo if(proc0) call c2r(b,ri) end subroutine funnel2cr subroutine rfunnel2(a, b) real, dimension(m_low:,n_low:) :: a real, dimension(:,:) :: b integer :: m, n do n=m_lo%nl_world, m_lo%nu_world do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, n)) then b(m,n) = a(m,n) else call receive(b(m,n), proc_id(m_lo, m, n)) endif elseif (idx_local(m_lo, m, n)) then call send(a(m,n), 0) endif call barrier enddo enddo end subroutine rfunnel2 subroutine funnel1(a, b) real, dimension(m_low:) :: a real, dimension(:) :: b integer :: m do m=m_lo%ml_world, m_lo%mu_world if(proc0) then if(idx_local(m_lo, m, 1)) then b(m) = a(m) else call receive(b(m), proc_id(m_lo, m, 1)) endif elseif (idx_local(m_lo, m, 1)) then call send(a(m), 0) endif call barrier enddo end subroutine funnel1 end module funnel