! ! $Header$ ! SUBROUTINE advx(limit,dtx,pbaru,sm,s0, & sx,sy,sz,lati,latf) USE dimensions_mod, ONLY: iim, jjm, llm, ndm USE paramet_mod_h IMPLICIT NONE !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC ! C ! first-order moments (FOM) advection of tracer in X direction C ! C ! Source : Pascal Simon (Meteo,CNRM) C ! Adaptation : A.Armengaud (LGGE) juin 94 C ! C ! limit,dtx,pbaru,pbarv,sm,s0,sx,sy,sz C ! sont des arguments d'entree pour le s-pg... C ! C ! sm,s0,sx,sy,sz C ! sont les arguments de sortie pour le s-pg C ! C !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC ! ! parametres principaux du modele ! ! Arguments : ! ----------- ! dtx : frequence fictive d'appel du transport ! pbaru, pbarv : flux de masse en x et y en Pa.m2.s-1 INTEGER :: ntra PARAMETER (ntra = 1) ! ATTENTION partout ou on trouve ntra, insertion de boucle ! possible dans l'avenir. REAL :: dtx REAL :: pbaru ( iip1,jjp1,llm ) ! moments: SM total mass in each grid box ! S0 mass of tracer in each grid box ! Si 1rst order moment in i direction ! REAL :: SM(iip1,jjp1,llm),S0(iip1,jjp1,llm,ntra) REAL :: sx(iip1,jjp1,llm,ntra) & ,sy(iip1,jjp1,llm,ntra) REAL :: sz(iip1,jjp1,llm,ntra) ! Local : ! ------- ! mass fluxes across the boundaries (UGRI,VGRI,WGRI) ! mass fluxes in kg ! declaration : REAL :: UGRI(iip1,jjp1,llm) ! Rem : VGRI et WGRI ne sont pas utilises dans ! cette subroutine ( advection en x uniquement ) ! ! Ti are the moments for the current latitude and level ! REAL :: TM(iim) REAL :: T0(iim,ntra),TX(iim,ntra) REAL :: TY(iim,ntra),TZ(iim,ntra) REAL :: TEMPTM ! just a temporary variable ! ! the moments F are similarly defined and used as temporary ! storage for portions of the grid boxes in transit ! REAL :: FM(iim) REAL :: F0(iim,ntra),FX(iim,ntra) REAL :: FY(iim,ntra),FZ(iim,ntra) ! ! work arrays ! REAL :: ALF(iim),ALF1(iim),ALFQ(iim),ALF1Q(iim) ! REAL :: SMNEW(iim),UEXT(iim) ! REAL :: sqi,sqf LOGICAL :: LIMIT INTEGER :: NUM(jjp1),LONK,NUMK INTEGER :: lon,lati,latf,niv INTEGER :: i,i2,i3,j,jv,l,k,itrac lon = iim niv = llm ! *** Test de passage d'arguments ****** ! ------------------------------------- DO j = 1,jjp1 NUM(j) = 1 END DO sqi = 0. sqf = 0. DO l = 1,llm DO j = 1,jjp1 DO i = 1,iim !IM 240305 sqi = sqi + S0(i,j,l,9) sqi = sqi + S0(i,j,l,ntra) ENDDO ENDDO ENDDO PRINT*,'-------- DIAG DANS ADVX - ENTREE ---------' PRINT*,'sqi=',sqi ! Interface : adaptation nouveau modele ! ------------------------------------- ! ! --------------------------------------------------------- ! Conversion des flux de masses en kg/s ! pbaru est en N/s d'ou : ! ugri est en kg/s DO l = 1,llm DO j = 1,jjm+1 DO i = 1,iip1 ! ugri (i,j,llm+1-l) = pbaru (i,j,l) * ( dsig(l) / g ) ugri (i,j,llm+1-l) = pbaru (i,j,l) END DO END DO END DO ! --------------------------------------------------------- ! --------------------------------------------------------- ! --------------------------------------------------------- ! start here ! ! boucle principale sur les niveaux et les latitudes ! DO L=1,NIV DO K=lati,latf ! ! initialisation ! ! program assumes periodic boundaries in X ! DO I=2,LON SMNEW(I)=SM(I,K,L)+(UGRI(I-1,K,L)-UGRI(I,K,L))*DTX END DO SMNEW(1)=SM(1,K,L)+(UGRI(LON,K,L)-UGRI(1,K,L))*DTX ! ! modifications for extended polar zones ! NUMK=NUM(K) LONK=LON/NUMK ! IF(NUMK.GT.1) THEN ! DO I=1,LON TM(I)=0. END DO DO JV=1,NTRA DO I=1,LON T0(I,JV)=0. TX(I,JV)=0. TY(I,JV)=0. TZ(I,JV)=0. END DO END DO ! DO I2=1,NUMK ! DO I=1,LONK I3=(I-1)*NUMK+I2 TM(I)=TM(I)+SM(I3,K,L) ALF(I)=SM(I3,K,L)/TM(I) ALF1(I)=1.-ALF(I) END DO ! DO JV=1,NTRA DO I=1,LONK I3=(I-1)*NUMK+I2 TEMPTM=-ALF(I)*T0(I,JV)+ALF1(I) & *S0(I3,K,L,JV) T0(I,JV)=T0(I,JV)+S0(I3,K,L,JV) TX(I,JV)=ALF(I) *sx(I3,K,L,JV)+ & ALF1(I)*TX(I,JV) +3.*TEMPTM TY(I,JV)=TY(I,JV)+sy(I3,K,L,JV) TZ(I,JV)=TZ(I,JV)+sz(I3,K,L,JV) ENDDO ENDDO ! END DO ! ELSE ! DO I=1,LON TM(I)=SM(I,K,L) END DO DO JV=1,NTRA DO I=1,LON T0(I,JV)=S0(I,K,L,JV) TX(I,JV)=sx(I,K,L,JV) TY(I,JV)=sy(I,K,L,JV) TZ(I,JV)=sz(I,K,L,JV) END DO END DO ! ENDIF ! DO I=1,LONK UEXT(I)=UGRI(I*NUMK,K,L) END DO ! ! place limits on appropriate moments before transport ! (if flux-limiting is to be applied) ! IF(.NOT.LIMIT) GO TO 13 ! DO JV=1,NTRA DO I=1,LONK TX(I,JV)=SIGN(AMIN1(AMAX1(T0(I,JV),0.),ABS(TX(I,JV))),TX(I,JV)) END DO END DO ! 13 CONTINUE ! ! calculate flux and moments between adjacent boxes ! 1- create temporary moments/masses for partial boxes in transit ! 2- reajusts moments remaining in the box ! ! flux from IP to I if U(I).lt.0 ! DO I=1,LONK-1 IF(UEXT(I).LT.0.) THEN FM(I)=-UEXT(I)*DTX ALF(I)=FM(I)/TM(I+1) TM(I+1)=TM(I+1)-FM(I) ENDIF END DO ! I=LONK IF(UEXT(I).LT.0.) THEN FM(I)=-UEXT(I)*DTX ALF(I)=FM(I)/TM(1) TM(1)=TM(1)-FM(I) ENDIF ! ! flux from I to IP if U(I).gt.0 ! DO I=1,LONK IF(UEXT(I).GE.0.) THEN FM(I)=UEXT(I)*DTX ALF(I)=FM(I)/TM(I) TM(I)=TM(I)-FM(I) ENDIF END DO ! DO I=1,LONK ALFQ(I)=ALF(I)*ALF(I) ALF1(I)=1.-ALF(I) ALF1Q(I)=ALF1(I)*ALF1(I) END DO ! DO JV=1,NTRA DO I=1,LONK-1 ! IF(UEXT(I).LT.0.) THEN ! F0(I,JV)=ALF (I)* ( T0(I+1,JV)-ALF1(I)*TX(I+1,JV) ) FX(I,JV)=ALFQ(I)*TX(I+1,JV) FY(I,JV)=ALF (I)*TY(I+1,JV) FZ(I,JV)=ALF (I)*TZ(I+1,JV) ! T0(I+1,JV)=T0(I+1,JV)-F0(I,JV) TX(I+1,JV)=ALF1Q(I)*TX(I+1,JV) TY(I+1,JV)=TY(I+1,JV)-FY(I,JV) TZ(I+1,JV)=TZ(I+1,JV)-FZ(I,JV) ! ENDIF ! END DO END DO ! I=LONK IF(UEXT(I).LT.0.) THEN ! DO JV=1,NTRA ! F0 (I,JV)=ALF (I)* ( T0(1,JV)-ALF1(I)*TX(1,JV) ) FX (I,JV)=ALFQ(I)*TX(1,JV) FY (I,JV)=ALF (I)*TY(1,JV) FZ (I,JV)=ALF (I)*TZ(1,JV) ! T0(1,JV)=T0(1,JV)-F0(I,JV) TX(1,JV)=ALF1Q(I)*TX(1,JV) TY(1,JV)=TY(1,JV)-FY(I,JV) TZ(1,JV)=TZ(1,JV)-FZ(I,JV) ! END DO ! ENDIF ! DO JV=1,NTRA DO I=1,LONK ! IF(UEXT(I).GE.0.) THEN ! F0(I,JV)=ALF (I)* ( T0(I,JV)+ALF1(I)*TX(I,JV) ) FX(I,JV)=ALFQ(I)*TX(I,JV) FY(I,JV)=ALF (I)*TY(I,JV) FZ(I,JV)=ALF (I)*TZ(I,JV) ! T0(I,JV)=T0(I,JV)-F0(I,JV) TX(I,JV)=ALF1Q(I)*TX(I,JV) TY(I,JV)=TY(I,JV)-FY(I,JV) TZ(I,JV)=TZ(I,JV)-FZ(I,JV) ! ENDIF ! END DO END DO ! ! puts the temporary moments Fi into appropriate neighboring boxes ! DO I=1,LONK IF(UEXT(I).LT.0.) THEN TM(I)=TM(I)+FM(I) ALF(I)=FM(I)/TM(I) ENDIF END DO ! DO I=1,LONK-1 IF(UEXT(I).GE.0.) THEN TM(I+1)=TM(I+1)+FM(I) ALF(I)=FM(I)/TM(I+1) ENDIF END DO ! I=LONK IF(UEXT(I).GE.0.) THEN TM(1)=TM(1)+FM(I) ALF(I)=FM(I)/TM(1) ENDIF ! DO I=1,LONK ALF1(I)=1.-ALF(I) END DO ! DO JV=1,NTRA DO I=1,LONK ! IF(UEXT(I).LT.0.) THEN ! TEMPTM=-ALF(I)*T0(I,JV)+ALF1(I)*F0(I,JV) T0(I,JV)=T0(I,JV)+F0(I,JV) TX(I,JV)=ALF(I)*FX(I,JV)+ALF1(I)*TX(I,JV)+3.*TEMPTM TY(I,JV)=TY(I,JV)+FY(I,JV) TZ(I,JV)=TZ(I,JV)+FZ(I,JV) ! ENDIF ! END DO END DO ! DO JV=1,NTRA DO I=1,LONK-1 ! IF(UEXT(I).GE.0.) THEN ! TEMPTM=ALF(I)*T0(I+1,JV)-ALF1(I)*F0(I,JV) T0(I+1,JV)=T0(I+1,JV)+F0(I,JV) TX(I+1,JV)=ALF(I)*FX(I,JV)+ALF1(I)*TX(I+1,JV)+3.*TEMPTM TY(I+1,JV)=TY(I+1,JV)+FY(I,JV) TZ(I+1,JV)=TZ(I+1,JV)+FZ(I,JV) ! ENDIF ! END DO END DO ! I=LONK IF(UEXT(I).GE.0.) THEN DO JV=1,NTRA TEMPTM=ALF(I)*T0(1,JV)-ALF1(I)*F0(I,JV) T0(1,JV)=T0(1,JV)+F0(I,JV) TX(1,JV)=ALF(I)*FX(I,JV)+ALF1(I)*TX(1,JV)+3.*TEMPTM TY(1,JV)=TY(1,JV)+FY(I,JV) TZ(1,JV)=TZ(1,JV)+FZ(I,JV) END DO ENDIF ! ! retour aux mailles d'origine (passage des Tij aux Sij) ! IF(NUMK.GT.1) THEN ! DO I2=1,NUMK ! DO I=1,LONK ! I3=I2+(I-1)*NUMK SM(I3,K,L)=SMNEW(I3) ALF(I)=SMNEW(I3)/TM(I) TM(I)=TM(I)-SMNEW(I3) ! ALFQ(I)=ALF(I)*ALF(I) ALF1(I)=1.-ALF(I) ALF1Q(I)=ALF1(I)*ALF1(I) ! END DO END DO ! DO JV=1,NTRA DO I=1,LONK ! I3=I2+(I-1)*NUMK S0(I3,K,L,JV)=ALF (I) & * (T0(I,JV)-ALF1(I)*TX(I,JV)) sx(I3,K,L,JV)=ALFQ(I)*TX(I,JV) sy(I3,K,L,JV)=ALF (I)*TY(I,JV) sz(I3,K,L,JV)=ALF (I)*TZ(I,JV) ! ! reajusts moments remaining in the box ! T0(I,JV)=T0(I,JV)-S0(I3,K,L,JV) TX(I,JV)=ALF1Q(I)*TX(I,JV) TY(I,JV)=TY(I,JV)-sy(I3,K,L,JV) TZ(I,JV)=TZ(I,JV)-sz(I3,K,L,JV) ENDDO ENDDO ! ! ELSE ! DO I=1,LON SM(I,K,L)=TM(I) END DO DO JV=1,NTRA DO I=1,LON S0(I,K,L,JV)=T0(I,JV) sx(I,K,L,JV)=TX(I,JV) sy(I,K,L,JV)=TY(I,JV) sz(I,K,L,JV)=TZ(I,JV) END DO END DO ! ENDIF ! END DO END DO ! ! ----------- AA Test en fin de ADVX ------ Controle des S* ! OK ! DO 9998 l = 1, llm ! DO 9998 j = 1, jjp1 ! DO 9998 i = 1, iip1 ! IF (S0(i,j,l,ntra).lt.0..and.LIMIT) THEN ! PRINT*, '-------------------' ! PRINT*, 'En fin de ADVX' ! PRINT*,'SM(',i,j,l,')=',SM(i,j,l) ! PRINT*,'S0(',i,j,l,')=',S0(i,j,l,ntra) ! print*, 'sx(',i,j,l,')=',sx(i,j,l,ntra) ! print*, 'sy(',i,j,l,')=',sy(i,j,l,ntra) ! print*, 'sz(',i,j,l,')=',sz(i,j,l,ntra) ! WRITE (*,*) 'On arrete !! - pbl en fin de ADVX1' !c STOP ! ENDIF ! 9998 CONTINUE ! ! ---------- bouclage cyclique DO itrac=1,ntra DO l = 1,llm DO j = lati,latf SM(iip1,j,l) = SM(1,j,l) S0(iip1,j,l,itrac) = S0(1,j,l,itrac) sx(iip1,j,l,itrac) = sx(1,j,l,itrac) sy(iip1,j,l,itrac) = sy(1,j,l,itrac) sz(iip1,j,l,itrac) = sz(1,j,l,itrac) END DO END DO ENDDO ! ----------- qqtite totale de traceur dans tte l'atmosphere DO l = 1, llm DO j = 1, jjp1 DO i = 1, iim !IM 240405 sqf = sqf + S0(i,j,l,9) sqf = sqf + S0(i,j,l,ntra) END DO END DO END DO ! PRINT*,'------ DIAG DANS ADVX - SORTIE -----' PRINT*,'sqf=',sqf !------------- RETURN END SUBROUTINE advx !_________________________________________________________________ !_________________________________________________________________