/[lmdze]/trunk/Sources/dyn3d/Dissipation/inidissip.f
ViewVC logotype

Contents of /trunk/Sources/dyn3d/Dissipation/inidissip.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 55 - (show annotations)
Mon Dec 12 13:25:01 2011 UTC (12 years, 5 months ago) by guez
Original Path: trunk/libf/dyn3d/Dissipation/inidissip.f90
File size: 4110 byte(s)
-- In procedure "bilan_dyn", replaced average of "zvq" by integral of
"zvq", following a comment of Francis Codron :

Le calcul actuel donne des unités peu pratiques : transports de
chaleur en K m / s par exemple. C'est bien pour les sorties à 2
dimensions, latitude et pression, car alors le transport ne dépend pas
de l'espacement des niveaux, mieux pour comparer ou tracer en latitude
et pression. Par contre, quand on somme sur la verticale, on
préfèrerait avoir des transports d'énergie en watts, ou au moins an K
kg / s (à multiplier par "Cp" ou "L"). On doit pouvoir recalculer le
transport intégré à partir des fichiers de sortie, mais c'est embêtant
(calcul de "cv").

-- Gathered files in directory Dissipation.

1 module inidissip_m
2
3 use dimens_m, only: llm
4
5 IMPLICIT NONE
6
7 private llm
8
9 REAL dtdiss
10 integer idissip ! période de la dissipation (en pas de temps)
11 real tetaudiv(llm), tetaurot(llm), tetah(llm)
12 real cdivu, crot, cdivh
13
14 contains
15
16 SUBROUTINE inidissip
17
18 ! From dyn3d/inidissip.F, version 1.1.1.1 2004/05/19 12:53:06
19
20 ! Initialisation de la dissipation horizontale. Calcul des valeurs
21 ! propres des opérateurs par méthode itérative.
22
23 USE comconst, ONLY : dtvr
24 use comdissnew, only: lstardis, nitergdiv, nitergrot, niterh, tetagdiv, &
25 tetagrot, tetatemp
26 USE comvert, ONLY : preff, presnivs
27 USE conf_gcm_m, ONLY : iperiod
28 USE dimens_m, ONLY : iim, jjm, llm
29 USE paramet_m, ONLY : jjp1
30 use jumble, only: new_unit
31 use filtreg_m, only: filtreg
32 use gradiv2_m, only: gradiv2
33
34 ! Variables local to the procedure:
35 REAL zvert(llm), max_zvert
36 REAL, dimension(iim + 1, jjm + 1):: zh, zu
37 real zv(iim + 1, jjm), deltap(iim + 1, jjm + 1, llm)
38 REAL zllm
39 INTEGER l, seed_size, ii, unit
40 REAL tetamin ! in s
41
42 !-----------------------------------------------------------------------
43
44 PRINT *, 'Call sequence information: inidissip'
45 call random_seed(size=seed_size)
46 call random_seed(put=(/(0, ii = 1, seed_size)/))
47
48 PRINT *, 'Calcul des valeurs propres de divgrad'
49 deltap = 1.
50 call random_number(zh)
51 zh = zh - 0.5
52 CALL filtreg(zh, jjp1, 1, 2, 1, .TRUE., 1)
53
54 DO l = 1, 50
55 IF (lstardis) THEN
56 CALL divgrad2(1, zh, deltap, niterh, zh, -1.)
57 ELSE
58 CALL divgrad(1, zh, niterh, zh, -1.)
59 END IF
60
61 zllm = abs(maxval(zh))
62 zh = zh / zllm
63 END DO
64
65 IF (lstardis) THEN
66 cdivh = 1. / zllm
67 ELSE
68 cdivh = zllm**(- 1. / niterh)
69 END IF
70 PRINT *, 'cdivh = ', cdivh
71
72 PRINT *, 'Calcul des valeurs propres de gradiv'
73 call random_number(zu)
74 zu = zu - 0.5
75 CALL filtreg(zu, jjp1, 1, 2, 1, .TRUE., 1)
76 call random_number(zv)
77 zv = zv - 0.5
78 CALL filtreg(zv, jjm, 1, 2, 1, .FALSE., 1)
79
80 DO l = 1, 50
81 IF (lstardis) THEN
82 CALL gradiv2(1, zu, zv, nitergdiv, zu, zv, -1.)
83 ELSE
84 CALL gradiv(1, zu, zv, nitergdiv, zu, zv, -1.)
85 END IF
86
87 zllm = max(abs(maxval(zu)), abs(maxval(zv)))
88 zu = zu / zllm
89 zv = zv / zllm
90 end DO
91
92 IF (lstardis) THEN
93 cdivu = 1. / zllm
94 ELSE
95 cdivu = zllm**(- 1. / nitergdiv)
96 END IF
97 PRINT *, 'cdivu = ', cdivu
98
99 PRINT *, 'Calcul des valeurs propres de nxgrarot'
100 call random_number(zu)
101 zu = zu - 0.5
102 CALL filtreg(zu, jjp1, 1, 2, 1, .TRUE., 1)
103 call random_number(zv)
104 zv = zv - 0.5
105 CALL filtreg(zv, jjm, 1, 2, 1, .FALSE., 1)
106
107 DO l = 1, 50
108 IF (lstardis) THEN
109 CALL nxgraro2(1, zu, zv, nitergrot, zu, zv, -1.)
110 ELSE
111 CALL nxgrarot(1, zu, zv, nitergrot, zu, zv, -1.)
112 END IF
113
114 zllm = max(abs(maxval(zu)), abs(maxval(zv)))
115 zu = zu / zllm
116 zv = zv / zllm
117 end DO
118
119 IF (lstardis) THEN
120 crot = 1. / zllm
121 ELSE
122 crot = zllm**(-1. / nitergrot)
123 END IF
124 PRINT *, 'crot = ', crot
125
126 ! Variation verticale du coefficient de dissipation :
127 zvert = 2. - 1. / (1. + (preff / presnivs - 1.)**2)
128 ! (between 1 and 2)
129
130 tetaudiv = zvert / tetagdiv
131 tetaurot = zvert / tetagrot
132 tetah = zvert / tetatemp
133
134 call new_unit(unit)
135 open(unit, file="inidissip.csv", status="replace", action="write")
136 write(unit, fmt=*) "tetaudiv tetaurot tetah" ! title line
137 do l = 1, llm
138 write(unit, fmt=*) tetaudiv(l), tetaurot(l), tetah(l)
139 end do
140 close(unit)
141 print *, 'Created file "inidissip.csv".'
142
143 max_zvert = maxval(zvert)
144 tetamin = min(1e6, tetagdiv / max_zvert, tetagrot / max_zvert, &
145 tetatemp / max_zvert)
146 PRINT *, 'tetamin = ', tetamin
147 idissip = max(1, int(tetamin / (2 * dtvr * iperiod))) * iperiod
148 PRINT *, 'idissip = ', idissip
149 dtdiss = idissip * dtvr
150 PRINT *, 'dtdiss = ', dtdiss
151
152 END SUBROUTINE inidissip
153
154 end module inidissip_m

  ViewVC Help
Powered by ViewVC 1.1.21