Sources/dyn3d/sortvarc.f

SUBROUTINE sortvarc(itau,ucov,teta,ps,masse,pk,phis,vorpot,phi,bern,dp, &
     time)


  ! From dyn3d/sortvarc.F,v 1.1.1.1 2004/05/19 12:53:07

  !   Auteur:    P. Le Van
  !   Objet:
  !   sortie des variables de controle

  USE dimens_m, ONLY : iim, jjm, llm
  USE paramet_m, ONLY : iip1, iip2, ijp1llm, ip1jm, ip1jmp1, jjp1
  USE comconst, ONLY : daysec, dtvr, g, omeg, rad
  USE comgeom, ONLY : aire, cu, rlatu
  USE dynetat0_m, ONLY : day_ini
  USE ener, ONLY : ang, ang0, etot, etot0, ptot, ptot0, rmsdpdt, rmsv, &
       stot, stot0, ztot, ztot0
  use filtreg_m, only: filtreg

  IMPLICIT NONE

  !   Arguments:

  INTEGER, INTENT (IN) :: itau
  REAL :: ucov(ip1jmp1,llm), teta(ip1jmp1,llm), masse(ip1jmp1,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)

  !   Local:

  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

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

  dtvrs1j = dtvr/daysec
  rjour = float(int(itau*dtvrs1j))
  heure = (itau*dtvrs1j-rjour)*24.
  imjmp1 = iim*jjp1
  IF (abs(heure-24.)<=0.0001) heure = 0.

  CALL massbarxy(masse,massebxy)

  !   .....  Calcul  de  rmsdpdt  .....

  ge(:) = dp(:)*dp(:)

  rmsdpdt = ssum(ip1jmp1,ge,1) - ssum(jjp1,ge,iip1)

  rmsdpdt = daysec*1.E-2*sqrt(rmsdpdt/imjmp1)

  CALL scopy(ijp1llm,bern,1,bernf,1)
  CALL filtreg(bernf,jjp1,llm,-2,2,.TRUE.,1)

  !   .....  Calcul du moment  angulaire   .....

  radsg = rad/g
  radomeg = rad*omeg

  DO ij = iip2, ip1jm
     cosphi(ij) = cos(rlatu((ij-1)/iip1+1))
     omegcosp(ij) = radomeg*cosphi(ij)
  END DO

  !  ...  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)
     END DO
     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)+bernf(ij,l)-phi( &
             ij,l))
     END DO
     etotl(l) = ssum(ip1jmp1,ge,1) - ssum(jjp1,ge,iip1)

     DO ij = 1, ip1jmp1
        ge(ij) = masse(ij,l)*teta(ij,l)
     END DO
     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.)
     END DO
     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)
     END DO
     angl(l) = radsg*(ssum(ip1jm-iip1,ge(iip2),1)-ssum(jjm-1,ge(iip2),iip1) &
          )
  END DO

  DO ij = 1, ip1jmp1
     ge(ij) = ps(ij)*aire(ij)
  END DO
  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))

  IF (ptot0==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 SUBROUTINE sortvarc
1	SUBROUTINE sortvarc(itau,ucov,teta,ps,masse,pk,phis,vorpot,phi,bern,dp, &
2	time)
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 = (itaudtvrs1j-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 = daysec1.E-2sqrt(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