source: LMDZ6/trunk/libf/phylmd/ecrad/radiation/radiation_scheme.F90 @ 4787

Last change on this file since 4787 was 4787, checked in by idelkadi, 5 months ago

Continued work on implementing the Ecrad code in LMDZ (Integration of tropospheric aerosols when using the RRTMG and ECCKD gas modules):

  • Added aerosol optical properties files
  • Addition of output dignostics for tropospheric aerosols
File size: 27.3 KB
Line 
1! AI mars 2021
2! ====================== Interface between ECRAD and LMDZ ====================
3! radiation_scheme.F90 appelee dans radlwsw_m.F90 si iflag_rttm = 2
4! revoir toutes les parties avec "AI ATTENTION"
5! Mars 2021 :
6!             - Revoir toutes les parties commentees AI ATTENTION
7!             1. Traitement des aerosols
8!             2. Verifier les parametres times issus de LMDZ (calcul issed)
9!             3. Configuration a partir de namelist
10!             4. frac_std = 0.75     
11! Juillet 2023 :
12!             
13! ============================================================================
14module interface_lmdz_ecrad
15
16IMPLICIT NONE       
17
18contains
19
20SUBROUTINE RADIATION_SCHEME &
21! Inputs
22     & (KIDIA, KFDIA, KLON, KLEV, KAEROSOL, NSW, &
23     &  namelist_file, ok_3Deffect, IDAY, TIME, &
24     &  PSOLAR_IRRADIANCE, &
25     &  PMU0, PTEMPERATURE_SKIN, &
26     &  PALBEDO_DIF, PALBEDO_DIR, &
27     &  PEMIS, PEMIS_WINDOW, &
28     &  PGELAM, PGEMU, &
29     &  PPRESSURE_H, PTEMPERATURE_H, PQ, PQSAT, &
30     &  PCO2, PCH4, PN2O, PNO2, PCFC11, PCFC12, PHCFC22, &
31     &  PCCL4, PO3, PO2, &
32     &  PCLOUD_FRAC, PQ_LIQUID, PQ_ICE, PQ_SNOW, &
33     &  ZRE_LIQUID_UM, ZRE_ICE_UM, &
34     &  PAEROSOL_OLD, PAEROSOL, &
35! Outputs
36     &  PFLUX_SW, PFLUX_LW, PFLUX_SW_CLEAR, PFLUX_LW_CLEAR, &
37     &  PFLUX_SW_DN, PFLUX_LW_DN, PFLUX_SW_DN_CLEAR, PFLUX_LW_DN_CLEAR, &
38     &  PFLUX_SW_UP, PFLUX_LW_UP, PFLUX_SW_UP_CLEAR, PFLUX_LW_UP_CLEAR, &
39     &  PFLUX_DIR, PFLUX_DIR_CLEAR, PFLUX_DIR_INTO_SUN, &
40     &  PFLUX_UV, PFLUX_PAR, PFLUX_PAR_CLEAR, &
41     &  PEMIS_OUT, PLWDERIVATIVE, &
42     &  PSWDIFFUSEBAND, PSWDIRECTBAND, &
43     &  ecrad_cloud_cover_sw)
44
45! RADIATION_SCHEME - Interface to modular radiation scheme
46!
47! (C) Copyright 2015- ECMWF.
48!
49! This software is licensed under the terms of the Apache Licence Version 2.0
50! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
51!
52! In applying this licence, ECMWF does not waive the privileges and immunities
53! granted to it by virtue of its status as an intergovernmental organisation
54! nor does it submit to any jurisdiction.
55!
56! PURPOSE
57! -------
58!   The modular radiation scheme is contained in a separate
59!   library. This routine puts the the IFS arrays into appropriate
60!   objects, computing the additional data that is required, and sends
61!   it to the radiation scheme.  It returns net fluxes and surface
62!   flux components needed by the rest of the model.
63!
64!   Lower case is used for variables and types taken from the
65!   radiation library
66!
67! INTERFACE
68! ---------
69!    RADIATION_SCHEME is called from RADLSWR. The
70!    SETUP_RADIATION_SCHEME routine (in the RADIATION_SETUP module)
71!    should have been run first.
72!
73! AUTHOR
74! ------
75!   Robin Hogan, ECMWF
76!   Original: 2015-09-16
77!
78! MODIFICATIONS
79! -------------
80!
81! TO DO
82! -----
83!
84!-----------------------------------------------------------------------
85
86! Modules from ifs or ifsaux libraries
87USE PARKIND1 , ONLY : JPIM, JPRB
88USE YOMHOOK  , ONLY : LHOOK, DR_HOOK
89USE RADIATION_SETUP
90USE YOMCST   , ONLY : RSIGMA ! Stefan-Boltzmann constant
91!USE RADIATION_SETUP, ONLY : SETUP_RADIATION_SCHEME, &
92!                         &  config_type, driver_config_type, &
93!                         &  NWEIGHT_UV,  IBAND_UV,  WEIGHT_UV, &
94!                         &  NWEIGHT_PAR, IBAND_PAR, WEIGHT_PAR, &
95!                         &  ITYPE_TROP_BG_AER,  TROP_BG_AER_MASS_EXT, &
96!                         &  ITYPE_STRAT_BG_AER, STRAT_BG_AER_MASS_EXT, &
97!                         &  ISolverSpartacus
98
99! Modules from radiation library
100USE radiation_single_level,   ONLY : single_level_type
101USE radiation_thermodynamics, ONLY : thermodynamics_type
102USE radiation_gas
103USE radiation_cloud,          ONLY : cloud_type
104USE radiation_aerosol,        ONLY : aerosol_type
105USE radiation_flux,           ONLY : flux_type
106USE radiation_interface,      ONLY : radiation, set_gas_units
107USE radiation_save,           ONLY : save_inputs
108
109USE mod_phys_lmdz_para
110
111IMPLICIT NONE
112
113! INPUT ARGUMENTS
114! *** Array dimensions and ranges
115INTEGER(KIND=JPIM),INTENT(IN) :: KIDIA    ! Start column to process
116INTEGER(KIND=JPIM),INTENT(IN) :: KFDIA    ! End column to process
117!INTEGER, INTENT(IN) :: KIDIA, KFDIA
118INTEGER(KIND=JPIM),INTENT(IN) :: KLON     ! Number of columns
119INTEGER(KIND=JPIM),INTENT(IN) :: KLEV     ! Number of levels
120!INTEGER, INTENT(IN) :: KLON, KLEV
121!INTEGER(KIND=JPIM),INTENT(IN) :: KAEROLMDZ ! Number of aerosol types
122INTEGER(KIND=JPIM),INTENT(IN) :: KAEROSOL
123INTEGER(KIND=JPIM),INTENT(IN) :: NSW ! Numbe of bands
124
125! AI ATTENTION
126!INTEGER, PARAMETER :: KAEROSOL = 12
127
128! *** Single-level fields
129REAL(KIND=JPRB),   INTENT(IN) :: PSOLAR_IRRADIANCE ! (W m-2)
130REAL(KIND=JPRB),   INTENT(IN) :: PMU0(KLON) ! Cosine of solar zenith ang
131REAL(KIND=JPRB),   INTENT(IN) :: PTEMPERATURE_SKIN(KLON) ! (K)
132! Diffuse and direct components of surface shortwave albedo
133!REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIF(KLON,YRERAD%NSW)
134!REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIR(KLON,YRERAD%NSW)
135REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIF(KLON,NSW)
136REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIR(KLON,NSW)
137! Longwave emissivity outside and inside the window region
138REAL(KIND=JPRB),   INTENT(IN) :: PEMIS(KLON)
139REAL(KIND=JPRB),   INTENT(IN) :: PEMIS_WINDOW(KLON)
140! Longitude (radians), sine of latitude
141REAL(KIND=JPRB),   INTENT(IN) :: PGELAM(KLON)
142REAL(KIND=JPRB),   INTENT(IN) :: PGEMU(KLON)
143! Land-sea mask
144!REAL(KIND=JPRB),   INTENT(IN) :: PLAND_SEA_MASK(KLON)
145
146! *** Variables on half levels
147REAL(KIND=JPRB),   INTENT(IN) :: PPRESSURE_H(KLON,KLEV+1)    ! (Pa)
148REAL(KIND=JPRB),   INTENT(IN) :: PTEMPERATURE_H(KLON,KLEV+1) ! (K)
149
150! *** Gas mass mixing ratios on full levels
151REAL(KIND=JPRB),   INTENT(IN) :: PQ(KLON,KLEV)
152! AI
153REAL(KIND=JPRB),   INTENT(IN) :: PQSAT(KLON,KLEV)
154REAL(KIND=JPRB),   INTENT(IN) :: PCO2
155REAL(KIND=JPRB),   INTENT(IN) :: PCH4
156REAL(KIND=JPRB),   INTENT(IN) :: PN2O
157REAL(KIND=JPRB),   INTENT(IN) :: PNO2
158REAL(KIND=JPRB),   INTENT(IN) :: PCFC11
159REAL(KIND=JPRB),   INTENT(IN) :: PCFC12
160REAL(KIND=JPRB),   INTENT(IN) :: PHCFC22
161REAL(KIND=JPRB),   INTENT(IN) :: PCCL4
162REAL(KIND=JPRB),   INTENT(IN) :: PO3(KLON,KLEV) ! AI (kg/kg) ATTENTION (Pa*kg/kg)
163REAL(KIND=JPRB),   INTENT(IN) :: PO2
164
165! *** Cloud fraction and hydrometeor mass mixing ratios
166REAL(KIND=JPRB),   INTENT(IN) :: PCLOUD_FRAC(KLON,KLEV)
167REAL(KIND=JPRB),   INTENT(IN) :: PQ_LIQUID(KLON,KLEV)
168REAL(KIND=JPRB),   INTENT(IN) :: PQ_ICE(KLON,KLEV)
169!REAL(KIND=JPRB),   INTENT(IN) :: PQ_RAIN(KLON,KLEV)
170REAL(KIND=JPRB),   INTENT(IN) :: PQ_SNOW(KLON,KLEV)
171
172! *** Aerosol mass mixing ratios
173REAL(KIND=JPRB),   INTENT(IN) :: PAEROSOL_OLD(KLON,6,KLEV)
174REAL(KIND=JPRB),   INTENT(IN) :: PAEROSOL(KLON,KLEV,KAEROSOL)
175
176!REAL(KIND=JPRB),   INTENT(IN) :: PCCN_LAND(KLON)
177!REAL(KIND=JPRB),   INTENT(IN) :: PCCN_SEA(KLON)
178
179!AI mars 2021
180INTEGER(KIND=JPIM), INTENT(IN) :: IDAY
181REAL(KIND=JPRB), INTENT(IN)    :: TIME
182
183! Name of file names specified on command line
184character(len=512), INTENT(IN) :: namelist_file
185logical, INTENT(IN)            :: ok_3Deffect
186
187! OUTPUT ARGUMENTS
188
189! *** Net fluxes on half-levels (W m-2)
190REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW(KLON,KLEV+1)
191REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW(KLON,KLEV+1)
192REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_CLEAR(KLON,KLEV+1)
193REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_CLEAR(KLON,KLEV+1)
194
195!*** DN and UP flux on half-levels (W m-2)
196REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN(KLON,KLEV+1)
197REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN(KLON,KLEV+1)
198REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN_CLEAR(KLON,KLEV+1)
199REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN_CLEAR(KLON,KLEV+1)
200REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP(KLON,KLEV+1)
201REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP(KLON,KLEV+1)
202REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP_CLEAR(KLON,KLEV+1)
203REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP_CLEAR(KLON,KLEV+1)
204
205! Direct component of surface flux into horizontal plane
206REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR(KLON)
207REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_CLEAR(KLON)
208! As PFLUX_DIR but into a plane perpendicular to the sun
209REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_INTO_SUN(KLON)
210
211! *** Ultraviolet and photosynthetically active radiation (W m-2)
212REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_UV(KLON)
213REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_PAR(KLON)
214REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_PAR_CLEAR(KLON)
215
216! Diagnosed longwave surface emissivity across the whole spectrum
217REAL(KIND=JPRB),  INTENT(OUT) :: PEMIS_OUT(KLON)   
218
219! Partial derivative of total-sky longwave upward flux at each level
220! with respect to upward flux at surface, used to correct heating
221! rates at gridpoints/timesteps between calls to the full radiation
222! scheme.  Note that this version uses the convention of level index
223! increasing downwards, unlike the local variable ZLwDerivative that
224! is returned from the LW radiation scheme.
225REAL(KIND=JPRB),  INTENT(OUT) :: PLWDERIVATIVE(KLON,KLEV+1)
226
227! Surface diffuse and direct downwelling shortwave flux in each
228! shortwave albedo band, used in RADINTG to update the surface fluxes
229! accounting for high-resolution albedo information
230REAL(KIND=JPRB),  INTENT(OUT) :: PSWDIFFUSEBAND(KLON,NSW)
231REAL(KIND=JPRB),  INTENT(OUT) :: PSWDIRECTBAND (KLON,NSW)
232
233!AI Nov 2023
234REAL(KIND=JPRB),  INTENT(OUT) :: ecrad_cloud_cover_sw(KLON)
235
236! LOCAL VARIABLES
237! AI ATTENTION
238type(config_type),save         :: rad_config
239!!$OMP THREADPRIVATE(rad_config)
240type(driver_config_type),save  :: driver_config
241!!$OMP THREADPRIVATE(driver_config)
242!type(config_type)        :: rad_config
243!type(driver_config_type)  :: driver_config
244TYPE(single_level_type)   :: single_level
245TYPE(thermodynamics_type) :: thermodynamics
246TYPE(gas_type)            :: gas
247TYPE(cloud_type)          :: cloud
248TYPE(aerosol_type)        :: aerosol
249TYPE(flux_type)           :: flux
250
251! Mass mixing ratio of ozone (kg/kg)
252REAL(KIND=JPRB)           :: ZO3(KLON,KLEV)
253
254! Cloud effective radii in microns
255REAL(KIND=JPRB)           :: ZRE_LIQUID_UM(KLON,KLEV)
256REAL(KIND=JPRB)           :: ZRE_ICE_UM(KLON,KLEV)
257
258! Cloud overlap decorrelation length for cloud boundaries in km
259REAL(KIND=JPRB)           :: ZDECORR_LEN_KM(KLON)
260
261! Ratio of cloud overlap decorrelation length for cloud water
262! inhomogeneities to that for cloud boundaries (typically 0.5)
263!REAL(KIND=JPRB)           :: ZDECORR_LEN_RATIO = 0.5_jprb
264
265! The surface net longwave flux if the surface was a black body, used
266! to compute the effective broadband surface emissivity
267REAL(KIND=JPRB)           :: ZBLACK_BODY_NET_LW(KIDIA:KFDIA)
268
269! Layer mass in kg m-2
270REAL(KIND=JPRB)           :: ZLAYER_MASS(KIDIA:KFDIA,KLEV)
271
272! Time integers
273INTEGER :: ITIM
274
275! Loop indices
276INTEGER :: JLON, JLEV, JBAND, JB_ALBEDO, JAER
277
278REAL(KIND=JPRB) :: ZHOOK_HANDLE
279
280! AI ATTENTION traitement aerosols
281INTEGER, PARAMETER :: NAERMACC = 1
282
283logical :: loutput=.true.
284logical :: lprint_input=.false.
285logical :: lprint_config=.false.
286logical, save :: debut_ecrad=.true.
287!$OMP THREADPRIVATE(debut_ecrad)
288integer, save :: itap_ecrad=0
289!$OMP THREADPRIVATE(itap_ecrad)
290
291REAL(KIND=JPRB) ::  inv_cloud_effective_size(KLON,KLEV)
292REAL(KIND=JPRB) ::  inv_inhom_effective_size(KLON,KLEV)
293
294integer :: irang
295
296
297IF (LHOOK) CALL DR_HOOK('RADIATION_SCHEME',0,ZHOOK_HANDLE)
298
299! A.I juillet 2023 :
300! Initialisation dans radiation_setup au 1er passage dans Ecrad
301!$OMP MASTER
302if (.not.ok_3Deffect) then
303  if (debut_ecrad) then
304   call SETUP_RADIATION_SCHEME(loutput,namelist_file,rad_config,driver_config)
305   debut_ecrad=.false.
306  endif
307else
308   call SETUP_RADIATION_SCHEME(loutput,namelist_file,rad_config,driver_config)
309endif
310!$OMP END MASTER
311!$OMP BARRIER
312! Fin partie initialisation et configuration
313
314!AI print fichiers namelist utilise
315!if (is_omp_root) then
316! itap_ecrad=itap_ecrad+1
317! print*,'Dans radiation_scheme itap_ecrad, mpi_rank, omp_rank, namelist_file : ', &
318!         itap_ecrad, mpi_rank, omp_rank, namelist_file
319!else
320! print*,'mpi_rank omp_rank, namelist_file :', mpi_rank, omp_rank, namelist_file
321!endif
322
323! AI 11 23 Allocates depplaces au debut
324print*,'*********** ALLOCATES *******************************'
325! AI ATTENTION
326! Allocate memory in radiation objects
327! emissivite avec une seule bande
328CALL single_level%allocate(KLON, NSW, 1, &
329     &                     use_sw_albedo_direct=.TRUE.)
330CALL thermodynamics%allocate(KLON, KLEV, use_h2o_sat=.true.)
331CALL cloud%allocate(KLON, KLEV)
332CALL aerosol%allocate(KLON, 1, KLEV, KAEROSOL)
333CALL gas%allocate(KLON, KLEV)
334CALL flux%allocate(rad_config, 1, KLON, KLEV)
335
336print*,'************* THERMO (input) ************************************'
337! Set thermodynamic profiles: simply copy over the half-level
338! pressure and temperature
339! AI
340! pressure_hl > paprs
341! temperature_hl calculee dans radlsw de la meme facon que pour RRTM
342thermodynamics%pressure_hl   (KIDIA:KFDIA,:) = PPRESSURE_H   (KIDIA:KFDIA,:)
343thermodynamics%temperature_hl(KIDIA:KFDIA,:) = PTEMPERATURE_H(KIDIA:KFDIA,:)
344!print*,'Compute saturation specific humidity'
345! Compute saturation specific humidity, used to hydrate aerosols. The
346! "2" for the last argument indicates that the routine is not being
347! called from within the convection scheme.
348!CALL SATUR(KIDIA, KFDIA, KLON, 1, KLEV, &
349!     &  PPRESSURE, PTEMPERATURE, thermodynamics%h2o_sat_liq, 2)
350! Alternative approximate version using temperature and pressure from
351! the thermodynamics structure
352!CALL thermodynamics%calc_saturation_wrt_liquid(KIDIA, KFDIA)
353!AI ATTENTION
354thermodynamics%h2o_sat_liq = PQSAT
355
356print*,'********** SINGLE LEVEL VARS **********************************'
357!AI ATTENTION
358! Set single-level fileds
359single_level%solar_irradiance              = PSOLAR_IRRADIANCE
360single_level%cos_sza(KIDIA:KFDIA)          = PMU0(KIDIA:KFDIA)
361single_level%skin_temperature(KIDIA:KFDIA) = PTEMPERATURE_SKIN(KIDIA:KFDIA)
362single_level%sw_albedo(KIDIA:KFDIA,:)      = PALBEDO_DIF(KIDIA:KFDIA,:)
363single_level%sw_albedo_direct(KIDIA:KFDIA,:)=PALBEDO_DIR(KIDIA:KFDIA,:)
364single_level%lw_emissivity(KIDIA:KFDIA,1)  = PEMIS(KIDIA:KFDIA)
365!single_level%lw_emissivity(KIDIA:KFDIA,2)  = PEMIS_WINDOW(KIDIA:KFDIA)
366
367! Create the relevant seed from date and time get the starting day
368! and number of minutes since start
369!IDAY = NDD(NINDAT)
370!cur_day
371!ITIM = NINT(NSTEP * YRRIP%TSTEP / 60.0_JPRB)
372!ITIM = NINT(TIME / 60.0_JPRB)
373!current_time
374!allocate(single_level%iseed(KIDIA:KFDIA))
375!DO JLON = KIDIA, KFDIA
376  ! This method gives a unique value for roughly every 1-km square
377  ! on the globe and every minute.  ASIN(PGEMU)*60 gives rough
378  ! latitude in degrees, which we multiply by 100 to give a unique
379  ! value for roughly every km. PGELAM*60*100 gives a unique number
380  ! for roughly every km of longitude around the equator, which we
381  ! multiply by 180*100 so there is no overlap with the latitude
382  ! values.  The result can be contained in a 32-byte integer (but
383  ! since random numbers are generated with the help of integer
384  ! overflow, it should not matter if the number did overflow).
385!  single_level%iseed(JLON) = ITIM + IDAY &
386!       &  +  NINT(PGELAM(JLON)*108000000.0_JPRB &
387!       &          + ASIN(PGEMU(JLON))*6000.0_JPRB)
388!ENDDO
389!AI Nov 23
390! Simple initialization of the seeds for the Monte Carlo scheme
391call single_level%init_seed_simple(kidia, kfdia)
392
393print*,'********** CLOUDS (allocate + input) *******************************************'
394!print*,'Appel Allocate clouds'
395! Set cloud fields
396cloud%q_liq(KIDIA:KFDIA,:)    = PQ_LIQUID(KIDIA:KFDIA,:)
397cloud%q_ice(KIDIA:KFDIA,:)    = PQ_ICE(KIDIA:KFDIA,:) + PQ_SNOW(KIDIA:KFDIA,:)
398cloud%fraction(KIDIA:KFDIA,:) = PCLOUD_FRAC(KIDIA:KFDIA,:)
399!!! ok AI ATTENTION a voir avec JL
400! Compute effective radi and convert to metres
401! AI. : on passe directement les champs de LMDZ
402cloud%re_liq(KIDIA:KFDIA,:) = ZRE_LIQUID_UM(KIDIA:KFDIA,:)
403cloud%re_ice(KIDIA:KFDIA,:) = ZRE_ICE_UM(KIDIA:KFDIA,:)
404! Get the cloud overlap decorrelation length (for cloud boundaries),
405! in km, according to the parameterization specified by NDECOLAT,
406! and insert into the "cloud" object. Also get the ratio of
407! decorrelation lengths for cloud water content inhomogeneities and
408! cloud boundaries, and set it in the "rad_config" object.
409! IFS :
410!CALL CLOUD_OVERLAP_DECORR_LEN(KIDIA, KFDIA, KLON, PGEMU, YRERAD%NDECOLAT, &
411!     &    ZDECORR_LEN_KM, PDECORR_LEN_RATIO=ZDECORR_LEN_RATIO)
412! AI valeur dans namelist
413! rad_config%cloud_inhom_decorr_scaling = ZDECORR_LEN_RATIO
414!AI ATTENTION meme valeur que dans offline
415! A mettre dans namelist
416ZDECORR_LEN_KM = driver_config%overlap_decorr_length
417DO JLON = KIDIA,KFDIA
418  CALL cloud%set_overlap_param(thermodynamics, &
419       &                 ZDECORR_LEN_KM(JLON), &
420       &                 istartcol=JLON, iendcol=JLON)
421ENDDO
422! IFS :
423! Cloud water content fractional standard deviation is configurable
424! from namelist NAERAD but must be globally constant. Before it was
425! hard coded at 1.0.
426!CALL cloud%create_fractional_std(KLON, KLEV, YRERAD%RCLOUD_FRAC_STD)
427! AI ATTENTION frac_std=0.75 meme valeur que dans la version offline
428CALL cloud%create_fractional_std(KLON, KLEV, driver_config%frac_std)
429
430if (ok_3Deffect) then
431 if (driver_config%ok_effective_size) then
432     call cloud%create_inv_cloud_effective_size_eta(klon, klev, &
433               &  thermodynamics%pressure_hl, &
434               &  driver_config%low_inv_effective_size, &
435               &  driver_config%middle_inv_effective_size, &
436               &  driver_config%high_inv_effective_size, 0.8_jprb, 0.45_jprb, &
437               &  KIDIA, KFDIA)
438   else if (driver_config%ok_separation) then
439     call cloud%param_cloud_effective_separation_eta(klon, klev, &
440               &  thermodynamics%pressure_hl, &
441               &  driver_config%cloud_separation_scale_surface, &
442               &  driver_config%cloud_separation_scale_toa, &
443               &  driver_config%cloud_separation_scale_power, &
444               &  driver_config%cloud_inhom_separation_factor, &
445               &  KIDIA, KFDIA)
446   endif     
447 else 
448 if (rad_config%i_solver_sw == ISolverSPARTACUS &
449      & .or.   rad_config%i_solver_lw == ISolverSPARTACUS) then
450   ! AI ! Read cloud properties needed by SPARTACUS
451   if (driver_config%ok_effective_size) then
452     call cloud%create_inv_cloud_effective_size_eta(klon, klev, &
453               &  thermodynamics%pressure_hl, &
454               &  driver_config%low_inv_effective_size, &
455               &  driver_config%middle_inv_effective_size, &
456               &  driver_config%high_inv_effective_size, 0.8_jprb, 0.45_jprb, &
457               &  KIDIA, KFDIA)
458   else if (driver_config%ok_separation) then
459     call cloud%param_cloud_effective_separation_eta(klon, klev, &
460               &  thermodynamics%pressure_hl, &
461               &  driver_config%cloud_separation_scale_surface, &
462               &  driver_config%cloud_separation_scale_toa, &
463               &  driver_config%cloud_separation_scale_power, &
464               &  driver_config%cloud_inhom_separation_factor, &
465               &  KIDIA, KFDIA)
466   endif
467 endif 
468endif
469
470print*,'******** AEROSOLS (input) **************************************'
471!IF (NAERMACC > 0) THEN
472!ELSE
473!  CALL aerosol%allocate(KLON, 1, KLEV, 6) ! Tegen climatology
474!ENDIF
475! Compute the dry mass of each layer neglecting humidity effects, in
476! kg m-2, needed to scale some of the aerosol inputs
477! AI commente ATTENTION
478!CALL thermodynamics%get_layer_mass(ZLAYER_MASS)
479
480! Copy over aerosol mass mixing ratio
481!IF (NAERMACC > 0) THEN
482
483  ! MACC aerosol climatology - this is already in mass mixing ratio
484  ! units with the required array orientation so we can copy it over
485  ! directly
486  aerosol%mixing_ratio(KIDIA:KFDIA,:,:) = PAEROSOL(KIDIA:KFDIA,:,:)
487
488  ! Add the tropospheric and stratospheric backgrounds contained in the
489  ! old Tegen arrays - this is very ugly!
490! AI ATTENTION
491!  IF (TROP_BG_AER_MASS_EXT > 0.0_JPRB) THEN
492!    aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_TROP_BG_AER) &
493!         &  = aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_TROP_BG_AER) &
494!         &  + PAEROSOL_OLD(KIDIA:KFDIA,1,:) &
495!         &  / (ZLAYER_MASS * TROP_BG_AER_MASS_EXT)
496!  ENDIF
497!  IF (STRAT_BG_AER_MASS_EXT > 0.0_JPRB) THEN
498!    aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_STRAT_BG_AER) &
499!         &  = aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_STRAT_BG_AER) &
500!         &  + PAEROSOL_OLD(KIDIA:KFDIA,6,:) &
501!         &  / (ZLAYER_MASS * STRAT_BG_AER_MASS_EXT)
502!  ENDIF
503
504!ELSE
505
506  ! Tegen aerosol climatology - the array PAEROSOL_OLD contains the
507  ! 550-nm optical depth in each layer. The optics data file
508  ! aerosol_ifs_rrtm_tegen.nc does not contain mass extinction
509  ! coefficient, but a scaling factor that the 550-nm optical depth
510  ! should be multiplied by to obtain the optical depth in each
511  ! spectral band.  Therefore, in order for the units to work out, we
512  ! need to divide by the layer mass (in kg m-2) to obtain the 550-nm
513  ! cross-section per unit mass of dry air (so in m2 kg-1).  We also
514  ! need to permute the array.
515!  DO JLEV = 1,KLEV
516!    DO JAER = 1,6
517!      aerosol%mixing_ratio(KIDIA:KFDIA,JLEV,JAER) &
518!         &  = PAEROSOL_OLD(KIDIA:KFDIA,JAER,JLEV) &
519!         &  / ZLAYER_MASS(KIDIA:KFDIA,JLEV)
520!    ENDDO
521!  ENDDO
522!ENDIF
523
524print*,'********** GAS (input) ************************************************'
525!print*,'Appel Allocate gas'
526! Convert ozone Pa*kg/kg to kg/kg
527! AI ATTENTION
528!DO JLEV = 1,KLEV
529!  DO JLON = KIDIA,KFDIA
530!    ZO3(JLON,JLEV) = PO3_DP(JLON,JLEV) &
531!         &         / (PPRESSURE_H(JLON,JLEV+1)-PPRESSURE_H(JLON,JLEV))
532!  ENDDO
533!ENDDO
534!  Insert gas mixing ratios
535!print*,'Insert gas mixing ratios'
536CALL gas%put(IH2O,    IMassMixingRatio, PQ)
537CALL gas%put(IO3,     IMassMixingRatio, PO3)
538CALL gas%put_well_mixed(ICO2,    IMAssMixingRatio, PCO2)
539CALL gas%put_well_mixed(ICH4,    IMassMixingRatio, PCH4)
540CALL gas%put_well_mixed(IN2O,    IMassMixingRatio, PN2O)
541CALL gas%put_well_mixed(ICFC11,  IMassMixingRatio, PCFC11)
542CALL gas%put_well_mixed(ICFC12,  IMassMixingRatio, PCFC12)
543CALL gas%put_well_mixed(IHCFC22, IMassMixingRatio, PHCFC22)
544CALL gas%put_well_mixed(ICCL4,   IMassMixingRatio, PCCL4)
545CALL gas%put_well_mixed(IO2,     IMassMixingRatio, PO2)
546! Ensure the units of the gas mixing ratios are what is required by
547! the gas absorption model
548call set_gas_units(rad_config, gas)
549
550! Call radiation scheme
551!print*,'*** Appel radiation *** namelist **** omp_rank ****', &
552!         omp_rank, namelist_file
553!   if (rad_config%i_solver_sw == ISolverSPARTACUS) then
554!    if (driver_config%ok_separation) then
555!      print*,'Avant radiation, mpi_rank, omp_rank, size, chape inv_cloud = ',&
556!              mpi_rank, omp_rank, &
557!              shape(cloud%inv_cloud_effective_size), &
558!              size(cloud%inv_cloud_effective_size)     
559!      do jlon=KIDIA, KFDIA     
560!        do jlev=1,klev
561!         print*,' Avant radiation mpi_rank, omp_rank, jlon, jlev, &
562!            &     cloud%inv_cloud_effective_size =', mpi_rank, &
563!            &     omp_rank, jlon, jlev, &
564!            &     cloud%inv_cloud_effective_size(jlon,jlev)
565!        enddo
566!      enddo
567!       cloud%inv_cloud_effective_size=inv_cloud_effective_size
568!       cloud%inv_inhom_effective_size=inv_inhom_effective_size
569!    endif
570!   endif
571CALL radiation(KLON, KLEV, KIDIA, KFDIA, rad_config, &
572     &  single_level, thermodynamics, gas, cloud, aerosol, flux)
573
574 if (rad_config%use_aerosols) then
575    if (rad_config%i_gas_model == IGasModelIFSRRTMG) then     
576       CALL aeropt_5wv_ecrad(kidia, kfdia, 1, klev, &
577                             rad_config,thermodynamics,aerosol)
578    endif                 
579 endif               
580
581print*,'*********** Sortie flux ****************'
582! Cloud cover
583ecrad_cloud_cover_sw = flux%cloud_cover_sw
584! Compute required output fluxes
585! DN and UP flux
586PFLUX_SW_DN(KIDIA:KFDIA,:) = flux%sw_dn(KIDIA:KFDIA,:)
587PFLUX_SW_UP(KIDIA:KFDIA,:) = flux%sw_up(KIDIA:KFDIA,:)
588PFLUX_LW_DN(KIDIA:KFDIA,:) = flux%lw_dn(KIDIA:KFDIA,:)
589PFLUX_LW_UP(KIDIA:KFDIA,:) = flux%lw_up(KIDIA:KFDIA,:)
590PFLUX_SW_DN_CLEAR(KIDIA:KFDIA,:) = flux%sw_dn_clear(KIDIA:KFDIA,:)
591PFLUX_SW_UP_CLEAR(KIDIA:KFDIA,:) = flux%sw_up_clear(KIDIA:KFDIA,:)
592PFLUX_LW_DN_CLEAR(KIDIA:KFDIA,:) = flux%lw_dn_clear(KIDIA:KFDIA,:)
593PFLUX_LW_UP_CLEAR(KIDIA:KFDIA,:) = flux%lw_up_clear(KIDIA:KFDIA,:)
594! First the net fluxes
595PFLUX_SW(KIDIA:KFDIA,:) = flux%sw_dn(KIDIA:KFDIA,:) - flux%sw_up(KIDIA:KFDIA,:)
596PFLUX_LW(KIDIA:KFDIA,:) = flux%lw_dn(KIDIA:KFDIA,:) - flux%lw_up(KIDIA:KFDIA,:)
597PFLUX_SW_CLEAR(KIDIA:KFDIA,:) &
598     &  = flux%sw_dn_clear(KIDIA:KFDIA,:) - flux%sw_up_clear(KIDIA:KFDIA,:)
599PFLUX_LW_CLEAR(KIDIA:KFDIA,:) &
600     &  = flux%lw_dn_clear(KIDIA:KFDIA,:) - flux%lw_up_clear(KIDIA:KFDIA,:)
601! Now the surface fluxes
602!PFLUX_SW_DN_SURF(KIDIA:KFDIA) = flux%sw_dn(KIDIA:KFDIA,KLEV+1)
603!PFLUX_LW_DN_SURF(KIDIA:KFDIA) = flux%lw_dn(KIDIA:KFDIA,KLEV+1)
604!PFLUX_SW_UP_SURF(KIDIA:KFDIA) = flux%sw_up(KIDIA:KFDIA,KLEV+1)
605!PFLUX_LW_UP_SURF(KIDIA:KFDIA) = flux%lw_up(KIDIA:KFDIA,KLEV+1)
606!PFLUX_SW_DN_CLEAR_SURF(KIDIA:KFDIA) = flux%sw_dn_clear(KIDIA:KFDIA,KLEV+1)
607!PFLUX_LW_DN_CLEAR_SURF(KIDIA:KFDIA) = flux%lw_dn_clear(KIDIA:KFDIA,KLEV+1)
608!PFLUX_SW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%sw_up_clear(KIDIA:KFDIA,KLEV+1)
609!PFLUX_LW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%lw_up_clear(KIDIA:KFDIA,KLEV+1)
610PFLUX_DIR(KIDIA:KFDIA) = flux%sw_dn_direct(KIDIA:KFDIA,KLEV+1)
611PFLUX_DIR_CLEAR(KIDIA:KFDIA) = flux%sw_dn_direct_clear(KIDIA:KFDIA,KLEV+1)
612PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = 0.0_JPRB
613WHERE (PMU0(KIDIA:KFDIA) > EPSILON(1.0_JPRB))
614  PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = PFLUX_DIR(KIDIA:KFDIA) / PMU0(KIDIA:KFDIA)
615END WHERE
616! Top-of-atmosphere downwelling flux
617!PFLUX_SW_DN_TOA(KIDIA:KFDIA) = flux%sw_dn(KIDIA:KFDIA,1)
618!PFLUX_SW_UP_TOA(KIDIA:KFDIA) = flux%sw_up(KIDIA:KFDIA,1)
619!PFLUX_LW_DN_TOA(KIDIA:KFDIA) = flux%lw_dn(KIDIA:KFDIA,1)
620!PFLUX_LW_UP_TOA(KIDIA:KFDIA) = flux%lw_up(KIDIA:KFDIA,1)
621!AI ATTENTION
622if (0.eq.1) then
623PFLUX_UV       (KIDIA:KFDIA) = 0.0_JPRB
624DO JBAND = 1,NWEIGHT_UV
625  PFLUX_UV(KIDIA:KFDIA) = PFLUX_UV(KIDIA:KFDIA) + WEIGHT_UV(JBAND) &
626       &  * flux%sw_dn_surf_band(IBAND_UV(JBAND),KIDIA:KFDIA)
627ENDDO
628! Compute photosynthetically active radiation similarly
629PFLUX_PAR      (KIDIA:KFDIA) = 0.0_JPRB
630PFLUX_PAR_CLEAR(KIDIA:KFDIA) = 0.0_JPRB
631DO JBAND = 1,NWEIGHT_PAR
632  PFLUX_PAR(KIDIA:KFDIA) = PFLUX_PAR(KIDIA:KFDIA) + WEIGHT_PAR(JBAND) &
633       &  * flux%sw_dn_surf_band(IBAND_PAR(JBAND),KIDIA:KFDIA)
634  PFLUX_PAR_CLEAR(KIDIA:KFDIA) = PFLUX_PAR_CLEAR(KIDIA:KFDIA) &
635       &  + WEIGHT_PAR(JBAND) &
636       &  * flux%sw_dn_surf_clear_band(IBAND_PAR(JBAND),KIDIA:KFDIA)
637ENDDO
638endif
639! Compute effective broadband emissivity
640ZBLACK_BODY_NET_LW = flux%lw_dn(KIDIA:KFDIA,KLEV+1) &
641     &  - RSIGMA*PTEMPERATURE_SKIN(KIDIA:KFDIA)**4
642PEMIS_OUT(KIDIA:KFDIA) = PEMIS(KIDIA:KFDIA)
643WHERE (ABS(ZBLACK_BODY_NET_LW) > 1.0E-5)
644  PEMIS_OUT(KIDIA:KFDIA) = PFLUX_LW(KIDIA:KFDIA,KLEV+1) / ZBLACK_BODY_NET_LW
645END WHERE
646! Copy longwave derivatives
647! AI ATTENTION
648!IF (YRERAD%LAPPROXLWUPDATE) THEN
649IF (rad_config%do_lw_derivatives) THEN
650  PLWDERIVATIVE(KIDIA:KFDIA,:) = flux%lw_derivatives(KIDIA:KFDIA,:)
651END IF
652! Store the shortwave downwelling fluxes in each albedo band
653!AI ATTENTION
654!IF (YRERAD%LAPPROXSWUPDATE) THEN
655if (0.eq.1) then
656IF (rad_config%do_surface_sw_spectral_flux) THEN
657  PSWDIFFUSEBAND(KIDIA:KFDIA,:) = 0.0_JPRB
658  PSWDIRECTBAND (KIDIA:KFDIA,:) = 0.0_JPRB
659  DO JBAND = 1,rad_config%n_bands_sw
660    JB_ALBEDO = rad_config%i_albedo_from_band_sw(JBAND)
661    DO JLON = KIDIA,KFDIA
662      PSWDIFFUSEBAND(JLON,JB_ALBEDO) = PSWDIFFUSEBAND(JLON,JB_ALBEDO) &
663           &  + flux%sw_dn_surf_band(JBAND,JLON) &
664           &  - flux%sw_dn_direct_surf_band(JBAND,JLON)
665      PSWDIRECTBAND(JLON,JB_ALBEDO)  = PSWDIRECTBAND(JLON,JB_ALBEDO) &
666           &  + flux%sw_dn_direct_surf_band(JBAND,JLON)
667    ENDDO
668  ENDDO
669ENDIF
670endif
671
672print*,'********** DEALLOCATIONS ************************'
673CALL single_level%deallocate
674CALL thermodynamics%deallocate
675CALL gas%deallocate
676CALL cloud%deallocate
677CALL aerosol%deallocate
678CALL flux%deallocate
679
680IF (LHOOK) CALL DR_HOOK('RADIATION_SCHEME',1,ZHOOK_HANDLE)
681
682END SUBROUTINE RADIATION_SCHEME
683
684end module interface_lmdz_ecrad
Note: See TracBrowser for help on using the repository browser.