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

Diff of /trunk/dyn3d/sortvarc.f

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

trunk/libf/dyn3d/sortvarc.f revision 3 by guez, Wed Feb 27 13:16:39 2008 UTC trunk/libf/dyn3d/sortvarc.f90 revision 33 by guez, Fri Apr 9 10:56:14 2010 UTC
# Line 1  Line 1 
1  !  module sortvarc_m
 ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/sortvarc.F,v 1.1.1.1 2004/05/19 12:53:07 lmdzadmin Exp $  
 !  
       SUBROUTINE sortvarc  
      $(itau,ucov,teta,ps,masse,pk,phis,vorpot,phi,bern,dp,time ,  
      $ vcov )  
   
 c=======================================================================  
 c  
 c   Auteur:    P. Le Van  
 c   -------  
 c  
 c   Objet:  
 c   ------  
 c  
 c   sortie des variables de controle  
 c  
 c=======================================================================  
 c-----------------------------------------------------------------------  
 c   Declarations:  
 c   -------------  
   
       use dimens_m  
       use paramet_m  
       use comconst  
       use comvert  
       use logic  
       use comgeom  
       use temps  
       use ener  
       IMPLICIT NONE  
   
 c   Arguments:  
 c   ----------  
   
       INTEGER, intent(in):: itau  
       REAL ucov(ip1jmp1,llm),teta(ip1jmp1,llm),masse(ip1jmp1,llm)  
       REAL vcov(ip1jm,llm)  
       REAL ps(ip1jmp1),phis(ip1jmp1)  
       REAL vorpot(ip1jm,llm)  
       REAL phi(ip1jmp1,llm),bern(ip1jmp1,llm)  
       REAL dp(ip1jmp1)  
       REAL time  
       REAL, intent(in):: pk(ip1jmp1,llm)  
   
 c   Local:  
 c   ------  
   
       REAL vor(ip1jm),bernf(ip1jmp1,llm),ztotl(llm)  
       REAL etotl(llm),stotl(llm),rmsvl(llm),angl(llm),ge(ip1jmp1)  
       REAL cosphi(ip1jm),omegcosp(ip1jm)  
       REAL dtvrs1j,rjour,heure,radsg,radomeg  
       REAL rday, massebxy(ip1jm,llm)  
       INTEGER  l, ij, imjmp1  
   
       REAL       SSUM  
   
 c-----------------------------------------------------------------------  
   
        dtvrs1j   = dtvr/daysec  
        rjour     = FLOAT( INT( itau * dtvrs1j ))  
        heure     = ( itau*dtvrs1j-rjour ) * 24.  
        imjmp1    = iim * jjp1  
        IF(ABS(heure - 24.).LE.0.0001 ) heure = 0.  
 c  
        CALL massbarxy ( masse, massebxy )  
   
 c   .....  Calcul  de  rmsdpdt  .....  
   
        ge(:)=dp(:)*dp(:)  
   
        rmsdpdt = SSUM(ip1jmp1,ge,1) - SSUM(jjp1,ge,iip1)  
 c  
        rmsdpdt = daysec* 1.e-2 * SQRT(rmsdpdt/imjmp1)  
   
        CALL SCOPY( ijp1llm,bern,1,bernf,1 )  
        CALL filtreg(bernf,jjp1,llm,-2,2,.TRUE.,1)  
   
 c   .....  Calcul du moment  angulaire   .....  
   
        radsg    = rad /g  
        radomeg  = rad * omeg  
 c  
        DO ij=iip2,ip1jm  
           cosphi( ij ) = COS(rlatu((ij-1)/iip1+1))  
           omegcosp(ij) = radomeg   * cosphi(ij)  
        ENDDO  
   
 c  ...  Calcul  de l'energie,de l'enstrophie,de l'entropie et de rmsv  .  
   
        DO l=1,llm  
           DO ij = 1,ip1jm  
              vor(ij)=vorpot(ij,l)*vorpot(ij,l)*massebxy(ij,l)  
           ENDDO  
           ztotl(l)=(SSUM(ip1jm,vor,1)-SSUM(jjm,vor,iip1))  
   
           DO ij = 1,ip1jmp1  
              ge(ij)= masse(ij,l)*(phis(ij)+teta(ij,l)*pk(ij,l)  +  
      s        bernf(ij,l)-phi(ij,l))  
           ENDDO  
           etotl(l) = SSUM(ip1jmp1,ge,1) - SSUM(jjp1,ge,iip1)  
   
           DO   ij   = 1, ip1jmp1  
              ge(ij) = masse(ij,l)*teta(ij,l)  
           ENDDO  
           stotl(l)= SSUM(ip1jmp1,ge,1) - SSUM(jjp1,ge,iip1)  
   
           DO ij=1,ip1jmp1  
              ge(ij)=masse(ij,l)*AMAX1(bernf(ij,l)-phi(ij,l),0.)  
           ENDDO  
           rmsvl(l)=2.*(SSUM(ip1jmp1,ge,1)-SSUM(jjp1,ge,iip1))  
   
           DO ij =iip2,ip1jm  
              ge(ij)=(ucov(ij,l)/cu(ij)+omegcosp(ij))*masse(ij,l) *  
      *               cosphi(ij)  
           ENDDO  
           angl(l) = radsg *  
      s    (SSUM(ip1jm-iip1,ge(iip2),1)-SSUM(jjm-1,ge(iip2),iip1))  
       ENDDO  
   
           DO ij=1,ip1jmp1  
             ge(ij)= ps(ij)*aire(ij)  
           ENDDO  
       ptot  = SSUM(ip1jmp1,ge,1)-SSUM(jjp1,ge,iip1)  
       etot  = SSUM(     llm, etotl, 1 )  
       ztot  = SSUM(     llm, ztotl, 1 )  
       stot  = SSUM(     llm, stotl, 1 )  
       rmsv  = SSUM(     llm, rmsvl, 1 )  
       ang   = SSUM(     llm,  angl, 1 )  
   
       rday = FLOAT(INT ( day_ini + time ))  
 c  
       IF(ptot0.eq.0.)  THEN  
          PRINT 3500, itau, rday, heure,time  
          PRINT*,'WARNING!!! On recalcule les valeurs initiales de :'  
          PRINT*,'ptot,rmsdpdt,etot,ztot,stot,rmsv,ang'  
          PRINT *, ptot,rmsdpdt,etot,ztot,stot,rmsv,ang  
          etot0 = etot  
          ptot0 = ptot  
          ztot0 = ztot  
          stot0 = stot  
          ang0  = ang  
       END IF  
   
       etot= etot/etot0  
       rmsv= SQRT(rmsv/ptot)  
       ptot= ptot/ptot0  
       ztot= ztot/ztot0  
       stot= stot/stot0  
       ang = ang /ang0  
   
   
       PRINT 3500, itau, rday, heure, time  
       PRINT 4000, ptot,rmsdpdt,etot,ztot,stot,rmsv,ang  
   
       RETURN  
   
 3500   FORMAT(4x,'pas',i7,5x,'jour',f5.0,'heure',f5.1,4x  
      *   ,'date',f10.5)  
 4000   FORMAT(10x,'masse',4x,'rmsdpdt',7x,'energie',2x,'enstrophie'  
      * ,2x,'entropie',3x,'rmsv',4x,'mt.ang',/,'GLOB  '  
      .  ,f10.6,e13.6,5f10.3/  
      * )  
       END  
2    
3      IMPLICIT NONE
4    
5    contains
6    
7      SUBROUTINE sortvarc(itau, ucov, teta, ps, masse, pk, phis, vorpot, phi, &
8           bern, dp, time_0)
9    
10        ! From dyn3d/sortvarc.F, version 1.1.1.1 2004/05/19 12:53:07
11        ! Author: P. Le Van
12        ! Objet : sortie des variables de contrôle
13    
14        USE conf_gcm_m, ONLY: day_step
15        USE dimens_m, ONLY : iim, jjm, llm
16        USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, jjp1
17        USE comconst, ONLY : daysec, dtvr, g, omeg, rad
18        USE comgeom, ONLY : aire, cu, rlatu
19        USE dynetat0_m, ONLY : day_ini
20        USE ener, ONLY : ang, ang0, etot, etot0, ptot, ptot0, rmsdpdt, rmsv, &
21             stot, stot0, ztot, ztot0
22        use filtreg_m, only: filtreg
23    
24        ! Arguments:
25        INTEGER, INTENT (IN) :: itau
26        REAL :: ucov(ip1jmp1, llm), teta(ip1jmp1, llm), masse(ip1jmp1, llm)
27        REAL :: ps(ip1jmp1), phis(ip1jmp1)
28        REAL :: vorpot(ip1jm, llm)
29        REAL :: phi(ip1jmp1, llm), bern(ip1jmp1, llm)
30        REAL :: dp(ip1jmp1)
31        REAL, INTENT (IN):: time_0
32        REAL, INTENT (IN):: pk(ip1jmp1, llm)
33    
34        ! Local:
35        REAL :: vor(ip1jm), bernf(ip1jmp1, llm), ztotl(llm)
36        REAL :: etotl(llm), stotl(llm), rmsvl(llm), angl(llm), ge(ip1jmp1)
37        REAL :: cosphi(ip1jm), omegcosp(ip1jm)
38        REAL :: dtvrs1j, rjour, heure, radsg, radomeg
39        real massebxy(ip1jm, llm)
40        INTEGER :: l, ij, imjmp1
41        REAL :: ssum
42        real time
43    
44        !-----------------------------------------------------------------------
45    
46        time = real(itau) / day_step + time_0
47        dtvrs1j = dtvr/daysec
48        rjour = real(int(itau*dtvrs1j))
49        heure = (itau*dtvrs1j-rjour)*24.
50        imjmp1 = iim*jjp1
51        IF (abs(heure-24.)<=0.0001) heure = 0.
52    
53        CALL massbarxy(masse, massebxy)
54    
55        ! Calcul  de  rmsdpdt
56        ge(:) = dp(:)*dp(:)
57        rmsdpdt = ssum(ip1jmp1, ge, 1) - ssum(jjp1, ge, iip1)
58        rmsdpdt = daysec*1.E-2*sqrt(rmsdpdt/imjmp1)
59        CALL scopy(ijp1llm, bern, 1, bernf, 1)
60        CALL filtreg(bernf, jjp1, llm, -2, 2, .TRUE., 1)
61    
62        ! Calcul du moment  angulaire
63        radsg = rad/g
64        radomeg = rad*omeg
65        DO ij = iip2, ip1jm
66           cosphi(ij) = cos(rlatu((ij-1)/iip1+1))
67           omegcosp(ij) = radomeg*cosphi(ij)
68        END DO
69    
70        ! Calcul  de l'energie, de l'enstrophie, de l'entropie et de rmsv
71        DO l = 1, llm
72           DO ij = 1, ip1jm
73              vor(ij) = vorpot(ij, l)*vorpot(ij, l)*massebxy(ij, l)
74           END DO
75           ztotl(l) = (ssum(ip1jm, vor, 1)-ssum(jjm, vor, iip1))
76    
77           DO ij = 1, ip1jmp1
78              ge(ij) = masse(ij, l) * (phis(ij) + teta(ij, l) * pk(ij, l) &
79                   + bernf(ij, l)-phi(ij, l))
80           END DO
81           etotl(l) = ssum(ip1jmp1, ge, 1) - ssum(jjp1, ge, iip1)
82    
83           DO ij = 1, ip1jmp1
84              ge(ij) = masse(ij, l)*teta(ij, l)
85           END DO
86           stotl(l) = ssum(ip1jmp1, ge, 1) - ssum(jjp1, ge, iip1)
87    
88           DO ij = 1, ip1jmp1
89              ge(ij) = masse(ij, l)*amax1(bernf(ij, l)-phi(ij, l), 0.)
90           END DO
91           rmsvl(l) = 2.*(ssum(ip1jmp1, ge, 1)-ssum(jjp1, ge, iip1))
92    
93           DO ij = iip2, ip1jm
94              ge(ij) = (ucov(ij, l)/cu(ij)+omegcosp(ij))*masse(ij, l)*cosphi(ij)
95           END DO
96           angl(l) = radsg * (ssum(ip1jm-iip1, ge(iip2), 1) &
97                - ssum(jjm-1, ge(iip2), iip1))
98        END DO
99    
100        DO ij = 1, ip1jmp1
101           ge(ij) = ps(ij)*aire(ij)
102        END DO
103        ptot = ssum(ip1jmp1, ge, 1) - ssum(jjp1, ge, iip1)
104        etot = ssum(llm, etotl, 1)
105        ztot = ssum(llm, ztotl, 1)
106        stot = ssum(llm, stotl, 1)
107        rmsv = ssum(llm, rmsvl, 1)
108        ang = ssum(llm, angl, 1)
109    
110        IF (ptot0 == 0.) THEN
111           PRINT *, 'WARNING!!! On recalcule les valeurs initiales de :'
112           PRINT *, 'ptot, rmsdpdt, etot, ztot, stot, rmsv, ang'
113           PRINT *, ptot, rmsdpdt, etot, ztot, stot, rmsv, ang
114           etot0 = etot
115           ptot0 = ptot
116           ztot0 = ztot
117           stot0 = stot
118           ang0 = ang
119        END IF
120    
121        etot = etot/etot0
122        rmsv = sqrt(rmsv/ptot)
123        ptot = ptot/ptot0
124        ztot = ztot/ztot0
125        stot = stot/stot0
126        ang = ang/ang0
127    
128        PRINT 3500, itau, int(day_ini + time), heure, time
129        PRINT 4000, ptot, rmsdpdt, etot, ztot, stot, rmsv, ang
130    
131    3500 FORMAT (4X, 'pas', I7, 5X, 'jour', i5, 1X, 'heure', F5.1, 4X, 'date', &
132              F10.5)
133    4000 FORMAT (10X, 'masse', 4X, 'rmsdpdt', 7X, 'energie', 2X, 'enstrophie', &
134              2X, 'entropie', 3X, 'rmsv', 4X, 'mt.ang', /, 'GLOB  ', F10.6, &
135              E13.6, 5F10.3/)
136    
137      END SUBROUTINE sortvarc
138    
139    end module sortvarc_m

Legend:
Removed from v.3  
changed lines
  Added in v.33

  ViewVC Help
Powered by ViewVC 1.1.21