GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: dyn3d/integrd.F Lines: 65 72 90.3 %
Date: 2023-06-30 12:56:34 Branches: 42 44 95.5 %

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