11 INTEGER,
DIMENSION(366) ::
cblday = 0
13 CHARACTER(len=200),
DIMENSION(366) ::
filesonde =
""
15 REAL(kind(1d0)) ::
wsb
16 REAL(kind(1d0)),
DIMENSION(1:10) ::
cbldata
17 REAL(kind(1d0)),
DIMENSION(:, :),
ALLOCATABLE ::
inicbldata
24 REAL(kind(1d0)) ::
c2k = 273.16
52 REAL(kind(1d0)),
DIMENSION(6) ::
y
69 SUBROUTINE cbl(iy, id, it, imin, ir, Gridiv, qh_choice, dectime, Temp_C, Press_hPa, avkdn, avu1, avrh, &
70 avcp, avdens, es_hPa, lv_J_kg, nsh_real, tstep, UStar, psih, is, NumberOfGrids, qhforCBL, qeforCBL, &
71 ReadLinesMetdata, dataOutBL)
73 INTEGER,
PARAMETER :: ncolumnsdataOutBL = 22
75 INTEGER,
INTENT(IN) :: tstep, is, NumberOfGrids, Gridiv, ReadLinesMetdata, ir
76 REAL(KIND(1D0)),
INTENT(IN),
DIMENSION(NumberOfGrids) :: qhforCBL, qeforCBL
77 REAL(KIND(1D0)),
INTENT(IN) :: avkdn, nsh_real, UStar, psih
78 INTEGER,
INTENT(INOUT) :: qh_choice, iy, id, it, imin
79 REAL(KIND(1D0)),
INTENT(INOUT) :: dectime, Press_hPa, avu1, avrh, es_hPa, avcp, avdens, lv_J_kg
80 REAL(KIND(1D0)),
INTENT(OUT) :: Temp_C
81 REAL(KIND(1D0)),
INTENT(OUT),
DIMENSION(ReadLinesMetdata, ncolumnsdataOutBL, NumberOfGrids) :: dataOutBL
83 REAL(KIND(1D0)) :: gas_ct_dry = 8.31451/0.028965
85 REAL(KIND(1D0)) :: qh_use, qe_use, tm_K_zm, qm_gkg_zm
86 REAL(KIND(1D0)) :: Temp_C1, avrh1, es_hPa1
87 REAL(KIND(1D0)) :: secs0, secs1, Lv
88 REAL(KIND(1D0)) :: NAN = -999
89 INTEGER :: idoy, startflag
131 CALL nbl(iy, id, it, imin, dectime, ir, qh_choice, qh_use, qe_use, tm_k_zm, qm_gkg_zm, startflag, gridiv, &
132 psih, ustar, temp_c, &
133 numberofgrids, qhforcbl, qeforcbl, &
134 press_hpa, avu1, avrh, &
135 readlinesmetdata, dataoutbl, &
136 avcp, avdens, es_hpa, lv_j_kg)
141 IF (startflag == 0)
THEN
143 dataoutbl(ir, 1:ncolumnsdataoutbl, gridiv) &
144 = (/real(iy, 8), real(id, 8), real(it, 8), real(imin, 8), dectime,
blh_m,
tm_k, &
149 qh_use = qhforcbl(gridiv)
150 qe_use = qeforcbl(gridiv)
151 IF (qh_use < -900 .OR. qe_use < -900)
THEN
152 CALL errorhint(22,
'Unrealistic qh or qe_value for CBL in CBL.', qh_use, qe_use, qh_choice)
176 cbldata(1) = float(it) + float(imin)/60.
191 secs1 = secs0 + float(tstep)
250 IF (it == 0 .AND. imin == (nsh_real - 1)/nsh_real*60) idoy = id - 1
253 IF ((qh_choice == 1) .OR. (qh_choice == 2))
THEN
258 lv = (2500.25 - 2.365*temp_c)*1000
262 CALL errorhint(34,
'subroutine CBL dectime, relative humidity', idoy +
cbldata(1)/24.0, avrh, 100)
267 dataoutbl(ir, 1:ncolumnsdataoutbl, gridiv) &
268 = (/real(iy, 8), real(id, 8), real(it, 8), real(imin, 8), dectime,
blh_m,
tm_k, &
270 temp_c, avrh,
cbldata([2, 3, 9, 7, 8, 4, 5, 6]), &
272 ELSEIF (qh_choice == 3)
THEN
276 lv = (2500.25 - 2.365*temp_c1)*1000
279 IF (avrh1 > 100)
THEN
280 CALL errorhint(34,
'subroutine CBL dectime, relative humidity', idoy +
cbldata(1)/24.0, avrh1, 100)
285 dataoutbl(ir, 1:ncolumnsdataoutbl, gridiv) &
286 = (/real(iy, 8), real(id, 8), real(it, 8), real(imin, 8), dectime,
blh_m,
tm_k, &
288 temp_c1, avrh1,
cbldata([2, 3, 9, 7, 8, 4, 5, 6]), &
303 INTEGER,
INTENT(INOUT) :: qh_choice
304 CHARACTER(len=150),
INTENT(IN) :: FileInputPath
320 OPEN (51, file=trim(fileinputpath)//
'CBLInput.nml', status=
'old', err=24)
321 READ (51, nml=cblinput, err=24)
330 READ (52, *, iostat=ios) l
331 IF (ios < 0 .OR. l == -9)
EXIT
35424
CALL errorhint(24,
'CBLInput.nml', 0.00d0, 0.000d0, 0)
361 SUBROUTINE cbl_initial(qh_use, qe_use, tm_K_zm, qm_gkg_zm, startflag, ir, Gridiv)
379 REAL(KIND(1D0)) :: qh_use, qe_use, tm_K_zm, qm_gkg_zm
380 REAL(KIND(1D0)) :: lv
381 INTEGER :: i, nLineDay, ir, Gridiv, startflag
385 IF (qh_use < -900 .OR. qe_use < -900)
THEN
386 CALL errorhint(22,
'Unrealistic qh or qe_value for CBL in CBL_initial.', qh_use, qe_use,
qh_choice)
408 dataoutbl(ir, 1:
ncolumnsdataoutbl, gridiv) = (/real(
iy, 8), real(
id, 8), real(
it, 8), real(
imin, 8),
dectime, &
414 nlineday = nlineday + 1
434 lv = (2500.25 - 2.365*
temp_c)*1000
446 lv = (2500.25 - 2.365*
temp_c)*1000
482 SUBROUTINE nbl(iy, id, it, imin, dectime, ir, qh_choice, qh_use, qe_use, &
483 tm_K_zm, qm_gkg_zm, startflag, Gridiv, &
484 psih, UStar, Temp_C, NumberOfGrids, qhforCBL, qeforCBL, Press_hPa, avu1, avrh, &
485 ReadLinesMetdata, dataOutBL, &
486 avcp, avdens, es_hPa, lv_J_kg)
489 INTEGER,
PARAMETER :: ncolumnsdataOutBL = 22
491 INTEGER,
INTENT(IN) :: qh_choice, iy, id, it, imin, NumberOfGrids, ReadLinesMetdata, ir
492 REAL(KIND(1D0)),
INTENT(IN) :: Press_hPa, psih, UStar
493 REAL(KIND(1D0)),
INTENT(OUT) :: Temp_C, tm_K_zm, qm_gkg_zm
494 REAL(KIND(1D0)),
INTENT(INOUT) :: dectime, avu1, avRH, avcp, avdens, es_hPa, lv_J_kg
495 REAL(KIND(1D0)),
INTENT(IN),
DIMENSION(NumberOfGrids) :: qhforCBL, qeforCBL
496 REAL(KIND(1D0)),
INTENT(OUT),
DIMENSION(ReadLinesMetdata, ncolumnsdataOutBL, NumberOfGrids) :: dataOutBL
500 gas_ct_dry = 8.31451/0.028965
501 REAL(KIND(1D0)) :: qh_use, qe_use
502 REAL(KIND(1D0)) :: lv
503 INTEGER :: i, nLineDay, Gridiv, startflag
505 qh_use = qhforcbl(gridiv)
506 qe_use = qeforcbl(gridiv)
507 IF (qh_use < -900 .OR. qe_use < -900)
THEN
508 CALL errorhint(22,
'Unrealistic qh or qe value for CBL in NBL.', qh_use, qe_use, qh_choice)
512 DO i = 1, nlineindata
513 IF (int(inicbldata(i, 1)) <= id)
THEN
514 nlineday = nlineday + 1
524 tm_k = inicbldata(nlineday, 7)
525 qm_gkg = inicbldata(nlineday, 8)
531 gamt_km = inicbldata(nlineday, 3)
532 gamq_gkgm = inicbldata(nlineday, 4)
533 tp_k = inicbldata(nlineday, 5)
534 qp_gkg = inicbldata(nlineday, 6)
535 tm_k = inicbldata(nlineday, 7)
536 qm_gkg = inicbldata(nlineday, 8)
542 temp_c = tm_k/((1000/press_hpa)**(gas_ct_dry/avcp)) - c2k
543 es_hpa = sat_vap_press_x(temp_c, press_hpa, 1, dectime)
544 lv = (2500.25 - 2.365*temp_c)*1000
546 avrh = 100*((qm_gkg*press_hpa/(622 + qm_gkg))/es_hpa)
548 CALL errorhint(34,
'subroutine CBL dectime, relative humidity', dectime, avrh, 100)
555 dataoutbl(ir, 1:ncolumnsdataoutbl, gridiv) = &
556 (/real(iy, 8), real(id, 8), real(it, 8), real(imin, 8), dectime, &
557 blh_m, tm_k, qm_gkg, &
559 temp_c, avrh, qh_use, qe_use, press_hpa, avu1, ustar, avdens, lv_j_kg, avcp, &
560 gamt_km, gamq_kgkgm/)
562 IF (initialdata_use == 2)
THEN
563 blh_m = inicbldata(nlineday, 2)
564 gamt_km = inicbldata(nlineday, 3)
565 gamq_gkgm = inicbldata(nlineday, 4)
566 tp_k = inicbldata(nlineday, 5)
567 qp_gkg = inicbldata(nlineday, 6)
568 tm_k = inicbldata(nlineday, 7)
569 qm_gkg = inicbldata(nlineday, 8)
570 ELSEIF (initialdata_use == 1 .AND. inicbldata(nlineday, 1) == id)
THEN
571 blh_m = inicbldata(nlineday, 2)
572 gamt_km = inicbldata(nlineday, 3)
573 gamq_gkgm = inicbldata(nlineday, 4)
574 tm_k_zm = (temp_c + c2k)*((1000/press_hpa)**(gas_ct_dry/avcp))
575 tm_k = tm_k_zm - psih*qh_use/(k*ustar*avcp*avdens)
576 es_hpa = sat_vap_press_x(temp_c, press_hpa, 1, dectime)
577 qm_gkg_zm = 622*avrh/(100*press_hpa/es_hpa - avrh)
578 lv = (2500.25 - 2.365*temp_c)*1000
579 qm_gkg = qm_gkg_zm - psih*qe_use/(k*ustar*avdens*lv)
582 ELSEIF (initialdata_use == 0)
THEN
586 tm_k_zm = (temp_c + c2k)*((1000/press_hpa)**(gas_ct_dry/avcp))
587 tm_k = tm_k_zm - psih*qh_use/(k*ustar*avcp*avdens)
588 es_hpa = sat_vap_press_x(temp_c, press_hpa, 1, dectime)
589 qm_gkg_zm = 622*avrh/(100*press_hpa/es_hpa - avrh)
590 lv = (2500.25 - 2.365*temp_c)*1000
591 qm_gkg = es_hpa - psih*qe_use/(k*ustar*avdens*lv)
596 gamq_kgkgm = gamq_gkgm/1000.
597 qp_kgkg = qp_gkg/1000
598 qm_kgkg = qm_gkg/1000
603 IF (sondeflag == 1 .AND. inicbldata(id, 2) /= -999)
THEN
612 IF (qp_kgkg > qsatf(tp_c, press_hpa) .OR. qp_kgkg < 0)
THEN
613 qp_kgkg = qsatf(tp_c, press_hpa)
615 IF (qm_kgkg > qsatf(tm_c, press_hpa) .OR. qm_kgkg < 0)
THEN
616 qm_kgkg = qsatf(tm_c, press_hpa)
631 SUBROUTINE rkutta(neqn_use, XA, XB, y_use, NSTEPS)
659 INTEGER :: ns, nsteps, nj, n, neqn_use
660 REAL(KIND(1D0)),
DIMENSION(neqn_use) :: y_use
661 REAL(KIND(1D0)),
DIMENSION(21) :: dydx, arg
662 REAL(KIND(1D0)),
DIMENSION(21, 5) :: rk
663 REAL(KIND(1D0)),
DIMENSION(4) :: coef
664 REAL(KIND(1D0)) :: XA, XB, step, X, xx
671 step = (xb - xa)/nsteps
677 x = xa + (ns - 1)*step
682 xx = x + coef(n)*step
686 arg(nj) = y_use(nj) + coef(n)*rk(nj, n)
689 CALL diff(xx, arg, dydx)
692 rk(nj, n + 1) = step*dydx(nj)
698 y_use(nj) = y_use(nj) + rk(nj, n + 1)/(6*coef(n))
726 REAL(KIND(1D0)),
DIMENSION(neqn) :: dyds, y1
727 REAL(KIND(1D0)) :: zero = 0.0
728 REAL(KIND(1D0)) :: h1, t_K, q_kgkg, c, cp, ws, s, foo
730 REAL(KIND(1D0)) :: delt_K, delq_kgkg, delc
731 REAL(KIND(1D0)) :: gamtv_Km, deltv_K, ftv_Kms
732 REAL(KIND(1D0)) :: ftva_Kms, delb, qs2, qs3
733 REAL(KIND(1D0)) :: dhds, dtds, dqds, dcds, dtpds, dqpds
734 REAL(KIND(1D0)) :: conk, conn, cona, conc, cont
759 deltv_k = delt_k + 0.61*
tm_k*delq_kgkg
762 ftva_kms = max(ftv_kms, zero)
763 ws = (h1*ftva_kms*
grav/
tm_k)**0.3333333333
768 dhds = ftva_kms/(h1*gamtv_km)
772 IF (deltv_k <= 0.01)
THEN
773 dhds = ftva_kms/(h1*gamtv_km)
774 CALL errorhint(30,
"subroutine diff [CBL: Deltv_K<0.01 EntrainmentType=1], deltv_K,delt_K,", deltv_k, delt_k,
notusedi)
781 dhds = (conc*ws**3 + cona*
cbldata(8)**3)/(h1*delb)
787 IF (deltv_k <= 0.01)
THEN
788 dhds = ftva_kms/(h1*gamtv_km)
789 CALL errorhint(31,
'subroutine difflfnout: [CBL: deltv_K<0.01 EntrainmentType=4],deltv_K', &
796 dhds =
alpha3*ftva_kms/deltv_k
807 qs3 = ws**3 + (conn*
cbldata(8))**3
808 qs2 = qs3**(0.6666666667)
810 IF (deltv_k <= 0.01)
THEN
811 dhds = ftva_kms/(h1*gamtv_km)
812 CALL errorhint(31,
'subroutine difflfnout: [CBL: deltv_K<0.01 EntrainmentType=3],deltv_K', &
817 dhds = (conk*qs3)/(cont*qs2 + h1*delb)
828 dcds =
fcbl/h1 + delc*(dhds -
wsb)/h1
833 dtds =
fhbl_kms/h1 + delt_k*(dhds)/h1
835 dcds =
fcbl/h1 + delc*(dhds)/h1
859 INTEGER :: i, fn = 101, izm = 500, notusedi = -9999, id
860 CHARACTER(len=200) :: FileN
861 REAL(KIND(1D0)) :: dxx
862 REAL(KIND(1D0)),
PARAMETER :: notUsed = -9999.99
865 OPEN (fn, file=filen, status=
"old", err=24)
872 READ (fn, *,
END=900, err=25) gtheta(i, 1), dxx, gtheta(i, 2), ghum(i, 1), dxx, ghum(i, 2)
873 ghum(i, 2) = ghum(i, 2)
877 CALL errorhint(23, filen, real(zmax, kind(1d0)), notused, izm)
88124
CALL errorhint(24, filen, notused, notused, notusedi)
88225
CALL errorhint(25, filen, notused, notused, i)
892 REAL(KIND(1D0)) :: gamtt, gamqq
907 gamqq =
ghum(j - 1, 2)
integer, parameter ncolumnsdataoutbl
real(kind(1d0)), dimension(:, :, :), allocatable dataoutbl
real(kind(1d0)), parameter k
subroutine cbl_initial(qh_use, qe_use, tm_K_zm, qm_gkg_zm, startflag, ir, Gridiv)
subroutine cbl_readinputdata(FileInputPath, qh_choice)
subroutine cbl(iy, id, it, imin, ir, Gridiv, qh_choice, dectime, Temp_C, Press_hPa, avkdn, avu1, avrh, avcp, avdens, es_hPa, lv_J_kg, nsh_real, tstep, UStar, psih, is, NumberOfGrids, qhforCBL, qeforCBL, ReadLinesMetdata, dataOutBL)
subroutine rkutta(neqn_use, XA, XB, y_use, NSTEPS)
subroutine diff(s, y1, dyds)
subroutine nbl(iy, id, it, imin, dectime, ir, qh_choice, qh_use, qe_use, tm_K_zm, qm_gkg_zm, startflag, Gridiv, psih, UStar, Temp_C, NumberOfGrids, qhforCBL, qeforCBL, Press_hPa, avu1, avrh, ReadLinesMetdata, dataOutBL, avcp, avdens, es_hPa, lv_J_kg)
real(kind(1d0)) gamq_kgkgm
real(kind(1d0)), dimension(6) y
real(kind(1d0)), dimension(:, :), allocatable inicbldata
integer, dimension(366) cblday
real(kind(1d0)) febl_kgkgms
real(kind(1d0)) gamq_gkgm
real(kind(1d0)), dimension(0:500, 2) gtheta
real(kind(1d0)), dimension(0:500, 2) ghum
character(len=200) initialdatafilename
real(kind(1d0)), dimension(1:10) cbldata
character(len=200), dimension(366) filesonde
real(kind(1d0)) press_hpa
character(len=150) fileinputpath
real(kind(1d0)) gas_ct_dry
real(kind(1d0)) function qsatf(T, PMB)
real(kind(1d0)) function sat_vap_press_x(Temp_c, PRESS_hPa, from, dectime)
real(kind(1d0)), dimension(:), allocatable qeforcbl
real(kind(1d0)), dimension(:), allocatable qhforcbl
subroutine errorhint(errh, ProblemFile, VALUE, value2, valueI)