/[lmdze]/trunk/Sources/phylmd/clesphys2.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/clesphys2.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 209 - (show annotations)
Wed Dec 7 17:37:21 2016 UTC (7 years, 4 months ago) by guez
File size: 1310 byte(s)
The program did not work with cycle_diurne set to false. mu0 in
physiq, which is supposed to be a cosine, was set to -999.999. So prmu
in swu had a value of the order of 1e3. So zrmum1 in sw2s had a value
of the order of 1e3. So zrayl in sw2s had a value of the order of
1e15. So ztray and ptauaz in swclr also had a large value. So zcorae
at line 138 in swclr had a large negative value, which resulted in
overflow at line 138.

This assignment of -999.999 to mu0 dates from somewhere between
revisions 348 and 524 of LMDZ. It was corrected in revision 1068 of
LMDZ with a call to angle which was present in revision 348. However,
procedure angle was removed from LMDZE in revision 22 because it was
not used. Hesitated to bring back angle but, finally, just removed the
option of having no diurnal cycle.

1 module clesphys2
2
3 ! From version 1.3 2005/06/06 13:16:33
4
5 implicit none
6
7 LOGICAL:: soil_model = .TRUE.
8 ! Choix du modele de sol (Thermique ?)
9
10 LOGICAL:: new_oliq = .TRUE.
11 ! Permet de mettre en route la nouvelle parametrisation de l'eau liquide
12
13 ! Pour l'orographie:
14 LOGICAL:: ok_orodr = .TRUE.
15 LOGICAL:: ok_orolf = .TRUE.
16
17 LOGICAL:: ok_limitvrai = .FALSE.
18 ! On peut forcer le modele a lire le fichier SST de la bonne
19 ! annee.
20
21 INTEGER:: nbapp_rad = 12
22 ! nombre d'appels des routines de rayonnements par jour
23
24 logical:: conv_emanuel = .true. ! convection scheme of Emanuel, else Tiedtke
25
26 contains
27
28 subroutine read_clesphys2
29
30 use unit_nml_m, only: unit_nml
31 use nr_util, only: assert
32 use conf_gcm_m, only: day_step, iphysiq
33
34 namelist /clesphys2_nml/soil_model, new_oliq, ok_orodr, ok_orolf, &
35 ok_limitvrai, nbapp_rad, conv_emanuel
36
37 !------------------------------------
38
39 print *, "Enter namelist 'clesphys2_nml'."
40 read(unit=*, nml=clesphys2_nml)
41 write(unit_nml, nml=clesphys2_nml)
42 call assert(mod(day_step / iphysiq, nbapp_rad) == 0, &
43 "read_clesphys2 nbapp_rad")
44 call assert(nbapp_rad >= 4, &
45 "read_clesphys2: minimum 4 calls to radiative transfer per day")
46
47 end subroutine read_clesphys2
48
49 end module clesphys2

  ViewVC Help
Powered by ViewVC 1.1.21