1 |
|
|
! |
2 |
|
|
! $Id: phyaqua_mod.F90 4523 2023-05-03 16:21:08Z evignon $ |
3 |
|
|
! |
4 |
|
|
MODULE phyaqua_mod |
5 |
|
|
! Routines complementaires pour la physique planetaire. |
6 |
|
|
IMPLICIT NONE |
7 |
|
|
|
8 |
|
|
CONTAINS |
9 |
|
|
|
10 |
|
|
SUBROUTINE iniaqua(nlon,year_len,iflag_phys) |
11 |
|
|
|
12 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
13 |
|
|
! Creation d'un etat initial et de conditions aux limites |
14 |
|
|
! (resp startphy.nc et limit.nc) pour des configurations idealisees |
15 |
|
|
! du modele LMDZ dans sa version terrestre. |
16 |
|
|
! iflag_phys est un parametre qui controle |
17 |
|
|
! iflag_phys = N |
18 |
|
|
! de 100 a 199 : aqua planetes avec SST forcees |
19 |
|
|
! N-100 determine le type de SSTs |
20 |
|
|
! de 200 a 299 : terra planetes avec Ts calcule |
21 |
|
|
|
22 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
23 |
|
|
|
24 |
|
|
USE dimphy, ONLY: klon |
25 |
|
|
USE geometry_mod, ONLY : latitude |
26 |
|
|
USE surface_data, ONLY: type_ocean, ok_veget |
27 |
|
|
USE pbl_surface_mod, ONLY: pbl_surface_init |
28 |
|
|
USE fonte_neige_mod, ONLY: fonte_neige_init |
29 |
|
|
USE phys_state_var_mod |
30 |
|
|
USE time_phylmdz_mod, ONLY: day_ref, ndays, pdtphys, & |
31 |
|
|
day_ini,day_end |
32 |
|
|
USE indice_sol_mod |
33 |
|
|
USE nrtype, ONLY: pi |
34 |
|
|
! USE ioipsl |
35 |
|
|
USE mod_phys_lmdz_para, ONLY: is_master |
36 |
|
|
USE mod_phys_lmdz_transfert_para, ONLY: bcast |
37 |
|
|
USE mod_grid_phy_lmdz |
38 |
|
|
USE ioipsl_getin_p_mod, ONLY : getin_p |
39 |
|
|
USE phys_cal_mod , ONLY: calend, year_len_phy => year_len |
40 |
|
|
IMPLICIT NONE |
41 |
|
|
|
42 |
|
|
include "YOMCST.h" |
43 |
|
|
include "clesphys.h" |
44 |
|
|
include "dimsoil.h" |
45 |
|
|
|
46 |
|
|
INTEGER, INTENT (IN) :: nlon, year_len, iflag_phys |
47 |
|
|
! IM ajout latfi, lonfi |
48 |
|
|
! REAL, INTENT (IN) :: lonfi(nlon), latfi(nlon) |
49 |
|
|
|
50 |
|
|
INTEGER type_profil, type_aqua |
51 |
|
|
|
52 |
|
|
! Ajouts initialisation des surfaces |
53 |
|
|
REAL :: run_off_lic_0(nlon) |
54 |
|
|
REAL :: qsolsrf(nlon, nbsrf), snsrf(nlon, nbsrf) |
55 |
|
|
REAL :: tsoil(nlon, nsoilmx, nbsrf) |
56 |
|
|
REAL :: tslab(nlon), seaice(nlon) |
57 |
|
|
REAL fder(nlon) |
58 |
|
|
|
59 |
|
|
|
60 |
|
|
|
61 |
|
|
! Arguments : |
62 |
|
|
! ----------- |
63 |
|
|
|
64 |
|
|
! integer radpas |
65 |
|
|
INTEGER it, unit, i, k, itap |
66 |
|
|
|
67 |
|
|
REAL rugos, albedo |
68 |
|
|
REAL tsurf |
69 |
|
|
REAL time, timestep, day, day0 |
70 |
|
|
REAL qsol_f |
71 |
|
|
REAL rugsrel(nlon) |
72 |
|
|
LOGICAL alb_ocean |
73 |
|
|
|
74 |
|
|
CHARACTER *80 ans, file_forctl, file_fordat, file_start |
75 |
|
|
CHARACTER *100 file, var |
76 |
|
|
CHARACTER *2 cnbl |
77 |
|
|
|
78 |
|
|
REAL phy_nat(nlon, year_len) |
79 |
|
|
REAL phy_alb(nlon, year_len) |
80 |
|
|
REAL phy_sst(nlon, year_len) |
81 |
|
|
REAL phy_bil(nlon, year_len) |
82 |
|
|
REAL phy_rug(nlon, year_len) |
83 |
|
|
REAL phy_ice(nlon, year_len) |
84 |
|
|
REAL phy_fter(nlon, year_len) |
85 |
|
|
REAL phy_foce(nlon, year_len) |
86 |
|
|
REAL phy_fsic(nlon, year_len) |
87 |
|
|
REAL phy_flic(nlon, year_len) |
88 |
|
|
|
89 |
|
|
INTEGER, SAVE :: read_climoz = 0 ! read ozone climatology |
90 |
|
|
!$OMP THREADPRIVATE(read_climoz) |
91 |
|
|
|
92 |
|
|
! ------------------------------------------------------------------------- |
93 |
|
|
! declaration pour l'appel a phyredem |
94 |
|
|
! ------------------------------------------------------------------------- |
95 |
|
|
|
96 |
|
|
! real pctsrf(nlon,nbsrf),ftsol(nlon,nbsrf) |
97 |
|
|
REAL falbe(nlon, nbsrf), falblw(nlon, nbsrf) |
98 |
|
|
! real pbl_tke(nlon,llm,nbsrf) |
99 |
|
|
! real rain_fall(nlon),snow_fall(nlon) |
100 |
|
|
! real solsw(nlon), sollw(nlon),radsol(nlon) |
101 |
|
|
! real t_ancien(nlon,llm),q_ancien(nlon,llm),rnebcon(nlon,llm) |
102 |
|
|
! real ratqs(nlon,llm) |
103 |
|
|
! real clwcon(nlon,llm) |
104 |
|
|
|
105 |
|
|
INTEGER longcles |
106 |
|
|
PARAMETER (longcles=20) |
107 |
|
|
REAL clesphy0(longcles) |
108 |
|
|
|
109 |
|
|
|
110 |
|
|
! ----------------------------------------------------------------------- |
111 |
|
|
! dynamial tendencies : |
112 |
|
|
! --------------------- |
113 |
|
|
|
114 |
|
|
INTEGER l, ierr, aslun |
115 |
|
|
|
116 |
|
|
REAL paire |
117 |
|
|
|
118 |
|
|
! Local |
119 |
|
|
CHARACTER (LEN=20) :: modname='phyaqua' |
120 |
|
|
CHARACTER (LEN=80) :: abort_message |
121 |
|
|
|
122 |
|
|
|
123 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
124 |
|
|
! INITIALISATIONS |
125 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
126 |
|
|
|
127 |
|
|
! ----------------------------------------------------------------------- |
128 |
|
|
! Initialisations des constantes |
129 |
|
|
! ------------------------------- |
130 |
|
|
|
131 |
|
|
!IF (calend .EQ. "earth_360d") Then |
132 |
|
|
year_len_phy = year_len |
133 |
|
|
!END IF |
134 |
|
|
|
135 |
|
|
if (year_len.ne.360) then |
136 |
|
|
write (*,*) year_len |
137 |
|
|
call abort_physic("iniaqua", 'iniaqua: 360 day calendar is required !', 1) |
138 |
|
|
endif |
139 |
|
|
|
140 |
|
|
type_aqua = iflag_phys/100 |
141 |
|
|
type_profil = iflag_phys - type_aqua*100 |
142 |
|
|
PRINT *, 'iniaqua:type_aqua, type_profil', type_aqua, type_profil |
143 |
|
|
|
144 |
|
|
IF (klon/=nlon) THEN |
145 |
|
|
WRITE (*, *) 'iniaqua: klon=', klon, ' nlon=', nlon |
146 |
|
|
abort_message= 'probleme de dimensions dans iniaqua' |
147 |
|
|
CALL abort_physic(modname,abort_message,1) |
148 |
|
|
END IF |
149 |
|
|
CALL phys_state_var_init(read_climoz) |
150 |
|
|
|
151 |
|
|
|
152 |
|
|
read_climoz = 0 |
153 |
|
|
day0 = 217. |
154 |
|
|
day = day0 |
155 |
|
|
it = 0 |
156 |
|
|
time = 0. |
157 |
|
|
|
158 |
|
|
! ----------------------------------------------------------------------- |
159 |
|
|
! initialisations de la physique |
160 |
|
|
! ----------------------------------------------------------------------- |
161 |
|
|
|
162 |
|
|
day_ini = day_ref |
163 |
|
|
day_end = day_ini + ndays |
164 |
|
|
|
165 |
|
|
nbapp_rad = 24 |
166 |
|
|
CALL getin_p('nbapp_rad', nbapp_rad) |
167 |
|
|
|
168 |
|
|
! --------------------------------------------------------------------- |
169 |
|
|
! Creation des conditions aux limites: |
170 |
|
|
! ------------------------------------ |
171 |
|
|
! Initialisations des constantes |
172 |
|
|
! Ajouter les manquants dans planete.def... (albedo etc) |
173 |
|
|
co2_ppm = 348. |
174 |
|
|
CALL getin_p('co2_ppm', co2_ppm) |
175 |
|
|
|
176 |
|
|
solaire = 1365. |
177 |
|
|
CALL getin_p('solaire', solaire) |
178 |
|
|
|
179 |
|
|
! CALL getin('albedo',albedo) ! albedo is set below, depending on |
180 |
|
|
! type_aqua |
181 |
|
|
alb_ocean = .TRUE. |
182 |
|
|
CALL getin_p('alb_ocean', alb_ocean) |
183 |
|
|
|
184 |
|
|
WRITE (*, *) 'iniaqua: co2_ppm=', co2_ppm |
185 |
|
|
WRITE (*, *) 'iniaqua: solaire=', solaire |
186 |
|
|
WRITE (*, *) 'iniaqua: alb_ocean=', alb_ocean |
187 |
|
|
|
188 |
|
|
radsol = 0. |
189 |
|
|
qsol_f = 10. |
190 |
|
|
|
191 |
|
|
! Conditions aux limites: |
192 |
|
|
! ----------------------- |
193 |
|
|
|
194 |
|
|
qsol(:) = qsol_f |
195 |
|
|
rugsrel = 0.0 ! (rugsrel = rugoro) |
196 |
|
|
rugoro = 0.0 |
197 |
|
|
u_ancien = 0.0 |
198 |
|
|
v_ancien = 0.0 |
199 |
|
|
agesno = 50.0 |
200 |
|
|
! Relief plat |
201 |
|
|
zmea = 0. |
202 |
|
|
zstd = 0. |
203 |
|
|
zsig = 0. |
204 |
|
|
zgam = 0. |
205 |
|
|
zthe = 0. |
206 |
|
|
zpic = 0. |
207 |
|
|
zval = 0. |
208 |
|
|
|
209 |
|
|
! Une seule surface |
210 |
|
|
pctsrf = 0. |
211 |
|
|
IF (type_aqua==1) THEN |
212 |
|
|
rugos = 1.E-4 |
213 |
|
|
albedo = 0.19 |
214 |
|
|
pctsrf(:, is_oce) = 1. |
215 |
|
|
ELSE IF (type_aqua==2) THEN |
216 |
|
|
rugos = 0.03 |
217 |
|
|
albedo = 0.1 |
218 |
|
|
pctsrf(:, is_ter) = 1. |
219 |
|
|
END IF |
220 |
|
|
|
221 |
|
|
CALL getin_p('rugos', rugos) |
222 |
|
|
|
223 |
|
|
WRITE (*, *) 'iniaqua: rugos=', rugos |
224 |
|
|
zmasq(:) = pctsrf(:, is_ter) |
225 |
|
|
|
226 |
|
|
! pctsrf_pot(:,is_oce) = 1. - zmasq(:) |
227 |
|
|
! pctsrf_pot(:,is_sic) = 1. - zmasq(:) |
228 |
|
|
|
229 |
|
|
! Si alb_ocean on calcule un albedo oceanique moyen |
230 |
|
|
! if (alb_ocean) then |
231 |
|
|
! Voir pourquoi on avait ca. |
232 |
|
|
! CALL ini_alb_oce(phy_alb) |
233 |
|
|
! else |
234 |
|
|
phy_alb(:, :) = albedo ! albedo land only (old value condsurf_jyg=0.3) |
235 |
|
|
! endif !alb_ocean |
236 |
|
|
|
237 |
|
|
DO i = 1, year_len |
238 |
|
|
! IM Terraplanete phy_sst(:,i) = 260.+50.*cos(rlatd(:))**2 |
239 |
|
|
! IM ajout calcul profil sst selon le cas considere (cf. FBr) |
240 |
|
|
|
241 |
|
|
phy_nat(:, i) = 1.0 ! 0=ocean libre, 1=land, 2=glacier, 3=banquise |
242 |
|
|
phy_bil(:, i) = 1.0 ! ne sert que pour les slab_ocean |
243 |
|
|
phy_rug(:, i) = rugos ! longueur rugosite utilisee sur land only |
244 |
|
|
phy_ice(:, i) = 0.0 ! fraction de glace (?) |
245 |
|
|
phy_fter(:, i) = pctsrf(:, is_ter) ! fraction de glace (?) |
246 |
|
|
phy_foce(:, i) = pctsrf(:, is_oce) ! fraction de glace (?) |
247 |
|
|
phy_fsic(:, i) = pctsrf(:, is_sic) ! fraction de glace (?) |
248 |
|
|
phy_flic(:, i) = pctsrf(:, is_lic) ! fraction de glace (?) |
249 |
|
|
END DO |
250 |
|
|
! IM calcul profil sst |
251 |
|
|
CALL profil_sst(nlon, latitude, type_profil, phy_sst) |
252 |
|
|
|
253 |
|
|
IF (grid_type==unstructured) THEN |
254 |
|
|
CALL writelim_unstruct(klon, phy_nat, phy_alb, phy_sst, phy_bil, phy_rug, phy_ice, & |
255 |
|
|
phy_fter, phy_foce, phy_flic, phy_fsic) |
256 |
|
|
ELSE |
257 |
|
|
|
258 |
|
|
CALL writelim(klon, phy_nat, phy_alb, phy_sst, phy_bil, phy_rug, phy_ice, & |
259 |
|
|
phy_fter, phy_foce, phy_flic, phy_fsic) |
260 |
|
|
ENDIF |
261 |
|
|
|
262 |
|
|
! --------------------------------------------------------------------- |
263 |
|
|
! Ecriture de l'etat initial: |
264 |
|
|
! --------------------------- |
265 |
|
|
|
266 |
|
|
|
267 |
|
|
! Ecriture etat initial physique |
268 |
|
|
|
269 |
|
|
timestep = pdtphys |
270 |
|
|
radpas = nint(rday/timestep/float(nbapp_rad)) |
271 |
|
|
|
272 |
|
|
DO i = 1, longcles |
273 |
|
|
clesphy0(i) = 0. |
274 |
|
|
END DO |
275 |
|
|
clesphy0(1) = float(iflag_con) |
276 |
|
|
clesphy0(2) = float(nbapp_rad) |
277 |
|
|
! IF( cycle_diurne ) clesphy0(3) = 1. |
278 |
|
|
clesphy0(3) = 1. ! cycle_diurne |
279 |
|
|
clesphy0(4) = 1. ! soil_model |
280 |
|
|
clesphy0(5) = 1. ! new_oliq |
281 |
|
|
clesphy0(6) = 0. ! ok_orodr |
282 |
|
|
clesphy0(7) = 0. ! ok_orolf |
283 |
|
|
clesphy0(8) = 0. ! ok_limitvrai |
284 |
|
|
|
285 |
|
|
|
286 |
|
|
! ======================================================================= |
287 |
|
|
! Profils initiaux |
288 |
|
|
! ======================================================================= |
289 |
|
|
|
290 |
|
|
! On initialise les temperatures de surfaces comme les sst |
291 |
|
|
DO i = 1, nlon |
292 |
|
|
ftsol(i, :) = phy_sst(i, 1) |
293 |
|
|
tsoil(i, :, :) = phy_sst(i, 1) |
294 |
|
|
tslab(i) = phy_sst(i, 1) |
295 |
|
|
END DO |
296 |
|
|
|
297 |
|
|
falbe(:, :) = albedo |
298 |
|
|
falblw(:, :) = albedo |
299 |
|
|
rain_fall(:) = 0. |
300 |
|
|
snow_fall(:) = 0. |
301 |
|
|
solsw(:) = 0. |
302 |
|
|
sollw(:) = 0. |
303 |
|
|
radsol(:) = 0. |
304 |
|
|
|
305 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
306 |
|
|
! intialisation bidon mais pas grave |
307 |
|
|
t_ancien(:, :) = 0. |
308 |
|
|
q_ancien(:, :) = 0. |
309 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
310 |
|
|
rnebcon = 0. |
311 |
|
|
ratqs = 0. |
312 |
|
|
clwcon = 0. |
313 |
|
|
pbl_tke = 1.E-8 |
314 |
|
|
|
315 |
|
|
! variables supplementaires pour appel a plb_surface_init |
316 |
|
|
fder(:) = 0. |
317 |
|
|
seaice(:) = 0. |
318 |
|
|
run_off_lic_0 = 0. |
319 |
|
|
fevap = 0. |
320 |
|
|
|
321 |
|
|
|
322 |
|
|
! Initialisations necessaires avant phyredem |
323 |
|
|
type_ocean = 'force' |
324 |
|
|
CALL fonte_neige_init(run_off_lic_0) |
325 |
|
|
qsolsrf(:, :) = qsol(1) ! humidite du sol des sous surface |
326 |
|
|
snsrf(:, :) = 0. ! couverture de neige des sous surface |
327 |
|
|
z0m(:, :) = rugos ! couverture de neige des sous surface |
328 |
|
|
z0h=z0m |
329 |
|
|
|
330 |
|
|
|
331 |
|
|
CALL pbl_surface_init(fder, snsrf, qsolsrf, tsoil) |
332 |
|
|
|
333 |
|
|
PRINT *, 'iniaqua: before phyredem' |
334 |
|
|
|
335 |
|
|
pbl_tke(:,:,:) = 1.e-8 |
336 |
|
|
falb1 = albedo |
337 |
|
|
falb2 = albedo |
338 |
|
|
zmax0 = 0. |
339 |
|
|
f0 = 0. |
340 |
|
|
sig1 = 0. |
341 |
|
|
w01 = 0. |
342 |
|
|
wake_deltat = 0. |
343 |
|
|
wake_deltaq = 0. |
344 |
|
|
wake_s = 0. |
345 |
|
|
wake_dens = 0. |
346 |
|
|
wake_cstar = 0. |
347 |
|
|
wake_pe = 0. |
348 |
|
|
wake_fip = 0. |
349 |
|
|
fm_therm = 0. |
350 |
|
|
entr_therm = 0. |
351 |
|
|
detr_therm = 0. |
352 |
|
|
ale_bl = 0. |
353 |
|
|
ale_bl_trig =0. |
354 |
|
|
alp_bl =0. |
355 |
|
|
treedrg(:,:,:)=0. |
356 |
|
|
|
357 |
|
|
u10m = 0. |
358 |
|
|
v10m = 0. |
359 |
|
|
|
360 |
|
|
ql_ancien = 0. |
361 |
|
|
qs_ancien = 0. |
362 |
|
|
qbs_ancien = 0. |
363 |
|
|
u_ancien = 0. |
364 |
|
|
v_ancien = 0. |
365 |
|
|
prw_ancien = 0. |
366 |
|
|
prlw_ancien = 0. |
367 |
|
|
prsw_ancien = 0. |
368 |
|
|
prbsw_ancien= 0. |
369 |
|
|
|
370 |
|
|
ale_wake = 0. |
371 |
|
|
ale_bl_stat = 0. |
372 |
|
|
|
373 |
|
|
|
374 |
|
|
!ym error : the sub surface dimension is the third not second : forgotten for iniaqua |
375 |
|
|
! falb_dir(:,is_ter,:)=0.08; falb_dir(:,is_lic,:)=0.6 |
376 |
|
|
! falb_dir(:,is_oce,:)=0.5; falb_dir(:,is_sic,:)=0.6 |
377 |
|
|
falb_dir(:,:,is_ter)=0.08; falb_dir(:,:,is_lic)=0.6 |
378 |
|
|
falb_dir(:,:,is_oce)=0.5; falb_dir(:,:,is_sic)=0.6 |
379 |
|
|
|
380 |
|
|
!ym falb_dif has been forgotten, initialize with defaukt value found in phyetat0 or 0 ? |
381 |
|
|
!ym probably the uninitialized value was 0 for standard (regular grid) case |
382 |
|
|
falb_dif(:,:,:)=0 |
383 |
|
|
|
384 |
|
|
|
385 |
|
|
CALL phyredem('startphy.nc') |
386 |
|
|
|
387 |
|
|
PRINT *, 'iniaqua: after phyredem' |
388 |
|
|
CALL phys_state_var_end |
389 |
|
|
|
390 |
|
|
RETURN |
391 |
|
|
END SUBROUTINE iniaqua |
392 |
|
|
|
393 |
|
|
|
394 |
|
|
! ==================================================================== |
395 |
|
|
! ==================================================================== |
396 |
|
|
SUBROUTINE zenang_an(cycle_diurne, gmtime, rlat, rlon, rmu0, fract) |
397 |
|
|
USE dimphy |
398 |
|
|
IMPLICIT NONE |
399 |
|
|
! ==================================================================== |
400 |
|
|
! ============================================================= |
401 |
|
|
! CALL zenang(cycle_diurne,gmtime,rlat,rlon,rmu0,fract) |
402 |
|
|
! Auteur : A. Campoy et F. Hourdin |
403 |
|
|
! Objet : calculer les valeurs moyennes du cos de l'angle zenithal |
404 |
|
|
! et l'ensoleillement moyen entre gmtime1 et gmtime2 |
405 |
|
|
! connaissant la declinaison, la latitude et la longitude. |
406 |
|
|
|
407 |
|
|
! Dans cette version particuliere, on calcule le rayonnement |
408 |
|
|
! moyen sur l'année à chaque latitude. |
409 |
|
|
! angle zenithal calculé pour obtenir un |
410 |
|
|
! Fit polynomial de l'ensoleillement moyen au sommet de l'atmosphere |
411 |
|
|
! en moyenne annuelle. |
412 |
|
|
! Spécifique de la terre. Utilisé pour les aqua planetes. |
413 |
|
|
|
414 |
|
|
! Rque : Different de la routine angle en ce sens que zenang |
415 |
|
|
! fournit des moyennes de pmu0 et non des valeurs |
416 |
|
|
! instantanees, du coup frac prend toutes les valeurs |
417 |
|
|
! entre 0 et 1. |
418 |
|
|
! Date : premiere version le 13 decembre 1994 |
419 |
|
|
! revu pour GCM le 30 septembre 1996 |
420 |
|
|
! =============================================================== |
421 |
|
|
! longi----INPUT : la longitude vraie de la terre dans son plan |
422 |
|
|
! solaire a partir de l'equinoxe de printemps (degre) |
423 |
|
|
! gmtime---INPUT : temps universel en fraction de jour |
424 |
|
|
! pdtrad---INPUT : pas de temps du rayonnement (secondes) |
425 |
|
|
! lat------INPUT : latitude en degres |
426 |
|
|
! long-----INPUT : longitude en degres |
427 |
|
|
! pmu0-----OUTPUT: angle zenithal moyen entre gmtime et gmtime+pdtrad |
428 |
|
|
! frac-----OUTPUT: ensoleillement moyen entre gmtime et gmtime+pdtrad |
429 |
|
|
! ================================================================ |
430 |
|
|
include "YOMCST.h" |
431 |
|
|
! ================================================================ |
432 |
|
|
LOGICAL cycle_diurne |
433 |
|
|
REAL gmtime |
434 |
|
|
REAL rlat(klon), rlon(klon), rmu0(klon), fract(klon) |
435 |
|
|
! ================================================================ |
436 |
|
|
INTEGER i |
437 |
|
|
REAL gmtime1, gmtime2 |
438 |
|
|
REAL pi_local |
439 |
|
|
|
440 |
|
|
|
441 |
|
|
REAL rmu0m(klon), rmu0a(klon) |
442 |
|
|
|
443 |
|
|
|
444 |
|
|
pi_local = 4.0*atan(1.0) |
445 |
|
|
|
446 |
|
|
! ================================================================ |
447 |
|
|
! Calcul de l'angle zenithal moyen sur la journee |
448 |
|
|
! ================================================================ |
449 |
|
|
|
450 |
|
|
DO i = 1, klon |
451 |
|
|
fract(i) = 1. |
452 |
|
|
! Calcule du flux moyen |
453 |
|
|
IF (abs(rlat(i))<=28.75) THEN |
454 |
|
|
rmu0m(i) = (210.1924+206.6059*cos(0.0174533*rlat(i))**2)/1365. |
455 |
|
|
ELSE IF (abs(rlat(i))<=43.75) THEN |
456 |
|
|
rmu0m(i) = (187.4562+236.1853*cos(0.0174533*rlat(i))**2)/1365. |
457 |
|
|
ELSE IF (abs(rlat(i))<=71.25) THEN |
458 |
|
|
rmu0m(i) = (162.4439+284.1192*cos(0.0174533*rlat(i))**2)/1365. |
459 |
|
|
ELSE |
460 |
|
|
rmu0m(i) = (172.8125+183.7673*cos(0.0174533*rlat(i))**2)/1365. |
461 |
|
|
END IF |
462 |
|
|
END DO |
463 |
|
|
|
464 |
|
|
! ================================================================ |
465 |
|
|
! Avec ou sans cycle diurne |
466 |
|
|
! ================================================================ |
467 |
|
|
|
468 |
|
|
IF (cycle_diurne) THEN |
469 |
|
|
|
470 |
|
|
! On redecompose flux au sommet suivant un cycle diurne idealise |
471 |
|
|
! identique a toutes les latitudes. |
472 |
|
|
|
473 |
|
|
DO i = 1, klon |
474 |
|
|
rmu0a(i) = 2.*rmu0m(i)*sqrt(2.)*pi_local/(4.-pi_local) |
475 |
|
|
rmu0(i) = rmu0a(i)*abs(sin(pi_local*gmtime+pi_local*rlon(i)/360.)) - & |
476 |
|
|
rmu0a(i)/sqrt(2.) |
477 |
|
|
END DO |
478 |
|
|
|
479 |
|
|
DO i = 1, klon |
480 |
|
|
IF (rmu0(i)<=0.) THEN |
481 |
|
|
rmu0(i) = 0. |
482 |
|
|
fract(i) = 0. |
483 |
|
|
ELSE |
484 |
|
|
fract(i) = 1. |
485 |
|
|
END IF |
486 |
|
|
END DO |
487 |
|
|
|
488 |
|
|
! Affichage de l'angel zenitale |
489 |
|
|
! print*,'************************************' |
490 |
|
|
! print*,'************************************' |
491 |
|
|
! print*,'************************************' |
492 |
|
|
! print*,'latitude=',rlat(i),'longitude=',rlon(i) |
493 |
|
|
! print*,'rmu0m=',rmu0m(i) |
494 |
|
|
! print*,'rmu0a=',rmu0a(i) |
495 |
|
|
! print*,'rmu0=',rmu0(i) |
496 |
|
|
|
497 |
|
|
ELSE |
498 |
|
|
|
499 |
|
|
DO i = 1, klon |
500 |
|
|
fract(i) = 0.5 |
501 |
|
|
rmu0(i) = rmu0m(i)*2. |
502 |
|
|
END DO |
503 |
|
|
|
504 |
|
|
END IF |
505 |
|
|
|
506 |
|
|
RETURN |
507 |
|
|
END SUBROUTINE zenang_an |
508 |
|
|
|
509 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
510 |
|
|
|
511 |
|
|
SUBROUTINE writelim_unstruct(klon, phy_nat, phy_alb, phy_sst, phy_bil, phy_rug, & |
512 |
|
|
phy_ice, phy_fter, phy_foce, phy_flic, phy_fsic) |
513 |
|
|
|
514 |
|
|
USE mod_phys_lmdz_para, ONLY: is_omp_master, klon_mpi |
515 |
|
|
USE mod_phys_lmdz_transfert_para, ONLY: gather_omp |
516 |
|
|
#ifdef CPP_XIOS |
517 |
|
|
USE xios |
518 |
|
|
#endif |
519 |
|
|
IMPLICIT NONE |
520 |
|
|
|
521 |
|
|
include "netcdf.inc" |
522 |
|
|
|
523 |
|
|
INTEGER, INTENT (IN) :: klon |
524 |
|
|
REAL, INTENT (IN) :: phy_nat(klon, 360) |
525 |
|
|
REAL, INTENT (IN) :: phy_alb(klon, 360) |
526 |
|
|
REAL, INTENT (IN) :: phy_sst(klon, 360) |
527 |
|
|
REAL, INTENT (IN) :: phy_bil(klon, 360) |
528 |
|
|
REAL, INTENT (IN) :: phy_rug(klon, 360) |
529 |
|
|
REAL, INTENT (IN) :: phy_ice(klon, 360) |
530 |
|
|
REAL, INTENT (IN) :: phy_fter(klon, 360) |
531 |
|
|
REAL, INTENT (IN) :: phy_foce(klon, 360) |
532 |
|
|
REAL, INTENT (IN) :: phy_flic(klon, 360) |
533 |
|
|
REAL, INTENT (IN) :: phy_fsic(klon, 360) |
534 |
|
|
|
535 |
|
|
REAL :: phy_mpi(klon_mpi, 360) ! temporary variable, to store phy_***(:) |
536 |
|
|
! on the whole physics grid |
537 |
|
|
|
538 |
|
|
#ifdef CPP_XIOS |
539 |
|
|
PRINT *, 'writelim: Ecriture du fichier limit' |
540 |
|
|
|
541 |
|
|
CALL gather_omp(phy_foce, phy_mpi) |
542 |
|
|
IF (is_omp_master) CALL xios_send_field('foce_limout',phy_mpi) |
543 |
|
|
|
544 |
|
|
CALL gather_omp(phy_fsic, phy_mpi) |
545 |
|
|
IF (is_omp_master) CALL xios_send_field('fsic_limout',phy_mpi) |
546 |
|
|
|
547 |
|
|
CALL gather_omp(phy_fter, phy_mpi) |
548 |
|
|
IF (is_omp_master) CALL xios_send_field('fter_limout',phy_mpi) |
549 |
|
|
|
550 |
|
|
CALL gather_omp(phy_flic, phy_mpi) |
551 |
|
|
IF (is_omp_master) CALL xios_send_field('flic_limout',phy_mpi) |
552 |
|
|
|
553 |
|
|
CALL gather_omp(phy_sst, phy_mpi) |
554 |
|
|
IF (is_omp_master) CALL xios_send_field('sst_limout',phy_mpi) |
555 |
|
|
|
556 |
|
|
CALL gather_omp(phy_bil, phy_mpi) |
557 |
|
|
IF (is_omp_master) CALL xios_send_field('bils_limout',phy_mpi) |
558 |
|
|
|
559 |
|
|
CALL gather_omp(phy_alb, phy_mpi) |
560 |
|
|
IF (is_omp_master) CALL xios_send_field('alb_limout',phy_mpi) |
561 |
|
|
|
562 |
|
|
CALL gather_omp(phy_rug, phy_mpi) |
563 |
|
|
IF (is_omp_master) CALL xios_send_field('rug_limout',phy_mpi) |
564 |
|
|
#endif |
565 |
|
|
END SUBROUTINE writelim_unstruct |
566 |
|
|
|
567 |
|
|
|
568 |
|
|
|
569 |
|
|
SUBROUTINE writelim(klon, phy_nat, phy_alb, phy_sst, phy_bil, phy_rug, & |
570 |
|
|
phy_ice, phy_fter, phy_foce, phy_flic, phy_fsic) |
571 |
|
|
|
572 |
|
|
USE mod_phys_lmdz_para, ONLY: is_master |
573 |
|
|
USE mod_grid_phy_lmdz, ONLY: klon_glo |
574 |
|
|
USE mod_phys_lmdz_transfert_para, ONLY: gather |
575 |
|
|
USE phys_cal_mod, ONLY: year_len |
576 |
|
|
use netcdf, only: nf90_def_var, nf90_double, nf90_float |
577 |
|
|
IMPLICIT NONE |
578 |
|
|
include "netcdf.inc" |
579 |
|
|
|
580 |
|
|
INTEGER, INTENT (IN) :: klon |
581 |
|
|
REAL, INTENT (IN) :: phy_nat(klon, year_len) |
582 |
|
|
REAL, INTENT (IN) :: phy_alb(klon, year_len) |
583 |
|
|
REAL, INTENT (IN) :: phy_sst(klon, year_len) |
584 |
|
|
REAL, INTENT (IN) :: phy_bil(klon, year_len) |
585 |
|
|
REAL, INTENT (IN) :: phy_rug(klon, year_len) |
586 |
|
|
REAL, INTENT (IN) :: phy_ice(klon, year_len) |
587 |
|
|
REAL, INTENT (IN) :: phy_fter(klon, year_len) |
588 |
|
|
REAL, INTENT (IN) :: phy_foce(klon, year_len) |
589 |
|
|
REAL, INTENT (IN) :: phy_flic(klon, year_len) |
590 |
|
|
REAL, INTENT (IN) :: phy_fsic(klon, year_len) |
591 |
|
|
|
592 |
|
|
REAL :: phy_glo(klon_glo, year_len) ! temporary variable, to store phy_***(:) |
593 |
|
|
! on the whole physics grid |
594 |
|
|
INTEGER :: k |
595 |
|
|
INTEGER ierr |
596 |
|
|
INTEGER dimfirst(3) |
597 |
|
|
INTEGER dimlast(3) |
598 |
|
|
|
599 |
|
|
INTEGER nid, ndim, ntim |
600 |
|
|
INTEGER dims(2), debut(2), epais(2) |
601 |
|
|
INTEGER id_tim |
602 |
|
|
INTEGER id_nat, id_sst, id_bils, id_rug, id_alb |
603 |
|
|
INTEGER id_fter, id_foce, id_fsic, id_flic |
604 |
|
|
|
605 |
|
|
IF (is_master) THEN |
606 |
|
|
|
607 |
|
|
PRINT *, 'writelim: Ecriture du fichier limit' |
608 |
|
|
|
609 |
|
|
ierr = nf_create('limit.nc', IOR(NF_CLOBBER,NF_64BIT_OFFSET), nid) |
610 |
|
|
|
611 |
|
|
ierr = nf_put_att_text(nid, nf_global, 'title', 30, & |
612 |
|
|
'Fichier conditions aux limites') |
613 |
|
|
! ! ierr = NF_DEF_DIM (nid, "points_physiques", klon, ndim) |
614 |
|
|
ierr = nf_def_dim(nid, 'points_physiques', klon_glo, ndim) |
615 |
|
|
ierr = nf_def_dim(nid, 'time', nf_unlimited, ntim) |
616 |
|
|
|
617 |
|
|
dims(1) = ndim |
618 |
|
|
dims(2) = ntim |
619 |
|
|
|
620 |
|
|
#ifdef NC_DOUBLE |
621 |
|
|
ierr = nf90_def_var(nid, 'TEMPS', nf90_double, [ntim], id_tim) |
622 |
|
|
#else |
623 |
|
|
ierr = nf90_def_var(nid, 'TEMPS', nf90_float, [ntim], id_tim) |
624 |
|
|
#endif |
625 |
|
|
ierr = nf_put_att_text(nid, id_tim, 'title', 17, 'Jour dans l annee') |
626 |
|
|
|
627 |
|
|
#ifdef NC_DOUBLE |
628 |
|
|
ierr = nf90_def_var(nid, 'NAT', nf90_double, dims, id_nat) |
629 |
|
|
#else |
630 |
|
|
ierr = nf90_def_var(nid, 'NAT', nf90_float, dims, id_nat) |
631 |
|
|
#endif |
632 |
|
|
ierr = nf_put_att_text(nid, id_nat, 'title', 23, & |
633 |
|
|
'Nature du sol (0,1,2,3)') |
634 |
|
|
|
635 |
|
|
#ifdef NC_DOUBLE |
636 |
|
|
ierr = nf90_def_var(nid, 'SST', nf90_double, dims, id_sst) |
637 |
|
|
#else |
638 |
|
|
ierr = nf90_def_var(nid, 'SST', nf90_float, dims, id_sst) |
639 |
|
|
#endif |
640 |
|
|
ierr = nf_put_att_text(nid, id_sst, 'title', 35, & |
641 |
|
|
'Temperature superficielle de la mer') |
642 |
|
|
|
643 |
|
|
#ifdef NC_DOUBLE |
644 |
|
|
ierr = nf90_def_var(nid, 'BILS', nf90_double, dims, id_bils) |
645 |
|
|
#else |
646 |
|
|
ierr = nf90_def_var(nid, 'BILS', nf90_float, dims, id_bils) |
647 |
|
|
#endif |
648 |
|
|
ierr = nf_put_att_text(nid, id_bils, 'title', 32, & |
649 |
|
|
'Reference flux de chaleur au sol') |
650 |
|
|
|
651 |
|
|
#ifdef NC_DOUBLE |
652 |
|
|
ierr = nf90_def_var(nid, 'ALB', nf90_double, dims, id_alb) |
653 |
|
|
#else |
654 |
|
|
ierr = nf90_def_var(nid, 'ALB', nf90_float, dims, id_alb) |
655 |
|
|
#endif |
656 |
|
|
ierr = nf_put_att_text(nid, id_alb, 'title', 19, 'Albedo a la surface') |
657 |
|
|
|
658 |
|
|
#ifdef NC_DOUBLE |
659 |
|
|
ierr = nf90_def_var(nid, 'RUG', nf90_double, dims, id_rug) |
660 |
|
|
#else |
661 |
|
|
ierr = nf90_def_var(nid, 'RUG', nf90_float, dims, id_rug) |
662 |
|
|
#endif |
663 |
|
|
ierr = nf_put_att_text(nid, id_rug, 'title', 8, 'Rugosite') |
664 |
|
|
|
665 |
|
|
#ifdef NC_DOUBLE |
666 |
|
|
ierr = nf90_def_var(nid, 'FTER', nf90_double, dims, id_fter) |
667 |
|
|
#else |
668 |
|
|
ierr = nf90_def_var(nid, 'FTER', nf90_float, dims, id_fter) |
669 |
|
|
#endif |
670 |
|
|
ierr = nf_put_att_text(nid, id_fter, 'title',10,'Frac. Land') |
671 |
|
|
#ifdef NC_DOUBLE |
672 |
|
|
ierr = nf90_def_var(nid, 'FOCE', nf90_double, dims, id_foce) |
673 |
|
|
#else |
674 |
|
|
ierr = nf90_def_var(nid, 'FOCE', nf90_float, dims, id_foce) |
675 |
|
|
#endif |
676 |
|
|
ierr = nf_put_att_text(nid, id_foce, 'title',11,'Frac. Ocean') |
677 |
|
|
#ifdef NC_DOUBLE |
678 |
|
|
ierr = nf90_def_var(nid, 'FSIC', nf90_double, dims, id_fsic) |
679 |
|
|
#else |
680 |
|
|
ierr = nf90_def_var(nid, 'FSIC', nf90_float, dims, id_fsic) |
681 |
|
|
#endif |
682 |
|
|
ierr = nf_put_att_text(nid, id_fsic, 'title',13,'Frac. Sea Ice') |
683 |
|
|
#ifdef NC_DOUBLE |
684 |
|
|
ierr = nf90_def_var(nid, 'FLIC', nf90_double, dims, id_flic) |
685 |
|
|
#else |
686 |
|
|
ierr = nf90_def_var(nid, 'FLIC', nf90_float, dims, id_flic) |
687 |
|
|
#endif |
688 |
|
|
ierr = nf_put_att_text(nid, id_flic, 'title',14,'Frac. Land Ice') |
689 |
|
|
|
690 |
|
|
ierr = nf_enddef(nid) |
691 |
|
|
IF (ierr/=nf_noerr) THEN |
692 |
|
|
WRITE (*, *) 'writelim error: failed to end define mode' |
693 |
|
|
WRITE (*, *) nf_strerror(ierr) |
694 |
|
|
END IF |
695 |
|
|
|
696 |
|
|
|
697 |
|
|
! write the 'times' |
698 |
|
|
DO k = 1, year_len |
699 |
|
|
#ifdef NC_DOUBLE |
700 |
|
|
ierr = nf_put_var1_double(nid, id_tim, k, dble(k)) |
701 |
|
|
#else |
702 |
|
|
ierr = nf_put_var1_real(nid, id_tim, k, float(k)) |
703 |
|
|
#endif |
704 |
|
|
IF (ierr/=nf_noerr) THEN |
705 |
|
|
WRITE (*, *) 'writelim error with temps(k),k=', k |
706 |
|
|
WRITE (*, *) nf_strerror(ierr) |
707 |
|
|
END IF |
708 |
|
|
END DO |
709 |
|
|
|
710 |
|
|
END IF ! of if (is_master) |
711 |
|
|
|
712 |
|
|
! write the fields, after having collected them on master |
713 |
|
|
|
714 |
|
|
CALL gather(phy_nat, phy_glo) |
715 |
|
|
IF (is_master) THEN |
716 |
|
|
#ifdef NC_DOUBLE |
717 |
|
|
ierr = nf_put_var_double(nid, id_nat, phy_glo) |
718 |
|
|
#else |
719 |
|
|
ierr = nf_put_var_real(nid, id_nat, phy_glo) |
720 |
|
|
#endif |
721 |
|
|
IF (ierr/=nf_noerr) THEN |
722 |
|
|
WRITE (*, *) 'writelim error with phy_nat' |
723 |
|
|
WRITE (*, *) nf_strerror(ierr) |
724 |
|
|
END IF |
725 |
|
|
END IF |
726 |
|
|
|
727 |
|
|
CALL gather(phy_sst, phy_glo) |
728 |
|
|
IF (is_master) THEN |
729 |
|
|
#ifdef NC_DOUBLE |
730 |
|
|
ierr = nf_put_var_double(nid, id_sst, phy_glo) |
731 |
|
|
#else |
732 |
|
|
ierr = nf_put_var_real(nid, id_sst, phy_glo) |
733 |
|
|
#endif |
734 |
|
|
IF (ierr/=nf_noerr) THEN |
735 |
|
|
WRITE (*, *) 'writelim error with phy_sst' |
736 |
|
|
WRITE (*, *) nf_strerror(ierr) |
737 |
|
|
END IF |
738 |
|
|
END IF |
739 |
|
|
|
740 |
|
|
CALL gather(phy_bil, phy_glo) |
741 |
|
|
IF (is_master) THEN |
742 |
|
|
#ifdef NC_DOUBLE |
743 |
|
|
ierr = nf_put_var_double(nid, id_bils, phy_glo) |
744 |
|
|
#else |
745 |
|
|
ierr = nf_put_var_real(nid, id_bils, phy_glo) |
746 |
|
|
#endif |
747 |
|
|
IF (ierr/=nf_noerr) THEN |
748 |
|
|
WRITE (*, *) 'writelim error with phy_bil' |
749 |
|
|
WRITE (*, *) nf_strerror(ierr) |
750 |
|
|
END IF |
751 |
|
|
END IF |
752 |
|
|
|
753 |
|
|
CALL gather(phy_alb, phy_glo) |
754 |
|
|
IF (is_master) THEN |
755 |
|
|
#ifdef NC_DOUBLE |
756 |
|
|
ierr = nf_put_var_double(nid, id_alb, phy_glo) |
757 |
|
|
#else |
758 |
|
|
ierr = nf_put_var_real(nid, id_alb, phy_glo) |
759 |
|
|
#endif |
760 |
|
|
IF (ierr/=nf_noerr) THEN |
761 |
|
|
WRITE (*, *) 'writelim error with phy_alb' |
762 |
|
|
WRITE (*, *) nf_strerror(ierr) |
763 |
|
|
END IF |
764 |
|
|
END IF |
765 |
|
|
|
766 |
|
|
CALL gather(phy_rug, phy_glo) |
767 |
|
|
IF (is_master) THEN |
768 |
|
|
#ifdef NC_DOUBLE |
769 |
|
|
ierr = nf_put_var_double(nid, id_rug, phy_glo) |
770 |
|
|
#else |
771 |
|
|
ierr = nf_put_var_real(nid, id_rug, phy_glo) |
772 |
|
|
#endif |
773 |
|
|
IF (ierr/=nf_noerr) THEN |
774 |
|
|
WRITE (*, *) 'writelim error with phy_rug' |
775 |
|
|
WRITE (*, *) nf_strerror(ierr) |
776 |
|
|
END IF |
777 |
|
|
END IF |
778 |
|
|
|
779 |
|
|
CALL gather(phy_fter, phy_glo) |
780 |
|
|
IF (is_master) THEN |
781 |
|
|
#ifdef NC_DOUBLE |
782 |
|
|
ierr = nf_put_var_double(nid, id_fter, phy_glo) |
783 |
|
|
#else |
784 |
|
|
ierr = nf_put_var_real(nid, id_fter, phy_glo) |
785 |
|
|
#endif |
786 |
|
|
IF (ierr/=nf_noerr) THEN |
787 |
|
|
WRITE (*, *) 'writelim error with phy_fter' |
788 |
|
|
WRITE (*, *) nf_strerror(ierr) |
789 |
|
|
END IF |
790 |
|
|
END IF |
791 |
|
|
|
792 |
|
|
CALL gather(phy_foce, phy_glo) |
793 |
|
|
IF (is_master) THEN |
794 |
|
|
#ifdef NC_DOUBLE |
795 |
|
|
ierr = nf_put_var_double(nid, id_foce, phy_glo) |
796 |
|
|
#else |
797 |
|
|
ierr = nf_put_var_real(nid, id_foce, phy_glo) |
798 |
|
|
#endif |
799 |
|
|
IF (ierr/=nf_noerr) THEN |
800 |
|
|
WRITE (*, *) 'writelim error with phy_foce' |
801 |
|
|
WRITE (*, *) nf_strerror(ierr) |
802 |
|
|
END IF |
803 |
|
|
END IF |
804 |
|
|
|
805 |
|
|
CALL gather(phy_fsic, phy_glo) |
806 |
|
|
IF (is_master) THEN |
807 |
|
|
#ifdef NC_DOUBLE |
808 |
|
|
ierr = nf_put_var_double(nid, id_fsic, phy_glo) |
809 |
|
|
#else |
810 |
|
|
ierr = nf_put_var_real(nid, id_fsic, phy_glo) |
811 |
|
|
#endif |
812 |
|
|
IF (ierr/=nf_noerr) THEN |
813 |
|
|
WRITE (*, *) 'writelim error with phy_fsic' |
814 |
|
|
WRITE (*, *) nf_strerror(ierr) |
815 |
|
|
END IF |
816 |
|
|
END IF |
817 |
|
|
|
818 |
|
|
CALL gather(phy_flic, phy_glo) |
819 |
|
|
IF (is_master) THEN |
820 |
|
|
#ifdef NC_DOUBLE |
821 |
|
|
ierr = nf_put_var_double(nid, id_flic, phy_glo) |
822 |
|
|
#else |
823 |
|
|
ierr = nf_put_var_real(nid, id_flic, phy_glo) |
824 |
|
|
#endif |
825 |
|
|
IF (ierr/=nf_noerr) THEN |
826 |
|
|
WRITE (*, *) 'writelim error with phy_flic' |
827 |
|
|
WRITE (*, *) nf_strerror(ierr) |
828 |
|
|
END IF |
829 |
|
|
END IF |
830 |
|
|
|
831 |
|
|
! close file: |
832 |
|
|
IF (is_master) THEN |
833 |
|
|
ierr = nf_close(nid) |
834 |
|
|
END IF |
835 |
|
|
|
836 |
|
|
END SUBROUTINE writelim |
837 |
|
|
|
838 |
|
|
! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
839 |
|
|
|
840 |
|
|
SUBROUTINE profil_sst(nlon, rlatd, type_profil, phy_sst) |
841 |
|
|
USE dimphy |
842 |
|
|
USE phys_cal_mod , ONLY: year_len |
843 |
|
|
IMPLICIT NONE |
844 |
|
|
|
845 |
|
|
INTEGER nlon, type_profil, i, k, j |
846 |
|
|
REAL :: rlatd(nlon), phy_sst(nlon, year_len) |
847 |
|
|
INTEGER imn, imx, amn, amx, kmn, kmx |
848 |
|
|
INTEGER p, pplus, nlat_max |
849 |
|
|
PARAMETER (nlat_max=72) |
850 |
|
|
REAL x_anom_sst(nlat_max) |
851 |
|
|
CHARACTER (LEN=20) :: modname='profil_sst' |
852 |
|
|
CHARACTER (LEN=80) :: abort_message |
853 |
|
|
|
854 |
|
|
IF (klon/=nlon) THEN |
855 |
|
|
abort_message='probleme de dimensions dans profil_sst' |
856 |
|
|
CALL abort_physic(modname,abort_message,1) |
857 |
|
|
ENDIF |
858 |
|
|
WRITE (*, *) ' profil_sst: type_profil=', type_profil |
859 |
|
|
DO i = 1, year_len |
860 |
|
|
! phy_sst(:,i) = 260.+50.*cos(rlatd(:))**2 |
861 |
|
|
|
862 |
|
|
! Rajout fbrlmd |
863 |
|
|
|
864 |
|
|
IF (type_profil==1) THEN |
865 |
|
|
! Méthode 1 "Control" faible plateau à l'Equateur |
866 |
|
|
DO j = 1, klon |
867 |
|
|
phy_sst(j, i) = 273. + 27.*(1-sin(1.5*rlatd(j))**2) |
868 |
|
|
! PI/3=1.047197551 |
869 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
870 |
|
|
phy_sst(j, i) = 273. |
871 |
|
|
END IF |
872 |
|
|
END DO |
873 |
|
|
END IF |
874 |
|
|
IF (type_profil==2) THEN |
875 |
|
|
! Méthode 2 "Flat" fort plateau à l'Equateur |
876 |
|
|
DO j = 1, klon |
877 |
|
|
phy_sst(j, i) = 273. + 27.*(1-sin(1.5*rlatd(j))**4) |
878 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
879 |
|
|
phy_sst(j, i) = 273. |
880 |
|
|
END IF |
881 |
|
|
END DO |
882 |
|
|
END IF |
883 |
|
|
|
884 |
|
|
|
885 |
|
|
IF (type_profil==3) THEN |
886 |
|
|
! Méthode 3 "Qobs" plateau réel à l'Equateur |
887 |
|
|
DO j = 1, klon |
888 |
|
|
phy_sst(j, i) = 273. + 0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5* & |
889 |
|
|
rlatd(j))**4) |
890 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
891 |
|
|
phy_sst(j, i) = 273. |
892 |
|
|
END IF |
893 |
|
|
END DO |
894 |
|
|
END IF |
895 |
|
|
|
896 |
|
|
IF (type_profil==4) THEN |
897 |
|
|
! Méthode 4 : Méthode 3 + SST+2 "Qobs" plateau réel à l'Equateur |
898 |
|
|
DO j = 1, klon |
899 |
|
|
phy_sst(j, i) = 273. + 0.5*29.*(2-sin(1.5*rlatd(j))**2-sin(1.5* & |
900 |
|
|
rlatd(j))**4) |
901 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
902 |
|
|
phy_sst(j, i) = 273. |
903 |
|
|
END IF |
904 |
|
|
END DO |
905 |
|
|
END IF |
906 |
|
|
|
907 |
|
|
IF (type_profil==5) THEN |
908 |
|
|
! Méthode 5 : Méthode 3 + +2K "Qobs" plateau réel à l'Equateur |
909 |
|
|
DO j = 1, klon |
910 |
|
|
phy_sst(j, i) = 273. + 2. + 0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5 & |
911 |
|
|
*rlatd(j))**4) |
912 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
913 |
|
|
phy_sst(j, i) = 273. + 2. |
914 |
|
|
END IF |
915 |
|
|
|
916 |
|
|
END DO |
917 |
|
|
END IF |
918 |
|
|
|
919 |
|
|
IF (type_profil==6) THEN |
920 |
|
|
! Méthode 6 "cst" valeur constante de SST |
921 |
|
|
DO j = 1, klon |
922 |
|
|
phy_sst(j, i) = 288. |
923 |
|
|
END DO |
924 |
|
|
END IF |
925 |
|
|
|
926 |
|
|
|
927 |
|
|
IF (type_profil==7) THEN |
928 |
|
|
! Méthode 7 "cst" valeur constante de SST +2 |
929 |
|
|
DO j = 1, klon |
930 |
|
|
phy_sst(j, i) = 288. + 2. |
931 |
|
|
END DO |
932 |
|
|
END IF |
933 |
|
|
|
934 |
|
|
p = 0 |
935 |
|
|
IF (type_profil==8) THEN |
936 |
|
|
! Méthode 8 profil anomalies SST du modèle couplé AR4 |
937 |
|
|
DO j = 1, klon |
938 |
|
|
IF (rlatd(j)==rlatd(j-1)) THEN |
939 |
|
|
phy_sst(j, i) = 273. + x_anom_sst(pplus) + & |
940 |
|
|
0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5*rlatd(j))**4) |
941 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
942 |
|
|
phy_sst(j, i) = 273. + x_anom_sst(pplus) |
943 |
|
|
END IF |
944 |
|
|
ELSE |
945 |
|
|
p = p + 1 |
946 |
|
|
pplus = 73 - p |
947 |
|
|
phy_sst(j, i) = 273. + x_anom_sst(pplus) + & |
948 |
|
|
0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5*rlatd(j))**4) |
949 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
950 |
|
|
phy_sst(j, i) = 273. + x_anom_sst(pplus) |
951 |
|
|
END IF |
952 |
|
|
WRITE (*, *) rlatd(j), x_anom_sst(pplus), phy_sst(j, i) |
953 |
|
|
END IF |
954 |
|
|
END DO |
955 |
|
|
END IF |
956 |
|
|
|
957 |
|
|
IF (type_profil==9) THEN |
958 |
|
|
! Méthode 5 : Méthode 3 + -2K "Qobs" plateau réel à l'Equateur |
959 |
|
|
DO j = 1, klon |
960 |
|
|
phy_sst(j, i) = 273. - 2. + 0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5 & |
961 |
|
|
*rlatd(j))**4) |
962 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
963 |
|
|
phy_sst(j, i) = 273. - 2. |
964 |
|
|
END IF |
965 |
|
|
END DO |
966 |
|
|
END IF |
967 |
|
|
|
968 |
|
|
|
969 |
|
|
IF (type_profil==10) THEN |
970 |
|
|
! Méthode 10 : Méthode 3 + +4K "Qobs" plateau réel à l'Equateur |
971 |
|
|
DO j = 1, klon |
972 |
|
|
phy_sst(j, i) = 273. + 4. + 0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5 & |
973 |
|
|
*rlatd(j))**4) |
974 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
975 |
|
|
phy_sst(j, i) = 273. + 4. |
976 |
|
|
END IF |
977 |
|
|
END DO |
978 |
|
|
END IF |
979 |
|
|
|
980 |
|
|
IF (type_profil==11) THEN |
981 |
|
|
! Méthode 11 : Méthode 3 + 4CO2 "Qobs" plateau réel à l'Equateur |
982 |
|
|
DO j = 1, klon |
983 |
|
|
phy_sst(j, i) = 273. + 0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5* & |
984 |
|
|
rlatd(j))**4) |
985 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
986 |
|
|
phy_sst(j, i) = 273. |
987 |
|
|
END IF |
988 |
|
|
END DO |
989 |
|
|
END IF |
990 |
|
|
|
991 |
|
|
IF (type_profil==12) THEN |
992 |
|
|
! Méthode 12 : Méthode 10 + 4CO2 "Qobs" plateau réel à l'Equateur |
993 |
|
|
DO j = 1, klon |
994 |
|
|
phy_sst(j, i) = 273. + 4. + 0.5*27.*(2-sin(1.5*rlatd(j))**2-sin(1.5 & |
995 |
|
|
*rlatd(j))**4) |
996 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
997 |
|
|
phy_sst(j, i) = 273. + 4. |
998 |
|
|
END IF |
999 |
|
|
END DO |
1000 |
|
|
END IF |
1001 |
|
|
|
1002 |
|
|
IF (type_profil==13) THEN |
1003 |
|
|
! Méthode 13 "Qmax" plateau réel à l'Equateur augmenté ! |
1004 |
|
|
DO j = 1, klon |
1005 |
|
|
phy_sst(j, i) = 273. + 0.5*29.*(2-sin(1.5*rlatd(j))**2-sin(1.5* & |
1006 |
|
|
rlatd(j))**4) |
1007 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
1008 |
|
|
phy_sst(j, i) = 273. |
1009 |
|
|
END IF |
1010 |
|
|
END DO |
1011 |
|
|
END IF |
1012 |
|
|
|
1013 |
|
|
IF (type_profil==14) THEN |
1014 |
|
|
! Méthode 13 "Qmax2K" plateau réel à l'Equateur augmenté +2K ! |
1015 |
|
|
DO j = 1, klon |
1016 |
|
|
phy_sst(j, i) = 273. + 2. + 0.5*29.*(2-sin(1.5*rlatd(j))**2-sin(1.5 & |
1017 |
|
|
*rlatd(j))**4) |
1018 |
|
|
IF ((rlatd(j)>1.0471975) .OR. (rlatd(j)<-1.0471975)) THEN |
1019 |
|
|
phy_sst(j, i) = 273. |
1020 |
|
|
END IF |
1021 |
|
|
END DO |
1022 |
|
|
END IF |
1023 |
|
|
|
1024 |
|
|
if (type_profil.EQ.20) then |
1025 |
|
|
print*,'Profile SST 20' |
1026 |
|
|
! Méthode 13 "Qmax2K" plateau réel �| l'Equateur augmenté +2K |
1027 |
|
|
|
1028 |
|
|
do j=1,klon |
1029 |
|
|
phy_sst(j,i)=248.+55.*(1-sin(rlatd(j))**2) |
1030 |
|
|
enddo |
1031 |
|
|
endif |
1032 |
|
|
|
1033 |
|
|
if (type_profil.EQ.21) then |
1034 |
|
|
print*,'Profile SST 21' |
1035 |
|
|
! Méthode 13 "Qmax2K" plateau réel �| l'Equateur augmenté +2K |
1036 |
|
|
do j=1,klon |
1037 |
|
|
phy_sst(j,i)=252.+55.*(1-sin(rlatd(j))**2) |
1038 |
|
|
enddo |
1039 |
|
|
endif |
1040 |
|
|
|
1041 |
|
|
|
1042 |
|
|
|
1043 |
|
|
END DO |
1044 |
|
|
|
1045 |
|
|
! IM beg : verif profil SST: phy_sst |
1046 |
|
|
amn = min(phy_sst(1,1), 1000.) |
1047 |
|
|
amx = max(phy_sst(1,1), -1000.) |
1048 |
|
|
imn = 1 |
1049 |
|
|
kmn = 1 |
1050 |
|
|
imx = 1 |
1051 |
|
|
kmx = 1 |
1052 |
|
|
DO k = 1, year_len |
1053 |
|
|
DO i = 2, nlon |
1054 |
|
|
IF (phy_sst(i,k)<amn) THEN |
1055 |
|
|
amn = phy_sst(i, k) |
1056 |
|
|
imn = i |
1057 |
|
|
kmn = k |
1058 |
|
|
END IF |
1059 |
|
|
IF (phy_sst(i,k)>amx) THEN |
1060 |
|
|
amx = phy_sst(i, k) |
1061 |
|
|
imx = i |
1062 |
|
|
kmx = k |
1063 |
|
|
END IF |
1064 |
|
|
END DO |
1065 |
|
|
END DO |
1066 |
|
|
|
1067 |
|
|
PRINT *, 'profil_sst: imn, kmn, phy_sst(imn,kmn) ', imn, kmn, amn |
1068 |
|
|
PRINT *, 'profil_sst: imx, kmx, phy_sst(imx,kmx) ', imx, kmx, amx |
1069 |
|
|
! IM end : verif profil SST: phy_sst |
1070 |
|
|
|
1071 |
|
|
RETURN |
1072 |
|
|
END SUBROUTINE profil_sst |
1073 |
|
|
|
1074 |
|
|
END MODULE phyaqua_mod |