%% General parameters
rho     =   0.998e6;            %Water density, g/m3
V_      =   0.2155/1000;        %Molal volume sucrose, m3 mol-1
R       =   8.3143*1e-6;        %Gas constant, m3 MPa K-1 mol-1
T       =   293;                %Absolute temperature, K
mu_x    =   1.0019e-9;          %Dynamic visc water, MPa s
Lr      =   5e-8;               %Radial hydraulic conductivity, m MPa-1 s-1
t_switch = 3600*48;

%Phloem
NoP     =   100;                %Number of sieve elements
r_p0    =   8.4e-6;             %Sieve element radius, m
A_ax_p0 =   NoP*(pi*r_p0^2)   ; %Initial axial surface area phloem, m2
ap_p    =   3.82e-12;           %Axial permeability phloem, m2
e_phloem=   30;                 %Elastic modulus, MPa
    %Sucrose
v_load    = 1e-8 ;            %sucrose loading rate, mol/s
v_unload1 = 0.148 ; K_unload1 = 100;
v_unload2 = 0.148 ; K_unload2 = 100;

%Xylem
NoX     =   150;             %Number of xylem conduits
r_x     =   30e-6;           %Xylem radius, m
e_xylem =   1000;            %Elastic modulus, MPa
A_ax_x  =   NoX*(pi*r_x^2);  %Axial surface area xylem, m2
%v_unload1 = 0.01 * 0.5e-10 / (N_un1); K_unload1 = 1000;
%v_unload2 = 0.5e-10   / (N_un2);      K_unload2 = 1000;

%Apoplast
Rmax  = 0.4e-8;    %Removal rate, mol/s
KmRes = 5;                  %Removal rate affinity, mol/m3

%% Stem architecture;
l       =   1;                  %Total plant length, m

%Loading zone (mature leaves)
N_load   = 10; l_load = 0.1;  x_load   = l_load*l/N_load;
A_rad0_load  =   NoP*2*pi*r_p0*x_load ;  %Radial surface area, m2
M_p_load     =   A_ax_p0*x_load*rho;     %Mass phloem, g
M_x_load     =   A_ax_x*x_load*rho;      %Mass xylem, g

N_1      = 10; l_1    = sc.l1;  x_1    = l_1   *l /N_1;
A_rad0_1  =   NoP*2*pi*r_p0*x_1;
M_p_1     =   A_ax_p0*x_1*rho;
M_x_1     =   A_ax_x*x_1*rho;

N_1_un   = 1; l_1_un = 0.1;  x_1_un = l_1_un*l /N_1_un;
A_rad0_1_un  =   NoP*2*pi*r_p0*x_1_un ;
M_p_1_un     =   A_ax_p0*x_1_un*rho;
M_x_1_un     =   A_ax_x*x_1_un*rho;

N_2      = 10; l_2    = sc.l2;  x_2    = l_2   *l /N_2; 
A_rad0_2  =   NoP*2*pi*r_p0*x_2;
M_p_2     =   A_ax_p0*x_2*rho;
M_x_2     =   A_ax_x*x_2*rho;

N_2_un   = 1; l_2_un = 0.1;  x_2_un = l_2_un*l /N_2_un;
A_rad0_2_un  =   NoP*2*pi*r_p0*x_2_un ;
M_p_2_un     =   A_ax_p0*x_2_un*rho;
M_x_2_un     =   A_ax_x*x_2_un*rho;
