/[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 40 by guez, Tue Feb 22 13:49:36 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), teta(ip1jmp1, llm)
24        REAL q(:, :, :, :) ! (iim + 1, jjm + 1, llm, nq)
25        REAL ps(ip1jmp1), masse(ip1jmp1, llm)
26    
27        REAL vcovm1(ip1jm, llm), ucovm1(ip1jmp1, llm)
28        REAL tetam1(ip1jmp1, llm), psm1(ip1jmp1), massem1(ip1jmp1, llm)
29    
30        REAL dv(ip1jm, llm), du(ip1jmp1, llm)
31        REAL dteta(ip1jmp1, llm), dp(ip1jmp1)
32        REAL finvmaold(ip1jmp1, llm)
33        LOGICAL, INTENT (IN) :: leapf
34        real, intent(in):: dt
35    
36        ! Local:
37    
38        INTEGER nq
39        REAL vscr(ip1jm), uscr(ip1jmp1), hscr(ip1jmp1), pscr(ip1jmp1)
40        REAL massescr(ip1jmp1, llm), finvmasse(ip1jmp1, llm)
41        REAL p(ip1jmp1, llmp1)
42        REAL tpn, tps, tppn(iim), tpps(iim)
43        REAL qpn, qps, qppn(iim), qpps(iim)
44        REAL deltap(ip1jmp1, llm)
45    
46        INTEGER l, ij, iq
47    
48        REAL ssum
49    
50        !-----------------------------------------------------------------------
51    
52        call assert(size(q, 1) == iim + 1, size(q, 2) == jjm + 1, &
53             size(q, 3) == llm, "integrd")
54        nq = size(q, 4)
55    
56        DO l = 1, llm
57           DO ij = 1, iip1
58              ucov(ij, l) = 0.
59              ucov(ij+ip1jm, l) = 0.
60              uscr(ij) = 0.
61              uscr(ij+ip1jm) = 0.
62           END DO
63        END DO
64    
65        ! integration de ps
66    
67        CALL scopy(ip1jmp1*llm, masse, 1, massescr, 1)
68    
69        DO ij = 1, ip1jmp1
70           pscr(ij) = ps(ij)
71           ps(ij) = psm1(ij) + dt*dp(ij)
72        END DO
73    
74        DO ij = 1, ip1jmp1
75           IF (ps(ij)<0.) THEN
76              PRINT *, ' Au point ij = ', ij, ' , pression sol neg. ', ps(ij)
77              STOP 'integrd'
78           END IF
79        END DO
80    
81        DO ij = 1, iim
82           tppn(ij) = aire(ij)*ps(ij)
83           tpps(ij) = aire(ij+ip1jm)*ps(ij+ip1jm)
84        END DO
85        tpn = ssum(iim, tppn, 1)/apoln
86        tps = ssum(iim, tpps, 1)/apols
87        DO ij = 1, iip1
88           ps(ij) = tpn
89           ps(ij+ip1jm) = tps
90        END DO
91    
92        ! Calcul de la nouvelle masse d'air au dernier temps integre t+1
93    
94        forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
95        CALL massdair(p, masse)
96    
97        CALL scopy(ijp1llm, masse, 1, finvmasse, 1)
98        CALL filtreg(finvmasse, jjp1, llm, -2, 2, .TRUE., 1)
99    
100        ! integration de ucov, vcov, h
101    
102        DO l = 1, llm
103           DO ij = iip2, ip1jm
104              uscr(ij) = ucov(ij, l)
105              ucov(ij, l) = ucovm1(ij, l) + dt*du(ij, l)
106           END DO
107    
108           DO ij = 1, ip1jm
109              vscr(ij) = vcov(ij, l)
110              vcov(ij, l) = vcovm1(ij, l) + dt*dv(ij, l)
111           END DO
112    
113           DO ij = 1, ip1jmp1
114              hscr(ij) = teta(ij, l)
115              teta(ij, l) = tetam1(ij, l)*massem1(ij, l)/masse(ij, l) + &
116                   dt*dteta(ij, l)/masse(ij, l)
117           END DO
118    
119           ! Calcul de la valeur moyenne, unique aux poles pour teta
120    
121           DO ij = 1, iim
122              tppn(ij) = aire(ij)*teta(ij, l)
123              tpps(ij) = aire(ij+ip1jm)*teta(ij+ip1jm, l)
124           END DO
125           tpn = ssum(iim, tppn, 1)/apoln
126           tps = ssum(iim, tpps, 1)/apols
127    
128           DO ij = 1, iip1
129              teta(ij, l) = tpn
130              teta(ij+ip1jm, l) = tps
131           END DO
132    
133           IF (leapf) THEN
134              CALL scopy(ip1jmp1, uscr(1), 1, ucovm1(1, l), 1)
135              CALL scopy(ip1jm, vscr(1), 1, vcovm1(1, l), 1)
136              CALL scopy(ip1jmp1, hscr(1), 1, tetam1(1, l), 1)
137           END IF
138        END DO
139    
140        DO l = 1, llm
141           DO ij = 1, ip1jmp1
142              deltap(ij, l) = p(ij, l) - p(ij, l+1)
143           END DO
144        END DO
145    
146        CALL qminimum(q, nq, deltap)
147    
148        ! Calcul de la valeur moyenne, unique aux poles pour q
149    
150        DO iq = 1, nq
151           DO l = 1, llm
152              DO ij = 1, iim
153                 qppn(ij) = aire(ij)*q(ij, 1, l, iq)
154                 qpps(ij) = aire(ij+ip1jm)*q(ij, jjm + 1, l, iq)
155              END DO
156              qpn = ssum(iim, qppn, 1)/apoln
157              qps = ssum(iim, qpps, 1)/apols
158    
159              DO ij = 1, iip1
160                 q(ij, 1, l, iq) = qpn
161                 q(ij, jjm + 1, l, iq) = qps
162              END DO
163           END DO
164        END DO
165    
166        CALL scopy(ijp1llm, finvmasse, 1, finvmaold, 1)
167    
168        ! Fin de l'integration de q
169    
170        IF (leapf) THEN
171           CALL scopy(ip1jmp1, pscr, 1, psm1, 1)
172           CALL scopy(ip1jmp1*llm, massescr, 1, massem1, 1)
173        END IF
174    
175      END SUBROUTINE integrd
176    
177    end module integrd_m

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

  ViewVC Help
Powered by ViewVC 1.1.21