trunk/dyn3d/integrd.f

module integrd_m

  IMPLICIT NONE

contains

  SUBROUTINE integrd(nq,vcovm1,ucovm1,tetam1,psm1,massem1,dv,du,dteta,dp, &
       vcov,ucov,teta,q,ps,masse,finvmaold,leapf, dt)

    ! From dyn3d/integrd.F,v 1.1.1.1 2004/05/19 12:53:05
    !   Auteur:  P. Le Van                                                  
    !   objet:                                                              
    !   Incrementation des tendances dynamiques                             

    USE dimens_m, ONLY : iim, llm
    USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, jjp1, llmp1
    USE comvert, ONLY : ap, bp
    USE comgeom, ONLY : aire, apoln, apols
    USE filtreg_m, ONLY : filtreg

    !   Arguments:                                                          

    INTEGER, intent(in):: nq

    REAL vcov(ip1jm,llm), ucov(ip1jmp1,llm), teta(ip1jmp1,llm)
    REAL q(ip1jmp1,llm,nq)
    REAL ps(ip1jmp1), masse(ip1jmp1,llm)

    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 finvmaold(ip1jmp1,llm)
    LOGICAL, INTENT (IN) :: leapf
    real, intent(in):: dt

    !   Local:                                                              

    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

    !-----------------------------------------------------------------------

    DO l = 1, llm
       DO ij = 1, iip1
          ucov(ij,l) = 0.
          ucov(ij+ip1jm,l) = 0.
          uscr(ij) = 0.
          uscr(ij+ip1jm) = 0.
       END DO
    END DO


    !    ............    integration  de       ps         ..............    

    CALL scopy(ip1jmp1*llm,masse,1,massescr,1)

    DO ij = 1, ip1jmp1
       pscr(ij) = ps(ij)
       ps(ij) = psm1(ij) + dt*dp(ij)
    END DO

    DO ij = 1, ip1jmp1
       IF (ps(ij)<0.) THEN
          PRINT *, ' Au point ij = ', ij, ' , pression sol neg. ', ps(ij)
          STOP 'integrd'
       END IF
    END DO

    DO ij = 1, iim
       tppn(ij) = aire(ij)*ps(ij)
       tpps(ij) = aire(ij+ip1jm)*ps(ij+ip1jm)
    END DO
    tpn = ssum(iim,tppn,1)/apoln
    tps = ssum(iim,tpps,1)/apols
    DO ij = 1, iip1
       ps(ij) = tpn
       ps(ij+ip1jm) = tps
    END DO

    !  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 .

    forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
    CALL massdair(p,masse)

    CALL scopy(ijp1llm,masse,1,finvmasse,1)
    CALL filtreg(finvmasse,jjp1,llm,-2,2,.TRUE.,1)


    !    ............   integration  de  ucov, vcov,  h     ..............  

    DO  l = 1, llm

       DO ij = iip2, ip1jm
          uscr(ij) = ucov(ij,l)
          ucov(ij,l) = ucovm1(ij,l) + dt*du(ij,l)
       END DO

       DO ij = 1, ip1jm
          vscr(ij) = vcov(ij,l)
          vcov(ij,l) = vcovm1(ij,l) + dt*dv(ij,l)
       END DO

       DO 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)
       END DO

       !   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    .


       DO ij = 1, iim
          tppn(ij) = aire(ij)*teta(ij,l)
          tpps(ij) = aire(ij+ip1jm)*teta(ij+ip1jm,l)
       END DO
       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
       END DO


       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

    END DO

    DO l = 1, llm
       DO ij = 1, ip1jmp1
          deltap(ij,l) = p(ij,l) - p(ij,l+1)
       END DO
    END DO

    CALL qminimum(q,nq,deltap)

    !    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....


    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)
          END DO
          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
          END DO

       END DO
    END DO


    CALL scopy(ijp1llm,finvmasse,1,finvmaold,1)


    !     .....   FIN  de l'integration  de   q    .......                  

    IF (leapf) THEN
       CALL scopy(ip1jmp1,pscr,1,psm1,1)
       CALL scopy(ip1jmp1*llm,massescr,1,massem1,1)
    END IF

  END SUBROUTINE integrd

end module integrd_m
1	guez	32	module integrd_m
2	guez	3
3	guez	32	IMPLICIT NONE
4	guez	3
5	guez	32	contains
6	guez	3
7	guez	32	SUBROUTINE integrd(nq,vcovm1,ucovm1,tetam1,psm1,massem1,dv,du,dteta,dp, &
8			vcov,ucov,teta,q,ps,masse,finvmaold,leapf, dt)
9	guez	3
10	guez	32	! From dyn3d/integrd.F,v 1.1.1.1 2004/05/19 12:53:05
11			! Auteur: P. Le Van
12			! objet:
13			! Incrementation des tendances dynamiques
14	guez	3
15	guez	32	USE dimens_m, ONLY : iim, llm
16			USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, jjp1, llmp1
17			USE comvert, ONLY : ap, bp
18			USE comgeom, ONLY : aire, apoln, apols
19			USE filtreg_m, ONLY : filtreg
20	guez	3
21	guez	32	! Arguments:
22	guez	3
23	guez	32	INTEGER, intent(in):: nq
24	guez	3
25	guez	32	REAL vcov(ip1jm,llm), ucov(ip1jmp1,llm), teta(ip1jmp1,llm)
26			REAL q(ip1jmp1,llm,nq)
27			REAL ps(ip1jmp1), masse(ip1jmp1,llm)
28	guez	3
29	guez	32	REAL vcovm1(ip1jm,llm), ucovm1(ip1jmp1,llm)
30			REAL tetam1(ip1jmp1,llm), psm1(ip1jmp1), massem1(ip1jmp1,llm)
31	guez	3
32	guez	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	guez	3
38	guez	32	! Local:
39	guez	3
40	guez	32	REAL vscr(ip1jm), uscr(ip1jmp1), hscr(ip1jmp1), pscr(ip1jmp1)
41			REAL massescr(ip1jmp1,llm), finvmasse(ip1jmp1,llm)
42			REAL p(ip1jmp1,llmp1)
43			REAL tpn, tps, tppn(iim), tpps(iim)
44			REAL qpn, qps, qppn(iim), qpps(iim)
45			REAL deltap(ip1jmp1,llm)
46	guez	3
47	guez	32	INTEGER l, ij, iq
48	guez	3
49	guez	32	REAL ssum
50	guez	3
51	guez	32	!-----------------------------------------------------------------------
52	guez	3
53	guez	32	DO l = 1, llm
54			DO ij = 1, iip1
55			ucov(ij,l) = 0.
56			ucov(ij+ip1jm,l) = 0.
57			uscr(ij) = 0.
58			uscr(ij+ip1jm) = 0.
59			END DO
60			END DO
61	guez	3
62
63	guez	32	! ............ integration de ps ..............
64	guez	3
65	guez	32	CALL scopy(ip1jmp1*llm,masse,1,massescr,1)
66	guez	3
67	guez	32	DO ij = 1, ip1jmp1
68			pscr(ij) = ps(ij)
69			ps(ij) = psm1(ij) + dt*dp(ij)
70			END DO
71	guez	3
72	guez	32	DO ij = 1, ip1jmp1
73			IF (ps(ij)<0.) THEN
74			PRINT *, ' Au point ij = ', ij, ' , pression sol neg. ', ps(ij)
75			STOP 'integrd'
76			END IF
77			END DO
78	guez	3
79	guez	32	DO ij = 1, iim
80			tppn(ij) = aire(ij)*ps(ij)
81			tpps(ij) = aire(ij+ip1jm)*ps(ij+ip1jm)
82			END DO
83			tpn = ssum(iim,tppn,1)/apoln
84			tps = ssum(iim,tpps,1)/apols
85			DO ij = 1, iip1
86			ps(ij) = tpn
87			ps(ij+ip1jm) = tps
88			END DO
89	guez	3
90	guez	32	! ... Calcul de la nouvelle masse d'air au dernier temps integre t+1 .
91	guez	3
92	guez	37	forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
93	guez	32	CALL massdair(p,masse)
94	guez	3
95	guez	32	CALL scopy(ijp1llm,masse,1,finvmasse,1)
96			CALL filtreg(finvmasse,jjp1,llm,-2,2,.TRUE.,1)
97	guez	3
98
99	guez	32	! ............ integration de ucov, vcov, h ..............
100	guez	3
101	guez	32	DO l = 1, llm
102	guez	3
103	guez	32	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	guez	3
108	guez	32	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	guez	3
113	guez	32	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	guez	3
119	guez	32	! .... Calcul de la valeur moyenne, unique aux poles pour teta .
120	guez	3
121
122	guez	32	DO ij = 1, iim
123			tppn(ij) = aire(ij)*teta(ij,l)
124			tpps(ij) = aire(ij+ip1jm)*teta(ij+ip1jm,l)
125			END DO
126			tpn = ssum(iim,tppn,1)/apoln
127			tps = ssum(iim,tpps,1)/apols
128	guez	3
129	guez	32	DO ij = 1, iip1
130			teta(ij,l) = tpn
131			teta(ij+ip1jm,l) = tps
132			END DO
133	guez	3
134
135	guez	32	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	guez	3
141	guez	32	END DO
142	guez	3
143	guez	32	DO l = 1, llm
144			DO ij = 1, ip1jmp1
145			deltap(ij,l) = p(ij,l) - p(ij,l+1)
146			END DO
147			END DO
148	guez	3
149	guez	32	CALL qminimum(q,nq,deltap)
150	guez	3
151	guez	32	! ..... Calcul de la valeur moyenne, unique aux poles pour q .....
152	guez	3
153
154	guez	32	DO iq = 1, nq
155			DO l = 1, llm
156	guez	28
157	guez	32	DO ij = 1, iim
158			qppn(ij) = aire(ij)*q(ij,l,iq)
159			qpps(ij) = aire(ij+ip1jm)*q(ij+ip1jm,l,iq)
160			END DO
161			qpn = ssum(iim,qppn,1)/apoln
162			qps = ssum(iim,qpps,1)/apols
163	guez	28
164	guez	32	DO ij = 1, iip1
165			q(ij,l,iq) = qpn
166			q(ij+ip1jm,l,iq) = qps
167			END DO
168	guez	28
169	guez	32	END DO
170			END DO
171	guez	28
172
173	guez	32	CALL scopy(ijp1llm,finvmasse,1,finvmaold,1)
174	guez	28
175
176	guez	32	! ..... FIN de l'integration de q .......
177
178			IF (leapf) THEN
179			CALL scopy(ip1jmp1,pscr,1,psm1,1)
180			CALL scopy(ip1jmp1*llm,massescr,1,massem1,1)
181			END IF
182
183			END SUBROUTINE integrd
184
185			end module integrd_m