Commit 37847675 authored by Koren, Gerbrand's avatar Koren, Gerbrand
Browse files

Upload New File

parent 3996610c
def leaf_model_D47(key1,var1,Cin,d13Cin,d18Oin,D47in,An,Aleaf,key2,var2,gm13C,d18O_leaf):
"""
Name: leaf_model_D47
Author: Gerbrand Koren (Wageningen University)
Date: June 29th, 2017
----------------------------------------------------------------------------------
Description:
This model calculates the outgoing D47 signature for a leaf in the cuvette experiment as described by Hofmann et al.
List of input variables:
- Experimental conditions
*airflow rate of air entering cuvette [L/min]
*drawdown CO2 drawdown [ppm]
- Ingoing air characteristics
Cin ingoing CO2 mixing ratio [ppm]
d13Cin ingoing d13C signature [permil VPDB]
d18Oin ingoing d18O signature [permil VSMOW]
D47in ingoing D47 signature [permil]
- Leaf properties
An net assimilation rate [umol/(m2s)]
Aleaf leaf area [cm2]
**gs stomatal conductance [mol/(m2s)]
**Cm_Ca CO2-in-mesophyll to CO2-in-amtosphere ratio [-]
gm13C conductance between intercellular space and chloroplast [mol/(m2s)]
- Leaf water characteristics
d18O_leaf leaf water d18O signature [permil VSMOW]
Notes:
* The user can call the function with the value of aiflow or drawdown for var1. The user is
required to specifiy which input parameter is given by using the key1 variable.
** The user can call the function with the value of gs or Cm_Ca for var2. The user is required
to specifiy which input parameter is given by using the key2 variable.
----------------------------------------------------------------------------------
"""
#################################################
# #
# INITIALIZATION #
# #
#################################################
# -- Import statement
import numpy as np
# -- Output settings
print_diag_An_13C = False
print_diag_d13Ca = False
print_diag_d18Oa = False
print_diag_D47a = False
print_table = False
#################################################
# #
# MODEL SETTINGS #
# #
#################################################
# -- Test key1
if (key1=='airflow'):
airflow = var1 # L/min
elif (key1=='drawdown'):
drawdown = var1 # ppm
else:
print 'Error: key1 not recognized'
return
# -- Test key2
if (key2=='gs'):
gs = var2 # mol/(m2s)
elif (key2=='Cm_Ca'):
Cm_Ca = var2 # -
else:
print 'Error: key2 not recognized'
return
#################################################
# #
# MODEL PARAMETERS #
# #
#################################################
# -- Isotope reference
R13_VPDB = 0.0112372 # -
R18_VSMOW = 0.0020052 # -
# -- Physical constants
Vm = 24.5 # L/mol
# -- Leaf water signature
D47_leaf = 0.95 # permil
# -- 13C fractionations
a_13C = 4.4 # permil
es_13C = 1.1 # permil
al_13C = 0.7 # permil
b_13C = 28.2 # permil
# -- 18O fractionations
a18O_diff = 8.8 # permil (Farquhar et al., 1993)
a18O_liq = 0.8 # permil (Farquhar et al., 1993)
a18O_eq = 42.2 # permil (Brenninkmeijer et al., 1983)
# -- 13C18OO fractionations
a47_diff = 12.6 # permil
a47_liq = 1.0 # permil
#################################################
# #
# CALCULATION #
# #
#################################################
# ==== CO2 MIXING RATIOS
# -- CO2 assimilation
CO2_assimilation = An*1.0e-6*Aleaf*1.0e-4 # molCO2/s
# -- Test key1
if (key1=='airflow'):
# -- CO2 fluxes
CO2_inflow = airflow/60.0/Vm*Cin*1.0e-6 # molCO2/s
CO2_outflow = CO2_inflow-CO2_assimilation # molCO2/s
# -- Cuvette condition
Ca = CO2_outflow/(airflow/60.0)*Vm*1.0e6 # ppm
# -- Test key1
if (key1=='drawdown'):
# -- Airflow
airflow = CO2_assimilation*60.0*Vm/(drawdown*1.0e-6) # L/min
# -- Cuvette condition
Ca = Cin-drawdown # ppm
# -- Mesophyll conductance
gm18O = 3.0*gm13C # mol/(m2s)
# -- Test key2
if (key2=='Cm_Ca'):
# -- Mesophyll CO2 mixing ratio
Cm = Cm_Ca*Ca # ppm
# -- Combined stomatal-mesophyll conductance
gtot = An/(Ca-Cm) # mol/(m2s)
# -- Stomatal conductance
gs = 1.0/(1.0/gtot-1.0/gm18O) # mol/(m2s)
# -- Intercellular CO2 mixing ratio
Ci = Ca-An/gs # ppm
# -- Test key2
if (key2=='gs'):
# -- Mesophyll CO2 mixing ratio
Cm = Ci-An/gm18O # ppm
# -- Chloroplast CO2 mixing ratio
Cc = Ci-An/gm13C # ppm
# ==== 13CO2 MIXING RATIOS
# -- Assimilation factor
k = An/Cc # mol/(m2s)
k_13C = k*(1.0-b_13C/1000.0) # mol/(m2s)
# -- 13CO2 conductances
gs_13C = gs*(1.0-a_13C/1000.0) # mol/(m2s)
gm18O_13C = gm18O*(1.0-es_13C/1000.0) # mol/(m2s)
gm13C_13C = gm13C*(1.0-(es_13C+al_13C)/1000.0) # mol/(m2s)
# -- 13CO2 mixing ratio
R13in = (d13Cin/1000.0+1.0)*R13_VPDB # -
Cin_13C = R13in*Cin # ppm
# -- 13CO2 signature
d13Ca = d13Cin # permil VPDB
max_jj = 30 # -
relax_jj = 0.5 # -
for jj in np.arange(max_jj):
# -- 13CO2 mixing ratio
R13a = (d13Ca/1000.0+1.0)*R13_VPDB # -
Ca_13C = R13a*Ca # ppm
# -- 13CO2 assimilation rate
An_13C = An*R13a # umol/(m2s)
max_ii = 30 # -
relax_ii = 0.3 # -
for ii in np.arange(max_ii):
# -- Calculate 13CO2 mixing ratios
Ci_13C = Ca_13C-An_13C/gs_13C # ppm
Cm_13C = Ci_13C-An_13C/gm18O_13C # ppm
Cc_13C = Ci_13C-An_13C/gm13C_13C # ppm
# -- Calculate new 13CO2 assimliation rate
An_13C_diff = k_13C*Cc_13C-An_13C # umol/(m2s)
An_13C = An_13C+relax_ii*An_13C_diff # umol/(m2s)
# -- Calculate 13CO2 ratios
R13i = Ci_13C/Ci # -
R13m = Cm_13C/Cm # -
R13c = Cc_13C/Cc # -
# -- Calcualte 13CO2 signatures
d13Ci = 1000.0*(R13i/R13_VPDB-1.0) # permil VPDB
d13Cm = 1000.0*(R13m/R13_VPDB-1.0) # permil VPDB
d13Cc = 1000.0*(R13c/R13_VPDB-1.0) # permil VPDB
# -- Print diagnostics
print_1 = print_diag_d13Ca and ii==(max_ii-1)
print_2 = print_diag_An_13C and jj==(max_jj-1)
print_3 = print_diag_An_13C and print_diag_d13Ca
if (print_1 or print_2 or print_3):
print '---------------------------------------'
print 'd13Ca iteration', jj+1
print 'An_13C iteration', ii+1
print '---------------------------------------'
print 'An_13C ', An_13C
print 'Ca_13C ', Ca_13C
print 'Ci_13C ', Ci_13C
print 'Cm_13C ', Cm_13C
print 'Cc_13C ', Cc_13C
print 'R13a ', R13a
print 'R13i ', R13i
print 'R13m ', R13m
print 'R13c ', R13c
print 'd13Ca ', d13Ca
print 'd13Ci ', d13Ci
print 'd13Cm ', d13Cm
print 'd13Cc ', d13Cc
print '---------------------------------------'
print
# -- Clear loop index
del ii
# -- Calculate 13CO2 fluxes
inflow_13CO2 = airflow/60.0/Vm*Cin_13C*1.0e-6 # mol13CO2/s
assimilation_13CO2 = An_13C*1.0e-6*Aleaf*1.0e-4 # mol13CO2/s
outflow_13CO2 = inflow_13CO2-assimilation_13CO2 # mol13CO2/s
# -- Copy previous d13Ca signature
d13Ca_prev = d13Ca # permil VPDB
# -- Calculate new 13CO2 signature
Ca_13C = outflow_13CO2/(airflow/60.0)*Vm*1.0e6 # ppm
R13a = Ca_13C/Ca # -
d13Ca_diff = 1000.0*(R13a/R13_VPDB-1.0)-d13Ca # permil VPDB
d13Ca = d13Ca+relax_jj*d13Ca_diff # permil VPDB
# -- Clear loop index
del jj
# -- Determine convergence
conv_An_13C = np.abs(An_13C_diff/An_13C) # -
conv_d13Ca = np.abs((d13Ca-d13Ca_prev)/d13Ca_prev) # -
# -- Test convergence
if (conv_An_13C > 1.0e-6):
print 'WARNING: Convergence in 13CO2 assimilation less than 1.0e-6'
print 'conv ', conv_An_13C
print
# -- Test convergence
if (conv_d13Ca > 1.0e-6):
print 'WARNING: Convergence in outgoing d13C signature less than 1.0e-6'
print 'conv ', conv_d13Ca
print
# ==== C18OO MIXING RATIOS
# -- Mesophyll C18OO mixing ratio
d18Om = d18O_leaf+a18O_eq # permil VSMOW
R18m = (d18Om/1000.0+1.0)*R18_VSMOW # -
Cm_18O = 2*R18m*Cm # ppm
# -- C18OO conductances
gs_18O = gs*(1.0-a18O_diff/1000.0) # mol/(m2s)
gm18O_18O = gm18O*(1.0-a18O_liq/1000.0) # mol/(m2s)
# -- Ingoing C18OO mixing ratio
R18in = (d18Oin/1000.0+1.0)*R18_VSMOW # -
Cin_18O = 2*R18in*Cin # ppm
# -- C18OO signature
d18Oa = d18Oin # permil VSMOW
max_pp = 100 # -
relax_pp = 0.5 # -
for pp in np.arange(max_pp):
# -- Atmospheric C18OO mixing ratio
R18a = (d18Oa/1000.0+1.0)*R18_VSMOW # -
Ca_18O = 2*R18a*Ca # ppm
# -- Intercellular d18O signature
Ci_18O = (Ca_18O*gs_18O+Cm_18O*gm18O_18O)/(gs_18O+gm18O_18O) # ppm
R18i = 0.5*Ci_18O/Ci # -
d18Oi = 1000.0*(R18i/R18_VSMOW-1.0) # permil VSMOW
# -- Calculate C18OO fluxes
inflow_C18OO = airflow/60.0/Vm*Cin_18O*1.0e-6 # molC18OO/s
leafuptake_C18OO = (Ca_18O-Ci_18O)*gs_18O*1.0e-6*Aleaf*1.0e-4 # molC18OO/s
outflow_C18OO = inflow_C18OO-leafuptake_C18OO # molC18OO/s
# -- Copy previous d18Oa signature
d18Oa_prev = d18Oa # permil VSMOW
# -- Calculate new d18Oa signature
Ca_18O = outflow_C18OO/(airflow/60.0)*Vm*1.0e6 # ppm
R18a = 0.5*Ca_18O/Ca # -
d18Oa_diff = 1000.0*(R18a/R18_VSMOW-1.0)-d18Oa # permil VSMOW
d18Oa = d18Oa+relax_pp*d18Oa_diff # permil VSMOW
# -- Print diagnostics
if (print_diag_d18Oa):
print '---------------------------------------'
print 'd18Oa iteration ', pp+1
print '---------------------------------------'
print 'Ca_18O ', Ca_18O
print 'Ci_18O ', Ci_18O
print 'Cm_18O ', Cm_18O
print 'R18a ', R18a
print 'R18i ', R18i
print 'R18m ', R18m
print 'd18Oa ', d18Oa
print 'd18Oi ', d18Oi
print 'd18Om ', d18Om
print '---------------------------------------'
print
# -- Clear loop index
del pp
# -- Determine convergence
conv_d18Oa = np.abs((d18Oa-d18Oa_prev)/d18Oa_prev) # -
# -- Test convergence
if (conv_d18Oa > 1.0e-6):
print 'WARNING: Convergence in outgoing d18O signature less than 1.0e-6'
print 'conv ', conv_d18Oa
print
# ==== 13C18OO MIXING RATIOS
# -- Mesophyll 13C18OO mixing ratio
D47m = D47_leaf # permil
R47m_stoch = 2*R13m*R18m # -
R47m = (D47m/1000.0+1.0)*R47m_stoch # -
Cm_13C18OO = R47m*Cm # ppm
# -- 13C18OO conductances
gs_13C18OO = gs*(1.0-a47_diff/1000.0) # mol/(m2s)
gm18O_13C18OO = gm18O*(1.0-a47_liq/1000.0) # mol/(m2s)
# -- Ingoing 13C18OO mixing ratio
R47in_stoch = 2*R13in*R18in # -
R47in = (D47in/1000.0+1.0)*R47in_stoch # -
Cin_13C18OO = R47in*Cin # ppm
# -- Atmospheric 13C18OO signature
D47a = D47in # permil VSMOW
max_kk = 100 # -
relax_kk = 0.5 # -
for kk in np.arange(max_kk):
# -- Atmospheric 13C18OO mixing ratio
R47a_stoch = 2*R13a*R18a # -
R47a = (D47a/1000.0+1.0)*R47a_stoch # -
Ca_13C18OO = R47a*Ca # ppm
# -- Intercellular D47 signature
Ci_13C18OO = (Ca_13C18OO*gs_13C18OO+Cm_13C18OO*gm18O_13C18OO)/(gs_13C18OO+gm18O_13C18OO) # ppm
R47i = Ci_13C18OO/Ci # -
R47i_stoch = 2*R13i*R18i # -
D47i = 1000.0*(R47i/R47i_stoch-1.0) # permil
# -- Calculate 13C18OO fluxes
inflow_13C18OO = airflow/60.0/Vm*Cin_13C18OO*1.0e-6 # mol13C18OO/s
leafuptake_13C18OO = (Ca_13C18OO-Ci_13C18OO)*gs_13C18OO*1.0e-6*Aleaf*1.0e-4 # mol13C18OO/s
outflow_13C18OO = inflow_13C18OO-leafuptake_13C18OO # mol13C18OO/s
# -- Copy previous D47a signature
D47a_prev = D47a # permil
# -- Calculate new D47a signature
Ca_13C18OO = outflow_13C18OO/(airflow/60.0)*Vm*1.0e6 # ppm
R47a = Ca_13C18OO/Ca # -
D47a_diff = 1000.0*(R47a/R47a_stoch-1.0)-D47a # permil
D47a = D47a+relax_kk*D47a_diff # permil VSMOW
# -- Print diagnostics
if (print_diag_D47a):
print '---------------------------------------'
print 'D47a iteration ', kk+1
print '---------------------------------------'
print 'Ca_13C18OO ', Ca_13C18OO
print 'Ci_13C18OO ', Ci_13C18OO
print 'Cm_13C18OO ', Cm_13C18OO
print 'R47a ', R47a
print 'R47i ', R47i
print 'R47m ', R47m
print 'D47a ', D47a
print 'D47i ', D47i
print 'D47m ', D47m
print '---------------------------------------'
print
# -- Clear loop index
del kk
# -- Determine convergence
conv_D47a = np.abs((D47a-D47a_prev)/D47a_prev) # -
# -- Test convergence
if (conv_D47a > 1.0e-6):
print 'WARNING: Convergence in outgoing D47 signature less than 1.0e-6'
print 'conv ', conv_D47a
print
#################################################
# #
# RESULTS #
# #
#################################################
# -- Test print settings
if print_table:
# -- Print paramters
print '================================================================'
print 'CO2 mixing ratios'
print ' atmosphere \t\t\t Ca \t %.1f \t ppm' %(Ca)
print ' intercellular \t\t Ci \t %.1f \t ppm' %(Ci)
print ' mesophyll \t\t\t Cm \t %.1f \t ppm' %(Cm)
print ' chloroplast \t\t\t Cc \t %.1f \t ppm' %(Cc)
print
print 'd13C signatures'
print ' atmosphere \t\t\t d13Ca \t %.1f \t permil VPDB' %(d13Ca)
print ' intercellular \t\t d13Ci \t %.1f \t permil VPDB' %(d13Ci)
print ' mesophyll \t\t\t d13Cm \t %.1f \t permil VPDB' %(d13Cm)
print ' chloroplast \t\t\t d13Cc \t %.1f \t permil VPDB' %(d13Cc)
print
print 'd18O signatures'
print ' atmosphere \t\t\t d18Oa \t %.1f \t permil VSMOW' %(d18Oa)
print ' intercellular \t\t d18Oi \t %.1f \t permil VSMOW' %(d18Oi)
print ' mesophyll \t\t\t d18Om \t %.1f \t permil VSMOW' %(d18Om)
print
print 'D47 signature'
print ' atmosphere \t\t\t D47a \t %.2f \t permil' %(D47a)
print ' intercellular \t\t D47i \t %.2f \t permil' %(D47i)
print ' mesophyll \t\t\t D47m \t %.2f \t permil' %(D47m)
print '================================================================'
print
# -- Return outgoing D47 signature
return D47a
\ No newline at end of file
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment