/[lmdze]/trunk/libf/dyn3d/sortvarc.f90
ViewVC logotype

Diff of /trunk/libf/dyn3d/sortvarc.f90

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

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

  ViewVC Help
Powered by ViewVC 1.1.21