Skip to content
Snippets Groups Projects
Commit f9cf182e authored by a_medvedev's avatar a_medvedev
Browse files

тест машинной точности

parent 3200361e
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
CCC.f90
+ 22
− 20
View file @ f9cf182e
......@@ -175,11 +175,6 @@ do t = 1, ntime
! вывод данных:
if (hour /= hour_mem) then
hour_mem = hour
endif
if (month /= month_mem) then
print'(" ",i4.4,"-",i2.2,"-",i2.2," ",i2.2,":",i2.2,":",i2.2)', year, month, day, hour, minute, second
print*, adc6(i0,j0,nv), adc8(i0,j0,nv), adc8b(i0,j0,nv)
print*, Cveg(i0,j0), Csoil(i0,j0), Csoilb(i0,j0), Catm(i0,j0)
......@@ -188,18 +183,25 @@ do t = 1, ntime
if (pool_mem(i0,j0,n_Csoilb) > 0) print*, (pool(i0,j0,n_Csoilb) - pool_mem(i0,j0,n_Csoilb))/pool_mem(i0,j0,n_Csoilb)
print*, 'eps = ', (adc6(i0,j0,nv) - Cveg(i0,j0)) / adc6(i0,j0,nv) * 100, '%'
print*, '---'
print*, 'fpsn', fpsn(nv), fpsn_new(i0,j0)/umol2kg
print*, 'frg', frg(nv), frg_new(i0,j0)/umol2kg
print*, 'frmf', frmf(nv), frmf_new(i0,j0)/umol2kg
print*, 'frms', frms(nv), frms_new(i0,j0)/umol2kg
print*, 'frmr', frmr(nv), frmr_new(i0,j0)/umol2kg
print*, 'fmicr', fmicr(nv), fmicr_new(i0,j0)/umol2kg
print*, 'fmicrb', fmicrb(nv), fmicrb_new(i0,j0)/umol2kg
print*, 'fpsn', fpsn(nv), fpsn_new(i0,j0)!/umol2kg
print*, 'frg', frg(nv), frg_new(i0,j0)!/umol2kg
print*, 'frmf', frmf(nv), frmf_new(i0,j0)!/umol2kg
print*, 'frms', frms(nv), frms_new(i0,j0)!/umol2kg
print*, 'frmr', frmr(nv), frmr_new(i0,j0)!/umol2kg
print*, 'fmicr', fmicr(nv), fmicr_new(i0,j0)!/umol2kg
print*, 'fmicrb', fmicrb(nv), fmicrb_new(i0,j0)!/umol2kg
print*, 'dv68', dv68(nv)/umol2kg, dv68_new(i0,j0)/umol2kg
print*, 'ddc6a', ddc6a(nv)/umol2kg/dt, ddc6a_new(i0,j0)/umol2kg
print*, 'ddc6b', ddc6b(nv)/umol2kg/dt, ddc6b_new(i0,j0)/umol2kg
print*, 'ddc8', ddc8(nv)/umol2kg/dt, ddc8_new(i0,j0)/umol2kg
print*, 'ddc6a', ddc6a(nv)/dt/umol2kg, ddc6a_new(i0,j0)/umol2kg
print*, 'ddc6b', ddc6b(nv)/dt/umol2kg, ddc6b_new(i0,j0)/umol2kg
print*, 'ddc8', ddc8(nv)/dt/umol2kg, ddc8_new(i0,j0)/umol2kg
print*
pause
if (hour /= hour_mem) then
hour_mem = hour
endif
if (month /= month_mem) then
month_mem = month
endif
......
core.f90
+ 3
− 2
View file @ f9cf182e
......@@ -5,7 +5,7 @@ implicit none
! ----------------------------------------------- Constant values ------------------------------------------------
! физические константы:
integer, parameter :: ti = 60*60*24 !< number of second per day
real, parameter :: pi = 4.*atan(1.) !< pi
real, parameter :: pi = 3.1415!4.*atan(1.) !< pi
real, parameter :: hPa2Pa = 100. !< гектопаскали в паскали
real, parameter :: umol2kg = 12.*1.e-9 !< микромоли в килограммы углерода
real, parameter :: Tz = 0. !< точка замерзания воды, С
......@@ -227,7 +227,8 @@ contains
laisha(k) = elai(k)*fsha(k)
fpsn(k)=amndf(k)*(psnsun(k)*laisun(k)+psnsha(k)*laisha(k))
f = fpsn(k) * veget_ich(k)*amsq*(1.0-sq(4)) * 12.e-9
! f = fpsn(k) * veget_ich(k)*amsq*(1.0-sq(4)) * umol2kg
f = fpsn(k)
end function
......
inmcm.f90
+ 32
− 17
View file @ f9cf182e
......@@ -42,7 +42,7 @@ integer, parameter :: n3 = 1 !< Number of respiration function
integer, parameter :: n4 = 1 !< Number of litterfall function
! ---
real, parameter :: pi_loc = 4.*atan(1.)
real, parameter :: pi_loc = 3.1415!4.*atan(1.)
real, parameter :: r_earth_loc = 6371000.
integer, parameter :: nv1 = 13 !< число типов поверхности
......@@ -420,10 +420,13 @@ select case (n1)
DO J = j0, j1
DO I = i0, i1
! IF(OLIM(I,J).EQ.1) THEN
DO NVV=1,NV2
SUMS(1)=SUMS(1)+ADC6(I,J,NVV)*DEFR(NVV)*AREA*1.E-12
SUMS(2)=SUMS(2)+ADC8(I,J,NVV)*SOER(NVV)*AREA*1.E-12
END DO
! DO NVV=1,NV2
NVV = 3
! SUMS(1)=SUMS(1)+ADC6(I,J,NVV)*DEFR(NVV)*AREA*1.E-12
! SUMS(2)=SUMS(2)+ADC8(I,J,NVV)*SOER(NVV)*AREA*1.E-12
SUMS(1)=SUMS(1)+DEFR(NVV)*AREA!*1.E-12
SUMS(2)=SUMS(2)+SOER(NVV)*AREA!*1.E-12
! END DO
! END IF
END DO
END DO
......@@ -599,22 +602,28 @@ select case (n1)
defor0b = defor0-defor0a
if (sc6 /= 0.) then
ddc6a(i) = dc6(i) * defr(i) * defor0a * adefr * dt / yrs / sc6
ddc6b(i) = dc6(i) * defr(i) * defor0b * adefr * dt / yrs / sc6
! ddc6a(i) = dc6(i) * defr(i) * defor0a * adefr * dt / yrs / sc6
! ddc6b(i) = dc6(i) * defr(i) * defor0b * adefr * dt / yrs / sc6
ddc6a(i) = defr(i) / sc6 * defor0a * adefr * dt / yrs
ddc6b(i) = defr(i) / sc6 * defor0b * adefr * dt / yrs
else
ddc6a(i) = 0.
ddc6b(i) = 0.
endif
if (sc8 /= 0.) then
ddc8(i) = dc8(i) * soer(i) * defor0 * (1.0 - adefr) * dt / yrs / sc8
! ddc8(i) = dc8(i) * soer(i) * defor0 * (1.0 - adefr) * dt / yrs / sc8
ddc8(i) = soer(i) / sc8 * defor0 * (1.0 - adefr) * dt / yrs
else
ddc8(i) = 0.
endif
ddc6a(i) = ddc6a(i) * area / (4.*pi_loc*r_earth_loc**2)
ddc6b(i) = ddc6b(i) * area / (4.*pi_loc*r_earth_loc**2)
ddc8(i) = ddc8(i) * area / (4.*pi_loc*r_earth_loc**2)
! ddc6a(i) = ddc6a(i) * area / (4.*pi_loc*r_earth_loc**2)
! ddc6b(i) = ddc6b(i) * area / (4.*pi_loc*r_earth_loc**2)
! ddc8(i) = ddc8(i) * area / (4.*pi_loc*r_earth_loc**2)
ddc6a(i) = ddc6a(i) * area / (4.*pi_loc*r_earth_loc**2) * 1.e+12
ddc6b(i) = ddc6b(i) * area / (4.*pi_loc*r_earth_loc**2) * 1.e+12
ddc8(i) = ddc8(i) * area / (4.*pi_loc*r_earth_loc**2) * 1.e+12
case default
ddc6a(i) = 1.
......@@ -649,9 +658,12 @@ select case (n3)
amns = stemb(i)/znam !< Tech variable
amnr = rootb(i)/znam !< Tech variable
frmf(i) = 1.00 * amndf(i) * rmf25(i) * tlai(i) * fnf * tf * rf * btran(i)
frms(i) = 0.35 * amndf(i) * rms25(i) * dc6(i) * amns * tf * rf
frmr(i) = 0.35 * amndf(i) * rmr25(i) * dc6(i) * amnr * tf * rf
! frmf(i) = 1.00 * amndf(i) * rmf25(i) * tlai(i) * fnf * tf * rf * btran(i)
! frms(i) = 0.35 * amndf(i) * rms25(i) * dc6(i) * amns * tf * rf
! frmr(i) = 0.35 * amndf(i) * rmr25(i) * dc6(i) * amnr * tf * rf
frmf(i) = 1.00 * amndf(i) * rmf25(i) * tlai(i) * fnf * rf * btran(i)
frms(i) = 0.35 * amndf(i) * rms25(i) * amns * rf
frmr(i) = 0.35 * amndf(i) * rmr25(i) * amnr * rf
frm(i) = frmf(i) + frms(i) + frmr(i)
! Growth respiration
......@@ -671,8 +683,10 @@ select case (n3)
end if
! fswfst = fsw * fst
fmicr(i) = 1.5 * fsw * fst * amrp(i) * dc8(i)
fmicrb(i) = 1.5 * fsw * fst * amrp(i) * dc8b(i) * cv81b
! fmicr(i) = 1.5 * fsw * fst * amrp(i) * dc8(i)
! fmicrb(i) = 1.5 * fsw * fst * amrp(i) * dc8b(i) * cv81b
fmicr(i) = 1.5 * fsw * amrp(i)
fmicrb(i) = 1.5 * fsw * amrp(i) * cv81b
!end do
case default
......@@ -687,7 +701,8 @@ end select
! ----------------------------------------------------- Litterfall ---------------------------------------------------
select case (n4)
case (1) !Type 1 (INMCM)
dv68(i) = dc6(i) / (al5(i) * yrs)
! dv68(i) = dc6(i) / (al5(i) * yrs)
dv68(i) = 1. / (al5(i) * yrs)
case default
dv68(i) = 1.
end select
......
my_carbon_model.f90
+ 36
− 19
View file @ f9cf182e
......@@ -86,7 +86,8 @@ contains
! dv68(i) = dc6(i) / (al5(i) * yrs)
f = 1
nmult(n,m,f) = 1
call set_fun(n, m, f, 1, 'lin', x = Cveg, k = 1./(al5(nv)*yrs), b = 0.)
! call set_fun(n, m, f, 1, 'lin', x = Cveg, k = 1./(al5(nv)*yrs), b = 0.)
call set_fun(n, m, f, 1, 'const', c = 1./(al5(nv)*yrs))
route(n,m,f) = +1
dv68_new(i0:i1,j0:j1) => flux(:,:,n,m,f)
......@@ -116,10 +117,12 @@ contains
call set_fun(n, m, f, 3, 'const', c = rmf25(nv))
call set_fun(n, m, f, 4, 'const_fict', c = tlai(nv))
call set_fun(n, m, f, 5, 'const', c = min(foln(nv)/max(folnmx(nv),1.e-6),1.))
call set_fun(n, m, f, 6, 'exp', x = Temp, a = 1., k = 0.1*log(arm(nv)), b = Kelvin0-298.16)
! call set_fun(n, m, f, 6, 'exp', x = Temp, a = 1., k = 0.1*log(arm(nv)), b = Kelvin0-298.16)
call set_fun(n, m, f, 6, 'const', c = 1.)
call set_fun(n, m, f, 7, 'step', x = igs_new, x0 = 0., y0 = 0.5, y1 = 1.)
call set_fun(n, m, f, 8, 'lin', x = btran_new, k = 1., b = 0.)
call set_fun(n, m, f, 9, 'const', c = umol2kg * veget_ich(nv)*amsq*(1.0-sq(4)))
! call set_fun(n, m, f, 9, 'const', c = umol2kg * veget_ich(nv)*amsq*(1.0-sq(4)))
call set_fun(n, m, f, 9, 'const', c = 1.)
route(n,m,f) = +1
frmf_new(i0:i1,j0:j1) => flux(:,:,n,m,f)
f = 4
......@@ -127,10 +130,13 @@ contains
call set_fun(n, m, f, 1, 'const', c = 0.35)
call set_fun(n, m, f, 2, 'lin', x = amndf3(:,:,nv), k = 1., b = 0.)
call set_fun(n, m, f, 3, 'const', c = rms25(nv))
call set_fun(n, m, f, 4, 'lin', x = Cveg, k = stemb(nv)/(stemb(nv)+rootb(nv)), b = 0.)
call set_fun(n, m, f, 5, 'exp', x = Temp, a = 1., k = 0.1*log(arm(nv)), b = Kelvin0-298.16)
! call set_fun(n, m, f, 4, 'lin', x = Cveg, k = stemb(nv)/(stemb(nv)+rootb(nv)), b = 0.)
call set_fun(n, m, f, 4, 'const', c = stemb(nv)/(stemb(nv)+rootb(nv)))
! call set_fun(n, m, f, 5, 'exp', x = Temp, a = 1., k = 0.1*log(arm(nv)), b = Kelvin0-298.16)
call set_fun(n, m, f, 5, 'const', c = 1.)
call set_fun(n, m, f, 6, 'step', x = igs_new, x0 = 0., y0 = 0.5, y1 = 1.)
call set_fun(n, m, f, 7, 'const', c = umol2kg)
! call set_fun(n, m, f, 7, 'const', c = umol2kg)
call set_fun(n, m, f, 7, 'const', c = 1.)
route(n,m,f) = +1
frms_new(i0:i1,j0:j1) => flux(:,:,n,m,f)
f = 5
......@@ -138,10 +144,13 @@ contains
call set_fun(n, m, f, 1, 'const', c = 0.35)
call set_fun(n, m, f, 2, 'lin', x = amndf3(:,:,nv), k = 1., b = 0.)
call set_fun(n, m, f, 3, 'const', c = rmr25(nv))
call set_fun(n, m, f, 4, 'lin', x = Cveg, k = rootb(nv)/(stemb(nv)+rootb(nv)), b = 0.)
call set_fun(n, m, f, 5, 'exp', x = Temp, a = 1., k = 0.1*log(arm(nv)), b = Kelvin0-298.16)
! call set_fun(n, m, f, 4, 'lin', x = Cveg, k = rootb(nv)/(stemb(nv)+rootb(nv)), b = 0.)
call set_fun(n, m, f, 4, 'const', c = rootb(nv)/(stemb(nv)+rootb(nv)))
! call set_fun(n, m, f, 5, 'exp', x = Temp, a = 1., k = 0.1*log(arm(nv)), b = Kelvin0-298.16)
call set_fun(n, m, f, 5, 'const', c = 1.)
call set_fun(n, m, f, 6, 'step', x = igs_new, x0 = 0., y0 = 0.5, y1 = 1.)
call set_fun(n, m, f, 7, 'const', c = umol2kg)
! call set_fun(n, m, f, 7, 'const', c = umol2kg)
call set_fun(n, m, f, 7, 'const', c = 1.)
route(n,m,f) = +1
frmr_new(i0:i1,j0:j1) => flux(:,:,n,m,f)
......@@ -155,11 +164,14 @@ contains
call set_fun(n, m, f, 1, 'const', c = 1.5)
call set_fun(n, m, f, 2, 'mm', x = thetaSoilTop, k = 1., b = 0.20)
call set_fun(n, m, f, 3, 'hyp', x = thetaSoilTop, k = 0.23, b = 0.23)
call set_fun(n, m, f, 4, 'exp', x = Tsoil(:,:,ms+1), a = 1., k = 0.1*log(2.), b = Kelvin0-283.16)
my_label(n,m,f,4) = 'fst'
! call set_fun(n, m, f, 4, 'exp', x = Tsoil(:,:,ms+1), a = 1., k = 0.1*log(2.), b = Kelvin0-283.16)
! my_label(n,m,f,4) = 'fst'
call set_fun(n, m, f, 4, 'const', c = 1.)
call set_fun(n, m, f, 5, 'const', c = amrp(nv))
call set_fun(n, m, f, 6, 'lin', x = Csoil, k = 1., b = 0.)
call set_fun(n, m, f, 7, 'const', c = umol2kg)
! call set_fun(n, m, f, 6, 'lin', x = Csoil, k = 1., b = 0.)
call set_fun(n, m, f, 6, 'const', c = 1.)
! call set_fun(n, m, f, 7, 'const', c = umol2kg)
call set_fun(n, m, f, 7, 'const', c = 1.)
route(n,m,f) = +1
fmicr_new(i0:i1,j0:j1) => flux(:,:,n,m,f)
......@@ -173,12 +185,15 @@ contains
call set_fun(n, m, f, 1, 'const', c = 1.5)
call set_fun(n, m, f, 2, 'mm', x = thetaSoilTop, k = 1., b = 0.20)
call set_fun(n, m, f, 3, 'hyp', x = thetaSoilTop, k = 0.23, b = 0.23)
call set_fun(n, m, f, 4, 'exp', x = Tsoil(:,:,ms+1), a = 1., k = 0.1*log(2.), b = Kelvin0-283.16)
my_label(n,m,f,4) = 'fst'
! call set_fun(n, m, f, 4, 'exp', x = Tsoil(:,:,ms+1), a = 1., k = 0.1*log(2.), b = Kelvin0-283.16)
! my_label(n,m,f,4) = 'fst'
call set_fun(n, m, f, 4, 'const', c = 1.)
call set_fun(n, m, f, 5, 'const', c = amrp(nv))
call set_fun(n, m, f, 6, 'lin', x = Csoilb, k = 1., b = 0.)
! call set_fun(n, m, f, 6, 'lin', x = Csoilb, k = 1., b = 0.)
call set_fun(n, m, f, 6, 'const', c = 1.)
call set_fun(n, m, f, 7, 'const', c = cv81b)
call set_fun(n, m, f, 8, 'const', c = umol2kg)
! call set_fun(n, m, f, 8, 'const', c = umol2kg)
call set_fun(n, m, f, 8, 'const', c = 1.)
route(n,m,f) = +1
fmicrb_new(i0:i1,j0:j1) => flux(:,:,n,m,f)
......@@ -288,8 +303,10 @@ contains
endif
! --- sc6, sc8 ---
dfr_weight(ii,jj) = pool(ii,jj,n_Cveg) * defr(nv) * area
ers_weight(ii,jj) = pool(ii,jj,n_Csoil) * soer(nv) * area
! dfr_weight(ii,jj) = pool(ii,jj,n_Cveg) * defr(nv) * area
! ers_weight(ii,jj) = pool(ii,jj,n_Csoil) * soer(nv) * area
dfr_weight(ii,jj) = defr(nv) * area
ers_weight(ii,jj) = soer(nv) * area
enddo
enddo
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment