/[lmdze]/trunk/dyn3d/etat0.f
ViewVC logotype

Diff of /trunk/dyn3d/etat0.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/dyn3d/etat0.f90 revision 15 by guez, Fri Aug 1 15:24:12 2008 UTC trunk/dyn3d/etat0.f revision 90 by guez, Wed Mar 12 21:16:36 2014 UTC
# Line 7  module etat0_mod Line 7  module etat0_mod
7    
8    REAL pctsrf(klon, nbsrf)    REAL pctsrf(klon, nbsrf)
9    ! ("pctsrf(i, :)" is the composition of the surface at horizontal    ! ("pctsrf(i, :)" is the composition of the surface at horizontal
10    !  position "i")    ! position "i")
11    
12    private nbsrf, klon    private nbsrf, klon
13    
# Line 15  contains Line 15  contains
15    
16    SUBROUTINE etat0    SUBROUTINE etat0
17    
18      ! From "etat0_netcdf.F", version 1.3 2005/05/25 13:10:09      ! From "etat0_netcdf.F", version 1.3, 2005/05/25 13:10:09
19    
20      ! This subroutine creates "mask".      use caldyn0_m, only: caldyn0
21        use comconst, only: cpp, kappa, iniconst
22      USE ioipsl, only: flinget, flinclo, flinopen_nozoom, flininfo, histclo      use comgeom, only: rlatu, rlonv, rlonu, rlatv, aire_2d, apoln, apols, &
23             cu_2d, cv_2d, inigeom
24      USE start_init_orog_m, only: start_init_orog, mask, phis      use conf_gcm_m, only: dayref, anneeref
     use start_init_phys_m, only: qsol_2d  
     use startdyn, only: start_inter_3d, start_init_dyn  
25      use dimens_m, only: iim, jjm, llm, nqmx      use dimens_m, only: iim, jjm, llm, nqmx
     use paramet_m, only: ip1jm, ip1jmp1  
     use comconst, only: dtvr, daysec, cpp, kappa, pi  
     use comdissnew, only: lstardis, nitergdiv, nitergrot, niterh, &  
          tetagdiv, tetagrot, tetatemp  
     use indicesol, only: is_oce, is_sic, is_ter, is_lic, epsfra  
     use comvert, only: ap, bp, preff, pa  
26      use dimphy, only: zmasq      use dimphy, only: zmasq
     use conf_gcm_m, only: day_step, iphysiq, dayref, anneeref  
     use comgeom, only: rlatu, rlonv, rlonu, rlatv, aire_2d, apoln, apols, &  
          cu_2d, cv_2d  
     use serre, only: alphax  
27      use dimsoil, only: nsoilmx      use dimsoil, only: nsoilmx
28      use temps, only: itau_dyn, itau_phy, annee_ref, day_ref, dt      use disvert_m, only: ap, bp, preff, pa, disvert
29        use dynredem0_m, only: dynredem0
30        use dynredem1_m, only: dynredem1
31        use exner_hyb_m, only: exner_hyb
32        use geopot_m, only: geopot
33      use grid_atob, only: grille_m      use grid_atob, only: grille_m
34      use grid_change, only: init_dyn_phy, dyn_phy      use grid_change, only: init_dyn_phy, dyn_phy
35      use q_sat_m, only: q_sat      use histclo_m, only: histclo
36      use exner_hyb_m, only: exner_hyb      use indicesol, only: is_oce, is_sic, is_ter, is_lic, epsfra
37      use advtrac_m, only: iniadvtrac      use iniadvtrac_m, only: iniadvtrac
38        use inifilr_m, only: inifilr
39        use massdair_m, only: massdair
40        use netcdf, only: nf90_nowrite
41        use netcdf95, only: nf95_close, nf95_get_var, nf95_gw_var, &
42             nf95_inq_varid, nf95_open
43        use nr_util, only: pi, assert
44        use paramet_m, only: ip1jm, ip1jmp1
45        use phyredem_m, only: phyredem
46      use pressure_var, only: pls, p3d      use pressure_var, only: pls, p3d
47      use dynredem0_m, only: dynredem0      use q_sat_m, only: q_sat
48      use regr_lat_time_coefoz_m, only: regr_lat_time_coefoz      use regr_lat_time_coefoz_m, only: regr_lat_time_coefoz
49      use regr_pr_o3_m, only: regr_pr_o3      use regr_pr_o3_m, only: regr_pr_o3
50      use phyredem_m, only: phyredem      use serre, only: alphax
51        use startdyn, only: start_init_dyn
52        USE start_init_orog_m, only: start_init_orog, mask
53        use start_init_phys_m, only: start_init_phys
54        use start_inter_3d_m, only: start_inter_3d
55        use temps, only: itau_phy, annee_ref, day_ref
56    
57      ! Variables local to the procedure:      ! Variables local to the procedure:
58    
59      REAL latfi(klon), lonfi(klon)      REAL latfi(klon), lonfi(klon)
60      ! (latitude and longitude of a point of the scalar grid identified      ! (latitude and longitude of a point of the scalar grid identified
61      ! by a simple index, in °)      ! by a simple index, in degrees)
62    
63      REAL, dimension(iim + 1, jjm + 1, llm):: uvent, t3d, tpot      REAL, dimension(iim + 1, jjm + 1, llm):: ucov, t3d, teta
64      REAL vvent(iim + 1, jjm, llm)      REAL vcov(iim + 1, jjm, llm)
65    
66      REAL q3d(iim + 1, jjm + 1, llm, nqmx)      REAL q(iim + 1, jjm + 1, llm, nqmx)
67      ! (mass fractions of trace species      ! (mass fractions of trace species
68      ! "q3d(i, j, l)" is at longitude "rlonv(i)", latitude "rlatu(j)"      ! "q(i, j, l)" is at longitude "rlonv(i)", latitude "rlatu(j)"
69      ! and pressure level "pls(i, j, l)".)      ! and pressure level "pls(i, j, l)".)
70    
71      real qsat(iim + 1, jjm + 1, llm) ! mass fraction of saturating water vapor      real qsat(iim + 1, jjm + 1, llm) ! mass fraction of saturating water vapor
# Line 76  contains Line 81  contains
81      real seaice(klon) ! kg m-2      real seaice(klon) ! kg m-2
82      REAL frugs(klon, nbsrf), agesno(klon, nbsrf)      REAL frugs(klon, nbsrf), agesno(klon, nbsrf)
83      REAL rugmer(klon)      REAL rugmer(klon)
84      real, dimension(iim + 1, jjm + 1):: relief, zstd_2d, zsig_2d, zgam_2d      REAL phis(iim + 1, jjm + 1) ! surface geopotential, in m2 s-2
85        real, dimension(iim + 1, jjm + 1):: zmea_2d, zstd_2d, zsig_2d, zgam_2d
86      real, dimension(iim + 1, jjm + 1):: zthe_2d, zpic_2d, zval_2d      real, dimension(iim + 1, jjm + 1):: zthe_2d, zpic_2d, zval_2d
87      real, dimension(iim + 1, jjm + 1):: tsol_2d, psol      real, dimension(iim + 1, jjm + 1):: tsol_2d, qsol_2d, ps
88      REAL zmea(klon), zstd(klon)      REAL zmea(klon), zstd(klon)
89      REAL zsig(klon), zgam(klon)      REAL zsig(klon), zgam(klon)
90      REAL zthe(klon)      REAL zthe(klon)
91      REAL zpic(klon), zval(klon)      REAL zpic(klon), zval(klon)
92      REAL t_ancien(klon, llm), q_ancien(klon, llm)      !      REAL t_ancien(klon, llm), q_ancien(klon, llm)
93      REAL run_off_lic_0(klon)      REAL run_off_lic_0(klon)
94      real clwcon(klon, llm), rnebcon(klon, llm), ratqs(klon, llm)      real clwcon(klon, llm), rnebcon(klon, llm), ratqs(klon, llm)
95      ! déclarations pour lecture glace de mer  
96      INTEGER iml_lic, jml_lic, llm_tmp, ttm_tmp      ! D\'eclarations pour lecture glace de mer :
97      INTEGER itaul(1), fid      INTEGER iml_lic, jml_lic
98      REAL lev(1), date      INTEGER ncid, varid
99      REAL, ALLOCATABLE:: lon_lic(:, :), lat_lic(:, :)      REAL, pointer:: dlon_lic(:), dlat_lic(:)
     REAL, ALLOCATABLE:: dlon_lic(:), dlat_lic(:)  
100      REAL, ALLOCATABLE:: fraclic(:, :) ! fraction land ice      REAL, ALLOCATABLE:: fraclic(:, :) ! fraction land ice
101      REAL flic_tmp(iim + 1, jjm + 1) !fraction land ice temporary      REAL flic_tmp(iim + 1, jjm + 1) ! fraction land ice temporary
102    
103      INTEGER l, ji      INTEGER l, ji
104    
# Line 101  contains Line 106  contains
106      real pks(iim + 1, jjm + 1)      real pks(iim + 1, jjm + 1)
107    
108      REAL masse(iim + 1, jjm + 1, llm)      REAL masse(iim + 1, jjm + 1, llm)
109      REAL phi(ip1jmp1, llm)      REAL phi(iim + 1, jjm + 1, llm)
110      REAL pbaru(ip1jmp1, llm), pbarv(ip1jm, llm)      REAL pbaru(ip1jmp1, llm), pbarv(ip1jm, llm)
111      REAL w(ip1jmp1, llm)      REAL w(ip1jmp1, llm)
112      REAL phystep  
113        real sig1(klon, llm) ! section adiabatic updraft
114        real w01(klon, llm) ! vertical velocity within adiabatic updraft
115    
116      !---------------------------------      !---------------------------------
117    
118      print *, "Call sequence information: etat0"      print *, "Call sequence information: etat0"
119    
120      ! Construct a grid:      CALL iniconst
121    
122      dtvr = daysec / real(day_step)      ! Construct a grid:
     print *, 'dtvr = ', dtvr  
123    
124      pa = 5e4      pa = 5e4
125      CALL iniconst      CALL disvert
126      CALL inigeom      CALL inigeom
127      CALL inifilr      CALL inifilr
128    
# Line 130  contains Line 136  contains
136      ! (with conversion to degrees)      ! (with conversion to degrees)
137      lonfi(klon) = 0.      lonfi(klon) = 0.
138    
139      call start_init_orog(relief, zstd_2d, zsig_2d, zgam_2d, zthe_2d, zpic_2d, &      call start_init_orog(phis, zmea_2d, zstd_2d, zsig_2d, zgam_2d, zthe_2d, &
140           zval_2d) ! also compute "mask" and "phis"           zpic_2d, zval_2d) ! also compute "mask"
141      call init_dyn_phy ! define the mask "dyn_phy" for distinct grid points      call init_dyn_phy ! define the mask "dyn_phy" for distinct grid points
142      zmasq = pack(mask, dyn_phy)      zmasq = pack(mask, dyn_phy)
143      PRINT *, 'Masque construit'      PRINT *, 'Masque construit'
144    
145      CALL start_init_dyn(tsol_2d, psol) ! also compute "qsol_2d"      call start_init_phys(tsol_2d, qsol_2d)
146        CALL start_init_dyn(tsol_2d, phis, ps)
147    
148      ! Compute pressure on intermediate levels:      ! Compute pressure on intermediate levels:
149      forall(l = 1: llm + 1) p3d(:, :, l) = ap(l) + bp(l) * psol      forall(l = 1: llm + 1) p3d(:, :, l) = ap(l) + bp(l) * ps
150      CALL exner_hyb(psol, p3d, pks, pk)      CALL exner_hyb(ps, p3d, pks, pk)
151      IF (MINVAL(pk) == MAXVAL(pk)) stop '"pk" should not be constant'      call assert(MINVAL(pk) /= MAXVAL(pk), '"pk" should not be constant')
152    
153      pls(:, :, :) = preff * (pk(:, :, :) / cpp)**(1. / kappa)      pls = preff * (pk / cpp)**(1. / kappa)
154      PRINT *, "minval(pls(:, :, :)) = ", minval(pls(:, :, :))      PRINT *, "minval(pls) = ", minval(pls)
155      print *, "maxval(pls(:, :, :)) = ", maxval(pls(:, :, :))      print *, "maxval(pls) = ", maxval(pls)
156    
157      uvent(:, :, :) = start_inter_3d('U', rlonv, rlatv, pls)      call start_inter_3d('U', rlonv, rlatv, pls, ucov)
158      forall (l = 1: llm) uvent(:iim, :, l) = uvent(:iim, :, l) * cu_2d(:iim, :)      forall (l = 1: llm) ucov(:iim, :, l) = ucov(:iim, :, l) * cu_2d(:iim, :)
159      uvent(iim+1, :, :) = uvent(1, :, :)      ucov(iim+1, :, :) = ucov(1, :, :)
160    
161      vvent(:, :, :) = start_inter_3d('V', rlonu, rlatu(:jjm), pls(:, :jjm, :))      call start_inter_3d('V', rlonu, rlatu(:jjm), pls(:, :jjm, :), vcov)
162      forall (l = 1: llm) vvent(:iim, :, l) = vvent(:iim, :, l) * cv_2d(:iim, :)      forall (l = 1: llm) vcov(:iim, :, l) = vcov(:iim, :, l) * cv_2d(:iim, :)
163      vvent(iim + 1, :, :) = vvent(1, :, :)      vcov(iim + 1, :, :) = vcov(1, :, :)
164    
165      t3d(:, :, :) = start_inter_3d('TEMP', rlonu, rlatv, pls)      call start_inter_3d('TEMP', rlonu, rlatv, pls, t3d)
166      PRINT *,  'minval(t3d(:, :, :)) = ', minval(t3d(:, :, :))      PRINT *, 'minval(t3d) = ', minval(t3d)
167      print *, "maxval(t3d(:, :, :)) = ", maxval(t3d(:, :, :))      print *, "maxval(t3d) = ", maxval(t3d)
168    
169      tpot(:iim, :, :) = t3d(:iim, :, :) * cpp / pk(:iim, :, :)      teta(:iim, :, :) = t3d(:iim, :, :) * cpp / pk(:iim, :, :)
170      tpot(iim + 1, :, :) = tpot(1, :, :)      teta(iim + 1, :, :) = teta(1, :, :)
171      DO l=1, llm      DO l = 1, llm
172         tpot(:, 1, l) = SUM(aire_2d(:, 1) * tpot(:, 1, l)) / apoln         teta(:, 1, l) = SUM(aire_2d(:, 1) * teta(:, 1, l)) / apoln
173         tpot(:, jjm + 1, l) = SUM(aire_2d(:, jjm + 1) * tpot(:, jjm + 1, l)) &         teta(:, jjm + 1, l) = SUM(aire_2d(:, jjm + 1) * teta(:, jjm + 1, l)) &
174              / apols              / apols
175      ENDDO      ENDDO
176    
177      ! Calcul de l'humidité à saturation :      ! Calcul de l'humidit\'e \`a saturation :
178      qsat(:, :, :) = q_sat(t3d, pls)      qsat = q_sat(t3d, pls)
179      PRINT *, "minval(qsat(:, :, :)) = ", minval(qsat(:, :, :))      PRINT *, "minval(qsat) = ", minval(qsat)
180      print *, "maxval(qsat(:, :, :)) = ", maxval(qsat(:, :, :))      print *, "maxval(qsat) = ", maxval(qsat)
181      IF (MINVAL(qsat) == MAXVAL(qsat)) stop '"qsat" should not be constant'      IF (MINVAL(qsat) == MAXVAL(qsat)) stop '"qsat" should not be constant'
182    
183      ! Water vapor:      ! Water vapor:
184      q3d(:, :, :, 1) = 0.01 * start_inter_3d('R', rlonu, rlatv, pls) * qsat      call start_inter_3d('R', rlonu, rlatv, pls, q(:, :, :, 1))
185      WHERE (q3d(:, :, :, 1) < 0.) q3d(:, :, :, 1) = 1E-10      q(:, :, :, 1) = 0.01 * q(:, :, :, 1) * qsat
186        WHERE (q(:, :, :, 1) < 0.) q(:, :, :, 1) = 1E-10
187      DO l = 1, llm      DO l = 1, llm
188         q3d(:, 1, l, 1) = SUM(aire_2d(:, 1) * q3d(:, 1, l, 1)) / apoln         q(:, 1, l, 1) = SUM(aire_2d(:, 1) * q(:, 1, l, 1)) / apoln
189         q3d(:, jjm + 1, l, 1) &         q(:, jjm + 1, l, 1) &
190              = SUM(aire_2d(:, jjm + 1) * q3d(:, jjm + 1, l, 1)) / apols              = SUM(aire_2d(:, jjm + 1) * q(:, jjm + 1, l, 1)) / apols
191      ENDDO      ENDDO
192    
193      q3d(:, :, :, 2:4) = 0. ! liquid water, radon and lead      q(:, :, :, 2:4) = 0. ! liquid water, radon and lead
194    
195      if (nqmx >= 5) then      if (nqmx >= 5) then
196         ! Ozone:         ! Ozone:
197         call regr_lat_time_coefoz         call regr_lat_time_coefoz
198         call regr_pr_o3(q3d(:, :, :, 5))         call regr_pr_o3(q(:, :, :, 5))
199         ! Convert from mole fraction to mass fraction:         ! Convert from mole fraction to mass fraction:
200         q3d(:, :, :, 5) = q3d(:, :, :, 5)  * 48. / 29.         q(:, :, :, 5) = q(:, :, :, 5) * 48. / 29.
201      end if      end if
202    
203      tsol = pack(tsol_2d, dyn_phy)      tsol = pack(tsol_2d, dyn_phy)
# Line 199  contains Line 207  contains
207      tslab = 0. ! IM "slab" ocean      tslab = 0. ! IM "slab" ocean
208      seaice = 0.      seaice = 0.
209      rugmer = 0.001      rugmer = 0.001
210      zmea = pack(relief, dyn_phy)      zmea = pack(zmea_2d, dyn_phy)
211      zstd = pack(zstd_2d, dyn_phy)      zstd = pack(zstd_2d, dyn_phy)
212      zsig = pack(zsig_2d, dyn_phy)      zsig = pack(zsig_2d, dyn_phy)
213      zgam = pack(zgam_2d, dyn_phy)      zgam = pack(zgam_2d, dyn_phy)
# Line 207  contains Line 215  contains
215      zpic = pack(zpic_2d, dyn_phy)      zpic = pack(zpic_2d, dyn_phy)
216      zval = pack(zval_2d, dyn_phy)      zval = pack(zval_2d, dyn_phy)
217    
218      ! On initialise les sous-surfaces:      ! On initialise les sous-surfaces.
219      ! Lecture du fichier glace de terre pour fixer la fraction de terre      ! Lecture du fichier glace de terre pour fixer la fraction de terre
220      ! et de glace de terre:      ! et de glace de terre :
221      CALL flininfo("landiceref.nc", iml_lic, jml_lic, llm_tmp, &  
222           ttm_tmp, fid)      call nf95_open("landiceref.nc", nf90_nowrite, ncid)
223      ALLOCATE(lat_lic(iml_lic, jml_lic))  
224      ALLOCATE(lon_lic(iml_lic, jml_lic))      call nf95_inq_varid(ncid, 'longitude', varid)
225      ALLOCATE(dlon_lic(iml_lic))      call nf95_gw_var(ncid, varid, dlon_lic)
226      ALLOCATE(dlat_lic(jml_lic))      iml_lic = size(dlon_lic)
227    
228        call nf95_inq_varid(ncid, 'latitude', varid)
229        call nf95_gw_var(ncid, varid, dlat_lic)
230        jml_lic = size(dlat_lic)
231    
232        call nf95_inq_varid(ncid, 'landice', varid)
233      ALLOCATE(fraclic(iml_lic, jml_lic))      ALLOCATE(fraclic(iml_lic, jml_lic))
234      CALL flinopen_nozoom("landiceref.nc", iml_lic, jml_lic, &      call nf95_get_var(ncid, varid, fraclic)
          llm_tmp, lon_lic, lat_lic, lev, ttm_tmp, itaul, date, dt,  &  
          fid)  
     CALL flinget(fid, 'landice', iml_lic, jml_lic, llm_tmp, ttm_tmp &  
          , 1, 1, fraclic)  
     CALL flinclo(fid)  
235    
236      ! Interpolation sur la grille T du modèle :      call nf95_close(ncid)
237      PRINT *, 'Dimensions de "landice"'  
238        ! Interpolation sur la grille T du mod\`ele :
239        PRINT *, 'Dimensions de "landiceref.nc"'
240      print *, "iml_lic = ", iml_lic      print *, "iml_lic = ", iml_lic
241      print *, "jml_lic = ", jml_lic      print *, "jml_lic = ", jml_lic
242    
243      ! Si les coordonnées sont en degrés, on les transforme :      ! Si les coordonn\'ees sont en degr\'es, on les transforme :
244      IF (MAXVAL( lon_lic ) > pi)  THEN      IF (MAXVAL(dlon_lic) > pi) THEN
245         lon_lic = lon_lic * pi / 180.         dlon_lic = dlon_lic * pi / 180.
246      ENDIF      ENDIF
247      IF (maxval( lat_lic ) > pi) THEN      IF (maxval(dlat_lic) > pi) THEN
248         lat_lic = lat_lic * pi/ 180.         dlat_lic = dlat_lic * pi/ 180.
249      ENDIF      ENDIF
250    
     dlon_lic = lon_lic(:, 1)  
     dlat_lic = lat_lic(1, :)  
   
251      flic_tmp(:iim, :) = grille_m(dlon_lic, dlat_lic, fraclic, rlonv(:iim), &      flic_tmp(:iim, :) = grille_m(dlon_lic, dlat_lic, fraclic, rlonv(:iim), &
252           rlatu)           rlatu)
253      flic_tmp(iim + 1, :) = flic_tmp(1, :)      flic_tmp(iim + 1, :) = flic_tmp(1, :)
254    
255        deallocate(dlon_lic, dlat_lic) ! pointers
256    
257      ! Passage sur la grille physique      ! Passage sur la grille physique
258      pctsrf = 0.      pctsrf = 0.
259      pctsrf(:, is_lic) = pack(flic_tmp, dyn_phy)      pctsrf(:, is_lic) = pack(flic_tmp, dyn_phy)
260      ! Adéquation avec le maque terre/mer      ! Ad\'equation avec le maque terre/mer
261      WHERE (pctsrf(:, is_lic) < EPSFRA ) pctsrf(:, is_lic) = 0.      WHERE (pctsrf(:, is_lic) < EPSFRA) pctsrf(:, is_lic) = 0.
262      WHERE (zmasq < EPSFRA) pctsrf(:, is_lic) = 0.      WHERE (zmasq < EPSFRA) pctsrf(:, is_lic) = 0.
263      pctsrf(:, is_ter) = zmasq      pctsrf(:, is_ter) = zmasq
264      where (zmasq > EPSFRA)      where (zmasq > EPSFRA)
# Line 264  contains Line 274  contains
274         end where         end where
275      end where      end where
276    
277      ! Sous-surface océan et glace de mer (pour démarrer on met glace      ! Sous-surface oc\'ean et glace de mer (pour d\'emarrer on met glace
278      ! de mer à 0) :      ! de mer \`a 0) :
279      pctsrf(:, is_oce) = 1. - zmasq      pctsrf(:, is_oce) = 1. - zmasq
280      WHERE (pctsrf(:, is_oce) < EPSFRA) pctsrf(:, is_oce) = 0.      WHERE (pctsrf(:, is_oce) < EPSFRA) pctsrf(:, is_oce) = 0.
281    
282      ! Vérification que somme des sous-surfaces vaut 1:      ! V\'erification que somme des sous-surfaces vaut 1 :
283      ji = count(abs(sum(pctsrf, dim = 2) - 1.) > EPSFRA)      ji = count(abs(sum(pctsrf, dim = 2) - 1.) > EPSFRA)
284      IF (ji /= 0) then      IF (ji /= 0) then
285         PRINT *, 'Problème répartition sous maille pour ', ji, 'points'         PRINT *, 'Probl\`eme r\'epartition sous maille pour ', ji, 'points'
286      end IF      end IF
287    
288      ! Calcul intermédiaire:      ! Calcul interm\'ediaire :
289      CALL massdair(p3d, masse)      CALL massdair(p3d, masse)
290    
291      print *, 'ALPHAX = ', alphax      print *, 'ALPHAX = ', alphax
292    
293      forall  (l = 1:llm)      forall (l = 1:llm)
294         masse(:, 1, l) = SUM(aire_2d(:iim, 1) * masse(:iim, 1, l)) / apoln         masse(:, 1, l) = SUM(aire_2d(:iim, 1) * masse(:iim, 1, l)) / apoln
295         masse(:, jjm + 1, l) = &         masse(:, jjm + 1, l) = &
296              SUM(aire_2d(:iim, jjm + 1) * masse(:iim, jjm + 1, l)) / apols              SUM(aire_2d(:iim, jjm + 1) * masse(:iim, jjm + 1, l)) / apols
# Line 288  contains Line 298  contains
298    
299      ! Initialisation pour traceurs:      ! Initialisation pour traceurs:
300      call iniadvtrac      call iniadvtrac
     ! Ecriture:  
     CALL inidissip(lstardis, nitergdiv, nitergrot, niterh, tetagdiv, &  
          tetagrot, tetatemp)  
     itau_dyn = 0  
301      itau_phy = 0      itau_phy = 0
302      day_ref = dayref      day_ref = dayref
303      annee_ref = anneeref      annee_ref = anneeref
304    
305      CALL geopot(ip1jmp1, tpot, pk , pks,  phis  , phi)      CALL geopot(teta, pk , pks, phis, phi)
306      CALL caldyn0(0, uvent, vvent, tpot, psol, masse, pk, phis, phi, w, &      CALL caldyn0(ucov, vcov, teta, ps, masse, pk, phis, phi, w, pbaru, &
307           pbaru, pbarv, 0)           pbarv)
308      CALL dynredem0("start.nc", dayref, phis)      CALL dynredem0("start.nc", dayref, phis)
309      CALL dynredem1("start.nc", 0., vvent, uvent, tpot, q3d, masse, psol)      CALL dynredem1("start.nc", vcov, ucov, teta, q, masse, ps, itau=0)
   
     ! Ecriture état initial physique:  
     print *, 'dtvr = ', dtvr  
     print *, "iphysiq = ", iphysiq  
     phystep   = dtvr * REAL(iphysiq)  
     print *, 'phystep = ', phystep  
310    
311      ! Initialisations :      ! Initialisations :
312      tsolsrf(:, is_ter) = tsol      tsolsrf(:, is_ter) = tsol
# Line 340  contains Line 340  contains
340      seaice = 0.      seaice = 0.
341    
342      frugs(:, is_oce) = rugmer      frugs(:, is_oce) = rugmer
343      frugs(:, is_ter) = MAX(1.e-05, zstd * zsig / 2)      frugs(:, is_ter) = MAX(1e-5, zstd * zsig / 2)
344      frugs(:, is_lic) = MAX(1.e-05, zstd * zsig / 2)      frugs(:, is_lic) = MAX(1e-5, zstd * zsig / 2)
345      frugs(:, is_sic) = 0.001      frugs(:, is_sic) = 0.001
346      fder = 0.      fder = 0.
347      clwcon = 0.      clwcon = 0.
348      rnebcon = 0.      rnebcon = 0.
349      ratqs = 0.      ratqs = 0.
350      run_off_lic_0 = 0.      run_off_lic_0 = 0.
351        sig1 = 0.
352        w01 = 0.
353    
354      call phyredem("startphy.nc", latfi, lonfi, pctsrf, &      call phyredem("startphy.nc", latfi, lonfi, pctsrf, &
355           tsolsrf, tsoil, tslab, seaice, qsolsrf, qsol, snsrf, albe, alblw, &           tsolsrf, tsoil, tslab, seaice, qsolsrf, qsol, snsrf, albe, alblw, &
356           evap, rain_fall, snow_fall, solsw, sollw, fder, radsol, frugs, &           evap, rain_fall, snow_fall, solsw, sollw, fder, radsol, frugs, &
357           agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &           agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
358           t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0)           t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, w01)
359      CALL histclo      CALL histclo
360    
361    END SUBROUTINE etat0    END SUBROUTINE etat0
Legend:
Removed from v.15  
changed lines
  Added in v.90
  ViewVC Help
Powered by ViewVC 1.1.21