Skip to content
Snippets Groups Projects
Commit c607d60a authored by Georgiy Faikin's avatar Georgiy Faikin
Browse files

Model Rothc calculation without evapotranspiration and precipitation

parent 53bf0588
No related branches found
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
source/carbon/carbon_models/carbon_model_ROTHC.f90
+ 122
− 24
View file @ c607d60a
......@@ -43,12 +43,12 @@ module carbon_model_rothc
! VLDR 1 0.00777 1.38303 0.03275 0.14961 0.119
! ROST 1 0.04936 0.799 0.11873 7.07025 0.814
! litterfall
! litterfall ---------------------
call set_flux(fid = n_F, pid_out = n_Cveg, pid_in = n_CDPM, name = 'fdpm*organic')
call set_mult(n_F, 'lin', x_ij = organic, k = fdpm)
call set_flux(fid = n_F, pid_out = n_Cveg, pid_in = n_CRPM, name = '(1-fdpm)*organic')
call set_mult(n_F, 'lin', x_ij = organic, k = (1 - fdpm))
call set_mult(n_F, 'lin', x_ij = organic, k = (1. - fdpm))
call set_flux(fid = n_F, pid_out = n_Cveg, pid_in = n_CRPM, name = '0.49*FYM')
call set_mult(n_F, 'lin', x_ij = FYM, k = fdpm2)
......@@ -57,106 +57,204 @@ module carbon_model_rothc
call set_mult(n_F, 'lin', x_ij = FYM, k = fdpm2)
call set_flux(fid = n_F, pid_out = n_Cveg, pid_in = n_CHUM, name = '0.02*FYM')
call set_mult(n_F, 'lin', x_ij = FYM, k = 1. - 2*fdpm2)
call set_mult(n_F, 'lin', x_ij = FYM, k = 1. - 2.*fdpm2)
! --------------------
! to atmosphere
! to atmosphere -----------------------
call set_flux(fid = n_F, pid_out = n_CDPM, pid_in = n_Catm, name = 'RDPM')
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(1))
call set_mult(n_F, 'lin', x_ijn = CDPM)
call set_flux(fid = n_F, pid_out = n_CRPM, pid_in = n_Catm, name = 'RRPM')
!call set_mult(n_F, 'const', c_n = 47.91)
!call set_mult(n_F, 'exp', x_ij = Tsoil, a = exp, k = 106.06, x0 = -18.27, y0 = 1.)
!call set_mult(n_F, 'step', x_ij = Tsoil, x0 = -4.99)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_flux(fid = n_F, pid_out = n_CRPM, pid_in = n_Catm, name = 'RRPM')
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(2))
call set_mult(n_F, 'lin', x_ijn = CRPM)
call set_flux(fid = n_F, pid_out = n_CBIO, pid_in = n_Catm, name = 'RBIO')
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(3))
call set_mult(n_F, 'lin', x_ijn = CBIO)
call set_flux(fid = n_F, pid_out = n_CHUM, pid_in = n_Catm, name = 'RHUM')
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(4))
call set_mult(n_F, 'lin', x_ijn = CHUM)
! atmosphere fall
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CBIO, name = 'BIO_prop*RDPM')
call set_mult(n_F, 'lin', x_ij = BIO_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(1))
call set_mult(n_F, 'lin', x_ijn = CDPM)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CBIO, name = 'BIO_prop*RRPM')
call set_mult(n_F, 'lin', x_ij = BIO_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(2))
call set_mult(n_F, 'lin', x_ijn = CRPM)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CBIO, name = 'BIO_prop*RBIO')
call set_mult(n_F, 'lin', x_ij = BIO_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(3))
call set_mult(n_F, 'lin', x_ijn = CBIO)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CBIO, name = 'BIO_prop*RHUM')
call set_mult(n_F, 'lin', x_ij = BIO_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(4))
call set_mult(n_F, 'lin', x_ijn = CHUM)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CHUM, name = 'HUM_prop*RDPM')
call set_mult(n_F, 'lin', x_ij = HUM_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(1))
call set_mult(n_F, 'lin', x_ijn = CDPM)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CHUM, name = 'HUM_prop*RRPM')
call set_mult(n_F, 'lin', x_ij = HUM_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(2))
call set_mult(n_F, 'lin', x_ijn = CRPM)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CHUM, name = 'HUM_prop*RBIO')
call set_mult(n_F, 'lin', x_ij = HUM_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(3))
call set_mult(n_F, 'lin', x_ijn = CBIO)
call set_flux(fid = n_F, pid_out = n_Catm, pid_in = n_CHUM, name = 'HUM_prop*RHUM')
call set_mult(n_F, 'lin', x_ij = HUM_prop)
! piecewise of Ft
call set_mult(n_F, 'lin', x = Ft)
call set_mult(n_F, 'step', x_ij = Tsoil(:,:,1), x0 = -5., k = 1.0)
! piecewise of Fs
call set_mult(n_F, 'PieWise2', x_ij = Wsoil(:,:,1), x0_ij = smin(:,:,1), y0_ij = s0(:,:,1), k = 0., &
& a_ij = 0.8/(s0(:,:,1) - smin(:,:,1)), amp = -0.8, c01 = 0.2, &
& c02_ij = 0.2 - 0.8*smin(:,:,1)/(s0(:,:,1) - smin(:,:,1)), c_ij = 1. + 0.8*s0(:,:,1))
call set_mult(n_F, 'lin', y0 = 0.6, k = -0.4, x_ij = veg, x0 = 1.)
call set_mult(n_F, 'lin', x = Fs)
call set_mult(n_F, 'lin', x = ks(4))
call set_mult(n_F, 'lin', x_ijn = CHUM)
......
source/carbon/carbon_models/carbon_model_ROTHC_aux.f90
+ 44
− 70
View file @ c607d60a
......@@ -61,7 +61,7 @@ module carbon_model_rothc_aux
! ------ Initialization ------
real :: cs = 0. !< Сумма пулов
! ----- The topsoil moisture deficit -------
real :: deep = 20. ! глубина слоя, см (@todo сопоставить с переменной Ml)
real :: deep = 30. ! глубина слоя, см (@todo сопоставить с переменной Ml, внедрить в environment_data_station)
real, allocatable :: TSMD (:,:) !< the topsoil moisture deficit
real, allocatable :: SMD (:,:) !< TSMD на актуальной глубине
real, allocatable :: bSMD (:,:) !< Случай TSMD для голой почвы (bare soil)
......@@ -69,9 +69,6 @@ module carbon_model_rothc_aux
real, allocatable :: aTSMD (:,:) !< accumulated topsoil mositure deficit
real, allocatable :: evpar_part(:,:) !< 75% from evaporation
real :: a,b
contains
subroutine carbon_model_init()
......@@ -86,32 +83,33 @@ contains
allocate(BIO_prop(i0:i1,j0:j1))
allocate(HUM_prop(i0:i1,j0:j1))
!! расчёт через осадки и эвотранспирацию
! init of TPMD
allocate(TSMD(i0:i1,j0:j1))
allocate(SMD(i0:i1,j0:j1))
allocate(bSMD(i0:i1,j0:j1))
allocate(moi_bal(i0:i1,j0:j1))
allocate(aTSMD(i0:i1,j0:j1))
allocate(evpar_part(i0:i1,j0:j1))
TSMD = 0.
SMD = 0.
bSMD = 0.
moi_bal = 0.
aTSMD = 0.
evpar_part = 0.
SA = 0.
CO2_prop = 0.
BIO_prop = 0.
HUM_prop = 0.
s0 = 0.
smin = 0.
TSMD(:,:) = -(20. + 1.3*bettar(:,:) - 0.01*bettar(:,:)**2)
SMD(:,:) = (deep/23.)*TSMD(:,:)
bSMD(:,:) = SMD(:,:)/1.8
! allocate(TSMD(i0:i1,j0:j1))
! allocate(SMD(i0:i1,j0:j1))
! allocate(bSMD(i0:i1,j0:j1))
! allocate(moi_bal(i0:i1,j0:j1))
! allocate(aTSMD(i0:i1,j0:j1))
! allocate(evpar_part(i0:i1,j0:j1))
! TSMD = 0.
! SMD = 0.
! bSMD = 0.
! moi_bal = 0.
! aTSMD = 0.
! evpar_part = 0.
! SA = 0.
! CO2_prop = 0.
! BIO_prop = 0.
! HUM_prop = 0.
! s0 = 0.
! smin = 0.
! TSMD(:,:) = -(20. + 1.3*bettar(:,:) - 0.01*bettar(:,:)**2)
! SMD(:,:) = (deep/23.)*TSMD(:,:)
! bSMD(:,:) = SMD(:,:)/1.8
! main init part of Rothc
do i = 1, ml
......@@ -143,50 +141,28 @@ contains
! ---- TPMD_calc_at_cell ----
evpar_part(ii,jj) = 0.75*Evpar(ii,jj)
moi_bal(ii,jj) = pr(ii,jj) - evpar_part(ii,jj)
if (veg(ii,jj) == 1) then
aTSMD(ii,jj) = max(SMD(ii,jj), min(0.,aTSMD(ii,jj) + moi_bal(ii,jj)))
else
aTSMD(ii,jj) = max(min(bSMD(ii,jj),aTSMD(ii,jj)),min(0.,aTSMD(ii,jj) + moi_bal(ii,jj)))
end if
! evpar_part(ii,jj) = 0.75*Evpar(ii,jj)
! moi_bal(ii,jj) = pr(ii,jj) - evpar_part(ii,jj)
! main part of rothc_calc_at_cell
!if (Wsoil(ii,jj,1) >= s0(ii,jj,1)) then
! Fs = 1. - 0.8*(Wsoil(ii,jj,1) - s0(ii,jj,1))
! else if (smin(ii,jj,1) < Wsoil(ii,jj,1) .and. Wsoil(ii,jj,1) < s0(ii,jj,1)) then
! Fs = 0.2 + 0.8*(Wsoil(ii,jj,1) - smin(ii,jj,1))/(s0(ii,jj,1) - smin(ii,jj,1))
! else if (Wsoil(ii,jj,1) <= smin(ii,jj,1)) then
! Fs = 0.2
! if (veg(ii,jj) == 1) then
! aTSMD(ii,jj) = max(SMD(ii,jj), min(0.,aTSMD(ii,jj) + moi_bal(ii,jj)))
! else
! aTSMD(ii,jj) = max(min(bSMD(ii,jj),aTSMD(ii,jj)),min(0.,aTSMD(ii,jj) + moi_bal(ii,jj)))
! end if
!a = Fs
if (aTSMD(ii,jj) > 0.444*SMD(ii,jj)) then
Fs = 1.
else
Fs = 0.2 + (1. - 0.2)*(SMD(ii,jj) - aTSMD(ii,jj))/(SMD(ii,jj) - 0.444*SMD(ii,jj))
end if
! if (aTSMD(ii,jj) > 0.444*SMD(ii,jj)) then
! Fs = 1.
! else
! Fs = 0.2 + (1. - 0.2)*(SMD(ii,jj) - aTSMD(ii,jj))/(SMD(ii,jj) - 0.444*SMD(ii,jj))
! end if
!b = Fs
! main part of rothc_calc_at_cell
if (Tsoil(ii,jj,1) >= -5.) then
Ft = 47.91/(1. + exp(106.06/(Tsoil(ii,jj,1) + 18.27)))
else
Ft = 0
end if
! ---- Initialization ------
Cs = CDPM(ii,jj,1) + CRPM(ii,jj,1) + CBIO(ii,jj,1) + CHUM(ii,jj,1)
print*, 'moi_bal:', moi_bal(ii,jj)
print*, 'pr', pr(ii,jj)
print*, 'Evpar:', Evpar(ii,jj)
print*, 'evpar_part:', evpar_part(ii,jj)
end subroutine
......@@ -208,8 +184,6 @@ contains
open(unit=500, file='results/'//trim(carbon_model_type)//'/'//trim(station_name)//'_'//trim(station_opt)//'.txt', status='unknown')
write(500,*) date_c%timestamp,';',CDPM(:,:,1),';',CRPM(:,:,1),';',CBIO(:,:,1),';',CHUM(:,:,1),';', CIOM(:,:,1)
write(513,*) date_c%timestamp,';',a,';',b
! ---- Initialization ------
!print*, 'finish thing 1, ', CDPM(:,:,1)/cs*100,'%'
!print*, 'finish thing 2, ', CRPM(:,:,1)/cs*100,'%'
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment