source: LMDZ5/branches/LMDZ6_rc0/libf/dyn3d/integrd.F @ 5394

Last change on this file since 5394 was 2160, checked in by Laurent Fairhead, 10 years ago

Merged trunk changes -r2070:2158 into testing branch. Compilation problems introduced by revision r2155 have been corrected by hand

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.0 KB
Line 
1!
2! $Id: integrd.F 2160 2014-11-28 15:36:29Z evignon $
3!
4      SUBROUTINE integrd
5     $  (  nq,vcovm1,ucovm1,tetam1,psm1,massem1,
6     $     dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis !,finvmaold
7     &  )
8
9      use control_mod, only : planet_type
10
11      IMPLICIT NONE
12
13
14c=======================================================================
15c
16c   Auteur:  P. Le Van
17c   -------
18c
19c   objet:
20c   ------
21c
22c   Incrementation des tendances dynamiques
23c
24c=======================================================================
25c-----------------------------------------------------------------------
26c   Declarations:
27c   -------------
28
29#include "dimensions.h"
30#include "paramet.h"
31#include "comconst.h"
32#include "comgeom.h"
33#include "comvert.h"
34#include "logic.h"
35#include "temps.h"
36#include "serre.h"
37#include "iniprint.h"
38
39c   Arguments:
40c   ----------
41
42      integer,intent(in) :: nq ! number of tracers to handle in this routine
43      real,intent(inout) :: vcov(ip1jm,llm) ! covariant meridional wind
44      real,intent(inout) :: ucov(ip1jmp1,llm) ! covariant zonal wind
45      real,intent(inout) :: teta(ip1jmp1,llm) ! potential temperature
46      real,intent(inout) :: q(ip1jmp1,llm,nq) ! advected tracers
47      real,intent(inout) :: ps(ip1jmp1) ! surface pressure
48      real,intent(inout) :: masse(ip1jmp1,llm) ! atmospheric mass
49      real,intent(in) :: phis(ip1jmp1) ! ground geopotential !!! unused
50      ! values at previous time step
51      real,intent(inout) :: vcovm1(ip1jm,llm)
52      real,intent(inout) :: ucovm1(ip1jmp1,llm)
53      real,intent(inout) :: tetam1(ip1jmp1,llm)
54      real,intent(inout) :: psm1(ip1jmp1)
55      real,intent(inout) :: massem1(ip1jmp1,llm)
56      ! the tendencies to add
57      real,intent(in) :: dv(ip1jm,llm)
58      real,intent(in) :: du(ip1jmp1,llm)
59      real,intent(in) :: dteta(ip1jmp1,llm)
60      real,intent(in) :: dp(ip1jmp1)
61      real,intent(in) :: dq(ip1jmp1,llm,nq) !!! unused
62!      real,intent(out) :: finvmaold(ip1jmp1,llm) !!! unused
63
64c   Local:
65c   ------
66
67      REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
68      REAL massescr( ip1jmp1,llm )
69!      REAL finvmasse(ip1jmp1,llm)
70      REAL p(ip1jmp1,llmp1)
71      REAL tpn,tps,tppn(iim),tpps(iim)
72      REAL qpn,qps,qppn(iim),qpps(iim)
73      REAL deltap( ip1jmp1,llm )
74
75      INTEGER  l,ij,iq,i,j
76
77      REAL SSUM
78
79c-----------------------------------------------------------------------
80
81      DO  l = 1,llm
82        DO  ij = 1,iip1
83         ucov(    ij    , l) = 0.
84         ucov( ij +ip1jm, l) = 0.
85         uscr(     ij      ) = 0.
86         uscr( ij +ip1jm   ) = 0.
87        ENDDO
88      ENDDO
89
90
91c    ............    integration  de       ps         ..............
92
93      CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
94
95      DO ij = 1,ip1jmp1
96       pscr (ij)    = ps(ij)
97       ps (ij)      = psm1(ij) + dt * dp(ij)
98      ENDDO
99c
100      DO ij = 1,ip1jmp1
101        IF( ps(ij).LT.0. ) THEN
102         write(lunout,*) "integrd: negative surface pressure ",ps(ij)
103         write(lunout,*) " at node ij =", ij
104         ! since ij=j+(i-1)*jjp1 , we have
105         j=modulo(ij,jjp1)
106         i=1+(ij-j)/jjp1
107         write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg",
108     &                   " lat = ",rlatu(j)*180./pi, " deg"
109         call abort_gcm("integrd", "", 1)
110        ENDIF
111      ENDDO
112c
113      DO  ij    = 1, iim
114       tppn(ij) = aire(   ij   ) * ps(  ij    )
115       tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
116      ENDDO
117       tpn      = SSUM(iim,tppn,1)/apoln
118       tps      = SSUM(iim,tpps,1)/apols
119      DO ij   = 1, iip1
120       ps(   ij   )  = tpn
121       ps(ij+ip1jm)  = tps
122      ENDDO
123c
124c  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...
125c
126      CALL pression ( ip1jmp1, ap, bp, ps, p )
127      CALL massdair (     p  , masse         )
128
129! Ehouarn : we don't use/need finvmaold and finvmasse,
130!           so might as well not compute them
131!      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
132!      CALL filtreg( finvmasse, jjp1, llm, -2, 2, .TRUE., 1  )
133c
134
135c    ............   integration  de  ucov, vcov,  h     ..............
136
137      DO l = 1,llm
138
139       DO ij = iip2,ip1jm
140        uscr( ij )   =  ucov( ij,l )
141        ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
142       ENDDO
143
144       DO ij = 1,ip1jm
145        vscr( ij )   =  vcov( ij,l )
146        vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
147       ENDDO
148
149       DO ij = 1,ip1jmp1
150        hscr( ij )    =  teta(ij,l)
151        teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l)
152     &                + dt * dteta(ij,l) / masse(ij,l)
153       ENDDO
154
155c   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......
156c
157c
158       DO  ij   = 1, iim
159        tppn(ij) = aire(   ij   ) * teta(  ij    ,l)
160        tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
161       ENDDO
162        tpn      = SSUM(iim,tppn,1)/apoln
163        tps      = SSUM(iim,tpps,1)/apols
164
165       DO ij   = 1, iip1
166        teta(   ij   ,l)  = tpn
167        teta(ij+ip1jm,l)  = tps
168       ENDDO
169c
170
171       IF(leapf)  THEN
172         CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
173         CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
174         CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
175       END IF
176
177      ENDDO ! of DO l = 1,llm
178
179
180c
181c   .......  integration de   q   ......
182c
183c$$$      IF( iadv(1).NE.3.AND.iadv(2).NE.3 )    THEN
184c$$$c
185c$$$       IF( forward. OR . leapf )  THEN
186c$$$        DO iq = 1,2
187c$$$        DO  l = 1,llm
188c$$$        DO ij = 1,ip1jmp1
189c$$$        q(ij,l,iq) = ( q(ij,l,iq)*finvmaold(ij,l) + dtvr *dq(ij,l,iq) )/
190c$$$     $                            finvmasse(ij,l)
191c$$$        ENDDO
192c$$$        ENDDO
193c$$$        ENDDO
194c$$$       ELSE
195c$$$         DO iq = 1,2
196c$$$         DO  l = 1,llm
197c$$$         DO ij = 1,ip1jmp1
198c$$$         q( ij,l,iq ) = q( ij,l,iq ) * finvmaold(ij,l) / finvmasse(ij,l)
199c$$$         ENDDO
200c$$$         ENDDO
201c$$$         ENDDO
202c$$$
203c$$$       END IF
204c$$$c
205c$$$      ENDIF
206
207      if (planet_type.eq."earth") then
208! Earth-specific treatment of first 2 tracers (water)
209        DO l = 1, llm
210          DO ij = 1, ip1jmp1
211            deltap(ij,l) =  p(ij,l) - p(ij,l+1)
212          ENDDO
213        ENDDO
214
215        CALL qminimum( q, nq, deltap )
216
217c
218c    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
219c
220
221       DO iq = 1, nq
222        DO l = 1, llm
223
224           DO ij = 1, iim
225             qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)
226             qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
227           ENDDO
228             qpn  =  SSUM(iim,qppn,1)/apoln
229             qps  =  SSUM(iim,qpps,1)/apols
230
231           DO ij = 1, iip1
232             q(   ij   ,l,iq)  = qpn
233             q(ij+ip1jm,l,iq)  = qps
234           ENDDO
235
236        ENDDO
237       ENDDO
238
239! Ehouarn: forget about finvmaold
240!      CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
241
242      endif ! of if (planet_type.eq."earth")
243c
244c
245c     .....   FIN  de l'integration  de   q    .......
246
247c    .................................................................
248
249
250      IF( leapf )  THEN
251         CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )
252         CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )
253      END IF
254
255      RETURN
256      END
Note: See TracBrowser for help on using the repository browser.