/[lmdze]/trunk/phylmd/Interface_surf/interfsur_lim.f
ViewVC logotype

Annotation of /trunk/phylmd/Interface_surf/interfsur_lim.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 191 - (hide annotations)
Mon May 9 19:56:28 2016 UTC (8 years, 1 month ago) by guez
Original Path: trunk/Sources/phylmd/Interface_surf/interfsur_lim.f
File size: 2435 byte(s)
Extracted the call to read_comdissnew out of conf_gcm.

Made ok_instan a variable of module clesphys, itau_phy a variable of
module phyetat0_m, nid_ins a variable of module ini_histins_m, itap a
variable of new module time_phylmdz, so that histwrite_phy can be
called from any procedure without the need to cascade those variables
into that procedure. Made itau_w a variable of module time_phylmdz so
that it is computed only once per time step of physics.

Extracted variables of module clesphys which were in namelist
conf_phys_nml into their own namelist, clesphys_nml, and created
procedure read_clesphys reading clesphys_nml, to avoid side effect.

No need for double precision in procedure getso4fromfile. Assume there
is a single variable for the whole year in the NetCDF file instead of
one variable per month.

Created generic procedure histwrite_phy and removed procedure
write_histins, following LMDZ. histwrite_phy has only two arguments,
can be called from anywhere, and should manage the logic of writing or
not writing into various history files with various operations. So the
test on ok_instan goes inside histwrite_phy.

Test for raz_date in phyetat0 instead of physiq to avoid side effect.

Created procedure increment_itap to avoid side effect.

Removed unnecessary differences between procedures readsulfate and
readsulfate_pi.

1 guez 54 module interfsur_lim_m
2    
3     implicit none
4    
5     contains
6    
7 guez 191 SUBROUTINE interfsur_lim(dtime, jour, knindex, debut, albedo, z0_new)
8 guez 54
9 guez 150 ! Cette routine sert d'interface entre le mod\`ele atmosph\'erique et
10 guez 54 ! un fichier de conditions aux limites.
11    
12 guez 104 ! Laurent FAIRHEAD, February 2000
13 guez 54
14 guez 98 USE dimphy, ONLY: klon
15     use netcdf, only: NF90_NOWRITE
16     use netcdf95, only: NF95_close, NF95_GET_VAR, NF95_INQ_VARID, NF95_OPEN
17 guez 191 use time_phylmdz, only: itap
18 guez 54
19 guez 98 real, intent(IN):: dtime ! pas de temps de la physique (en s)
20     integer, intent(IN):: jour ! jour a lire dans l'annee
21 guez 54
22 guez 106 integer, intent(in):: knindex(:) ! (knon)
23 guez 150 ! index des points de la surface \`a traiter
24 guez 98
25 guez 150 logical, intent(IN):: debut ! premier appel \`a la physique (initialisation)
26 guez 156 real, intent(out):: albedo(:) ! (knon) albedo lu
27 guez 150 real, intent(out):: z0_new(:) ! (klon) longueur de rugosit\'e lue
28 guez 98
29     ! Local:
30    
31 guez 106 integer knon ! nombre de points dans le domaine a traiter
32    
33 guez 98 integer, save:: lmt_pas ! frequence de lecture des conditions limites
34 guez 54 ! (en pas de physique)
35    
36 guez 98 logical, save:: deja_lu_sur
37 guez 150 ! jour \`a lire d\'ej\`a lu pour une surface pr\'ec\'edente
38 guez 54
39 guez 98 integer, save:: jour_lu_sur
40 guez 54
41 guez 98 ! Champs lus dans le fichier de conditions aux limites :
42     real, allocatable, save:: alb_lu(:), rug_lu(:)
43    
44     integer ncid, varid
45    
46 guez 54 !------------------------------------------------------------
47    
48 guez 106 knon = size(knindex)
49    
50 guez 54 if (debut) then
51 guez 191 lmt_pas = nint(86400. / dtime) ! pour une lecture une fois par jour
52 guez 54 jour_lu_sur = jour - 1
53     allocate(alb_lu(klon))
54     allocate(rug_lu(klon))
55     endif
56    
57 guez 98 if (jour - jour_lu_sur /= 0) deja_lu_sur = .false.
58 guez 54
59     ! Tester d'abord si c'est le moment de lire le fichier
60 guez 191 if (mod(itap - 1, lmt_pas) == 0 .and. .not. deja_lu_sur) then
61 guez 98 call NF95_OPEN('limit.nc', NF90_NOWRITE, ncid)
62 guez 54
63     ! Lecture Albedo
64 guez 98 call NF95_INQ_VARID(ncid, 'ALB', varid)
65     call NF95_GET_VAR(ncid, varid, alb_lu, start=(/1, jour/))
66 guez 54
67 guez 150 ! Lecture rugosit\'e
68 guez 98 call NF95_INQ_VARID(ncid, 'RUG', varid)
69     call NF95_GET_VAR(ncid, varid, rug_lu, start=(/1, jour/))
70 guez 54
71 guez 98 call NF95_CLOSE(ncid)
72 guez 54 deja_lu_sur = .true.
73     jour_lu_sur = jour
74     endif
75    
76     ! Recopie des variables dans les champs de sortie
77 guez 156 albedo = alb_lu(knindex)
78 guez 155 z0_new(:knon) = rug_lu(knindex)
79 guez 98 z0_new(knon + 1:) = 999999.
80 guez 54
81     END SUBROUTINE interfsur_lim
82    
83     end module interfsur_lim_m

  ViewVC Help
Powered by ViewVC 1.1.21