Parent Directory | Revision Log
Rename phisinit to phis in restart.nc: clearer, same name as Fortran variable. In aaam_bud, use rlat and rlon from phyetat0_m instead of having these module variables associated to actual arguments in physiq. In clmain, too many wind variables make the procedure hard to understand. Use yu(:knon, 1) and yv(:knon, 1) instead of u1lay(:knon) and v1lay(:knon). Note that when yu(:knon, 1) and yv(:knon, 1) are used as actual arguments, they are probably copied to new arrays since the elements are not contiguous. Rename yu10m to wind10m because this is the norm of wind vector, not its zonal component. Rename yustar to ustar. Rename uzon and vmer to u1 and v1 since these are wind components at first layer and u1 and v1 are the names of corresponding dummy arguments in stdlevvar. In clmain, rename yzlev to zlev. In clmain, screenc, stdlevvar and coefcdrag, remove the code corresponding to zxli true (not used in LMDZ either). Subroutine ustarhb becomes a function. Simplifications using the fact that zx_alf2 = 0 and zx_alf1 = 1 (discarding the possibility to change this). In procedure vdif_kcay, remove unused dummy argument plev. Remove useless computations of sss and sssq. In clouds_gno, exp(100.) would overflow in single precision. Set maximum to exp(80.) instead. In physiq, use u(:, 1) and v(:, 1) as arguments to phytrac instead of creating ad hoc variables yu1 and yv1. In stdlevvar, rename dummy argument u_10m to wind10m, following the corresponding modification in clmain. Simplifications using the fact that ok_pred = 0 and ok_corr = 1 (discarding the possibility to change this).
1 | module conf_phys_m |
2 | |
3 | implicit none |
4 | |
5 | integer:: iflag_pbl = 1 ! for the planetary boundary layer |
6 | REAL:: rad_chau1 = 13., rad_chau2 = 9. |
7 | real:: epmax = 0.993 ! \'efficacit\'e de pr\'ecipitation |
8 | integer:: iflag_clw = 0 |
9 | |
10 | contains |
11 | |
12 | subroutine conf_phys |
13 | |
14 | ! From phylmd/conf_phys.F90, version 1.7 2005/07/05 07:21:23 |
15 | |
16 | ! Configuration de la "physique" de LMDZ. |
17 | |
18 | USE clesphys, ONLY: read_clesphys |
19 | use clesphys2, only: read_clesphys2 |
20 | USE comfisrtilp, ONLY: cld_lc_con, cld_lc_lsc, cld_tau_con, & |
21 | cld_tau_lsc, coef_eva, ffallv_con, ffallv_lsc, iflag_pdf, reevap_ice |
22 | use nr_util, only: assert |
23 | use unit_nml_m, only: unit_nml |
24 | USE yomcst, ONLY: read_YOMCST |
25 | |
26 | namelist /conf_phys_nml/ epmax, iflag_clw, cld_lc_lsc, cld_lc_con, & |
27 | cld_tau_lsc, cld_tau_con, ffallv_lsc, ffallv_con, coef_eva, & |
28 | reevap_ice, iflag_pdf, iflag_pbl |
29 | |
30 | namelist /nuagecom/ rad_chau1, rad_chau2 |
31 | |
32 | !----------------------------------------------------------- |
33 | |
34 | print *, "Call sequence information: conf_phys" |
35 | call read_clesphys2 |
36 | call read_YOMCST |
37 | |
38 | cld_lc_lsc = 2.6e-4 |
39 | cld_lc_con = 2.6e-4 |
40 | cld_tau_lsc = 3600. |
41 | cld_tau_con = 3600. |
42 | ffallv_lsc = 1. |
43 | ffallv_con = 1. |
44 | coef_eva = 2.e-5 |
45 | reevap_ice = .false. |
46 | iflag_pdf = 0 |
47 | |
48 | print *, "Enter namelist 'conf_phys_nml'." |
49 | read(unit=*, nml=conf_phys_nml) |
50 | write(unit_nml, nml=conf_phys_nml) |
51 | |
52 | call assert(iflag_pbl <=2 .or. iflag_pbl >= 6, & |
53 | "conf_phys: bad value for iflag_pbl") |
54 | call read_clesphys |
55 | |
56 | print *, "Enter namelist 'nuagecom'." |
57 | read(unit=*, nml=nuagecom) |
58 | write(unit_nml, nml=nuagecom) |
59 | |
60 | end subroutine conf_phys |
61 | |
62 | end module conf_phys_m |
ViewVC Help | |
Powered by ViewVC 1.1.21 |