source: LMDZ.3.3/trunk/libf/dyn3d/advect.F @ 343

Last change on this file since 343 was 2, checked in by lmdz, 25 years ago

Initial revision
  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 4.5 KB
Line 
1      SUBROUTINE advect(ucov,vcov,teta,w,massebx,masseby,du,dv,dteta)
2
3      IMPLICIT NONE
4c=======================================================================
5c
6c   Auteurs:  P. Le Van , Fr. Hourdin  .
7c   -------
8c
9c   Objet:
10c   ------
11c
12c   *************************************************************
13c   .... calcul des termes d'advection vertic.pour u,v,teta,q ...
14c   *************************************************************
15c        ces termes sont ajoutes a du,dv,dteta et dq .
16c  Modif F.Forget 03/94 : on retire q de advect
17c
18c=======================================================================
19c-----------------------------------------------------------------------
20c   Declarations:
21c   -------------
22
23#include "dimensions.h"
24#include "paramet.h"
25#include "comconst.h"
26#include "comvert.h"
27#include "comgeom.h"
28#include "logic.h"
29#include "ener.h"
30
31c   Arguments:
32c   ----------
33
34      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)
35      REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm),w(ip1jmp1,llm)
36      REAL dv(ip1jm,llm),du(ip1jmp1,llm),dteta(ip1jmp1,llm)
37
38c   Local:
39c   ------
40
41      REAL uav(ip1jmp1,llm),vav(ip1jm,llm),wsur2(ip1jmp1)
42      REAL unsaire2(ip1jmp1), ge(ip1jmp1)
43      REAL deuxjour, ww, gt, uu, vv
44
45      INTEGER  ij,l
46
47      EXTERNAL  SSUM
48      REAL      SSUM
49
50c-----------------------------------------------------------------------
51c   2. Calculs preliminaires:
52c   -------------------------
53
54      IF (conser)  THEN
55         deuxjour = 2. * daysec
56
57         DO   1  ij   = 1, ip1jmp1
58         unsaire2(ij) = unsaire(ij) * unsaire(ij)
59   1     CONTINUE
60      END IF
61
62
63c------------------  -yy ----------------------------------------------
64c   .  Calcul de     u
65
66      DO  l=1,llm
67         DO    ij     = iip2, ip1jmp1
68            uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )
69         ENDDO
70         DO    ij     = iip2, ip1jm
71            uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)
72         ENDDO
73         DO      ij         = 1, iip1
74            uav(ij      ,l) = 0.
75            uav(ip1jm+ij,l) = 0.
76         ENDDO
77      ENDDO
78
79c------------------  -xx ----------------------------------------------
80c   .  Calcul de     v
81
82      DO  l=1,llm
83         DO    ij   = 2, ip1jm
84          vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )
85         ENDDO
86         DO    ij   = 1,ip1jm,iip1
87          vav(ij,l) = vav(ij+iim,l)
88         ENDDO
89         DO    ij   = 1, ip1jm-1
90          vav(ij,l) = vav(ij,l) + vav(ij+1,l)
91         ENDDO
92         DO    ij       = 1, ip1jm, iip1
93          vav(ij+iim,l) = vav(ij,l)
94         ENDDO
95      ENDDO
96
97c-----------------------------------------------------------------------
98
99c
100      DO 20 l = 1, llmm1
101
102
103c       ......   calcul de  - w/2.    au niveau  l+1   .......
104
105      DO 5   ij   = 1, ip1jmp1
106      wsur2( ij ) = - 0.5 * w( ij,l+1 )
107   5  CONTINUE
108
109
110c     .....................     calcul pour  du     ..................
111
112      DO 6 ij = iip2 ,ip1jm-1
113      ww        = wsur2 (  ij  )     + wsur2( ij+1 )
114      uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
115      du(ij,l)  = du(ij,l)   - ww * ( uu - uav(ij, l ) )/massebx(ij, l )
116      du(ij,l+1)= du(ij,l+1) + ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
117   6  CONTINUE
118
119c     .....  correction pour  du(iip1,j,l)  ........
120c     .....     du(iip1,j,l)= du(1,j,l)   .....
121
122CDIR$ IVDEP
123      DO   7  ij   = iip1 +iip1, ip1jm, iip1
124      du( ij, l  ) = du( ij -iim, l  )
125      du( ij,l+1 ) = du( ij -iim,l+1 )
126   7  CONTINUE
127
128c     .................    calcul pour   dv      .....................
129
130      DO 8 ij = 1, ip1jm
131      ww        = wsur2( ij+iip1 )   + wsur2( ij )
132      vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
133      dv(ij,l)  = dv(ij, l ) - ww * (vv - vav(ij, l ) )/masseby(ij, l )
134      dv(ij,l+1)= dv(ij,l+1) + ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
135   8  CONTINUE
136
137c
138
139c     ............................................................
140c     ...............    calcul pour   dh      ...................
141c     ............................................................
142
143c                       ---z
144c       calcul de  - d( teta  * w )      qu'on ajoute a   dh
145c                   ...............
146
147        DO 15 ij = 1, ip1jmp1
148         ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
149         dteta(ij, l ) = dteta(ij, l )  -  ww
150         dteta(ij,l+1) = dteta(ij,l+1)  +  ww
151  15    CONTINUE
152
153      IF( conser)  THEN
154        DO 17 ij = 1,ip1jmp1
155        ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)
156  17    CONTINUE
157        gt       = SSUM( ip1jmp1,ge,1 )
158        gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )
159      END IF
160
161  20  CONTINUE
162 
163      RETURN
164      END
Note: See TracBrowser for help on using the repository browser.