/[lmdze]/trunk/dyn3d/advect.f
ViewVC logotype

Diff of /trunk/dyn3d/advect.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/dyn3d/advect.f revision 3 by guez, Wed Feb 27 13:16:39 2008 UTC trunk/libf/dyn3d/advect.f90 revision 43 by guez, Fri Apr 8 12:43:31 2011 UTC
# Line 1  Line 1 
1  !  SUBROUTINE advect(ucov, vcov, teta, w, massebx, masseby, du, dv, dteta, conser)
2  ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/advect.F,v 1.1.1.1 2004/05/19 12:53:06 lmdzadmin Exp $  
3  !    ! From dyn3d/advect.F, v 1.1.1.1 2004/05/19 12:53:06
4        SUBROUTINE advect(ucov,vcov,teta,w,massebx,masseby,du,dv,dteta)    !   Auteurs:  P. Le Van , Fr. Hourdin  .                                
5      !   Objet:                                                              
6        use dimens_m    !   .... calcul des termes d'advection vertic.pour u, v, teta, q ...      
7        use paramet_m    !        ces termes sont ajoutes a du, dv, dteta et dq .                  
8        use comconst    !  Modif F.Forget 03/94 : on retire q de advect                        
9        use comvert  
10        use logic    USE dimens_m
11        use comgeom    USE paramet_m
12        use ener    USE comconst
13        IMPLICIT NONE    USE comvert
14  c=======================================================================    USE comgeom
15  c    USE ener
16  c   Auteurs:  P. Le Van , Fr. Hourdin  .  
17  c   -------    IMPLICIT NONE
18  c  
19  c   Objet:    !   Arguments:                                                          
20  c   ------  
21  c    REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm)
22  c   *************************************************************    real, intent(in):: teta(ip1jmp1, llm)
23  c   .... calcul des termes d'advection vertic.pour u,v,teta,q ...    REAL massebx(ip1jmp1, llm), masseby(ip1jm, llm)
24  c   *************************************************************    real, INTENT (IN):: w(ip1jmp1, llm)
25  c        ces termes sont ajoutes a du,dv,dteta et dq .    REAL dv(ip1jm, llm), du(ip1jmp1, llm), dteta(ip1jmp1, llm)
26  c  Modif F.Forget 03/94 : on retire q de advect    LOGICAL, INTENT (IN):: conser
27  c  
28  c=======================================================================    !   Local:                                                              
29  c-----------------------------------------------------------------------  
30  c   Declarations:    REAL uav(ip1jmp1, llm), vav(ip1jm, llm), wsur2(ip1jmp1)
31  c   -------------    REAL unsaire2(ip1jmp1), ge(ip1jmp1)
32      REAL deuxjour, ww, gt, uu, vv
33    
34  c   Arguments:    INTEGER ij, l
35  c   ----------  
36      REAL ssum
37        REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)  
38        REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm),w(ip1jmp1,llm)    !-----------------------------------------------------------------------
39        REAL dv(ip1jm,llm),du(ip1jmp1,llm),dteta(ip1jmp1,llm)    !   2. Calculs preliminaires:                                          
40      !   -------------------------                                          
41  c   Local:  
42  c   ------    IF (conser) THEN
43         deuxjour = 2.*daysec
44        REAL uav(ip1jmp1,llm),vav(ip1jm,llm),wsur2(ip1jmp1)  
45        REAL unsaire2(ip1jmp1), ge(ip1jmp1)       DO ij = 1, ip1jmp1
46        REAL deuxjour, ww, gt, uu, vv          unsaire2(ij) = unsaire(ij)*unsaire(ij)
47         end DO
48        INTEGER  ij,l    END IF
49    
50        REAL      SSUM  
51      !------------------  -yy ----------------------------------------------
52  c-----------------------------------------------------------------------    !   .  Calcul de     u                                                  
53  c   2. Calculs preliminaires:  
54  c   -------------------------    DO l = 1, llm
55         DO ij = iip2, ip1jmp1
56        IF (conser)  THEN          uav(ij, l) = 0.25*(ucov(ij, l)+ucov(ij-iip1, l))
57           deuxjour = 2. * daysec       END DO
58         DO ij = iip2, ip1jm
59           DO   1  ij   = 1, ip1jmp1          uav(ij, l) = uav(ij, l) + uav(ij+iip1, l)
60           unsaire2(ij) = unsaire(ij) * unsaire(ij)       END DO
61     1     CONTINUE       DO ij = 1, iip1
62        END IF          uav(ij, l) = 0.
63            uav(ip1jm+ij, l) = 0.
64         END DO
65  c------------------  -yy ----------------------------------------------    END DO
66  c   .  Calcul de     u  
67      !------------------  -xx ----------------------------------------------
68        DO  l=1,llm    !   .  Calcul de     v                                                  
69           DO    ij     = iip2, ip1jmp1  
70              uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )    DO l = 1, llm
71           ENDDO       DO ij = 2, ip1jm
72           DO    ij     = iip2, ip1jm          vav(ij, l) = 0.25*(vcov(ij, l)+vcov(ij-1, l))
73              uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)       END DO
74           ENDDO       DO ij = 1, ip1jm, iip1
75           DO      ij         = 1, iip1          vav(ij, l) = vav(ij+iim, l)
76              uav(ij      ,l) = 0.       END DO
77              uav(ip1jm+ij,l) = 0.       DO ij = 1, ip1jm - 1
78           ENDDO          vav(ij, l) = vav(ij, l) + vav(ij+1, l)
79        ENDDO       END DO
80         DO ij = 1, ip1jm, iip1
81  c------------------  -xx ----------------------------------------------          vav(ij+iim, l) = vav(ij, l)
82  c   .  Calcul de     v       END DO
83      END DO
84        DO  l=1,llm  
85           DO    ij   = 2, ip1jm    !-----------------------------------------------------------------------
86            vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )  
87           ENDDO  
88           DO    ij   = 1,ip1jm,iip1    DO l = 1, llmm1
89            vav(ij,l) = vav(ij+iim,l)  
90           ENDDO  
91           DO    ij   = 1, ip1jm-1       !       ......   calcul de  - w/2.    au niveau  l+1   .......          
92            vav(ij,l) = vav(ij,l) + vav(ij+1,l)  
93           ENDDO       DO ij = 1, ip1jmp1
94           DO    ij       = 1, ip1jm, iip1          wsur2(ij) = -0.5*w(ij, l+1)
95            vav(ij+iim,l) = vav(ij,l)       END DO
96           ENDDO  
97        ENDDO  
98         !     .....................     calcul pour  du     ..................  
99  c-----------------------------------------------------------------------  
100         DO ij = iip2, ip1jm - 1
101  c          ww = wsur2(ij) + wsur2(ij+1)
102        DO 20 l = 1, llmm1          uu = 0.5*(ucov(ij, l)+ucov(ij, l+1))
103            du(ij, l) = du(ij, l) - ww*(uu-uav(ij, l))/massebx(ij, l)
104            du(ij, l+1) = du(ij, l+1) + ww*(uu-uav(ij, l+1))/massebx(ij, l+1)
105  c       ......   calcul de  - w/2.    au niveau  l+1   .......       END DO
106    
107        DO 5   ij   = 1, ip1jmp1       !     .....  correction pour  du(iip1, j, l)  ........                    
108        wsur2( ij ) = - 0.5 * w( ij,l+1 )       !     .....     du(iip1, j, l)= du(1, j, l)   .....                        
109     5  CONTINUE  
110         !DIR$ IVDEP                                                            
111         DO ij = iip1 + iip1, ip1jm, iip1
112  c     .....................     calcul pour  du     ..................          du(ij, l) = du(ij-iim, l)
113            du(ij, l+1) = du(ij-iim, l+1)
114        DO 6 ij = iip2 ,ip1jm-1       END DO
115        ww        = wsur2 (  ij  )     + wsur2( ij+1 )  
116        uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )       !     .................    calcul pour   dv      .....................  
117        du(ij,l)  = du(ij,l)   - ww * ( uu - uav(ij, l ) )/massebx(ij, l )  
118        du(ij,l+1)= du(ij,l+1) + ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)       DO ij = 1, ip1jm
119     6  CONTINUE          ww = wsur2(ij+iip1) + wsur2(ij)
120            vv = 0.5*(vcov(ij, l)+vcov(ij, l+1))
121  c     .....  correction pour  du(iip1,j,l)  ........          dv(ij, l) = dv(ij, l) - ww*(vv-vav(ij, l))/masseby(ij, l)
122  c     .....     du(iip1,j,l)= du(1,j,l)   .....          dv(ij, l+1) = dv(ij, l+1) + ww*(vv-vav(ij, l+1))/masseby(ij, l+1)
123         END DO
124  CDIR$ IVDEP  
125        DO   7  ij   = iip1 +iip1, ip1jm, iip1  
126        du( ij, l  ) = du( ij -iim, l  )  
127        du( ij,l+1 ) = du( ij -iim,l+1 )       !     ............................................................      
128     7  CONTINUE       !     ...............    calcul pour   dh      ...................      
129         !     ............................................................      
130  c     .................    calcul pour   dv      .....................  
131         !                       ---z                                            
132        DO 8 ij = 1, ip1jm       !       calcul de  - d( teta  * w )      qu'on ajoute a   dh            
133        ww        = wsur2( ij+iip1 )   + wsur2( ij )       !                   ...............                                    
134        vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )  
135        dv(ij,l)  = dv(ij, l ) - ww * (vv - vav(ij, l ) )/masseby(ij, l )       DO ij = 1, ip1jmp1
136        dv(ij,l+1)= dv(ij,l+1) + ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)          ww = wsur2(ij)*(teta(ij, l)+teta(ij, l+1))
137     8  CONTINUE          dteta(ij, l) = dteta(ij, l) - ww
138            dteta(ij, l+1) = dteta(ij, l+1) + ww
139  c       end DO
140    
141  c     ............................................................       IF (conser) THEN
142  c     ...............    calcul pour   dh      ...................          DO ij = 1, ip1jmp1
143  c     ............................................................             ge(ij) = wsur2(ij)*wsur2(ij)*unsaire2(ij)
144            end DO
145  c                       ---z          gt = ssum(ip1jmp1, ge, 1)
146  c       calcul de  - d( teta  * w )      qu'on ajoute a   dh          gtot(l) = deuxjour*sqrt(gt/ip1jmp1)
147  c                   ...............       END IF
148    
149          DO 15 ij = 1, ip1jmp1    END DO
150           ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )  
151           dteta(ij, l ) = dteta(ij, l )  -  ww  END SUBROUTINE advect
          dteta(ij,l+1) = dteta(ij,l+1)  +  ww  
   15    CONTINUE  
   
       IF( conser)  THEN  
         DO 17 ij = 1,ip1jmp1  
         ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)  
   17    CONTINUE  
         gt       = SSUM( ip1jmp1,ge,1 )  
         gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )  
       END IF  
   
   20  CONTINUE  
   
       RETURN  
       END  

Legend:
Removed from v.3  
changed lines
  Added in v.43

  ViewVC Help
Powered by ViewVC 1.1.21