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

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

Parent Directory Parent Directory | Revision Log Revision Log


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

  ViewVC Help
Powered by ViewVC 1.1.21