* PACKAGE PGRIV !! river routing submodel * [HIS] 94/10 by kanae, 95/11 modified * 98/03/07(numaguti) AGCM5.5b * 99/08/11(hiro) parallel AGCM5.5b * 99/11/04(fuyuki) for more #PE than set when compiled * 06/11/24(kei) add OPT_TRIP ********************************************************************* SUBROUTINE RIVER !! river routing * #include "zcdim.F" /* # of grid point & wave */ #include "zpdim.F" /* physics etc. */ #include "zhdim.F" /* # of char. */ #include "zvdim.F" /* # of river grid point */ #include "zcorv.F" /* river coordinate */ #include "zccom.F" /* stand. physical const. */ * * [MODIFY] REAL * 8 TIME !! time * [INPUT] REAL * 8 TEND !! time for end LOGICAL OEND !! end of run? * * [INTERNAL SAVE] REAL*8 GDRIV ( IJVDIM, KRVMAX ) !! storage in river(kg/m**2) REAL*8 RDEST ( IJVDIM, KRVMAX ) !! destination speed/distance REAL*8 AREA ( IJVDIM ) !! grid area INTEGER IJDEST ( IJVDIM, KRVMAX ) !! Flow destination IJ INTEGER JDESTN ( IVDIM, KRVMAX ) !! Flow dest. at N-ridge INTEGER JDESTS ( IVDIM, KRVMAX ) !! Flow dest. at S-ridge INTEGER MDESTN ( KRVMAX ) !! No. flow dest. at N-ridge INTEGER MDESTS ( KRVMAX ) !! No. flow dest. at S-ridge INTEGER NDESTN ( KRVMAX ) !! No. flow dest. at N-ridge INTEGER NDESTS ( KRVMAX ) !! No. flow dest. at S-ridge REAL * 8 TSTART !! start time INTEGER NPROCS, MYRANK INTEGER IFPAR, JFPAR INTEGER ISTA , IEND INTEGER JSTA , JEND INTEGER ISTAL, IENDL LOGICAL OPCALC !! this rank run or not SAVE GDRIV, TSTART SAVE IJDEST, JDESTN, JDESTS SAVE MDESTN, MDESTS, NDESTN, NDESTS, RDEST, AREA SAVE NPROCS, MYRANK SAVE IFPAR, JFPAR SAVE ISTA , IEND SAVE ISTAL, IENDL SAVE OPCALC * * [INTERNAL WORK] #ifndef OPT_NOPHYSICS REAL*8 RUNOF ( IJVDIM, KRVMAX ) !! runoff REAL*8 RRUNOF ( IJVDIM, KRVMAX ) !! river runoff (kei) REAL*8 RFLOW ( IJVDIM, KRVMAX ) !! river flow (kg/sec) REAL*8 GSRIVN ( IVDIM, KRVMAX ) !! storage at N-ridge (send) REAL*8 GSRIVS ( IVDIM, KRVMAX ) !! storage at S-ridge (send) REAL*8 GRRIVN ( IVDIM, KRVMAX ) !! storage at N-ridge (recv) REAL*8 GRRIVS ( IVDIM, KRVMAX ) !! storage at S-ridge (recv) REAL*8 DELT !! time step REAL*8 DTIV INTEGER ID #endif INTEGER IJ, J, IJS, IJE, JN, JS, IERROR, IR INTEGER IMAP ( IVDIM*JVMAX ) !! flow direction (1 to 8) INTEGER IDESTN ( IVDIM ) !! Flow dest. at N-ridge INTEGER IDESTS ( IVDIM ) !! Flow dest. at S-ridge CHARACTER HTIME*20 REAL * 8 TINIT0 !! time of initial data REAL * 8 TINIT1 !! time of initial data REAL * 8 TINITD REAL * 8 TDURD REAL*8 PI CHARACTER HDFMT *(NCC) INTEGER IFLSTR, IEOD, KLEVS, JFLRST #ifdef OPT_TRIP REAL*8 GDRIVO ( IJVDIM, KRVMAX ) !! REAL*8 INFLW ( IJVDIM, KRVMAX ) !! river inflow (kg/sec) REAL*8 GDRIVALL,GDRIVALLO,RFLOWALL,INFLWALL,RIVALL,RUNOFALL REAL*8 AREAALL,RBUD #endif * * [INTERNAL PARM] REAL*8 VRIVER ( KRVMAX ) !! river flow speed (m/s) CHARACTER RIVMAP ( KRVMAX )*(NFILN) !! map file INTEGER IFILE LOGICAL OFIRST NAMELIST /NMRIVR/ VRIVER, RIVMAP, IFILE DATA VRIVER / KRVMAX*0.3 / DATA RIVMAP / KRVMAX*' ' / DATA IFILE / 96 / DATA OFIRST / .TRUE. / CHARACTER HCLAS*3 DATA HCLAS / 'RIV' / CHARACTER C3*3 * *========================================================================= ENTRY RIVSTP !! river step M ( TIME , I TEND ) * IF( .NOT. OPCALC ) THEN CALL CSS2CC ( HTIME, TIME ) WRITE (JFPAR,*) ' *** RIVER STEP TIME=', HTIME CALL SETTIM ( TIME ) TIME = TEND CALL SETTIM ( TIME ) RETURN ENDIF * #ifndef OPT_NOPHYSICS * CALL CSS2CC ( HTIME, TIME ) #ifdef OPT_MPE CALL MPE_Log_event(900, 0, 'start LND+RIV+OCN') #endif WRITE (JFPAR,*) ' *** RIVER STEP TIME=', HTIME CALL SETTIM ( TIME ) CALL HISTRT ( 4 ) DELT = TEND - TIME DTIV = 1.D0/DELT * CALL C2V( RUNOF, 'ROFF', KRVMAX, 0.D0 ) ** kei ** #ifdef OPT_TRIP CALL COPY (GDRIVO, GDRIV, IJVDIM*KRVMAX) CALL RESET (INFLW,IJVDIM*KRVMAX) #endif ** kei ** * DO 2000 IR = 1, KRVMAX #ifndef OPT_TRIP DO 2100 IJ = ISTA, IEND GDRIV( IJ,IR ) = GDRIV( IJ,IR ) + RUNOF( IJ,IR )*DELT RFLOW( IJ,IR ) = MAX( GDRIV( IJ,IR ),0.D0 )*AREA( IJ ) & * MIN( RDEST( IJ,IR ),DTIV ) 2100 CONTINUE JN = 0 JS = 0 DO 2200 IJ = ISTA, IEND ID = IJDEST( IJ,IR ) IF ( ID-ISTAL+1 .NE. IJ ) THEN GDRIV( IJ,IR ) = GDRIV( IJ,IR ) & - RFLOW( IJ,IR )/AREA( IJ )*DELT IF( ID .LT. ISTAL ) THEN JN = JN + 1 GSRIVN( JN, IR ) = RFLOW( IJ,IR ) ELSEIF( ID .GT. IENDL ) THEN JS = JS + 1 GSRIVS( JS, IR ) = RFLOW( IJ,IR ) ELSE GDRIV( ID-ISTAL+1,IR ) = GDRIV( ID-ISTAL+1,IR ) & + RFLOW( IJ,IR )/AREA( ID-ISTAL+1 )*DELT ENDIF RRUNOF( IJ,IR ) = RFLOW( IJ,IR )/AREA( IJ ) ENDIF 2200 CONTINUE * #else !! TRIP JN = 0 JS = 0 DO 2200 IJ = ISTA, IEND INFLW(IJ,IR)=INFLW(IJ,IR)+RUNOF(IJ,IR) ID = IJDEST( IJ,IR ) IF ( IMAP(IJ+ISTAL-1).GE.1 .AND. RDEST(IJ,IR).NE.0. ) THEN GDRIV(IJ,IR) = GDRIV(IJ,IR)*exp(-RDEST(IJ,IR)*DELT) $ + (1.D0 - exp(-(RDEST(IJ,IR)*DELT))) $ * RUNOF(IJ,IR)/RDEST(IJ,IR) RFLOW( IJ,IR ) = ((GDRIVO(IJ,IR)-GDRIV(IJ,IR))/DELT $ + INFLW(IJ,IR))*AREA(IJ) IF ( ID-ISTAL+1 .NE. IJ ) THEN IF( ID .LT. ISTAL ) THEN JN = JN + 1 GSRIVN( JN, IR ) = RFLOW( IJ,IR ) ELSEIF( ID .GT. IENDL ) THEN JS = JS + 1 GSRIVS( JS, IR ) = RFLOW( IJ,IR ) ELSE GDRIV( ID-ISTAL+1,IR ) = GDRIV( ID-ISTAL+1,IR ) & + RFLOW( IJ,IR )/AREA( ID-ISTAL+1 )*DELT INFLW(ID-ISTAL+1,IR)=INFLW(ID-ISTAL+1,IR) $ + RFLOW(IJ,IR)/AREA(ID-ISTAL+1) ENDIF ENDIF RRUNOF( IJ,IR ) = RFLOW( IJ,IR )/AREA( IJ ) ENDIF 2200 CONTINUE #endif * #ifndef OPT_NONPARALLEL * * communication at N/S boundary * CALL NSCOMMF O ( GRRIVS(1,IR) , GRRIVN(1,IR) , I GSRIVS(1,IR) , GSRIVN(1,IR) , I NDESTS(IR) , NDESTN(IR) , MDESTS(IR) , MDESTN(IR) , D IVDIM , IVDIM , IVDIM , IVDIM ) * IF( MYRANK .GT. 0 ) THEN DO 2400 IJ = 1, MDESTN(IR) ID = JDESTN( IJ, IR ) GDRIV( ID, IR ) = GDRIV( ID, IR ) & + GRRIVN( IJ,IR )/AREA( ID )*DELT #ifdef OPT_TRIP INFLW( ID, IR ) = INFLW( ID, IR ) & + GRRIVN( IJ,IR )/AREA( ID ) #endif 2400 CONTINUE ENDIF IF( MYRANK .LT. NPROCS-1 ) THEN DO 2500 IJ = 1, MDESTS(IR) ID = JDESTS( IJ, IR ) GDRIV( ID, IR ) = GDRIV( ID, IR ) & + GRRIVS( IJ,IR )/AREA( ID )*DELT #ifdef OPT_TRIP INFLW( ID, IR ) = INFLW( ID, IR ) & + GRRIVS( IJ,IR )/AREA( ID ) #endif 2500 CONTINUE ENDIF * #endif * DO 2300 IJ = ISTA, IEND ID = IJDEST( IJ,IR ) **kei IF ( ID-ISTAL+1 .EQ. IJ ) THEN IF ( ID-ISTAL+1 .EQ. IJ .AND. IMAP(IJ+ISTAL-1).EQ.0 ) THEN RRUNOF( IJ,IR ) = GDRIV( IJ,IR )/DELT RFLOW( IJ,IR ) = GDRIV( IJ,IR )*AREA( IJ )/DELT GDRIV( IJ,IR ) = 0.D0 ENDIF 2300 CONTINUE **kei #ifdef OPT_TRIP GDRIVALL=0.D0 GDRIVALLO=0.D0 RFLOWALL=0.D0 INFLWALL=0.D0 RIVALL=0.D0 RUNOFALL=0.D0 AREAALL=0.D0 DO 2350 IJ=ISTA,IEND IF ( IMAP(IJ+ISTAL-1).GE.1 ) THEN AREAALL=AREAALL+AREA(IJ) GDRIVALL=GDRIVALL+GDRIV(IJ,IR)*AREA(IJ) GDRIVALLO=GDRIVALLO+GDRIVO(IJ,IR)*AREA(IJ) RFLOWALL=RFLOWALL+RFLOW(IJ,IR)*DELT INFLWALL=INFLWALL+INFLW(IJ,IR)*AREA(IJ)*DELT RUNOFALL=RUNOFALL+RUNOF(IJ,IR)*AREA(IJ)*DELT IF ( IMAP(IJ+ISTAL-1).LE.8 ) THEN RIVALL=RIVALL+RFLOW(IJ,IR)*DELT ENDIF RBUD=GDRIV(IJ,IR)-GDRIVO(IJ,IR)- $ (INFLW(IJ,IR)-RRUNOF(IJ,IR))*DELT IF (ABS(RBUD).GT.0.1) THEN WRITE(JFPAR,*) 'Inbalance!',IJ,IMAP(IJ+ISTAL-1),RBUD ENDIF ENDIF 2350 CONTINUE WRITE(JFPAR,'(a20,2f20.7)') '*** PGRIV_TRIP BUD ', $ (GDRIVALL-GDRIVALLO-(INFLWALL-RFLOWALL))/AREAALL, $ (INFLWALL-RUNOFALL-RIVALL)/AREAALL WRITE(JFPAR,'(a20,e20.7)') '*** TRIP CHECK A ', $ AREAALL WRITE(JFPAR,'(a20,f20.7)') '*** TRIP CHECK S1 ', $ GDRIVALLO/AREAALL WRITE(JFPAR,'(a20,f20.7)') '*** TRIP CHECK S2 ', $ GDRIVALL/AREAALL WRITE(JFPAR,'(a20,f20.7)') '*** TRIP CHECK O ', $ RFLOWALL/AREAALL WRITE(JFPAR,'(a20,f20.7)') '*** TRIP CHECK I ', $ INFLWALL/AREAALL WRITE(JFPAR,'(a20,f20.7)') '*** TRIP CHECK R ', $ RUNOFALL/AREAALL #endif **kei * 2000 CONTINUE * **kei CALL V2C( RUNOF, 'ROF', KRVMAX, 0.D0 ) CALL V2C( RRUNOF, 'ROF', KRVMAX, DELT ) * CALL HISTIN(RFLOW ,'RFLOW' ,'river flow', 'kg/s', $ 'ASFC', HCLAS ) CALL HISTIN(GDRIV ,'GDRIV' ,'river water', 'kg/m**2', $ 'ASFC', HCLAS ) CALL HISTIN(RRUNOF ,'RRUNOF','river runoff','kg/m**2/s', $ 'ASFC', HCLAS ) **kei** #ifdef OPT_TRIP CALL HISTIN(INFLW ,'INFLW','river runoff','kg/m**2/s', $ 'ASFC', HCLAS ) #endif **kei** * #endif * TIME = TEND CALL SETTIM( TIME ) CALL HISTOU( .FALSE. ) * #ifdef OPT_MPE CALL MPE_Log_event(901, 0, 'end LND+RIV+OCN') #endif * RETURN *========================================================================= ENTRY RIVINI !! initialize river I ( TIME , TEND ) * CALL GETMPI O ( NPROCS, MYRANK, ISTA, IEND, IFPAR, JFPAR, I 'RIV' ) CALL GETMPI2 O ( NPROCS, MYRANK, JSTA, JEND, IFPAR, JFPAR, I 'RIV' ) * CALL GETOPCALC ( OPCALC ) IF( .NOT. OPCALC ) THEN GOTO 900 ENDIF * ISTAL = ( JSTA-1 ) * IVDIM + 1 IENDL = JEND * IVDIM OFIRST = .FALSE. WRITE ( JFPAR,* ) ' @@@ PGRIV: RIVER FLOW 99/04/05' CALL REWNML( IFPAR, JFPAR ) READ ( IFPAR, NMRIVR, END=190 ) 190 WRITE ( JFPAR, NMRIVR ) * CALL SETLOM O ( AVLON , DVLON, HVLON , O AVLAT , DVLAT, HVLAT , I IVDIV , JVDIV ) PI = 4.D0 * ATAN( 1.D0 ) DO 400 IJ = ISTA, IEND #ifndef OPT_LINGRD AREA( IJ ) = 4.D0*ER*ER*PI*DVLON( IJ )*DVLAT( IJ ) #else AREA( IJ ) = 2.D0*ER*ER*PI*PI*DVLON( IJ )*DVLAT( IJ ) & *COS( AVLAT( IJ ) ) #endif 400 CONTINUE * WRITE(C3,'(i3)') IVMAX DO 500 IR = 1, KRVMAX CALL IRESET( IMAP, IVDIM*JVMAX ) IF ( RIVMAP(IR) .NE. ' ' ) THEN CALL MKFILN( RIVMAP , JFPAR ) CALL FOPEN O ( IERROR, I IFILE , RIVMAP, 'READ', 'FORMATTED','SEQUENTIAL' ) IF ( IERROR .EQ. 0 ) THEN DO 200 J = 1, JVMAX #ifndef OPT_INVGRD IJS = (J-1)*IVDIM + 1 #else IJS = (JVMAX-J)*IVDIM + 1 #endif IJE = IJS+IVMAX-1 READ (IFILE,'('//C3//'I1)') (IMAP(IJ),IJ=IJS,IJE) 200 CONTINUE CLOSE ( IFILE ) ENDIF ELSE CALL RIVNER !! nearest route O ( IMAP , I ISTA , IEND , W IJDEST(1,IR) , RDEST(1,IR) ) ENDIF * CALL RIVDST !! destination O ( IJDEST(1,IR), IDESTN , IDESTS , O NDESTN(IR) , NDESTS(IR) , RDEST(1,IR) , I IMAP , ISTA , IEND , VRIVER(IR) ) * JN = NDESTN( IR ) JS = NDESTS( IR ) * #ifndef OPT_NONPARALLEL * CALL NSCOMMIF O ( MDESTS( IR ) , MDESTN( IR ) , I JS , JN , I 1 , 1 , 1 , 1 , D 1 , 1 , 1 , 1 ) * CALL NSCOMMIF O ( JDESTS(1,IR) , JDESTN(1,IR) , I IDESTS , IDESTN , I NDESTS(IR) , NDESTN(IR) , MDESTS(IR) , MDESTN(IR) , D IVDIM , IVDIM , IVDIM , IVDIM ) * #endif * DO 510 IJ = ISTA, IEND IF ( AREA( IJ ) .EQ. 0.D0 ) THEN IJDEST( IJ,IR ) = IJ+ISTAL-1 RDEST ( IJ,IR ) = 0.D0 ENDIF c IF (IMAP(IJ+ISTAL-1).GE.1) THEN c WRITE(JFPAR,'(i7,i3,6e15.7)') c $ IJ+ISTAL-1,IMAP(IJ+ISTAL-1),AREA(IJ),RDEST(IJ,IR), c $ AVLON(IJ),DVLON(IJ),AVLAT(IJ),DVLAT(IJ) c ENDIF 510 CONTINUE 500 CONTINUE * CALL OPNINI !! open initial data O ( IFLSTR, TINIT0, TINIT1, HDFMT , I TIME , 'RIV' ) CALL GDREAD O ( GDRIV , IEOD , TINITD, TDURD , KLEVS , I IFLSTR,'GRIVR', HDFMT , 'RIV' , TINIT0, TINIT1, I IVMAX , JVMAX , IVDIM , JVDIM , KRVMAX ) IF ( IEOD .NE. 0 ) THEN CALL REWNML ( IFPAR , JFPAR ) WRITE (JFPAR,*) ' ### RIVINI: Using zero as Init RIV Data.' CALL RESET( GDRIV, IJVDIM * KRVMAX ) ELSEIF( TINIT0 .LT. 0.D0 ) THEN TIME = TINITD ENDIF TSTART = TIME * 900 CALL SNDR8_COMM( TIME , 1 , 2 ) * RETURN *========================================================================= ENTRY RIVRST !! restart river I ( TIME , OEND ) * CALL GETOPCALC ( OPCALC ) IF( .NOT. OPCALC ) THEN RETURN ENDIF * CALL HISTOU ( OEND ) * CALL OPNRST !! open restart data O ( JFLRST, HDFMT , I TIME , TSTART, OEND , 'RIV' ) IF ( JFLRST .LE. 0 ) RETURN * CALL GDWRIT I ( GDRIV , 'GRIVR','river water storage','kg/m**2', I TIME , 1.D0 , JFLRST, HDFMT , 'RIV' , I IVMAX , JVMAX , IVDIM , JVDIM , KRVMAX, 'VSFC' ) * RETURN END ******************************************************************** SUBROUTINE RIVDST !! river destination O ( IJDEST, JDESTN , JDESTS , NDESTN , NDESTS , RDEST , I IMAP , ISTA , IEND , VRIVER ) * * [PARAM] #include "zxmpif.F" /* MPI comm. param. */ #include "zcdim.F" /* # of grid point & wave */ #include "zpdim.F" /* physics etc. */ #include "zhdim.F" /* # of char. */ #include "zvdim.F" /* # of river grid point */ #include "zcorv.F" /* river coordinate */ #include "zccom.F" /* stand. physical const. */ * * [OUTPUT] INTEGER IJDEST ( IJVDIM ) !! Flow destination IJ INTEGER JDESTN ( IVDIM ) !! Flow destination at N-ridge INTEGER JDESTS ( IVDIM ) !! Flow destination at S-ridge INTEGER NDESTN !! No. flow dest. at N-ridge INTEGER NDESTS !! No. flow dest. at S-ridge REAL*8 RDEST ( IJVDIM ) !! distination speed/distance * * [INPUT] INTEGER IMAP ( IVDIM*JVMAX ) !! flow direction (1 to 8) INTEGER ISTA, IEND !! start/end IJ index REAL*8 VRIVER !! river flow speed * * [WORK] INTEGER IDIST ( IVDIM*JVMAX ) !! distance from ocean INTEGER IDIST0 ( IJVDIM ) !! distance from ocean REAL*8 GIDX ( IJVDIM ) !! surface index * * [INTERNAL WORK] INTEGER IJ, IJL, IM, I, J, JX, ID, JN, JS, IJS, ISET INTEGER I1, J1 INTEGER IDMX, ITR, NP, IDM, IMX REAL*8 COSX REAL * 8 TIME INTEGER NPROCS, MYRANK INTEGER ISTAL, IENDL !! start/end IJ index (large) INTEGER JSTAS, JENDS !! start/end J index (small) INTEGER JSTAL, JENDL !! start/end J index (large) INTEGER JSTALN, JSTALS !! start J index (large) INTEGER IFPAR, JFPAR REAL*8 AVLATN( IVDIM ) REAL*8 AVLATS( IVDIM ) REAL*8 AVLONN( IVDIM ) REAL*8 AVLONS( IVDIM ) #ifndef OPT_NONPARALLEL REAL*8 WRK1( IVDIM ) REAL*8 WRK2( IVDIM ) INTEGER IDIR INTEGER IREQS( 0: MAXNP-1 ) INTEGER IREQR( 0: MAXNP-1) INTEGER ISTATUS(MPI_STATUS_SIZE) INTEGER MPI_INTNEW( 0: JVMAX-1) !! new data type INTEGER IRANK, IERR #endif LOGICAL OPFIRST DATA OPFIRST / .TRUE. / * * [INTERNAL SAVE] INTEGER IOFS ( 0:9 ) INTEGER JOFS ( 0:9 ) INTEGER IMPRF ( 8 ) * 0 1 2 3 4 5 6 7 8 9 DATA IOFS / 0, 0, 1, 1, 1, 0,-1,-1,-1, 0 / #ifndef OPT_INVGRD DATA JOFS / 0,-1,-1, 0, 1, 1, 1, 0,-1, 0 / #else DATA JOFS / 0, 1, 1, 0,-1,-1,-1, 0, 1, 0 / #endif DATA IMPRF / 3, 7, 1, 5, 2, 4, 6, 8 / * CHARACTER C3*3 !! kei * JSTALN = 1 JSTALS = 1 CALL GETMPI2 O ( NPROCS, MYRANK, JSTAL, JENDL, IFPAR, JFPAR, I 'RIV' ) #ifndef OPT_NONPARALLEL IF( MYRANK .GT. 0 ) CALL PARA_RANGE I ( 1 , JVMAX , NPROCS , MYRANK-1 , O JSTALN , IENDL ) IF( MYRANK .LT. NPROCS-1 ) CALL PARA_RANGE I ( 1 , JVMAX , NPROCS , MYRANK+1 , O JSTALS , IENDL ) * * communication at N/S boundary * CALL MCUT & ( WRK1, AVLAT, IJVDIM-IVDIM+1, IVDIM, IJVDIM, 1 ) CALL MCUT & ( WRK2, AVLAT, 1, IVDIM, IJVDIM, 1 ) CALL NSCOMM O ( AVLATN, AVLATS, I WRK1 , WRK2 , IVDIM , D IJVDIM ) CALL MCUT & ( WRK1, AVLON, IJVDIM-IVDIM+1, IVDIM, IJVDIM, 1 ) CALL MCUT & ( WRK2, AVLON, 1, IVDIM, IJVDIM, 1 ) CALL NSCOMM O ( AVLONN, AVLONS, I WRK1 , WRK2 , IVDIM , D IJVDIM ) #else #endif ISTAL = ( JSTAL-1 ) * IVDIM + 1 IENDL = JENDL * IVDIM * * identify river mouth DO 1050 J = 1, JVMAX DO 1050 I = 1, IVMAX IJ = ( J-1 )*IVDIM + I IM = MIN( MAX( IMAP( IJ ), 0 ), 9 ) I1 = I + IOFS( IM ) JX = J + JOFS( IM ) J1 = MIN( MAX( JX,1 ), JVMAX ) IF ( J1 .NE. JX ) I1 = I1 + IVMAX/2 IF ( I1 .LT. 1 ) I1 = I1 + IVMAX IF ( I1 .GT. IVMAX ) I1 = I1 - IVMAX ID = ( J1-1 )*IVDIM + I1 IF ( IMAP(ID) .EQ. 0 .AND. IJ .NE. ID ) IMAP( ID ) = 9 1050 CONTINUE * JN = 0 JS = 0 DO 1100 IJ = ISTA, IEND IJL = IJ + ( JSTAL-1 ) * IVDIM IM = MIN( MAX( IMAP( IJL ), 0 ), 9 ) I = MOD( IJ-1,IVDIM ) + 1 + IOFS( IM ) JX = (IJ-1)/IVDIM + JSTAL + JOFS( IM ) J = MIN( MAX( JX,1 ), JVMAX ) IF ( J .NE. JX ) I = I + IVMAX/2 IF ( I .LT. 1 ) I = I + IVMAX IF ( I .GT. IVMAX ) I = I - IVMAX ID = ( J-1 )*IVDIM + I IJDEST ( IJ ) = ID IF( ID .LT. ISTAL ) THEN JN = JN + 1 JDESTN( JN ) = ID - (JSTALN-1)*IVDIM ENDIF IF( ID .GT. IENDL ) THEN JS = JS + 1 JDESTS( JS ) = ID - (JSTALS-1)*IVDIM ENDIF ID = ID - ISTAL + 1 * RDEST ( IJ ) = 0.D0 IF ( IMAP(IJ+ISTAL-1).GE.1 ) THEN !! kei IF ( ID .NE. IJ ) THEN IF( ID .LT. ISTA ) THEN COSX = SIN( AVLAT( IJ ) )*SIN( AVLATS( ID+IVDIM ) ) & + COS( AVLAT( IJ ) )*COS( AVLATS( ID+IVDIM ) ) & * COS( AVLON( IJ ) - AVLONS( ID+IVDIM ) ) ELSEIF( ID .GT. IEND ) THEN COSX = SIN( AVLAT( IJ ) )*SIN( AVLATN( ID-IEND ) ) & + COS( AVLAT( IJ ) )*COS( AVLATN( ID-IEND ) ) & * COS( AVLON( IJ ) - AVLONN( ID-IEND ) ) ELSE COSX = SIN( AVLAT( IJ ) )*SIN( AVLAT( ID ) ) & + COS( AVLAT( IJ ) )*COS( AVLAT( ID ) ) & * COS( AVLON( IJ ) - AVLON( ID ) ) ENDIF IF ( COSX .LT. 1.D0 ) THEN RDEST( IJ ) = VRIVER/( ACOS( COSX )*ER ) ENDIF **kei ELSE IF ( MOD(IJ,IVDIM).EQ.IVMAX ) THEN COSX = SIN( AVLAT( IJ ) )*SIN( AVLAT( IJ ) ) & + COS( AVLAT( IJ ) )*COS( AVLAT( IJ ) ) & * COS( AVLON( IJ ) - AVLON( IJ-IVMAX+1 ) ) ELSE COSX = SIN( AVLAT( IJ ) )*SIN( AVLAT( IJ ) ) & + COS( AVLAT( IJ ) )*COS( AVLAT( IJ ) ) & * COS( AVLON( IJ ) - AVLON( IJ+1 ) ) ENDIF IF ( COSX .LT. 1.D0 ) THEN RDEST( IJ ) = VRIVER/( ACOS( COSX )*ER ) ENDIF ENDIF ENDIF 1100 CONTINUE NDESTN = JN NDESTS = JS * WRITE (JFPAR,*) ' *** RIVER FLOW MAP START ' WRITE(C3,'(i3)')IVMAX DO 1900 J = 1, JVMAX IJS = ( J-1 )*IVDIM + 1 WRITE(JFPAR,'('//C3//'I1)') (IMAP( IJ ), IJ=IJS, IJS+IVMAX-1) 1900 CONTINUE WRITE (JFPAR,*) ' *** RIVER FLOW MAP END ' * RETURN *==================================================================== ENTRY RIVNER !! river nearest route O ( IMAP , I ISTA , IEND , W IDIST0 , GIDX ) * CALL GETTIM( TIME ) CALL GETMPI O ( NPROCS, MYRANK, JSTAS, JENDS, IFPAR, JFPAR, I 'RIV' ) JSTAS = ( JSTAS - 1 ) / IVDIM + 1 JENDS = JENDS / IVDIM #ifndef OPT_NONPARALLEL CALL GETDTYPEi O ( MPI_INTNEW, JSTAL , JENDL , I IVDIM , JVMAX ) #endif CALL RDDAT ( GIDX, ISET, TIME, 'GRIDX', 'ASFC', .FALSE. ) * IDMX = 9999 DO 2100 IJ = 1, IVDIM*JVDIM IF ( GIDX( IJ ) .LT. 0.5D0 ) THEN IDIST0( IJ ) = 0 ELSE IDIST0( IJ ) = IDMX ENDIF 2100 CONTINUE * #ifndef OPT_NONPARALLEL IDIR = 0 CALL ISFTARRAY I ( IDIST0 , O IDIST , I IVDIM , JVDIM , JVMAX , 1 , I JSTAS , JENDS , JSTAL , JENDL , IDIR ) DO 2200 IRANK = 0, NPROCS - 1 IF( IRANK .NE. MYRANK ) THEN CALL MPI_ISEND( IDIST, 1, MPI_INTNEW(MYRANK), IRANK, 1, $ MPI_COMM_AGCM, IREQS(IRANK), IERR ) CALL MPI_IRECV( IDIST, 1, MPI_INTNEW(IRANK), IRANK, 1, $ MPI_COMM_AGCM, IREQR(IRANK), IERR ) ENDIF IF( OPFIRST ) WRITE( JFPAR, *) ' +++ IDIST Oobanburumai from', $ IRANK, ' +++++++ ' 2200 CONTINUE DO 2300 IRANK = 0, NPROCS - 1 IF( IRANK .NE. MYRANK ) THEN CALL MPI_WAIT( IREQS(IRANK), ISTATUS, IERR ) CALL MPI_WAIT( IREQR(IRANK), ISTATUS, IERR ) ENDIF 2300 CONTINUE #else CALL ICOPY ( IDIST , IDIST0 , IJVDIM ) #endif * DO 3000 ITR = 1, IDMX NP = 0 DO 3100 IJ = 1, IVDIM*JVMAX IMAP( IJ ) = IDIST( IJ ) 3100 CONTINUE DO 3200 IJ = 1, IVDIM*JVMAX IF ( IDIST( IJ ) .EQ. IDMX ) THEN DO 3210 IM = 1, 8 I = MOD( IJ-1,IVDIM ) + 1 + IOFS( IM ) JX = (IJ-1)/IVDIM + 1 + JOFS( IM ) J = MIN( MAX( JX,1 ), JVMAX ) IF ( J .NE. JX ) I = I + IVMAX/2 IF ( I .LT. 1 ) I = I + IVMAX IF ( I .GT. IVMAX ) I = I - IVMAX ID = ( J-1 )*IVDIM + I IDIST( IJ ) = MIN( IDIST( IJ ), IMAP( ID )+1 ) 3210 CONTINUE NP = NP + 1 ENDIF 3200 CONTINUE IF ( NP .EQ. 0 ) GOTO 3900 IF ( NP .EQ. IVDIM*JVMAX ) THEN WRITE (JFPAR,*) ' ### PGRIV: ALL GRID POINTS ARE LAND' GOTO 3900 ENDIF 3000 CONTINUE 3900 CONTINUE * DO 4100 IJ = 1, IVDIM*JVMAX IMAP( IJ ) = 0 IF ( IDIST( IJ ) .GT. 0 ) THEN IDM = IDMX DO 4110 IMX = 1, 8 IM = IMPRF( IMX ) I = MOD( IJ-1,IVDIM ) + 1 + IOFS( IM ) JX = (IJ-1)/IVDIM + 1 + JOFS( IM ) J = MIN( MAX( JX,1 ), JVMAX ) IF ( J .NE. JX ) I = I + IVMAX/2 IF ( I .LT. 1 ) I = I + IVMAX IF ( I .GT. IVMAX ) I = I - IVMAX ID = ( J-1 )*IVDIM + I IF ( IDIST( ID ) .LT. IDM ) THEN IDM = IDIST( ID ) IMAP( IJ ) = IM ENDIF 4110 CONTINUE ENDIF 4100 CONTINUE * OPFIRST = .FALSE. * RETURN END