** Warning **
Issuing rollback() due to DESTROY without explicit disconnect() of DBD::mysql::db handle dbname=MITgcm at /usr/local/share/lxr/lib/LXR/Common.pm line 1224.
Last-Modified: Wed, 5 Nov 2024 06:12:11 GMT
Content-Type: text/html; charset=utf-8
MITgcm/MITgcm/pkg/compon_communic/couprecv_r8tiles.F
File indexing completed on 2018-03-02 18:38:29 UTC
view on github raw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
bf4be02920 Jean* 0001
0002 subroutine couprecv_r8tiles ( component , dataname , Nx , Ny , arr )
0003 implicit none
0004
0005 character *(*) component
0006 character *(*) dataname
0007 integer Nx ,Ny
0008 real *8 arr (Nx ,Ny )
0009
0010 #include "CPLR_SIG.h "
0011
0012 #include "mpif.h"
0013 integer count ,dtype ,rank ,tag ,comm ,ierr
0014 integer stat (MPI_STATUS_SIZE)
0015
0016 integer mitcplr_match_comp
0017 integer generate_tag
0018
0019 integer compind ,numprocs
0020 integer i ,j ,ij ,n ,bibj
0021 integer Ni ,Io ,Nj ,Jo
0022 character *(MAXLEN_COMP_NAME ) recvdname
0023
0024
0025
0026 compind =mitcplr_match_comp ( component )
0027 if (compind .le. 0) stop 'couprecv_r8tiles: Bad component id'
0028 comm =MPI_COMM_compcplr ( compind )
0029 numprocs =num_component_procs (compind )
0030 if (numprocs .lt. 1) then
0031 write (LogUnit ,*) 'couprecv_r8tiles: compind = ' ,compind
0032 stop 'couprecv_r8tiles: numprocs < 1'
0033 endif
0034 if (VERB )
0035 & write (LogUnit ,*) 'couprecv_r8tiles: ' ,component_Name (compind )
0036 if (VERB )
0037 & write (LogUnit ,*) 'couprecv_r8tiles: dataname=' ,dataname
0038
0039
0040 do n =1,numprocs
0041
0042
0043 do bibj =1,component_num_tiles (n ,compind )
0044
0045
0046 count =MAX_R8_BUFLEN
0047 dtype =MPI_DOUBLE_PRECISION
0048 tag =generate_tag (103,bibj ,dataname )
0049 rank =rank_component_procs (n ,compind )
0050
0051 if (VERB ) then
0052 write (LogUnit ,*)
0053 & 'couprecv_r8tiles: calling MPI_Recv rank=' ,rank ,
0054 & ' proc=' ,n ,'/' ,numprocs ,' tile=' ,bibj
0055 call flush(LogUnit )
0056 endif
0057 call MPI_Recv(r8buf , count , dtype , rank , tag , comm , stat , ierr )
0058 if (VERB ) then
0059 write (LogUnit ,*) 'couprecv_r8tiles: returned ierr=' ,ierr
0060 call flush(LogUnit )
0061 endif
0062
0063 if (ierr .ne. 0) then
0064 write (LogUnit ,*) 'couprecv_r8tiles: rank(W,G)=' ,
0065 & my_rank_in_world ,my_rank_in_global ,
0066 & ' ierr=' ,ierr
0067 stop 'couprecv_r8tiles: MPI_Recv failed'
0068 endif
0069
0070
0071 Io =int(0.5+r8buf (1))
0072 Ni =int(0.5+r8buf (2))
0073 Jo =int(0.5+r8buf (3))
0074 Nj =int(0.5+r8buf (4))
0075
0076 if (Io +Ni -1.gt. Nx .or. Io .lt. 1) then
0077 write (LogUnit ,*) 'couprecv_r8tiles: Io,Ni=' ,Io ,Ni
0078 stop 'couprecv_r8tiles: Incompatible header/target array'
0079 endif
0080 if (Jo +Nj -1.gt. Ny .or. Jo .lt. 1) then
0081 write (LogUnit ,*) 'couprecv_r8tiles: Jo,Nj=' ,Jo ,Nj
0082 stop 'couprecv_r8tiles: Incompatible header/target array'
0083 endif
0084 if (Io .ne. component_tile_i0 (bibj ,n ,compind ))
0085 & stop 'couprecv_r8tiles: Io != component_tile_i0'
0086 if (Jo .ne. component_tile_j0 (bibj ,n ,compind ))
0087 & stop 'couprecv_r8tiles: Jo != component_tile_j0'
0088 if (Ni .ne. component_tile_nx (bibj ,n ,compind ))
0089 & stop 'couprecv_r8tiles: Ni != component_tile_nx'
0090 if (Nj .ne. component_tile_ny (bibj ,n ,compind ))
0091 & stop 'couprecv_r8tiles: Nj != component_tile_ny'
0092
e5266ce10c Jean* 0093 call mitcplr_dbl2char ( r8buf (9), recvdname )
bf4be02920 Jean* 0094 if (recvdname .ne. dataname ) then
0095 write (LogUnit ,*) 'couprecv_r8tiles: recvdname = ' ,recvdname
0096 write (LogUnit ,*) 'couprecv_r8tiles: dataname = ' ,dataname
0097 stop 'couprecv_r8tiles: recvdname != dataname'
0098 endif
0099
0100
0101 do j =1,Nj
0102 do i =1,Ni
0103 ij =HEADER_SIZE +i +Ni *(j -1)
0104 arr (Io +i -1,Jo +j -1)=r8buf (ij )
0105 enddo
0106 enddo
0107
0108 enddo
0109
0110 enddo
0111
0112
0113 return
0114 end
0115