Newer
Older
!MNH_LIC Copyright 2013-2019 CNRS, Meteo-France and Universite Paul Sabatier
!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
!MNH_LIC for details. version 1.
!-----------------------------------------------------------------
! ###############################
MODULE MODE_INI_LIMA_COLD_MIXED
! ###############################
!
!
! ###############################################
SUBROUTINE INI_LIMA_COLD_MIXED (PTSTEP, PDZMIN)
! ###############################################
!
!! PURPOSE
!! -------
!! The purpose of this routine is to initialize the constants used in the
!! microphysical scheme LIMA for the cold and mixed phase variables
!! and processes.
!!
!! AUTHOR
!! ------
!! J.-M. Cohard * Laboratoire d'Aerologie*
!! J.-P. Pinty * Laboratoire d'Aerologie*
!! S. Berthet * Laboratoire d'Aerologie*
!! B. Vié * Laboratoire d'Aerologie*
!!
!! MODIFICATIONS
!! -------------
!! Original ??/??/13
!! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
! C. Barthe 14/03/2022: add CIBU and RDSF
! J. Wurtz 03/2022: new snow characteristics
! M. Taufour 07/2022: add concentration for snow, graupel, hail
!
!-------------------------------------------------------------------------------
!
!* 0. DECLARATIONS
! ------------
!
USE MODD_CST
USE MODD_PARAMETERS
USE MODD_PARAM_LIMA
USE MODD_PARAM_LIMA_WARM
USE MODD_PARAM_LIMA_COLD
USE MODD_PARAM_LIMA_MIXED
!
use mode_msg
!
USE MODE_LIMA_FUNCTIONS, ONLY: MOMG, GAUHER
USE MODI_GAMMA
USE MODI_GAMMA_INC
USE MODE_RRCOLSS, ONLY: RRCOLSS
USE MODE_RZCOLX, ONLY: RZCOLX
USE MODE_RSCOLRG, ONLY: RSCOLRG
USE MODE_NRCOLSS, ONLY: NRCOLSS
USE MODE_NZCOLX, ONLY: NZCOLX
USE MODE_NSCOLRG, ONLY: NSCOLRG
USE MODE_LIMA_READ_XKER_RACCS, ONLY: LIMA_READ_XKER_RACCS
USE MODE_LIMA_READ_XKER_SDRYG, ONLY: LIMA_READ_XKER_SDRYG
USE MODE_LIMA_READ_XKER_RDRYG, ONLY: LIMA_READ_XKER_RDRYG
USE MODE_LIMA_READ_XKER_SWETH, ONLY: LIMA_READ_XKER_SWETH
USE MODE_LIMA_READ_XKER_GWETH, ONLY: LIMA_READ_XKER_GWETH
!
IMPLICIT NONE
!
!* 0.1 Declarations of dummy arguments :
!
REAL, INTENT(IN) :: PTSTEP ! Effective Time step
REAL, INTENT(IN) :: PDZMIN ! minimun vertical mesh size
!
!* 0.2 Declarations of local variables :
!
character(len=13) :: yval ! String for error message
INTEGER :: IKB ! Coordinates of the first physical
! points along z
INTEGER :: J1 ! Internal loop indexes
!
REAL, DIMENSION(8) :: ZGAMI ! parameters involving various moments
REAL, DIMENSION(2) :: ZGAMS ! of the generalized gamma law
!
REAL :: ZRHO00 ! Surface reference air density
REAL :: ZRATE ! Geometrical growth of Lbda in the tabulated
! functions and kernels
REAL :: ZBOUND ! XDCSLIM*Lbda_s: upper bound for the partial
! integration of the riming rate of the aggregates
REAL :: ZEGS, ZEGR, ZEHS, ZEHG! Bulk collection efficiencies
!
INTEGER :: IND ! Number of interval to integrate the kernels
REAL :: ZESR, ZESS ! Mean efficiency of rain-aggregate collection, aggregate-aggregate collection
REAL :: ZFDINFTY ! Factor used to define the "infinite" diameter
!
!
!INTEGER :: ILUOUT0 ! Logical unit number for output-listing
!LOGICAL :: GFLAG ! Logical flag for printing the constatnts on the output
! ! listing
REAL :: ZCONC_MAX ! Maximal concentration for snow
REAL :: ZFACT_NUCL! Amplification factor for the minimal ice concentration
!
INTEGER :: KND
INTEGER :: KACCLBDAS,KACCLBDAR,KDRYLBDAG,KDRYLBDAS,KDRYLBDAR
REAL :: PALPHAR,PALPHAS,PALPHAG,PALPHAH
REAL :: PNUR,PNUS,PNUG,PNUH
REAL :: PBR,PBS,PBG
REAL :: PCR,PCS,PCG,PCH
REAL :: PDR,PDS,PFVELOS,PDG,PDH
REAL :: PESR,PEGS,PEGR,PEHS,PEHG
REAL :: PFDINFTY
REAL :: PACCLBDAS_MAX,PACCLBDAR_MAX,PACCLBDAS_MIN,PACCLBDAR_MIN
REAL :: PDRYLBDAG_MAX,PDRYLBDAS_MAX,PDRYLBDAG_MIN,PDRYLBDAS_MIN
REAL :: PDRYLBDAR_MAX,PDRYLBDAR_MIN
REAL :: PWETLBDAS_MAX,PWETLBDAG_MAX,PWETLBDAS_MIN,PWETLBDAG_MIN
REAL :: PWETLBDAH_MAX,PWETLBDAH_MIN
INTEGER :: KWETLBDAS,KWETLBDAG,KWETLBDAH
!
REAL :: ZFAC_ZRNIC ! Zrnic factor used to decrease Long Kernels
REAL :: ZBOUND_CIBU_SMIN ! XDCSLIM*Lbda_s : lower bound used in the tabulated function
REAL :: ZBOUND_CIBU_GMIN ! XDCGLIM*Lbda_g : lower bound used in the tabulated function
REAL :: ZRATE_S ! Geometrical growth of Lbda_s in the tabulated function
REAL :: ZRATE_G ! Geometrical growth of Lbda_g in the tabulated function
!
REAL :: ZBOUND_RDSF_RMIN ! XDCRLIM*Lbda_r : lower bound used in the tabulated function
REAL :: ZRATE_R ! Geometrical growth of Lbda_r in the tabulated function
REAL :: ZKHI_LWM ! Coefficient of Lawson et al. (2015)
!
Marie TAUFOUR
committed
REAL :: ZRHOIW ! ice density
!
!-------------------------------------------------------------------------------
!
!
CALL PARAM_LIMA_COLD_ASSOCIATE()
CALL PARAM_LIMA_MIXED_ASSOCIATE()
!
!
!* 1. CHARACTERISTICS OF THE SPECIES
! ------------------------------
!
!CALL RAIN_ICE_PARAM_ASSOCIATE()
!
!* 1.2 Ice crystal characteristics
!
SELECT CASE (CPRISTINE_ICE_LIMA)
CASE('PLAT')
XAI = 0.82 ! Plates
XBI = 2.5 ! Plates
XC_I = 747. ! Plates
XDI = 1.0 ! Plates
Marie TAUFOUR
committed
XGAMMAI = 0.096
XDELTAI = 1.83
XC1I = 1./XPI ! Plates
CASE('COLU')
XAI = 2.14E-3 ! Columns
XBI = 1.7 ! Columns
XC_I = 1.96E5 ! Columns
XDI = 1.585 ! Columns
Marie TAUFOUR
committed
XGAMMAI = 0.659
XDELTAI = 2.0
XC1I = 0.8 ! Columns
CASE('BURO')
XAI = 44.0 ! Bullet rosettes
XBI = 3.0 ! Bullet rosettes
XC_I = 4.E5 ! Bullet rosettes
XDI = 1.663 ! Bullet rosettes
Marie TAUFOUR
committed
XGAMMAI = 0.062
XDELTAI = 1.81
XC1I = 0.5 ! Bullet rosettes
Marie TAUFOUR
committed
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
CASE('YPLA')
XAI = 0.745 ! Plates_from Yang et al (2013)
XBI = 2.47 ! Plates_from Yang et al (2013)
XC_I = 63. ! Plates_from Yang et al (2013)
XDI = 0.68 ! Plates_from Yang et al (2013)
XGAMMAI = 0.096
XDELTAI = 1.83
XC1I = 1./XPI ! Plates_from Yang et al (2013)
CASE('YCOL')
XAI = 261.102 ! Columns_from Yang et al (2013)
XBI = 2.99 ! Columns_from Yang et al (2013)
XC_I = 671 ! Columns_from Yang et al (2013)
XDI = 0.62 ! Columns_from Yang et al (2013)
XGAMMAI = 0.659
XDELTAI = 2.0
XC1I = 0.8 ! Columns_from Yang et al (2013)
CASE('YBUR')
XAI = 1.268 ! Bullet rosettes_from Yang et al (2013)
XBI = 2.59 ! Bullet rosettes_from Yang et al (2013)
XC_I = 128 ! Bullet rosettes_from Yang et al (2013)
XDI = 0.72 ! Bullet rosettes_from Yang et al (2013)
XGAMMAI = 0.062
XDELTAI = 1.81
XC1I = 0.5 ! Bullet rosettes_from Yang et al (2013)
CASE('YDRO')
XAI = 1.268 ! Droxtals_from Yang et al (2013)
XBI = 2.59 ! Droxtals_from Yang et al (2013)
XC_I = 128 ! Droxtals_from Yang et al (2013)
XDI = 0.72 ! Droxtals_from Yang et al (2013)
XGAMMAI = 0.673
XDELTAI = 2.0
XC1I = 0.5 ! Droxtals_from Yang et al (2013)
CASE('YHCO')
XAI = 217.586 ! Hollow_Columns_from Yang et al (2013)
XBI = 2.99 ! Hollow_Columns_from Yang et al (2013)
XC_I = 641 ! Hollow_Columns_from Yang et al (2013)
XDI = 0.63 ! Hollow_Columns_from Yang et al (2013)
XGAMMAI = 0.659
XDELTAI = 2.0
XC1I = 0.8 ! Hollow_Columns_from Yang et al (2013)
CASE('YHBU')
XAI = 1.258 ! Hollow_Bullet rosettes_from Yang et al (2013)
XBI = 2.61 ! Hollow_Bullet rosettes_from Yang et al (2013)
XC_I = 147 ! Hollow_Bullet rosettes_from Yang et al (2013)
XDI = 0.73 ! Hollow_Bullet rosettes_from Yang et al (2013)
XGAMMAI = 0.061
XDELTAI = 1.81
XC1I = 0.5 ! Hollow_Bullet rosettes_from Yang et al (2013)
END SELECT
!
! Note that XCCI=N_i (a locally predicted value) and XCXI=0.0, implicitly
!
XF0I = 1.00
! Correction BVIE XF2I from Pruppacher 1997 eq 13-88
!XF2I = 0.103
XF2I = 0.14
XF0IS = 0.86
XF1IS = 0.28
!
!* 1.3 Snowflakes/aggregates characteristics
!
XAS = 0.02
XBS = 1.9
IF (LSNOW_T) THEN
!Cas Gamma generalisee
XCS = 11.52
XDS = 0.39
XFVELOS =0.097
!Cas MP
!XCS = 13.2
!XDS = 0.423
!XFVELOS = 25.14
ELSE
XCS = 5.
XDS = 0.27
XFVELOS = 0.
END IF
XCCS = 5.0
XCXS = 1.0
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
XF0S = 0.86
XF1S = 0.28
!
XC1S = 1./XPI
!
!* 1.4 Graupel characteristics
!
XAG = 19.6 ! Lump graupel case
XBG = 2.8 ! Lump graupel case
XCG = 122. ! Lump graupel case
XDG = 0.66 ! Lump graupel case
!
XCCG = 5.E5
XCXG = -0.5
! XCCG = 4.E4 ! Test of Ziegler (1988)
! XCXG = -1.0 ! Test of Ziegler (1988)
!
XF0G = 0.86
XF1G = 0.28
!
XC1G = 1./2.
!
!* 2.5 Hailstone characteristics
!
!
XAH = 470.
XBH = 3.0
XCH = 201.
XDH = 0.64
!
!XCCH = 5.E-4
!XCXH = 2.0
!!!!!!!!!!!!
XCCH = 4.E4 ! Test of Ziegler (1988)
XCXH = -1.0 ! Test of Ziegler (1988)
!!! XCCH = 5.E5 ! Graupel_like
!!! XCXH = -0.5 ! Graupel_like
!!!!!!!!!!!!
!
XF0H = 0.86
XF1H = 0.28
!
XC1H = 1./2.
!
!-------------------------------------------------------------------------------
!
!
!* 2. DIMENSIONAL DISTRIBUTIONS OF THE SPECIES
! ----------------------------------------
!
!
!* 2.1 Ice, snow, graupel and hail distribution
!
!
XALPHAI = 3.0 ! Gamma law for the ice crystal volume
XNUI = 3.0 ! Gamma law with little dispersion
!
IF (LSNOW_T) THEN
!Cas GAMMAGEN
XALPHAS = .214 ! Generalized gamma law
XNUS = 43.7 ! Generalized gamma law
XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / &
( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) )
ELSE IF (NMOM_S.EQ.2) THEN
XALPHAS = 1.0 ! Gamma law
XNUS = 2.0 !
XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / &
( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) )
ELSE
XALPHAS = 1.0 ! Exponential law
XNUS = 1.0 ! Exponential law
XTRANS_MP_GAMMAS = 1.
END IF
if (NMOM_G.GE.2) then
XALPHAG = 1.0 !
XNUG = 2.0 !
else
XALPHAG = 1.0 ! Exponential law
XNUG = 1.0 ! Exponential law
end if
!
if (NMOM_H.GE.2) then
XALPHAH = 1.0 ! Gamma law
XNUH = 5.0 ! Gamma law with little dispersion
else
XALPHAH = 1.0 ! Gamma law
XNUH = 8.0 ! Gamma law with little dispersion
end if
!
!* 2.2 Constants for shape parameter
!
XLBEXI = 1.0/XBI
XLBI = XAI*MOMG(XALPHAI,XNUI,XBI)
!
XNS = 1.0/(XAS*MOMG(XALPHAS,XNUS,XBS))
IF (NMOM_S.EQ.1) THEN
XLBEXS = 1.0/(XCXS-XBS)
XLBS = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS)
ELSE
XLBEXS = 1./XBS
XLBS = XAS * MOMG(XALPHAS,XNUS,XBS)
END IF
XNG = 1.0/(XAG*MOMG(XALPHAG,XNUG,XBG))
IF (NMOM_G.EQ.1) THEN
XLBEXG = 1.0/(XCXG-XBG)
XLBG = ( XAG*XCCG*MOMG(XALPHAG,XNUG,XBG))**(-XLBEXG)
ELSE
XLBEXG = 1./XBG
XLBG = XAG * MOMG(XALPHAG,XNUG,XBG)
END IF
IF (NMOM_H.EQ.1) THEN
XLBEXH = 1.0/(XCXH-XBH)
XLBH = ( XAH*XCCH*MOMG(XALPHAH,XNUH,XBH) )**(-XLBEXH)
ELSE
XLBEXH = 1./XBH
XLBH = XAH * MOMG(XALPHAH,XNUH,XBH)
END IF
!!$GFLAG = .TRUE.
!!$IF (GFLAG) THEN
!!$ WRITE(UNIT=ILUOUT0,FMT='(" Shape Parameters")')
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XLBEXI =",E13.6," XLBI =",E13.6)') XLBEXI,XLBI
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XLBEXS =",E13.6," XLBS =",E13.6)') XLBEXS,XLBS
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XLBEXG =",E13.6," XLBG =",E13.6)') XLBEXG,XLBG
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XLBEXH =",E13.6," XLBH =",E13.6)') XLBEXH,XLBH
!!$END IF
XLBDAS_MAX = 1.E7 ! (eq to r~1E-7kg/kg) (for non MP PSD, use conversion XTRANS_MP_GAMMAS)
XLBDAS_MIN = 1. ! (eq to r~0.18kg/kg) (for non MP PSD, use conversion XTRANS_MP_GAMMAS)
XLBDAG_MAX = 100000.0
!
ZCONC_MAX = 1.E6 ! Maximal concentration for falling particules set to 1 per cc
!XLBDAS_MAX = ( ZCONC_MAX/XCCS )**(1./XCXS)
!XLBDAG_MAX = ( ZCONC_MAX/XCCG )**(1./XCXG)
!XLBDAH_MAX = ( ZCONC_MAX/XCCH )**(1./XCXH)
!
Marie TAUFOUR
committed
! constante for ecRad effective radius
ZRHOIW = 0.917
XREFFI = (3*XAI/(2*ZRHOIW*10**3*XGAMMAI)*MOMG(XALPHAI,XNUI,XBI)/MOMG(XALPHAI,XNUI,XDELTAI))*1E6
!
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
!-------------------------------------------------------------------------------
!
!
!* 3. CONSTANTS FOR THE SEDIMENTATION
! -------------------------------
!
!
!* 3.1 Exponent of the fall-speed air density correction
!
IKB = 1 + JPVEXT
! Correction
! ZRHO00 = XP00/(XRD*XTHVREFZ(IKB))
ZRHO00 = 1.2041 ! at P=1013.25hPa and T=20°C
!
!* 3.2 Constants for sedimentation
!
!! XEXRSEDI = (XBI+XDI)/XBI
!! XEXCSEDI = 1.0-XEXRSEDI
!! XFSEDI = (4.*XPI*900.)**(-XEXCSEDI) * &
!! XC_I*XAI*MOMG(XALPHAI,XNUI,XBI+XDI) * &
!! ((XAI*MOMG(XALPHAI,XNUI,XBI)))**(-XEXRSEDI) * &
!! (ZRHO00)**XCEXVT
!! !
!! ! Computations made for Columns
!! !
!! XEXRSEDI = 1.9324
!! XEXCSEDI =-0.9324
!! XFSEDI = 3.89745E11*MOMG(XALPHAI,XNUI,3.285)* &
!! MOMG(XALPHAI,XNUI,1.7)**(-XEXRSEDI)*(ZRHO00)**XCEXVT
!! XEXCSEDI =-0.9324*3.0
!! WRITE (ILUOUT0,FMT=*)' PRISTINE ICE SEDIMENTATION for columns XFSEDI=',XFSEDI
!
!
XFSEDRI = XC_I*GAMMA_X0D(XNUI+(XDI+XBI)/XALPHAI)/GAMMA_X0D(XNUI+XBI/XALPHAI)* &
(ZRHO00)**XCEXVT
XFSEDCI = XC_I*GAMMA_X0D(XNUI+XDI/XALPHAI)/GAMMA_X0D(XNUI)* &
(ZRHO00)**XCEXVT
!
IF (LSNOW_T) THEN
!HOUZE/HAIC
!XEXSEDS = -XDS !(2*XBS+XDS)
!XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) &
! *(ZRHO00)**XCEXVT
!LH_EXTENDED
XEXSEDS = -XDS-XBS
XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) &
*(ZRHO00)**XCEXVT
ELSE
XEXSEDS = (XBS+XDS-XCXS)/(XBS-XCXS)
XFSEDS = XCS*XAS*XCCS*MOMG(XALPHAS,XNUS,XBS+XDS)* &
(XAS*XCCS*MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS)*(ZRHO00)**XCEXVT
!
XEXSEDG = (XBG+XDG-XCXG)/(XBG-XCXG)
XFSEDG = XCG*XAG*XCCG*MOMG(XALPHAG,XNUG,XBG+XDG)* &
(XAG*XCCG*MOMG(XALPHAG,XNUG,XBG))**(-XEXSEDG)*(ZRHO00)**XCEXVT
!
XEXSEDH = (XBH+XDH-XCXH)/(XBH-XCXH)
XFSEDH = XCH*XAH*XCCH*MOMG(XALPHAH,XNUH,XBH+XDH)* &
(XAH*XCCH*MOMG(XALPHAH,XNUH,XBH))**(-XEXSEDH)*(ZRHO00)**XCEXVT
IF (NMOM_S.GE.2) THEN
XFSEDRS = XCS*GAMMA_X0D(XNUS+(XDS+XBS)/XALPHAS)/GAMMA_X0D(XNUS+XBS/XALPHAS)* &
(ZRHO00)**XCEXVT
XFSEDCS = XCS*GAMMA_X0D(XNUS+XDS/XALPHAS)/GAMMA_X0D(XNUS)* &
(ZRHO00)**XCEXVT
END IF
IF (NMOM_G.GE.2) THEN
XFSEDRG = XCG*GAMMA_X0D(XNUG+(XDG+XBG)/XALPHAG)/GAMMA_X0D(XNUG+XBG/XALPHAG)* &
(ZRHO00)**XCEXVT
XFSEDCG = XCG*GAMMA_X0D(XNUG+XDG/XALPHAG)/GAMMA_X0D(XNUG)* &
(ZRHO00)**XCEXVT
END IF
IF (NMOM_H.GE.2) THEN
XFSEDRH = XCH*GAMMA_X0D(XNUH+(XDH+XBH)/XALPHAH)/GAMMA_X0D(XNUH+XBH/XALPHAH)* &
(ZRHO00)**XCEXVT
XFSEDCH = XCH*GAMMA_X0D(XNUH+XDH/XALPHAH)/GAMMA_X0D(XNUH)* &
(ZRHO00)**XCEXVT
END IF
!
!
!
XLB(4) = XLBI
XLBEX(4) = XLBEXI
XD(4) = XDI
XFSEDR(4) = XFSEDRI
XFSEDC(4) = XFSEDCI
!
XLB(5) = XLBS
XLBEX(5) = XLBEXS
XD(5) = XDS
XFSEDR(5) = XCS*GAMMA_X0D(XNUS+(XDS+XBS)/XALPHAS)/GAMMA_X0D(XNUS+XBS/XALPHAS)* &
(ZRHO00)**XCEXVT
XFSEDC(5) = XCS*GAMMA_X0D(XNUS+XDS/XALPHAS)/GAMMA_X0D(XNUS)* &
(ZRHO00)**XCEXVT
!
XLB(6) = XLBG
XLBEX(6) = XLBEXG
XD(6) = XDG
XFSEDR(6) = XCG*GAMMA_X0D(XNUG+(XDG+XBG)/XALPHAG)/GAMMA_X0D(XNUG+XBG/XALPHAG)* &
(ZRHO00)**XCEXVT
XFSEDC(6) = XCG*GAMMA_X0D(XNUG+XDG/XALPHAG)/GAMMA_X0D(XNUG)* &
(ZRHO00)**XCEXVT
!
XLB(7) = XLBH
XLBEX(7) = XLBEXH
XD(7) = XDH
XFSEDR(7) = XCH*GAMMA_X0D(XNUH+(XDH+XBH)/XALPHAH)/GAMMA_X0D(XNUH+XBH/XALPHAH)* &
(ZRHO00)**XCEXVT
XFSEDC(7) = XCH*GAMMA_X0D(XNUH+XDH/XALPHAH)/GAMMA_X0D(XNUH)* &
(ZRHO00)**XCEXVT
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
!
!-------------------------------------------------------------------------------
!
!
!* 4. CONSTANTS FOR HETEROGENEOUS NUCLEATION
! --------------------------------------
!
!
! ***************
!* 4.1 LIMA_NUCLEATION
! ***************
!* 4.1.1 Constants for the computation of the number concentration
! of active IN
!
XRHO_CFDC = 0.76
!
XGAMMA = 2.
!
IF (NPHILLIPS == 13) THEN
XAREA1(1) = 2.0E-6 !DM1
XAREA1(2) = XAREA1(1) !DM2
XAREA1(3) = 1.0E-7 !BC
XAREA1(4) = 8.9E-7 !BIO
ELSE IF (NPHILLIPS == 8) THEN
XAREA1(1) = 2.0E-6 !DM1
XAREA1(2) = XAREA1(1) !DM2
XAREA1(3) = 2.7E-7 !BC
XAREA1(4) = 9.1E-7 !BIO
ELSE
call Print_msg( NVERB_FATAL, 'GEN', 'INI_LIMA_COLD_MIXED', 'NPHILLIPS should be equal to 8 or 13' )
END IF
!
!* 4.1.2 Constants for the computation of H_X (the fraction-redu-
! cing IN activity at low S_i and warm T) for X={DM1,DM2,BC,BIO}
!
!
IF (NPHILLIPS == 13) THEN
XDT0(1) = 5. !DM1
XDT0(2) = 5. !DM2
XDT0(3) = 10. !BC
XDT0(4) = 5. !BIOO
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
!
XT0(1) = -40. +273.15 !DM1
XT0(2) = XT0(1) !DM2
XT0(3) = -50. +273.15 !BC
XT0(4) = -20. +273.15 !BIO
!
XSW0 = 0.97
!
XDSI0(1) = 0.1 !DM1
XDSI0(2) = 0.1 !DM2
XDSI0(3) = 0.1 !BC
XDSI0(4) = 0.2 !BIO
!
XH(1) = 0.15 !DM1
XH(2) = 0.15 !DM2
XH(3) = 0. !BC
XH(4) = 0. !O
!
XTX1(1) = -30. +273.15 !DM1
XTX1(2) = XTX1(1) !DM2
XTX1(3) = -25. +273.15 !BC
XTX1(4) = -5. +273.15 !BIO
!
XTX2(1) = -10. +273.15 !DM1
XTX2(2) = XTX2(1) !DM2
XTX2(3) = -15. +273.15 !BC
XTX2(4) = -2. +273.15 !BIO
ELSE IF (NPHILLIPS == 8) THEN
XDT0(1) = 5. !DM1
XDT0(2) = 5. !DM2
XDT0(3) = 5. !BC
XDT0(4) = 5. !O
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
!
XT0(1) = -40. +273.15 !DM1
XT0(2) = XT0(1) !DM2
XT0(3) = -50. +273.15 !BC
XT0(4) = -50. +273.15 !BIO
!
XSW0 = 0.97
!
XDSI0(1) = 0.1 !DM1
XDSI0(2) = 0.1 !DM2
XDSI0(3) = 0.1 !BC
XDSI0(4) = 0.1 !BIO
!
XH(1) = 0.15 !DM1
XH(2) = 0.15 !DM2
XH(3) = 0. !BC
XH(4) = 0. !O
!
XTX1(1) = -5. +273.15 !DM1
XTX1(2) = XTX1(1) !DM2
XTX1(3) = -5. +273.15 !BC
XTX1(4) = -5. +273.15 !BIO
!
XTX2(1) = -2. +273.15 !DM1
XTX2(2) = XTX2(1) !DM2
XTX2(3) = -2. +273.15 !BC
XTX2(4) = -2. +273.15 !BIO
END IF
!
!* 4.1.3 Constants for the computation of the Gauss Hermitte
! quadrature method used for the integration of the total
! crystal number over T>-35°C
!
NDIAM = 70
!
CALL PARAM_LIMA_ALLOCATE('XABSCISS', NDIAM)
CALL PARAM_LIMA_ALLOCATE('XWEIGHT', NDIAM)
!
CALL GAUHER(XABSCISS, XWEIGHT, NDIAM)
!
! *****************
!* 4.2 MEYERS NUCLEATION
! *****************
!
ZFACT_NUCL = 1.0 ! Plates, Columns and Bullet rosettes
!
!* 5.2.1 Constants for nucleation from ice nuclei
!
XNUC_DEP = XFACTNUC_DEP*1000.*ZFACT_NUCL
XEXSI_DEP = 12.96E-2
XEX_DEP = -0.639
!
XNUC_CON = XFACTNUC_CON*1000.*ZFACT_NUCL
XEXTT_CON = -0.262
XEX_CON = -2.8
!
XMNU0 = 6.88E-13
!
!!$IF (LMEYERS) THEN
!!$ WRITE(UNIT=ILUOUT0,FMT='(" Heterogeneous nucleation")')
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XNUC_DEP=",E13.6," XEXSI=",E13.6," XEX=",E13.6)') &
!!$ XNUC_DEP,XEXSI_DEP,XEX_DEP
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XNUC_CON=",E13.6," XEXTT=",E13.6," XEX=",E13.6)') &
!!$ XNUC_CON,XEXTT_CON,XEX_CON
!!$ WRITE(UNIT=ILUOUT0,FMT='(" mass of embryo XMNU0=",E13.6)') XMNU0
!!$END IF
!
! *****************
!* 4.3 NUCLEATION for NMOM_I=1
! *****************
!
XNU10 = 50.*ZFACT_NUCL
XALPHA1 = 4.5
XBETA1 = 0.6
!
XNU20 = 1000.*ZFACT_NUCL
XALPHA2 = 12.96
XBETA2 = 0.639
!
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
!-------------------------------------------------------------------------------
!
!
!* 5. CONSTANTS FOR THE SLOW COLD PROCESSES
! -------------------------------------
!
!
!* 5.1.2 Constants for homogeneous nucleation from haze particules
!
XRHOI_HONH = 925.0
XCEXP_DIFVAP_HONH = 1.94
XCOEF_DIFVAP_HONH = (2.0*XPI)*0.211E-4*XP00/XTT**XCEXP_DIFVAP_HONH
XCRITSAT1_HONH = 2.583
XCRITSAT2_HONH = 207.83
XTMIN_HONH = 180.0
XTMAX_HONH = 240.0
XDLNJODT1_HONH = 4.37
XDLNJODT2_HONH = 0.03
XC1_HONH = 100.0
XC2_HONH = 22.6
XC3_HONH = 0.1
XRCOEF_HONH = (XPI/6.0)*XRHOI_HONH
!
!
!* 5.1.3 Constants for homogeneous nucleation from cloud droplets
!
XTEXP1_HONC = -606.3952*LOG(10.0)
XTEXP2_HONC = -52.6611*LOG(10.0)
XTEXP3_HONC = -1.7439*LOG(10.0)
XTEXP4_HONC = -0.0265*LOG(10.0)
XTEXP5_HONC = -1.536E-4*LOG(10.0)
IF (XALPHAC == 3.0) THEN
XC_HONC = XPI/6.0
XR_HONC = XPI/6.0
ELSE
write ( yval, '( E13.6 )' ) xalphac
call Print_msg( NVERB_FATAL, 'GEN', 'INI_LIMA_COLD_MIXED', 'homogeneous nucleation: XALPHAC='//trim(yval)// &
'/= 3. No algorithm developed for this case' )
END IF
!
!!$GFLAG = .TRUE.
!!$IF (GFLAG) THEN
!!$ WRITE(UNIT=ILUOUT0,FMT='(" Homogeneous nucleation")')
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XTEXP1_HONC=",E13.6)') XTEXP1_HONC
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XTEXP2_HONC=",E13.6)') XTEXP2_HONC
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XTEXP3_HONC=",E13.6)') XTEXP3_HONC
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XTEXP4_HONC=",E13.6)') XTEXP4_HONC
!!$ WRITE(UNIT=ILUOUT0,FMT='(" XTEXP5_HONC=",E13.6)') XTEXP5_HONC
!!$ WRITE(UNIT=ILUOUT0,FMT='("XC_HONC=",E13.6," XR_HONC=",E13.6)') XC_HONC,XR_HONC
!!$END IF
!
!
!* 5.2 Constants for vapor deposition on ice
!
XSCFAC = (0.63**(1./3.))*SQRT((ZRHO00)**XCEXVT) ! One assumes Sc=0.63
!
X0DEPI = (4.0*XPI)*XC1I*XF0I*MOMG(XALPHAI,XNUI,1.)
X2DEPI = (4.0*XPI)*XC1I*XF2I*XC_I*MOMG(XALPHAI,XNUI,XDI+2.0)
!
! Harrington parameterization for ice to snow conversion
!
XDICNVS_LIM = 125.E-6 ! size in microns
XLBDAICNVS_LIM = (50.0**(1.0/(XALPHAI)))/XDICNVS_LIM ! ZLBDAI Limitation
XC0DEPIS = ((4.0*XPI)/(XAI*XBI))*XC1I*XF0IS* &
(XALPHAI/GAMMA_X0D(XNUI))*XDICNVS_LIM**(1.0-XBI)
XC1DEPIS = ((4.0*XPI)/(XAI*XBI))*XC1I*XF1IS*SQRT(XC_I)* &
(XALPHAI/GAMMA_X0D(XNUI))*XDICNVS_LIM**(1.0-XBI+(XDI+1.0)/2.0)
XR0DEPIS = XC0DEPIS *(XAI*XDICNVS_LIM**XBI)
XR1DEPIS = XC1DEPIS *(XAI*XDICNVS_LIM**XBI)
!
! Harrington parameterization for snow to ice conversion
!
XLBDASCNVI_MAX = 6000.*XTRANS_MP_GAMMAS ! lbdas max after Field (1999)
!
XDSCNVI_LIM = 125.E-6 ! size in microns
XLBDASCNVI_LIM = (50.0**(1.0/(XALPHAS)))/XDSCNVI_LIM ! ZLBDAS Limitation
XC0DEPSI = ((4.0*XPI)/(XAS*XBS))*XC1S*XF0IS* &
(XALPHAS/GAMMA_X0D(XNUS))*XDSCNVI_LIM**(1.0-XBS)
XC1DEPSI = ((4.0*XPI)/(XAS*XBS))*XC1S*XF1IS*SQRT(XCS)* &
(XALPHAS/GAMMA_X0D(XNUS))*XDSCNVI_LIM**(1.0-XBS+(XDS+1.0)/2.0)
XR0DEPSI = XC0DEPSI *(XAS*XDSCNVI_LIM**XBS)
XR1DEPSI = XC1DEPSI *(XAS*XDSCNVI_LIM**XBS)
!
! Vapor deposition on snow and graupel and hail
!
X0DEPS = (4.0*XPI)*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.)
X1DEPS = (4.0*XPI)*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
XEX0DEPS = -1.0
XEX1DEPS = -0.5*(XDS+3.0)
X0DEPG = (4.0*XPI)*XC1G*XF0G*MOMG(XALPHAG,XNUG,1.)
X1DEPG = (4.0*XPI)*XC1G*XF1G*SQRT(XCG)*MOMG(XALPHAG,XNUG,0.5*XDG+1.5)
XEX0DEPG = -1.0
XEX1DEPG = -0.5*(XDG+3.0)
!
X0DEPH = (4.0*XPI)*XC1H*XF0H*MOMG(XALPHAH,XNUH,1.)
X1DEPH = (4.0*XPI)*XC1H*XF1H*SQRT(XCH)*MOMG(XALPHAH,XNUH,0.5*XDH+1.5)
XEX0DEPH = -1.0
XEX1DEPH = -0.5*(XDH+3.0)
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
!
!-------------------------------------------------------------------------------
!
!
!* 6. CONSTANTS FOR THE COALESCENCE PROCESSES
! ---------------------------------------
!
!
!* 6.0 Precalculation of the gamma function momentum
!
ZGAMI(1) = GAMMA_X0D(XNUI)
ZGAMI(2) = MOMG(XALPHAI,XNUI,3.)
ZGAMI(3) = MOMG(XALPHAI,XNUI,6.)
ZGAMI(4) = ZGAMI(3)-ZGAMI(2)**2 ! useful for Sig_I
ZGAMI(5) = MOMG(XALPHAI,XNUI,9.)
ZGAMI(6) = MOMG(XALPHAI,XNUI,3.+XBI)
ZGAMI(7) = MOMG(XALPHAI,XNUI,XBI)
ZGAMI(8) = MOMG(XALPHAI,XNUI,3.)/MOMG(XALPHAI,XNUI,2.)
!
ZGAMS(1) = GAMMA_X0D(XNUS)
ZGAMS(2) = MOMG(XALPHAS,XNUS,3.)
!
!
!* 6.1 Csts for the coalescence processes
!
ZFAC_ZRNIC = 0.1
XKER_ZRNIC_A1 = 2.59E15*ZFAC_ZRNIC**2! From Long a1=9.44E9 cm-3
! so XKERA1= 9.44E9*1E6*(PI/6)**2
XKER_ZRNIC_A2 = 3.03E3*ZFAC_ZRNIC ! From Long a2=5.78E3
! so XKERA2= 5.78E3* (PI/6)
!
!
!* 6.2 Csts for the pristine ice selfcollection process
!
XSELFI = XKER_ZRNIC_A1*ZGAMI(3)
XCOLEXII = 0.025 ! Temperature factor of the I+I collection efficiency
!
!
!* 6.3 Constants for pristine ice autoconversion
!
XTEXAUTI = 0.025 ! Temperature factor of the I+I collection efficiency
!
!!$GFLAG = .TRUE.
!!$IF (GFLAG) THEN
!!$ WRITE(UNIT=ILUOUT0,FMT='(" pristine ice autoconversion")')
!!$ WRITE(UNIT=ILUOUT0,FMT='(" Temp. factor XTEXAUTI=",E13.6)') XTEXAUTI
!!$END IF
!
XAUTO3 = 6.25E18*(ZGAMI(2))**(1./3.)*SQRT(ZGAMI(4))
XAUTO4 = 0.5E6*(ZGAMI(4))**(1./6.)
XLAUTS = 2.7E-2
XLAUTS_THRESHOLD = 0.4
XITAUTS= 0.27 ! (Notice that T2 of BR74 is uncorrect and that 0.27=1./3.7
XITAUTS_THRESHOLD = 7.5
!
!
!* 6.4 Constants for snow aggregation
!
XCOLEXIS = 0.05 ! Temperature factor of the I+S collection efficiency
XFIAGGS = (XPI/4.0)*0.25*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
XEXIAGGS = -XDS - 2.0
XAGGS_CLARGE1 = XKER_ZRNIC_A2*ZGAMI(2)
XAGGS_CLARGE2 = XKER_ZRNIC_A2*ZGAMS(2)
XAGGS_RLARGE1 = XKER_ZRNIC_A2*ZGAMI(6)*XAI
XAGGS_RLARGE2 = XKER_ZRNIC_A2*ZGAMI(7)*ZGAMS(2)*XAI
!
!!$GFLAG = .TRUE.
!!$IF (GFLAG) THEN
!!$ WRITE(UNIT=ILUOUT0,FMT='(" snow aggregation")')
!!$ WRITE(UNIT=ILUOUT0,FMT='(" Temp. factor XCOLEXIS=",E13.6)') XCOLEXIS
!!$END IF
!
!-------------------------------------------------------------------------------
!
!
!* 7. CONSTANTS FOR THE FAST COLD PROCESSES FOR THE AGGREGATES
! --------------------------------------------------------
!
!
!* 7.1 Constants for the riming of the aggregates
!
XDCSLIM = 0.007 ! D_cs^lim = 7 mm as suggested by Farley et al. (1989)
XCOLCS = 1.0
XEXCRIMSS= -XDS-2.0
XCRIMSS = (XPI/4.0)*XCOLCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
XEXCRIMSG= XEXCRIMSS
XCRIMSG = XCRIMSS
XSRIMCG = XAS*MOMG(XALPHAS,XNUS,XBS)
XEXSRIMCG= -XBS
!Pour Murakami 1990
XSRIMCG2 = XAG*MOMG(XALPHAS,XNUS,XBG)
XSRIMCG3 = 0.1
XEXSRIMCG2=XBG
!!$GFLAG = .TRUE.
!!$IF (GFLAG) THEN
!!$ WRITE(UNIT=ILUOUT0,FMT='(" riming of the aggregates")')
!!$ WRITE(UNIT=ILUOUT0,FMT='(" D_cs^lim (Farley et al.) XDCSLIM=",E13.6)') XDCSLIM
!!$ WRITE(UNIT=ILUOUT0,FMT='(" Coll. efficiency XCOLCS=",E13.6)') XCOLCS
!!$END IF
!!$!
PARAM_LIMA_MIXED%NGAMINC = 80
PARAM_LIMA_MIXED%XGAMINC_BOUND_MIN = (1000.*XTRANS_MP_GAMMAS*XDCSLIM)**XALPHAS !1.0E-1 ! Minimal value of (Lbda * D_cs^lim)**alpha
PARAM_LIMA_MIXED%XGAMINC_BOUND_MAX = (50000.*XTRANS_MP_GAMMAS*XDCSLIM)**XALPHAS !1.0E7 ! Maximal value of (Lbda * D_cs^lim)**alpha
ZRATE = EXP(LOG(PARAM_LIMA_MIXED%XGAMINC_BOUND_MAX/PARAM_LIMA_MIXED%XGAMINC_BOUND_MIN)/FLOAT(PARAM_LIMA_MIXED%NGAMINC-1))
CALL PARAM_LIMA_MIXED_ALLOCATE('XGAMINC_RIM1', NGAMINC)
CALL PARAM_LIMA_MIXED_ALLOCATE('XGAMINC_RIM2', NGAMINC)
CALL PARAM_LIMA_MIXED_ALLOCATE('XGAMINC_RIM4', NGAMINC)
!
DO J1=1,NGAMINC
ZBOUND = PARAM_LIMA_MIXED%XGAMINC_BOUND_MIN*ZRATE**(J1-1)
XGAMINC_RIM1(J1) = GAMMA_INC(XNUS+(2.0+XDS)/XALPHAS,ZBOUND)
XGAMINC_RIM2(J1) = GAMMA_INC(XNUS+XBS/XALPHAS ,ZBOUND)
XGAMINC_RIM4(J1) = GAMMA_INC(XNUS+XBG/XALPHAS ,ZBOUND) ! Pour Murakami 1990
END DO
!
XRIMINTP1 = XALPHAS / LOG(ZRATE)
XRIMINTP2 = 1.0 + XRIMINTP1*LOG( XDCSLIM/(XGAMINC_BOUND_MIN)**(1.0/XALPHAS) )
!
!* 7.1.1 Defining the constants for the Hallett-Mossop
! secondary ice nucleation process
!
XHMTMIN = XTT - 8.0
XHMTMAX = XTT - 3.0
XHM1 = 9.3E-3 ! Obsolete parameterization
XHM2 = 1.5E-3/LOG(10.0) ! from Ferrier (1995)
XHM_YIELD = 5.E-3 ! A splinter is produced after the riming of 200 droplets
XHM_COLLCS= 1.0 ! Collision efficiency snow/droplet (with Dc>25 microns)
XHM_FACTS = XHM_YIELD*(XHM_COLLCS/XCOLCS)
!
! Notice: One magnitude of lambda discretized over 10 points for the droplets
!
PARAM_LIMA_MIXED%XGAMINC_HMC_BOUND_MIN = 1.0E-3 ! Min value of (Lbda * (12,25) microns)**alpha
PARAM_LIMA_MIXED%XGAMINC_HMC_BOUND_MAX = 1.0E5 ! Max value of (Lbda * (12,25) microns)**alpha
ZRATE = EXP(LOG(PARAM_LIMA_MIXED%XGAMINC_HMC_BOUND_MAX/PARAM_LIMA_MIXED%XGAMINC_HMC_BOUND_MIN)/REAL(PARAM_LIMA_MIXED%NGAMINC-1))
CALL PARAM_LIMA_MIXED_ALLOCATE('XGAMINC_HMC', NGAMINC)
!
DO J1=1,NGAMINC
ZBOUND = PARAM_LIMA_MIXED%XGAMINC_HMC_BOUND_MIN*ZRATE**(J1-1)
XGAMINC_HMC(J1) = GAMMA_INC(XNUC,ZBOUND)
END DO
!
XHMSINTP1 = XALPHAC / LOG(ZRATE)
XHMSINTP2 = 1.0 + XHMSINTP1*LOG( 12.E-6/(XGAMINC_HMC_BOUND_MIN)**(1.0/XALPHAC) )
XHMLINTP1 = XALPHAC / LOG(ZRATE)
XHMLINTP2 = 1.0 + XHMLINTP1*LOG( 25.E-6/(XGAMINC_HMC_BOUND_MIN)**(1.0/XALPHAC) )
!
!
!* 7.2 Constants for the accretion of raindrops onto aggregates
!
XFRACCSS = XPI/4.0*XAR*(ZRHO00**XCEXVT)
XFNRACCSS = XPI/4.0*(ZRHO00**XCEXVT)
!
XLBRACCS1 = MOMG(XALPHAS,XNUS,2.)*MOMG(XALPHAR,XNUR,3.)
XLBRACCS2 = 2.*MOMG(XALPHAS,XNUS,1.)*MOMG(XALPHAR,XNUR,4.)
XLBRACCS3 = MOMG(XALPHAR,XNUR,5.)
XLBNRACCS1 = MOMG(XALPHAS,XNUS,2.)
XLBNRACCS2 = 2.*MOMG(XALPHAS,XNUS,1.)*MOMG(XALPHAR,XNUR,1.)
XLBNRACCS3 = MOMG(XALPHAR,XNUR,2.)
XFSACCRG = (XPI/4.0)*XAS*(ZRHO00**XCEXVT)
XFNSACCRG = (XPI/4.0)*(ZRHO00**XCEXVT)
!
XLBSACCR1 = MOMG(XALPHAR,XNUR,2.)*MOMG(XALPHAS,XNUS,XBS)
XLBSACCR2 = 2.*MOMG(XALPHAR,XNUR,1.)*MOMG(XALPHAS,XNUS,XBS+1.)
XLBSACCR3 = MOMG(XALPHAS,XNUS,XBS+2.)
XLBNSACCR1 = MOMG(XALPHAR,XNUR,2.)
XLBNSACCR2 = 2.*MOMG(XALPHAR,XNUR,1.)*MOMG(XALPHAS,XNUS,1.)
XLBNSACCR3 = MOMG(XALPHAS,XNUS,2.)
!
!* 7.2.1 Defining the ranges for the computation of the kernels
!
! Notice: One magnitude of lambda discretized over 10 points for rain
! Notice: One magnitude of lambda discretized over 10 points for snow
!
PARAM_LIMA_MIXED%NACCLBDAS = 40
PARAM_LIMA_MIXED%XACCLBDAS_MIN = 5.0E1*XTRANS_MP_GAMMAS !5.0E1*XTRANS_MP_GAMMAS ! Minimal value of Lbda_s to tabulate XKER_RACCS
PARAM_LIMA_MIXED%XACCLBDAS_MAX = 5.0E5*XTRANS_MP_GAMMAS !5.0E5*XTRANS_MP_GAMMAS ! Maximal value of Lbda_s to tabulate XKER_RACCS
ZRATE = LOG(XACCLBDAS_MAX/XACCLBDAS_MIN)/FLOAT(NACCLBDAS-1)
PARAM_LIMA_MIXED%XACCINTP1S = 1.0 / ZRATE
PARAM_LIMA_MIXED%XACCINTP2S = 1.0 - LOG( XACCLBDAS_MIN ) / ZRATE
PARAM_LIMA_MIXED%NACCLBDAR = 40
PARAM_LIMA_MIXED%XACCLBDAR_MIN = 1.0E3 ! Minimal value of Lbda_r to tabulate XKER_RACCS
PARAM_LIMA_MIXED%XACCLBDAR_MAX = 1.0E7 ! Maximal value of Lbda_r to tabulate XKER_RACCS
ZRATE = LOG(PARAM_LIMA_MIXED%XACCLBDAR_MAX/PARAM_LIMA_MIXED%XACCLBDAR_MIN)/REAL(NACCLBDAR-1)
XACCINTP1R = 1.0 / ZRATE
XACCINTP2R = 1.0 - LOG( PARAM_LIMA_MIXED%XACCLBDAR_MIN ) / ZRATE
!
!* 7.2.2 Computations of the tabulated normalized kernels
!
IND = 50 ! Interval number, collection efficiency and infinite diameter
ZESR = 1.0 ! factor used to integrate the dimensional distributions when
ZFDINFTY = 20.0 ! computing the kernels XKER_RACCSS, XKER_RACCS and XKER_SACCRG
!
CALL PARAM_LIMA_MIXED_ALLOCATE('XKER_RACCSS', NACCLBDAS,NACCLBDAR)
CALL PARAM_LIMA_MIXED_ALLOCATE('XKER_RACCS', NACCLBDAS,NACCLBDAR)
CALL PARAM_LIMA_MIXED_ALLOCATE('XKER_SACCRG', NACCLBDAR,NACCLBDAS)
CALL PARAM_LIMA_MIXED_ALLOCATE('XKER_N_RACCSS', NACCLBDAS,NACCLBDAR)
CALL PARAM_LIMA_MIXED_ALLOCATE('XKER_N_RACCS', NACCLBDAS,NACCLBDAR)
CALL PARAM_LIMA_MIXED_ALLOCATE('XKER_N_SACCRG', NACCLBDAR,NACCLBDAS)
CALL NRCOLSS ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
ZESR, XCS, XDS, XFVELOS, XCR, XDR, &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_N_RACCSS, XAG, XBS, XAS )
CALL NZCOLX ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
ZESR, XCS, XDS, XFVELOS, XCR, XDR, 0., & !
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_N_RACCS )
CALL NSCOLRG ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
ZESR, XCS, XDS, XFVELOS, XCR, XDR, &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_N_SACCRG,XAG, XBS, XAS )
!!$WRITE(UNIT=ILUOUT0,FMT='("!")')
!!$WRITE(UNIT=ILUOUT0,FMT='("IF( PRESENT(PKER_N_RACCSS) ) THEN")')
!!$DO J1 = 1 , NACCLBDAS
!!$ DO J2 = 1 , NACCLBDAR
!!$ WRITE(UNIT=ILUOUT0,FMT='(" PKER_N_RACCSS(",I3,",",I3,") = ",E13.6)') &
!!$ J1,J2,XKER_N_RACCSS(J1,J2)
!!$ END DO
!!$END DO
!!$WRITE(UNIT=ILUOUT0,FMT='("!")')
!!$WRITE(UNIT=ILUOUT0,FMT='("!")')
!!$WRITE(UNIT=ILUOUT0,FMT='("IF( PRESENT(PKER_N_RACCS) ) THEN")')
!!$DO J1 = 1 , NACCLBDAS
!!$ DO J2 = 1 , NACCLBDAR
!!$ WRITE(UNIT=ILUOUT0,FMT='(" PKER_N_RACCS(",I3,",",I3,") = ",E13.6)') &
!!$ J1,J2,XKER_N_RACCS(J1,J2)
!!$ END DO
!!$END DO
!!$WRITE(UNIT=ILUOUT0,FMT='("!")')
!!$WRITE(UNIT=ILUOUT0,FMT='("!")')
!!$WRITE(UNIT=ILUOUT0,FMT='("IF( PRESENT(PKER_N_SACCRG) ) THEN")')
!!$DO J1 = 1 , NACCLBDAR
!!$ DO J2 = 1 , NACCLBDAS
!!$ WRITE(UNIT=ILUOUT0,FMT='(" PKER_N_SACCRG(",I3,",",I3,") = ",E13.6)') &
!!$ J1,J2,XKER_N_SACCRG(J1,J2)