Changeset 4242 for dynamico_lmdz
- Timestamp:
- Jun 5, 2020, 1:59:19 PM (4 years ago)
- Location:
- dynamico_lmdz/simple_physics
- Files:
-
- 1 added
- 10 edited
-
config/DYNAMICO/TEST/field_def_dynamico.xml (modified) (4 diffs)
-
config/DYNAMICO/TEST/field_def_physics.xml (modified) (3 diffs)
-
config/DYNAMICO/TEST/file_def_dynamico.xml (modified) (3 diffs)
-
config/DYNAMICO/TEST/file_def_physics.xml (modified) (1 diff)
-
config/DYNAMICO/TEST/run.def (modified) (1 diff)
-
config/LMDZ/USEFUL (modified) (1 diff)
-
phyparam/DYNAMICO/icosa_phyparam_mod.F90 (modified) (3 diffs)
-
phyparam/physics/phyparam_mod.F90 (modified) (13 diffs)
-
phyparam/physics/radiative_mod.F90 (added)
-
phyparam/physics/radiative_sw.F90 (modified) (10 diffs)
-
phyparam/physics/surface.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
dynamico_lmdz/simple_physics/config/DYNAMICO/TEST/field_def_dynamico.xml
r4239 r4242 13 13 <field id="preff" grid_ref="scalar" unit="Pa"/> 14 14 <field id="ap" axis_ref="levp1" long_name="hybrid A coefficient at layer interface" /> 15 <field id="bp" axis_ref="levp1" long_name="hybrid B coefficient at layer interface" />16 <field id="mid_ap" axis_ref="lev" 17 <field id="mid_bp" axis_ref="lev" 15 <field id="bp" axis_ref="levp1" long_name="hybrid B coefficient at layer interface" /> 16 <field id="mid_ap" axis_ref="lev" long_name="hybrid A coefficient at midpoints" /> 17 <field id="mid_bp" axis_ref="lev" long_name="hybrid B coefficient at midpoints" /> 18 18 19 19 <field_group domain_ref="i"> 20 <field id="ps" /> 20 <field id="ps" standard_name="surface_pressure" long_name="Surface pressure" unit="Pa"/> 21 <field id="phis" standard_name="surface_geopotential" long_name="Surface geopotential" unit="m2/m2"/> 21 22 <field id="dps" /> 22 23 <field id="Ai" /> 23 <field id="phis" />24 24 <field id="phi" /> 25 25 <field id="precl" /> … … 27 27 <field id="Q2_col_int" /> 28 28 <field id="ps_init" /> 29 29 30 30 <field_group axis_ref="lev"> 31 31 <field id="mass" /> 32 <field id="p" />33 32 <field id="dmass"/> 34 <field id="theta" />35 33 <field id="dyn_q" /> 36 34 <field id="pk"/> 37 <field id="ulon"/> 38 <field id="ulat"/> 39 <field id="uz"/> 40 <field id="omega"/> 41 <field id="temp"/> 35 <field id="p" standard_name="air_pressure" long_name="Pressure" unit="Pa"/> 36 <field id="temp" standard_name="air_temperature" long_name="Temperature" unit="K"/> 37 <field id="theta" standard_name="pot_temp" long_name="Potential temperature" unit="m/s"/> 38 <field id="ulon" standard_name="eastward_wind" long_name="Zonal wind" unit="m/s"/> 39 <field id="ulat" standard_name="northward_wind" long_name="Meridional wind" unit="m/s"/> 40 <field id="uz" standard_name="upward_air_velocity" long_name="Vertical velocity" unit="m/s"/> 41 <field id="omega" standard_name="lagrangian_tendency_of_air_pressure" long_name="Vertical pressure velocity" unit="Pa/s"/> 42 42 <field id="kinetic_trisk"/> 43 43 <field id="kinetic"/> … … 93 93 94 94 </field_group> 95 <field id="u850" /> 96 <field id=" v850"/>97 <field id=" w850"/>98 <field id=" t850"/>99 <field id=" omega850"/>100 <field id=" u500"/>101 <field id=" v500"/>102 <field id=" w500"/>103 <field id="t500" />104 <field id="omega 500"/>105 <field id=" SST"/>106 95 96 <field id="u850" long_name="Zonal wind at 850 hPa" unit="m/s"/> 97 <field id="u500" long_name="Zonal wind at 500 hPa" unit="m/s"/> 98 <field id="v850" long_name="Meridional wind at 850 hPa" unit="m/s"/> 99 <field id="v500" long_name="Meridional wind at 850 hPa" unit="m/s"/> 100 <field id="w850" long_name="Vertical wind at 850 hPa" unit="m/s"/> 101 <field id="w500" long_name="Vertical wind at 500 hPa" unit="m/s"/> 102 <field id="t850" long_name="Temperature at 850 hPa" unit="K"/> 103 <field id="t500" long_name="Temperature at 500 hPa" unit="K"/> 104 <field id="omega850" long_name="Vertical pressure velocity at 850 hPa" unit="Pa/s"/> 105 <field id="omega500" long_name="Vertical pressure velocity at 500 hPa" unit="Pa/s"/> 106 <field id="SST"/> 107 107 <field id="ulon_850_500" field_ref="p" axis_ref="lev_pressure" /> 108 108 109 109 <field_group axis_ref="levp1"> 110 <field id="geopot" />110 <field id="geopot" standard_name="geopotential" long_name="geopotential" unit="m2/s2"/> 111 111 <field id="geopot_init"/> 112 112 </field_group> … … 124 124 125 125 </field_definition> 126 -
dynamico_lmdz/simple_physics/config/DYNAMICO/TEST/field_def_physics.xml
r4239 r4242 1 1 <!-- =========================================================================================================== --> 2 <!-- field_def_phy param.xml -->2 <!-- field_def_physics.xml --> 3 3 <!-- Definition of all existing variables that can be output from phyparam --> 4 4 <!-- =========================================================================================================== --> … … 8 8 <field_group id="phyparam_output"> 9 9 <field_group domain_ref="i"> 10 <field id="phyparam_mu0" /> 10 <field id="phyparam_ts" name="Tsurf" long_name="Surface temperature" unit="K"/> 11 <field id="phyparam_mu0" name="mu0" long_name="Cosine zenithal angle"/> 12 <field id="phyparam_swsurf" name="swsurf" long_name="Upward SW flux at surface"/> 13 <field id="phyparam_swtop" name="swtop" long_name="Downward SW fluw at TOA" unit="W/m2" /> 14 <field id="phyparam_lwsurf" name="lwsurf" long_name="LW flux at surface" unit="W/m2" /> 11 15 <field id="phyparam_alb" /> 12 <field id="phyparam_ts" />13 16 <field id="phyparam_ps" /> 14 17 <field id="phyparam_slp" /> 15 <field id="phyparam_swsurf" />16 <field id="phyparam_lwsurf" />17 18 18 19 <field id="phyparam_coslon" /> … … 20 21 <field id="phyparam_coslat" /> 21 22 <field id="phyparam_sinlat" /> 23 24 <field_group axis_ref="lev"> 25 <field id="phyparam_u" name="u" long_name="Zonal velocity" unit="m/s" /> 26 <field id="phyparam_v" name="v" long_name="Meridional velocity" unit="m/s" /> 27 <field id="phyparam_temp" name="T" long_name="Temperature" unit="K" /> 28 <field id="phyparam_geop" /> 29 <field id="phyparam_plev" /> 30 31 <field id="phyparam_du" name="du" long_name="Zonal velocity tendency" unit="m/s2" /> 32 <field id="phyparam_dv" name="dv" long_name="Meridional velocity tendency" unit="m/s2" /> 33 <field id="phyparam_dt" name="dT" long_name="Temperature tendency" unit="K/s" /> 34 <field id="phyparam_dtsw" name="dT_SW" long_name="Temperature tendency due to shortwave radiation" unit="K/s" /> 35 <field id="phyparam_dtlw" name="dT_LW" long_name="Temperature tendency due to longwave radiation" unit="K/s" /> 36 </field_group> 37 38 <field_group axis_ref="levp1"> 39 <field id="phyparam_swflux_down" name="swflux_down" long_name="Downward SW radiative flux" unit="W/m2" /> 40 </field_group> 41 22 42 </field_group> 23 43 </field_group> 24 44 25 45 </field_definition> 26 27 <!--28 call writefield('u','Vent zonal moy','m/s',pu)29 call writefield('v','Vent meridien moy','m/s',pv)30 call writefield('temp','Temperature','K',pt)31 call writefield('geop','Geopotential','m2/s2',pphi)32 call writefield('plev','plev','Pa',pplev(:,1:nlayer))33 34 call writefield('du','du',' ',pdu)35 call writefield('dv','du',' ',pdv)36 call writefield('dt','du',' ',pdt)37 call writefield('dtsw','dtsw',' ',zdtsw)38 call writefield('dtlw','dtlw',' ',zdtlw)39 40 call writefield('ts','Surface temper','K',tsurf)41 call writefield('coslon','coslon',' ',coslon)42 call writefield('sinlon','sinlon',' ',sinlon)43 call writefield('coslat','coslat',' ',coslat)44 call writefield('sinlat','sinlat',' ',sinlat)45 call writefield('mu0','mu0',' ',mu0)46 call writefield('alb','alb',' ',albedo)47 call writefield('ps','Surface pressure','Pa',pplev(:,1))48 call writefield('slp','Sea level pressure','Pa',zpmer)49 call writefield('swsurf','SW surf','Pa',zfluxsw)50 call writefield('lwsurf','LW surf','Pa',zfluxlw)51 --> -
dynamico_lmdz/simple_physics/config/DYNAMICO/TEST/file_def_dynamico.xml
r4239 r4242 7 7 <file id="dynamico_native" convention="CF" enabled="false" output_freq="3h" sync_freq="3h" output_level="10" timeseries="none" description="Simple physics" > 8 8 9 <field_group id=" dcmip2016_output_field_once" operation="once" freq_offset="0ts" ts_enabled="true">9 <field_group id="output_field_once" operation="once" freq_offset="0ts" ts_enabled="true"> 10 10 <field field_ref="timestep" name="mdt" enabled="false" /> 11 11 <field field_ref="preff" name="P0" long_name="reference pressure" /> … … 17 17 18 18 19 <field_group id="dcmip2016_output_field" ts_enabled="true"> 20 <field field_ref="ps" name="PS" standard_name="surface_pressure" long_name="Surface pressure" unit="Pa" /> 21 <field field_ref="phis" operation="once" freq_offset="0ts" name="PHIS" standard_name="surface_geopotential" long_name="Surface geopotential" unit="m2/m2"/> 19 <field_group id="output_field" ts_enabled="true"> 20 21 <field field_ref="phis" operation="once" freq_offset="0ts" /> 22 <field field_ref="ps" name="PS" /> 23 <field field_ref="ulon" /> 24 <field field_ref="ulat" /> 25 <field field_ref="temp" /> 26 <field field_ref="theta" /> 27 28 <field field_ref="u850" /> 29 <field field_ref="v850" /> 30 <field field_ref="t850" /> 31 <field field_ref="omega850" /> 32 <field field_ref="u500" /> 33 <field field_ref="v500" /> 34 <field field_ref="t500" /> 35 <field field_ref="omega500" /> 22 36 23 37 <!-- 24 <field field_ref="geopot" name="PHI" standard_name="geopotential" long_name="geopotential" unit="m2/s2"/>25 <field field_ref="uz" name="W" standard_name="upward_air_velocity" long_name="Vertical velocity" unit="m/s"/>26 <field field_ref="omega" name="OMEGA" standard_name="lagrangian_tendency_of_air_pressure" long_name="Vertical pressure velocity" unit="Pa/s"/>27 <field field_ref="p" name="P" standard_name="air_pressure" long_name="Pressure" unit="Pa"/>28 <field field_ref="w850" name="W850" long_name="Vertical velocity at 850 hPa" unit="m/s"/>29 <field field_ref="w500" name="W500" long_name="Vertical velocity at 500 hPa" unit="m/s"/>38 <field field_ref="geopot" /> 39 <field field_ref="uz" /> 40 <field field_ref="omega" /> 41 <field field_ref="p" /> 42 <field field_ref="w850" /> 43 <field field_ref="w500" /> 30 44 --> 31 <field field_ref="ulon" name="U" standard_name="eastward_wind" long_name="Zonal wind" unit="m/s"/>32 <field field_ref="ulat" name="V" standard_name="northward_wind" long_name="Meridional wind" unit="m/s"/>33 <field field_ref="temp" name="T" standard_name="air_temperature" long_name="Temperature" unit="K"/>34 45 35 <field field_ref="u850" name="U850" long_name="Zonal wind at 850 hPa" unit="m/s"/>36 <field field_ref="v850" name="V850" long_name="Meridional wind at 850 hPa" unit="m/s"/>37 <field field_ref="t850" name="T850" long_name="Temperature at 850 hPa" unit="K"/>38 <field field_ref="omega850" name="OMEGA850" long_name="Vertical pressure velocity at 850 hPa" unit="Pa/s"/>39 <field field_ref="u500" name="U500" long_name="Zonal wind at 500 hPa" unit="m/s"/>40 <field field_ref="v500" name="V500" long_name="Meridional wind at 850 hPa" unit="m/s"/>41 <field field_ref="t500" name="T500" long_name="Temperature at 850 hPa" unit="K"/>42 <field field_ref="omega500" name="OMEGA500" long_name="Vertical pressure velocity at 500 hPa" unit="Pa/s"/>43 46 44 47 </field_group> … … 53 56 <file id="dynamico_regular" enabled="true" output_freq="24h" sync_freq="24h" output_level="10" description="Simple physics" timeseries="none" > 54 57 55 <field_group group_ref="dcmip2016_output_field_once"/> 56 57 <field_group group_ref="dcmip2016_output_field" domain_ref="regular_one_degree" /> 58 58 <field_group group_ref="output_field_once"/> 59 <field_group group_ref="output_field" domain_ref="regular_one_degree" /> 59 60 <variable name="model" type="string" > dynamico </variable> 60 61 <variable name="grid" type="string" > hex </variable> -
dynamico_lmdz/simple_physics/config/DYNAMICO/TEST/file_def_physics.xml
r4239 r4242 4 4 5 5 <field_group id="phyparam_output_field" ts_enabled="true"> 6 <field field_ref="phyparam_ts" unit="K"/> 7 <field field_ref="phyparam_swsurf" unit="W/m2"/> 8 <field field_ref="phyparam_lwsurf" unit="W/m2"/> 6 7 <field field_ref="phyparam_mu0" /> 8 <field field_ref="phyparam_ts" /> 9 <field field_ref="phyparam_swsurf" /> 10 <field field_ref="phyparam_swtop" /> 11 <field field_ref="phyparam_lwsurf" /> 12 13 <field field_ref="phyparam_u" /> 14 <field field_ref="phyparam_v" /> 15 <field field_ref="phyparam_temp" /> 16 <field field_ref="phyparam_du" /> 17 <field field_ref="phyparam_dv" /> 18 <field field_ref="phyparam_dtsw" /> 19 <field field_ref="phyparam_dtlw" /> 20 <field field_ref="phyparam_swflux_down" /> 21 9 22 </field_group> 10 23 -
dynamico_lmdz/simple_physics/config/DYNAMICO/TEST/run.def
r4239 r4242 33 33 34 34 #-------------- Physics ------------- 35 itau_physics= 135 itau_physics=5 36 36 physics=plugin 37 #physics_log_level=DBG 38 37 39 # 0.25 day = 6h 38 40 period_sort=0.25 -
dynamico_lmdz/simple_physics/config/LMDZ/USEFUL
r4230 r4242 3 3 (cd TEST_PARAM/; time ./gcm.e | tee gcm.log | grep '* pas ' && ./check.sh IRENE) 4 4 grep WRN TEST_PARAM/gcm.log | sort -u 5 6 # IPSL / Camelot 7 module load gnu/7.2.0 -
dynamico_lmdz/simple_physics/phyparam/DYNAMICO/icosa_phyparam_mod.F90
r4240 r4242 199 199 CALL allocate_field(f_write_llmp1, field_t, type_real, llm+1, name='phyparam_write_llmp1') 200 200 writefield1_plugin => plugin_writefield1 201 writefield2_plugin => plugin_writefield2 201 202 END SUBROUTINE init_plugin_writefield 202 203 … … 206 207 CHARACTER(*), INTENT(IN) :: name, longname, unit 207 208 REAL, INTENT(IN) :: var(:) 208 WRITELOG(*,*) TRIM(name), ' : ', TRIM(longname), MINVAL(var), MAXVAL(var), inout%it 209 WRITELOG(*,*) TRIM(name), ' : ', TRIM(longname), SHAPE(var), inout%it 210 WRITELOG(*,*) TRIM(name), ' : ', MINVAL(var), MAXVAL(var) 209 211 LOG_INFO('writefield1') 210 212 CALL unpack_field(f_write2d, var) … … 213 215 214 216 SUBROUTINE plugin_writefield2(name,longname,unit, var) 217 USE physics_interface_mod, ONLY : unpack_field, inout => physics_inout 218 USE output_field_mod, ONLY : output_field 219 USE icosa, ONLY : llm 215 220 CHARACTER(*), INTENT(IN) :: name, longname, unit 216 221 REAL, INTENT(IN) :: var(:,:) 222 INTEGER :: nlev 223 WRITELOG(*,*) TRIM(name), ' : ', TRIM(longname), SHAPE(var), inout%it 224 WRITELOG(*,*) TRIM(name), ' : ', MINVAL(var), MAXVAL(var) 225 LOG_INFO('writefield2') 226 nlev = SIZE(var, 2) 227 IF(nlev==llm) THEN 228 CALL unpack_field(f_write_llm, var) 229 CALL output_field('phyparam_'//TRIM(name), f_write_llm) 230 ELSEIF(nlev==llm+1) THEN 231 CALL unpack_field(f_write_llmp1, var) 232 CALL output_field('phyparam_'//TRIM(name), f_write_llmp1) 233 END IF 217 234 END SUBROUTINE plugin_writefield2 218 235 -
dynamico_lmdz/simple_physics/phyparam/physics/phyparam_mod.F90
r4236 r4242 3 3 USE callkeys 4 4 USE comgeomfi 5 USE surface 5 6 IMPLICIT NONE 6 7 PRIVATE 7 8 SAVE 8 9 9 REAL :: nan ! initialized to NaN with adequate compiler options 10 11 REAL, PARAMETER :: pi=2*ASIN(1.), solarcst=1370., stephan=5.67e-08, & 12 ps_rad=1.e5, height_scale=10000., ref_temp=285., & 13 capcal_nosoil=1e5, tsoil_init=300. 14 15 ! internal state, written to / read from disk at checkpoint / restart 16 REAL, ALLOCATABLE :: tsurf(:), tsoil(:,:), capcal(:), fluxgrd(:) 17 !$OMP THREADPRIVATE(tsurf, tsoil, capcal, fluxgrd) 18 19 ! precomputed arrays 20 REAL, ALLOCATABLE :: rnatur(:), albedo(:),emissiv(:), z0(:), inertie(:) 21 !$OMP THREADPRIVATE( rnatur, albedo, emissiv, z0, inertie) 10 REAL, PARAMETER :: ref_temp=285., capcal_nosoil=1e5, tsoil_init=300. 22 11 23 12 INTEGER :: icount … … 36 25 & pdu,pdv,pdt,pdpsrf) BIND(C, name='phyparam_phyparam') 37 26 USE phys_const, ONLY : g, rcp, r, unjours 38 USE surface, ONLY : soil39 27 USE turbulence, ONLY : vdif 40 28 USE convection, ONLY : convadj 29 USE radiative_mod, ONLY : radiative_tendencies 41 30 USE writefield_mod, ONLY : writefield 42 31 … … 79 68 REAL, DIMENSION(ngrid) :: mu0 80 69 INTEGER :: j,l,ig,nlevel,igout 70 LOGICAL :: lwrite 81 71 ! 82 72 REAL :: zday, zdtime … … 90 80 REAL zdum2(ngrid,nlayer) 91 81 REAL zdum3(ngrid,nlayer) 92 REAL zdtlw(ngrid,nlayer),zdtsw(ngrid,nlayer) 93 REAL zfluxsw(ngrid),zfluxlw(ngrid), fluxrad(ngrid) 82 REAL fluxrad(ngrid) 94 83 95 84 WRITELOG(*,*) 'latitude0', ngrid, lati(1:2), lati(ngrid-1:ngrid) … … 134 123 zdtsrf(:) = 0. 135 124 125 IF(abs(zday-zday_last-period_sort)<=ptimestep/unjours/10.) THEN 126 WRITELOG(*,*) 'zday, zday_last SORTIE ', zday, zday_last 127 LOG_INFO('phyparam') 128 zday_last=zday 129 lwrite = .TRUE. 130 END IF 131 136 132 !----------------------------------------------------------------------- 137 133 ! calcul du geopotentiel aux niveaux intercouches … … 167 163 ! --------------------------------------- 168 164 169 IF(callrad) CALL radiative_tendencies( ngrid, igout, nlayer, &165 IF(callrad) CALL radiative_tendencies(lwrite, ngrid, igout, nlayer, & 170 166 gmtime, ptimestep*float(iradia), zday, pplev, pplay, pt, & 171 pdt, zdtlw, zfluxlw, zdtsw, zfluxsw, fluxrad, mu0)167 pdt, fluxrad) 172 168 173 169 !----------------------------------------------------------------------- … … 281 277 LOG_DBG('phyparam') 282 278 283 if(abs(zday-zday_last-period_sort)<=ptimestep/unjours/10.) then 284 WRITELOG(*,*) 'zday, zday_last SORTIE ', zday, zday_last 285 LOG_INFO('phyparam') 286 287 zday_last=zday 288 ! Ecriture/extension de la coordonnee temps 279 IF(lwrite) THEN 289 280 290 281 do ig=1,ngridmax … … 299 290 300 291 call writefield('du','du',' ',pdu) 301 call writefield('dv','du',' ',pdv) 302 call writefield('dt','du',' ',pdt) 303 call writefield('dtsw','dtsw',' ',zdtsw) 304 call writefield('dtlw','dtlw',' ',zdtlw) 292 call writefield('dv','dv',' ',pdv) 293 call writefield('dt','dt',' ',pdt) 305 294 306 295 call writefield('ts','Surface temper','K',tsurf) … … 309 298 call writefield('coslat','coslat',' ',coslat) 310 299 call writefield('sinlat','sinlat',' ',sinlat) 311 call writefield('mu0','mu0',' ',mu0)312 300 call writefield('alb','alb',' ',albedo) 313 301 call writefield('ps','Surface pressure','Pa',pplev(:,1)) 314 302 call writefield('slp','Sea level pressure','Pa',zpmer) 315 call writefield('swsurf','SW surf','Pa',zfluxsw) 316 call writefield('lwsurf','LW surf','Pa',zfluxlw) 317 318 endif 303 END IF 319 304 320 305 END SUBROUTINE phyparam 321 322 SUBROUTINE radiative_tendencies(ngrid, igout, nlayer, &323 gmtime, zdtime, zday, pplev, pplay, pt, &324 pdt, zdtlw, zfluxlw, zdtsw, zfluxsw, fluxrad, mu0)325 USE planet326 USE phys_const, ONLY : planet_rad, unjours327 USE astronomy, ONLY : orbite, solarlong328 USE solar, ONLY : solang, zenang, mucorr329 USE radiative_sw, ONLY : sw330 USE radiative_lw, ONLY : lw331 332 INTEGER, INTENT(IN) :: ngrid, igout, nlayer333 REAL, INTENT(IN) :: gmtime, zdtime, zday, &334 & pplev(ngrid,nlayer+1), pplay(ngrid, nlayer), &335 & pt(ngrid, nlayer+1)336 REAL, INTENT(INOUT) :: pdt(ngrid,nlayer)337 REAL, INTENT(OUT) :: zdtlw(ngrid,nlayer), & ! long-wave temperature tendency338 & zfluxlw(ngrid), & ! long-wave flux at surface339 & zdtsw(ngrid,nlayer), & ! short-wave temperature tendency340 & zfluxsw(ngrid), & ! short-wave flux at surface341 & fluxrad(ngrid), & ! surface flux342 mu0(ngrid) ! ??343 ! local variables344 REAL :: fract(ngrid), dtrad(ngrid, nlayer)345 REAL :: zls, zinsol, tsurf2, ztim1,ztim2,ztim3, dist_sol, declin346 REAL :: zplanck(ngrid)347 INTEGER :: ig, l348 349 ! 2.1 Insolation350 ! --------------------------------------------------351 352 CALL solarlong(zday,zls)353 CALL orbite(zls,dist_sol,declin)354 355 IF(diurnal) THEN356 IF ( .TRUE. ) then357 ztim1=SIN(declin)358 ztim2=COS(declin)*COS(2.*pi*(zday-.5))359 ztim3=-COS(declin)*SIN(2.*pi*(zday-.5))360 CALL solang(ngrid,sinlon,coslon,sinlat,coslat, &361 & ztim1,ztim2,ztim3, &362 & mu0,fract)363 ELSE364 CALL zenang(ngrid,zls,gmtime,zdtime,lati,long,mu0,fract)365 ENDIF366 367 IF(lverbose) THEN368 WRITELOG(*,*) 'day, declin, sinlon,coslon,sinlat,coslat'369 WRITELOG(*,*) zday, declin, &370 & sinlon(igout),coslon(igout), &371 & sinlat(igout),coslat(igout)372 LOG_DBG('radiative_tendencies')373 ENDIF374 ELSE375 WRITELOG(*,*) 'declin,ngrid,planet_rad', declin, ngrid, planet_rad376 LOG_DBG('radiative_tendencies')377 CALL mucorr(ngrid,declin,lati,mu0,fract,height_scale,planet_rad)378 ENDIF379 380 zinsol=solarcst/(dist_sol*dist_sol)381 382 ! 2.2 Radiative tendencies and fluxes:383 ! --------------------------------------------------384 385 CALL sw(ngrid,nlayer,diurnal,coefvis,albedo, &386 & pplev,ps_rad, &387 & mu0,fract,zinsol, &388 & zfluxsw,zdtsw, &389 & lverbose)390 391 CALL lw(ngrid,nlayer,coefir,emissiv, &392 & pplev,ps_rad,tsurf,pt, &393 & zfluxlw,zdtlw, &394 & lverbose)395 396 ! 2.4 surface fluxes397 ! ------------------------------398 399 DO ig=1,ngrid400 fluxrad(ig)=emissiv(ig)*zfluxlw(ig) &401 & +zfluxsw(ig)*(1.-albedo(ig))402 tsurf2 = tsurf(ig)*tsurf(ig)403 zplanck(ig)=emissiv(ig)*stephan*tsurf2*tsurf2404 fluxrad(ig)=fluxrad(ig)-zplanck(ig)405 ENDDO406 407 ! 2.5 Temperature tendencies408 ! --------------------------409 410 DO l=1,nlayer411 DO ig=1,ngrid412 dtrad(ig,l)=zdtsw(ig,l)+zdtlw(ig,l)413 pdt(ig,l)=pdt(ig,l)+dtrad(ig,l)414 ENDDO415 ENDDO416 417 IF(lverbose) THEN418 WRITELOG(*,*) 'Diagnostics for radiation'419 WRITELOG(*,*) 'albedo, emissiv, mu0,fract,Frad,Planck'420 WRITELOG(*,*) albedo(igout),emissiv(igout),mu0(igout), &421 & fract(igout), &422 & fluxrad(igout),zplanck(igout)423 WRITELOG(*,*) 'Tlay Play Plev dT/dt SW dT/dt LW (K/day)'424 WRITELOG(*,*) 'unjours',unjours425 DO l=1,nlayer426 WRITELOG(*,'(3f15.5,2e15.2)') pt(igout,l), &427 & pplay(igout,l),pplev(igout,l), &428 & zdtsw(igout,l),zdtlw(igout,l)429 ENDDO430 LOG_DBG('radiative_tendencies')431 ENDIF432 433 END SUBROUTINE radiative_tendencies434 306 435 307 SUBROUTINE alloc(ngrid, nlayer) BIND(C, name='phyparam_alloc') … … 437 309 !$cython wrapper def alloc(ngrid, nlayer) : phy.phyparam_alloc(ngrid, nlayer) 438 310 USE astronomy, ONLY : iniorbit 439 USE surface, ONLY : zc,zd440 311 INTEGER, INTENT(IN), VALUE :: ngrid, nlayer 441 LOGICAL, PARAMETER :: firstcall=.TRUE.442 312 ! allocate arrays for internal state 443 313 ALLOCATE(tsurf(ngrid)) … … 447 317 ALLOCATE(zc(ngrid,nsoilmx),zd(ngrid,nsoilmx)) 448 318 ! allocate precomputed arrays 449 ALLOCATE(rnatur(ngrid) )450 ALLOCATE( albedo(ngrid),emissiv(ngrid),z0(ngrid),inertie(ngrid))319 ALLOCATE(rnatur(ngrid), albedo(ngrid), emissiv(ngrid)) 320 ALLOCATE(z0(ngrid),inertie(ngrid)) 451 321 CALL iniorbit 452 322 END SUBROUTINE alloc … … 455 325 !$cython header void phyparam_precompute(); 456 326 !$cython wrapper def precompute() : phy.phyparam_precompute() 457 USE surface458 327 ! precompute time-independent arrays 459 328 rnatur(:) = 1. 329 inertie(:) = (1.-rnatur(:))*I_mer+rnatur(:)*I_ter 330 z0(:) = (1.-rnatur(:))*Cd_mer+rnatur(:)*Cd_ter 460 331 emissiv(:) = (1.-rnatur(:))*emi_mer+rnatur(:)*emi_ter 461 inertie(:) = (1.-rnatur(:))*I_mer+rnatur(:)*I_ter462 332 albedo(:) = (1.-rnatur(:))*alb_mer+rnatur(:)*alb_ter 463 z0(:) = (1.-rnatur(:))*Cd_mer+rnatur(:)*Cd_ter464 333 END SUBROUTINE precompute 465 334 … … 468 337 !$cython wrapper def coldstart (ngrid, timestep): phy.phyparam_coldstart(ngrid, timestep) 469 338 ! create internal state to start a run without a restart file 470 USE surface, ONLY : soil471 339 INTEGER, INTENT(IN), VALUE :: ngrid 472 340 REAL, INTENT(IN), VALUE :: ptimestep -
dynamico_lmdz/simple_physics/phyparam/physics/radiative_sw.F90
r4233 r4242 23 23 END SUBROUTINE monGATHER 24 24 25 PURE subroutine monscatter(n ,a,index,b)26 INTEGER, INTENT(IN) :: n ,index(n)27 REAL, INTENT(IN) :: b(n)28 REAL, INTENT(OUT) :: a(n)25 PURE subroutine monscatter(ngrid, n,a,index,b) 26 INTEGER, INTENT(IN) :: ngrid, n,index(n) 27 REAL, INTENT(IN) :: b(n) 28 REAL, INTENT(OUT) :: a(ngrid) 29 29 INTEGER :: i 30 DO i=1,n 31 a(index(i))=b(i) 32 END DO 30 IF(n<ngrid) THEN 31 a(:)=0. 32 DO i=1,n 33 a(index(i))=b(i) 34 END DO 35 ELSE 36 a(:)=b(:) 37 END IF 33 38 end subroutine monscatter 34 39 35 SUBROUTINE sw(ngrid,nlayer,ldiurn, & 36 coefvis,albedo, & 37 plevel,ps_rad,pmu,pfract,psolarf0, & 38 fsrfvis,dtsw, & 39 lwrite) 40 SUBROUTINE sw(ngrid,nlayer,ldiurn, coefvis,albedo, & 41 & plevel,ps_rad,pmu,pfract,psolarf0, & 42 & fsrfvis,dtsw, lverbose, lwrite) 40 43 USE phys_const, ONLY : cpp, g 44 USE writefield_mod, ONLY : writefield 45 41 46 !======================================================================= 42 47 ! … … 54 59 ! 55 60 INTEGER, INTENT(IN) :: ngrid,nlayer 56 LOGICAL, INTENT(IN) :: l write,ldiurn61 LOGICAL, INTENT(IN) :: ldiurn, lverbose, lwrite 57 62 REAL, INTENT(IN) :: albedo(ngrid), coefvis 58 63 REAL, INTENT(IN) :: pmu(ngrid), pfract(ngrid) … … 66 71 67 72 INTEGER ig,l,nlevel,index(ngrid),ncount,igout 68 REAL ztrdir(ngrid,nlayer+1),ztrref(ngrid,nlayer+1) 69 REAL zfsrfref(ngrid) 73 REAL ztrdir(ngrid,nlayer+1), & ! transmission coef for downward flux 74 ztrref(ngrid,nlayer+1) ! transmission coef for upward flux 75 REAL zfsrfref(ngrid) ! upward flux at surface 70 76 REAL z1(ngrid) 71 77 REAL zu(ngrid,nlayer+1) 72 78 REAL tau0 79 80 REAL :: flux_in(ngrid), & ! incoming solar flux 81 & flux_ref(ngrid), & ! flux reflected by surface 82 & flux_down(ngrid, nlayer+1), & ! downward flux 83 & flux_up(ngrid, nlayer+1) ! upward flux 84 REAL :: buf1(ngrid), buf2(ngrid, nlayer+1) ! buffers for output 73 85 74 86 !----------------------------------------------------------------------- … … 125 137 ! ----------------------------------------------------------- 126 138 139 DO ig=1,ncount 140 flux_in(ig) = psolarf0*zfract(ig)*zmu(ig) 141 ENDDO 142 127 143 DO l=1,nlevel 128 144 DO ig=1,ncount 129 ztrdir(ig,l)=exp(-zu(ig,l)/zmu(ig)) 130 ENDDO 131 ENDDO 132 133 IF (lwrite) THEN 145 ztrdir(ig,l)=exp(-zu(ig,l)/zmu(ig)) ! transmission coefficient 146 flux_down(ig,l) = flux_in(ig)*ztrdir(ig,l) 147 ENDDO 148 ENDDO 149 150 IF (lverbose) THEN 134 151 igout=ncount/2+1 135 152 WRITELOG(*,*) … … 149 166 DO l=1,nlayer 150 167 DO ig=1,ncount 151 zdtsw(ig,l)=g*psolarf0*zfract(ig)*zmu(ig)* & 152 (ztrdir(ig,l+1)-ztrdir(ig,l))/ & 153 (cpp*(zplev(ig,l)-zplev(ig,l+1))) 154 ENDDO 155 ENDDO 156 IF (lwrite) THEN 168 ! m.cp.dT = dflux/dz 169 ! m = -(dp/dz)/g 170 ! flux = ztrdir * psolarf0*zfract*zmu 171 IF(.FALSE.) THEN 172 zdtsw(ig,l)=g*psolarf0*zfract(ig)*zmu(ig) & 173 & *(ztrdir(ig,l+1)-ztrdir(ig,l)) & 174 & /(cpp*(zplev(ig,l)-zplev(ig,l+1))) 175 ELSE 176 zdtsw(ig,l)=(g/cpp) & 177 & * (flux_down(ig,l+1)-flux_down(ig,l)) & 178 & / (zplev(ig,l)-zplev(ig,l+1)) 179 END IF 180 ENDDO 181 ENDDO 182 IF (lverbose) THEN 157 183 WRITELOG(*,*) 158 184 WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires:' … … 169 195 170 196 DO ig=1,ncount 171 z1(ig)=zfract(ig)*zmu(ig)*psolarf0*ztrdir(ig,1) 172 zflux(ig)=(1.-zalb(ig))*z1(ig) 173 zfsrfref(ig)= zalb(ig)*z1(ig) 174 ENDDO 175 IF (l write) THEN197 z1(ig)=zfract(ig)*zmu(ig)*psolarf0*ztrdir(ig,1) ! total (down) 198 zflux(ig)=(1.-zalb(ig))*z1(ig) ! absorbed (net) 199 zfsrfref(ig)= zalb(ig)*z1(ig) ! reflected (up) 200 ENDDO 201 IF (lverbose) THEN 176 202 WRITELOG(*,*) 177 203 WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires:' … … 182 208 183 209 !----------------------------------------------------------------------- 184 ! 5.calcul des tra ansmissions depuis le sol, cas diffus:210 ! 5.calcul des transmissions depuis le sol, cas diffus: 185 211 ! ------------------------------------------------------ 186 212 … … 191 217 ENDDO 192 218 193 IF (l write) THEN219 IF (lverbose) THEN 194 220 WRITELOG(*,*) 195 221 WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires' … … 214 240 215 241 !----------------------------------------------------------------------- 216 ! 10. sor ites eventuelles:242 ! 10. sorties eventuelles: 217 243 ! ------------------------ 218 244 219 IF (l write) THEN245 IF (lverbose) THEN 220 246 WRITELOG(*,*) 221 247 WRITELOG(*,*) 'Diagnostique des taux de chauffage solaires:' … … 226 252 ENDIF 227 253 228 IF (ldiurn) THEN 229 fsrfvis(:)=0. 230 CALL monscatter(ncount,fsrfvis,index,zflux) 231 dtsw(:,:)=0. 232 DO l=1,nlayer 233 CALL monscatter(ncount,dtsw(1,l),index,zdtsw(1,l)) 234 ENDDO 235 ELSE 236 WRITELOG(*,*) 'NOT DIURNE' 237 fsrfvis(:)=zflux(:) 238 dtsw(:,:)=zdtsw(:,:) 239 ENDIF 240 ! call dump2d(iim,jjm-1,zflux(2),'ZFLUX ') 241 ! call dump2d(iim,jjm-1,fsrfvis(2),'FSRVIS ') 242 ! call dump2d(iim,jjm-1,ztrdir(2,1),'ztrdir ') 243 ! call dump2d(iim,jjm-1,pmu(2),'pmu ') 244 ! call dump2d(iim,jjm-1,pfract(2),'pfract ') 245 ! call dump2d(iim,jjm-1,albedo(2),'albedo ') 246 ! call dump2d(iim,jjm-1,ztrdir(2,1),'ztrdir ') 254 CALL monscatter(ngrid,ncount,fsrfvis,index,zflux) 255 DO l=1,nlayer 256 CALL monscatter(ngrid,ncount,dtsw(:,l),index,zdtsw(:,l)) 257 ENDDO 258 259 IF(lwrite) THEN 260 CALL monscatter(ngrid, ncount,buf1,index,flux_in) 261 CALL writefield('swtop','SW down TOA','W/m2',buf1) 262 263 DO l=1,nlayer+1 264 CALL monscatter(ngrid,ncount,buf2(:,l),index,flux_down(:,l)) 265 ENDDO 266 CALL writefield('swflux_down','Downwards SW flux','W/m2',buf2) 267 268 END IF 269 247 270 LOG_INFO('rad_sw') 248 271 -
dynamico_lmdz/simple_physics/phyparam/physics/surface.F90
r4233 r4242 15 15 ! precomputed variables 16 16 REAL :: lambda 17 REAL,ALLOCATABLE :: dz1(:),dz2(:) 18 !$OMP THREADPRIVATE(dz1,dz2,zc,lambda) 19 20 ! internal state variables needed for checkpoint/restart 17 REAL, ALLOCATABLE :: dz1(:),dz2(:) 18 !$OMP THREADPRIVATE(dz1,dz2,zc) 19 REAL, ALLOCATABLE :: rnatur(:), albedo(:),emissiv(:), z0(:), inertie(:) 20 !$OMP THREADPRIVATE( rnatur, albedo, emissiv, z0, inertie) 21 22 ! internal state, written to / read from disk at checkpoint / restart 21 23 REAL, ALLOCATABLE :: zc(:,:),zd(:,:) 22 24 !$OMP THREADPRIVATE(zc,zd) 23 24 PUBLIC :: soil, zc, zd 25 REAL, ALLOCATABLE :: tsurf(:), tsoil(:,:), capcal(:), fluxgrd(:) 26 !$OMP THREADPRIVATE(tsurf, tsoil, capcal, fluxgrd) 27 28 PUBLIC :: soil, zc, zd, & 29 rnatur, albedo, emissiv, z0, inertie, & 30 tsurf, tsoil, capcal, fluxgrd 25 31 26 32 CONTAINS
Note: See TracChangeset
for help on using the changeset viewer.
