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

Diff of /trunk/dyn3d/integrd.f

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

trunk/libf/dyn3d/integrd.f revision 18 by guez, Thu Aug 7 12:29:13 2008 UTC trunk/libf/dyn3d/integrd.f90 revision 44 by guez, Wed Apr 13 12:29:18 2011 UTC
# Line 1  Line 1 
1  !  module integrd_m
 ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/integrd.F,v 1.1.1.1 2004/05/19 12:53:05 lmdzadmin Exp $  
 !  
       SUBROUTINE integrd  
      $  (  nq,vcovm1,ucovm1,tetam1,psm1,massem1,  
      $     dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis,finvmaold,  
      $     leapf )  
   
       use dimens_m  
       use paramet_m  
       use comconst  
       use comvert  
       use logic  
       use comgeom  
       use serre  
       use temps  
       use iniadvtrac_m  
       use pression_m, only: pression  
   
       IMPLICIT NONE  
   
   
 c=======================================================================  
 c  
 c   Auteur:  P. Le Van  
 c   -------  
 c  
 c   objet:  
 c   ------  
 c  
 c   Incrementation des tendances dynamiques  
 c  
 c=======================================================================  
 c-----------------------------------------------------------------------  
 c   Declarations:  
 c   -------------  
   
   
 c   Arguments:  
 c   ----------  
   
       INTEGER nq  
   
       REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)  
       REAL q(ip1jmp1,llm,nq)  
       REAL ps(ip1jmp1),masse(ip1jmp1,llm),phis(ip1jmp1)  
   
       REAL vcovm1(ip1jm,llm),ucovm1(ip1jmp1,llm)  
       REAL tetam1(ip1jmp1,llm),psm1(ip1jmp1),massem1(ip1jmp1,llm)  
   
       REAL dv(ip1jm,llm),du(ip1jmp1,llm)  
       REAL dteta(ip1jmp1,llm),dp(ip1jmp1)  
       REAL dq(ip1jmp1,llm,nq), finvmaold(ip1jmp1,llm)  
       logical, intent(in):: leapf  
   
 c   Local:  
 c   ------  
   
       REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)  
       REAL massescr( ip1jmp1,llm ), finvmasse(ip1jmp1,llm)  
       REAL p(ip1jmp1,llmp1)  
       REAL tpn,tps,tppn(iim),tpps(iim)  
       REAL qpn,qps,qppn(iim),qpps(iim)  
       REAL deltap( ip1jmp1,llm )  
   
       INTEGER  l,ij,iq  
   
       REAL SSUM  
   
 c-----------------------------------------------------------------------  
   
       DO  l = 1,llm  
         DO  ij = 1,iip1  
          ucov(    ij    , l) = 0.  
          ucov( ij +ip1jm, l) = 0.  
          uscr(     ij      ) = 0.  
          uscr( ij +ip1jm   ) = 0.  
         ENDDO  
       ENDDO  
   
   
 c    ............    integration  de       ps         ..............  
   
       CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)  
   
       DO 2 ij = 1,ip1jmp1  
        pscr (ij)    = ps(ij)  
        ps (ij)      = psm1(ij) + dt * dp(ij)  
    2  CONTINUE  
 c  
       DO ij = 1,ip1jmp1  
         IF( ps(ij).LT.0. ) THEN  
          PRINT *,' Au point ij = ',ij, ' , pression sol neg. ', ps(ij)  
          STOP 'integrd'  
         ENDIF  
       ENDDO  
 c  
       DO  ij    = 1, iim  
        tppn(ij) = aire(   ij   ) * ps(  ij    )  
        tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)  
       ENDDO  
        tpn      = SSUM(iim,tppn,1)/apoln  
        tps      = SSUM(iim,tpps,1)/apols  
       DO ij   = 1, iip1  
        ps(   ij   )  = tpn  
        ps(ij+ip1jm)  = tps  
       ENDDO  
 c  
 c  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...  
 c  
       CALL pression ( ip1jmp1, ap, bp, ps, p )  
       CALL massdair (     p  , masse         )  
   
       CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )  
       CALL filtreg( finvmasse, jjp1, llm, -2, 2, .TRUE., 1  )  
 c  
   
 c    ............   integration  de  ucov, vcov,  h     ..............  
   
       DO 10 l = 1,llm  
   
       DO 4 ij = iip2,ip1jm  
       uscr( ij )   =  ucov( ij,l )  
       ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )  
    4  CONTINUE  
   
       DO 5 ij = 1,ip1jm  
       vscr( ij )   =  vcov( ij,l )  
       vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )  
    5  CONTINUE  
   
       DO 6 ij = 1,ip1jmp1  
       hscr( ij )    =  teta(ij,l)  
       teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l)  
      $                + dt * dteta(ij,l) / masse(ij,l)  
    6  CONTINUE  
   
 c   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......  
 c  
 c  
       DO  ij   = 1, iim  
         tppn(ij) = aire(   ij   ) * teta(  ij    ,l)  
         tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)  
       ENDDO  
         tpn      = SSUM(iim,tppn,1)/apoln  
         tps      = SSUM(iim,tpps,1)/apols  
   
       DO ij   = 1, iip1  
         teta(   ij   ,l)  = tpn  
         teta(ij+ip1jm,l)  = tps  
       ENDDO  
 c  
   
       IF(leapf)  THEN  
          CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )  
          CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )  
          CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )  
       END IF  
   
   10  CONTINUE  
   
          DO l = 1, llm  
           DO ij = 1, ip1jmp1  
            deltap(ij,l) =  p(ij,l) - p(ij,l+1)  
           ENDDO  
          ENDDO  
   
          CALL qminimum( q, nq, deltap )  
 c  
 c    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....  
 c  
   
       DO iq = 1, nq  
         DO l = 1, llm  
   
            DO ij = 1, iim  
              qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)  
              qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)  
            ENDDO  
              qpn  =  SSUM(iim,qppn,1)/apoln  
              qps  =  SSUM(iim,qpps,1)/apols  
   
            DO ij = 1, iip1  
              q(   ij   ,l,iq)  = qpn  
              q(ij+ip1jm,l,iq)  = qps  
            ENDDO  
   
         ENDDO  
       ENDDO  
   
   
          CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )  
 c  
 c  
 c     .....   FIN  de l'integration  de   q    .......  
   
 15    continue  
   
 c    .................................................................  
   
   
       IF( leapf )  THEN  
          CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )  
          CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )  
       END IF  
2    
3        RETURN    IMPLICIT NONE
4        END  
5    contains
6    
7      SUBROUTINE integrd(vcovm1, ucovm1, tetam1, psm1, massem1, dv, du, &
8           dteta, dp, vcov, ucov, teta, q, ps, masse, finvmaold, dt, leapf)
9    
10        ! From dyn3d/integrd.F, version 1.1.1.1 2004/05/19 12:53:05
11        ! Auteur: P. Le Van
12        ! Objet: incrémentation des tendances dynamiques
13    
14        USE comvert, ONLY : ap, bp
15        USE comgeom, ONLY : aire, apoln, apols
16        USE dimens_m, ONLY : iim, jjm, llm
17        USE filtreg_m, ONLY : filtreg
18        use nr_util, only: assert
19        USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, jjp1, llmp1
20    
21        ! Arguments:
22    
23        REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm)
24        real, intent(inout):: teta(ip1jmp1, llm)
25        REAL q(:, :, :, :) ! (iim + 1, jjm + 1, llm, nq)
26        REAL, intent(inout):: ps(ip1jmp1)
27        REAL masse(ip1jmp1, llm)
28    
29        REAL vcovm1(ip1jm, llm), ucovm1(ip1jmp1, llm)
30        REAL tetam1(ip1jmp1, llm), psm1(ip1jmp1), massem1(ip1jmp1, llm)
31    
32        REAL dv(ip1jm, llm), du(ip1jmp1, llm)
33        REAL dteta(ip1jmp1, llm), dp(ip1jmp1)
34        REAL finvmaold(ip1jmp1, llm)
35        LOGICAL, INTENT (IN) :: leapf
36        real, intent(in):: dt
37    
38        ! Local:
39    
40        INTEGER nq
41        REAL vscr(ip1jm), uscr(ip1jmp1), hscr(ip1jmp1), pscr(ip1jmp1)
42        REAL massescr(ip1jmp1, llm), finvmasse(ip1jmp1, llm)
43        REAL p(ip1jmp1, llmp1)
44        REAL tpn, tps, tppn(iim), tpps(iim)
45        REAL qpn, qps, qppn(iim), qpps(iim)
46        REAL deltap(ip1jmp1, llm)
47    
48        INTEGER l, ij, iq
49    
50        REAL ssum
51    
52        !-----------------------------------------------------------------------
53    
54        call assert(size(q, 1) == iim + 1, size(q, 2) == jjm + 1, &
55             size(q, 3) == llm, "integrd")
56        nq = size(q, 4)
57    
58        DO l = 1, llm
59           DO ij = 1, iip1
60              ucov(ij, l) = 0.
61              ucov(ij+ip1jm, l) = 0.
62              uscr(ij) = 0.
63              uscr(ij+ip1jm) = 0.
64           END DO
65        END DO
66    
67        ! integration de ps
68    
69        CALL scopy(ip1jmp1*llm, masse, 1, massescr, 1)
70    
71        DO ij = 1, ip1jmp1
72           pscr(ij) = ps(ij)
73           ps(ij) = psm1(ij) + dt*dp(ij)
74        END DO
75    
76        DO ij = 1, ip1jmp1
77           IF (ps(ij)<0.) THEN
78              PRINT *, ' Au point ij = ', ij, ' , pression sol neg. ', ps(ij)
79              STOP 'integrd'
80           END IF
81        END DO
82    
83        DO ij = 1, iim
84           tppn(ij) = aire(ij)*ps(ij)
85           tpps(ij) = aire(ij+ip1jm)*ps(ij+ip1jm)
86        END DO
87        tpn = ssum(iim, tppn, 1)/apoln
88        tps = ssum(iim, tpps, 1)/apols
89        DO ij = 1, iip1
90           ps(ij) = tpn
91           ps(ij+ip1jm) = tps
92        END DO
93    
94        ! Calcul de la nouvelle masse d'air au dernier temps integre t+1
95    
96        forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
97        CALL massdair(p, masse)
98    
99        CALL scopy(ijp1llm, masse, 1, finvmasse, 1)
100        CALL filtreg(finvmasse, jjp1, llm, -2, 2, .TRUE., 1)
101    
102        ! integration de ucov, vcov, h
103    
104        DO l = 1, llm
105           DO ij = iip2, ip1jm
106              uscr(ij) = ucov(ij, l)
107              ucov(ij, l) = ucovm1(ij, l) + dt*du(ij, l)
108           END DO
109    
110           DO ij = 1, ip1jm
111              vscr(ij) = vcov(ij, l)
112              vcov(ij, l) = vcovm1(ij, l) + dt*dv(ij, l)
113           END DO
114    
115           DO ij = 1, ip1jmp1
116              hscr(ij) = teta(ij, l)
117              teta(ij, l) = tetam1(ij, l) * massem1(ij, l) / masse(ij, l) &
118                   + dt * dteta(ij, l) / masse(ij, l)
119           END DO
120    
121           ! Calcul de la valeur moyenne, unique aux poles pour teta
122    
123           DO ij = 1, iim
124              tppn(ij) = aire(ij)*teta(ij, l)
125              tpps(ij) = aire(ij+ip1jm)*teta(ij+ip1jm, l)
126           END DO
127           tpn = ssum(iim, tppn, 1)/apoln
128           tps = ssum(iim, tpps, 1)/apols
129    
130           DO ij = 1, iip1
131              teta(ij, l) = tpn
132              teta(ij+ip1jm, l) = tps
133           END DO
134    
135           IF (leapf) THEN
136              CALL scopy(ip1jmp1, uscr(1), 1, ucovm1(1, l), 1)
137              CALL scopy(ip1jm, vscr(1), 1, vcovm1(1, l), 1)
138              CALL scopy(ip1jmp1, hscr(1), 1, tetam1(1, l), 1)
139           END IF
140        END DO
141    
142        DO l = 1, llm
143           DO ij = 1, ip1jmp1
144              deltap(ij, l) = p(ij, l) - p(ij, l+1)
145           END DO
146        END DO
147    
148        CALL qminimum(q, nq, deltap)
149    
150        ! Calcul de la valeur moyenne, unique aux poles pour q
151    
152        DO iq = 1, nq
153           DO l = 1, llm
154              DO ij = 1, iim
155                 qppn(ij) = aire(ij)*q(ij, 1, l, iq)
156                 qpps(ij) = aire(ij+ip1jm)*q(ij, jjm + 1, l, iq)
157              END DO
158              qpn = ssum(iim, qppn, 1)/apoln
159              qps = ssum(iim, qpps, 1)/apols
160    
161              DO ij = 1, iip1
162                 q(ij, 1, l, iq) = qpn
163                 q(ij, jjm + 1, l, iq) = qps
164              END DO
165           END DO
166        END DO
167    
168        CALL scopy(ijp1llm, finvmasse, 1, finvmaold, 1)
169    
170        ! Fin de l'integration de q
171    
172        IF (leapf) THEN
173           CALL scopy(ip1jmp1, pscr, 1, psm1, 1)
174           CALL scopy(ip1jmp1*llm, massescr, 1, massem1, 1)
175        END IF
176    
177      END SUBROUTINE integrd
178    
179    end module integrd_m

Legend:
Removed from v.18  
changed lines
  Added in v.44

  ViewVC Help
Powered by ViewVC 1.1.21