source: trunk/LMDZ.COMMON/libf/dyn3dpar/addfi_p.F @ 1189

Last change on this file since 1189 was 1189, checked in by emillour, 11 years ago

Common dynamics: a couple of bug fixes:

  • Correctly account for the change in pressure, mass, etc. after modifying surface pressure following a call to the physics.
  • Corrected tracer advection, which is computed using values at the beginning of the time step, so it is done at Matsuno forward step.

EM
File size: 6.5 KB
RevLine 
[1]1!
[7]2! $Id: addfi_p.F 1446 2010-10-22 09:27:25Z emillour $
[1]3!
4      SUBROUTINE addfi_p(pdt, leapf, forward,
5     S          pucov, pvcov, pteta, pq   , pps ,
6     S          pdufi, pdvfi, pdhfi,pdqfi, pdpfi  )
[1019]7      USE parallel_lmdz
[1]8      USE infotrac, ONLY : nqtot
[7]9      USE control_mod, ONLY : planet_type
[1]10      IMPLICIT NONE
11c
12c=======================================================================
13c
14c    Addition of the physical tendencies
15c
16c    Interface :
17c    -----------
18c
19c      Input :
20c      -------
21c      pdt                    time step of integration
22c      leapf                  logical
23c      forward                logical
24c      pucov(ip1jmp1,llm)     first component of the covariant velocity
25c      pvcov(ip1ip1jm,llm)    second component of the covariant velocity
26c      pteta(ip1jmp1,llm)     potential temperature
27c      pts(ip1jmp1,llm)       surface temperature
28c      pdufi(ip1jmp1,llm)     |
29c      pdvfi(ip1jm,llm)       |   respective
[108]30c      pdhfi(ip1jmp1)         |      tendencies  (unit/s)
[1]31c      pdtsfi(ip1jmp1)        |
32c
33c      Output :
34c      --------
35c      pucov
36c      pvcov
37c      ph
38c      pts
39c
40c
41c=======================================================================
42c
43c-----------------------------------------------------------------------
44c
45c    0.  Declarations :
46c    ------------------
47c
48#include "dimensions.h"
49#include "paramet.h"
50#include "comconst.h"
51#include "comgeom.h"
52#include "serre.h"
53c
54c    Arguments :
55c    -----------
56c
[1189]57      REAL,INTENT(IN) :: pdt ! time step for the integration (s)
[1]58c
[1189]59      REAL,INTENT(INOUT) :: pvcov(ip1jm,llm) ! covariant meridional wind
60      REAL,INTENT(INOUT) :: pucov(ip1jmp1,llm) ! covariant zonal wind
61      REAL,INTENT(INOUT) :: pteta(ip1jmp1,llm) ! potential temperature
62      REAL,INTENT(INOUT) :: pq(ip1jmp1,llm,nqtot) ! tracers
63      REAL,INTENT(INOUT) :: pps(ip1jmp1) ! surface pressure (Pa)
64c respective tendencies (.../s) to add
65      REAL,INTENT(IN) :: pdvfi(ip1jm,llm)
66      REAL,INTENT(IN) :: pdufi(ip1jmp1,llm)
67      REAL,INTENT(IN) :: pdqfi(ip1jmp1,llm,nqtot)
68      REAL,INTENT(IN) :: pdhfi(ip1jmp1,llm)
69      REAL,INTENT(IN) :: pdpfi(ip1jmp1)
[1]70c
[1189]71      LOGICAL,INTENT(IN) :: leapf,forward ! not used
[1]72c
73c
74c    Local variables :
75c    -----------------
76c
77      REAL xpn(iim),xps(iim),tpn,tps
78      INTEGER j,k,iq,ij
[1189]79      REAL,PARAMETER :: qtestw = 1.0e-15
80      REAL,PARAMETER :: qtestt = 1.0e-40
[1]81
82      REAL SSUM
83      EXTERNAL SSUM
84     
85      INTEGER :: ijb,ije
86c
87c-----------------------------------------------------------------------
88     
89      ijb=ij_begin
90      ije=ij_end
91     
92c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
93      DO k = 1,llm
94         DO j = ijb,ije
95            pteta(j,k)= pteta(j,k) + pdhfi(j,k) * pdt
96         ENDDO
97      ENDDO
98c$OMP END DO NOWAIT
99
100      if (pole_nord) then
101c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
102        DO  k    = 1, llm
103         DO  ij   = 1, iim
104           xpn(ij) = aire(   ij   ) * pteta(  ij    ,k)
105         ENDDO
106         tpn      = SSUM(iim,xpn,1)/ apoln
107
108         DO ij   = 1, iip1
109           pteta(   ij   ,k)  = tpn
110         ENDDO
111       ENDDO
112c$OMP END DO NOWAIT
113      endif
114
115      if (pole_sud) then
116c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
117        DO  k    = 1, llm
118         DO  ij   = 1, iim
119           xps(ij) = aire(ij+ip1jm) * pteta(ij+ip1jm,k)
120         ENDDO
121         tps      = SSUM(iim,xps,1)/ apols
122
123         DO ij   = 1, iip1
124           pteta(ij+ip1jm,k)  = tps
125         ENDDO
126       ENDDO
127c$OMP END DO NOWAIT
128      endif
129c
130
131      ijb=ij_begin
132      ije=ij_end
133      if (pole_nord) ijb=ij_begin+iip1
134      if (pole_sud)  ije=ij_end-iip1
135
136c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
137      DO k = 1,llm
138         DO j = ijb,ije
139            pucov(j,k)= pucov(j,k) + pdufi(j,k) * pdt
140         ENDDO
141      ENDDO
142c$OMP END DO NOWAIT
143
144      if (pole_nord) ijb=ij_begin
145
146c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
147      DO k = 1,llm
148         DO j = ijb,ije
149            pvcov(j,k)= pvcov(j,k) + pdvfi(j,k) * pdt
150         ENDDO
151      ENDDO
152c$OMP END DO NOWAIT
153
154c
155      if (pole_sud)  ije=ij_end
156c$OMP MASTER
157      DO j = ijb,ije
158         pps(j) = pps(j) + pdpfi(j) * pdt
159      ENDDO
160c$OMP END MASTER
161 
[7]162      if (planet_type=="earth") then
163      ! earth case, special treatment for first 2 tracers (water)
164       DO iq = 1, 2
[1]165c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
166         DO k = 1,llm
167            DO j = ijb,ije
168               pq(j,k,iq)= pq(j,k,iq) + pdqfi(j,k,iq) * pdt
169               pq(j,k,iq)= AMAX1( pq(j,k,iq), qtestw )
170            ENDDO
171         ENDDO
172c$OMP END DO NOWAIT
[7]173       ENDDO
[1]174
[7]175       DO iq = 3, nqtot
[1]176c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
177         DO k = 1,llm
178            DO j = ijb,ije
179               pq(j,k,iq)= pq(j,k,iq) + pdqfi(j,k,iq) * pdt
180               pq(j,k,iq)= AMAX1( pq(j,k,iq), qtestt )
181            ENDDO
182         ENDDO
183c$OMP END DO NOWAIT
[7]184       ENDDO
185      else
186      ! general case, treat all tracers equally)
187       DO iq = 1, nqtot
188c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
189         DO k = 1,llm
190            DO j = ijb,ije
191               pq(j,k,iq)= pq(j,k,iq) + pdqfi(j,k,iq) * pdt
192               pq(j,k,iq)= AMAX1( pq(j,k,iq), qtestt )
193            ENDDO
194         ENDDO
195c$OMP END DO NOWAIT
196       ENDDO
197      endif ! of if (planet_type=="earth")
[1]198
199c$OMP MASTER
200      if (pole_nord) then
201     
202        DO  ij   = 1, iim
203          xpn(ij) = aire(   ij   ) * pps(  ij     )
204        ENDDO
205
206        tpn      = SSUM(iim,xpn,1)/apoln
207
208        DO ij   = 1, iip1
209          pps (   ij     )  = tpn
210        ENDDO
211     
212      endif
213
214      if (pole_sud) then
215     
216        DO  ij   = 1, iim
217          xps(ij) = aire(ij+ip1jm) * pps(ij+ip1jm )
218        ENDDO
219
220        tps      = SSUM(iim,xps,1)/apols
221
222        DO ij   = 1, iip1
223          pps ( ij+ip1jm )  = tps
224        ENDDO
225     
226      endif
227c$OMP END MASTER
228
229      if (pole_nord) then
230        DO iq = 1, nqtot
231c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
232          DO  k    = 1, llm
233            DO  ij   = 1, iim
234              xpn(ij) = aire(   ij   ) * pq(  ij    ,k,iq)
235            ENDDO
236            tpn      = SSUM(iim,xpn,1)/apoln
237 
238            DO ij   = 1, iip1
239              pq (   ij   ,k,iq)  = tpn
240            ENDDO
241          ENDDO
242c$OMP END DO NOWAIT       
243        ENDDO
244      endif
245     
246      if (pole_sud) then
247        DO iq = 1, nqtot
248c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
249          DO  k    = 1, llm
250            DO  ij   = 1, iim
251              xps(ij) = aire(ij+ip1jm) * pq(ij+ip1jm,k,iq)
252            ENDDO
253            tps      = SSUM(iim,xps,1)/apols
254 
255            DO ij   = 1, iip1
256              pq (ij+ip1jm,k,iq)  = tps
257            ENDDO
258          ENDDO
259c$OMP END DO NOWAIT       
260        ENDDO
261      endif
262     
263     
264      RETURN
265      END
Note: See TracBrowser for help on using the repository browser.