Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:38:29 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
bf4be02920 Jean*0001 !=======================================================================
                0002       subroutine coupsend_r8( component, dataname, Nx, Ny, arr )
                0003       implicit none
fbcf25275b Jean*0004 ! Predefined constants/arrays
                0005 #include "CPLR_SIG.h"
                0006 ! MPI variables
                0007 #include "mpif.h"
bf4be02920 Jean*0008 ! Arguments
                0009       character*(*) component
                0010       character*(*) dataname
                0011       integer Nx,Ny
                0012       real*8 arr(Nx,Ny)
                0013 ! Functions
                0014       integer mitcplr_match_comp
                0015       integer generate_tag
fbcf25275b Jean*0016       external mitcplr_match_comp
                0017       external generate_tag
bf4be02920 Jean*0018 ! Local
fbcf25275b Jean*0019       integer count,dtype,dest,tag,comm,ierr
bf4be02920 Jean*0020       integer compind,numprocs
                0021       integer i,j,ij,n
                0022       integer Ni,Io,Nj,Jo
                0023 !     ------------------------------------------------------------------
                0024 
                0025 ! Establish who I am communicating with
                0026       compind=mitcplr_match_comp( component )
                0027       if (compind.le.0) stop 'coupsend_r8: Bad component id'
                0028       comm=MPI_COMM_compcplr( compind )
                0029       numprocs=num_component_procs(compind)
                0030       if (numprocs.lt.1) then
                0031        write(LogUnit,*) 'coupsend_r8: compind = ',compind
                0032        stop 'coupsend_r8: numprocs < 1'
                0033       endif
                0034       if (VERB)
                0035      &  write(LogUnit,*) 'coupsend_r8: ',component_Name(compind)
                0036       if (VERB)
                0037      &  write(LogUnit,*) 'coupsend_r8: dataname=',dataname
                0038 
                0039 ! Foreach component process
                0040       do n=1,numprocs
                0041 
                0042 ! Create header
                0043        Io=component_tile_i0(1,n,compind)
                0044        Jo=component_tile_j0(1,n,compind)
                0045        Ni=component_tile_nx(1,n,compind)
                0046        Nj=component_tile_ny(1,n,compind)
                0047        r8buf(1)=float( Io )
                0048        r8buf(2)=float( Jo )
                0049        r8buf(3)=float( Ni )
                0050        r8buf(4)=float( Nj )
e5266ce10c Jean*0051        call mitcplr_char2dbl( dataname, r8buf(9) )
bf4be02920 Jean*0052 
                0053 ! Pack data
                0054        do j=1,Nj
                0055         do i=1,Ni
                0056          ij=HEADER_SIZE+i+Ni*(j-1)
                0057          r8buf(ij)=arr(Io+i-1,Jo+j-1)
                0058         enddo
                0059        enddo
                0060 
                0061 ! Send message
                0062        count=HEADER_SIZE+Ni*Nj
                0063        dtype=MPI_DOUBLE_PRECISION
                0064        tag=generate_tag(121,n,dataname)
                0065        dest=rank_component_procs(n,compind)
                0066 
                0067        if (VERB) then
                0068         write(LogUnit,*) 'coupsend_r8: calling MPI_Send dest=',dest,
                0069      &    ' proc=',n,'/',numprocs
                0070         call flush(LogUnit)
                0071        endif
                0072        call MPI_Send( r8buf, count, dtype, dest, tag, comm, ierr )
                0073        if (VERB) then
                0074         write(LogUnit,*) 'coupsend_r8: returned ierr=',ierr
                0075         call flush(LogUnit)
                0076        endif
                0077 
                0078        if (ierr.ne.0) then
                0079         write(LogUnit,*) 'coupsend_r8tiles: rank(W,G)=',
                0080      &            my_rank_in_world,my_rank_in_global,
                0081      &            ' ierr=',ierr
fbcf25275b Jean*0082         stop 'coupsend_r8: MPI_Send failed'
bf4be02920 Jean*0083        endif
                0084 
                0085       enddo ! n
                0086 
                0087 !     ------------------------------------------------------------------
                0088       return
                0089       end
                0090 !=======================================================================