Merge pull request #174 from EcoExtreML/fix_163

fix issue #163

src/Root_properties.m

@@ -1,45 +1,36 @@
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-%   Subfunction - Root - Properties         %
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-%% REFERENCES
-function [Rl] = Root_properties(Rl, Ac, rroot, frac, bbx, KT)
-    %%% INPUTS
-    global DeltZ sfactor LAI_msr
-    % BR = 10:1:650; %% [gC /m^2 PFT]
-    % rroot =  0.5*1e-3 ; % 3.3*1e-4 ;%% [0.5-6 *10^-3] [m] root radius
-    %%% OUTPUTS
-    if KT < 2880
+function [Rl] = Root_properties(Rl, Ac, rroot, frac, bbx, KT, DeltZ, sfactor, LAI_msr)
+    %% REFERENCES
+    %{
+        This function is used to calculate the dynamic growth of root
+        This part can refer to Wang et al. (2021) "Intergrated modelling of photosynthesis and transfer of energy, mass, and momentum in the SPAC system"
+        Input:
+            Rl: root length index at last step
+            Ac: the net assimilation of canopy
+            rroot: radius of root
+        Output:
+            Rl: root length index at current step
+    %}
+
+    % this function is used to calculate the root fraction
+    if KT < 2880 % 2880 means the time step when the root stops growing
         fr = 0.3 * 3 * exp(-0.15 * LAI_msr(KT)) / (exp(-0.15 * LAI_msr(KT)) + 2 * sfactor);
         if fr < 0.15
-            fr = 0.15;
+           fr = 0.15;
         end
     else
         fr = 0.001;
     end
-    % if KT<3840
-    % fr=-0.091*LAI(KT)+0.40;
-    % else
-    % fr=0.001;
-    % end
-    % elseif KT<3200
-    % fr=0.05;
-    % else
-    % fr=0.02;
-    % end
+
     DeltZ0 = DeltZ' / 100;
     BR = Ac * fr * 1800 * 12 / 1000000;
-    root_den = 250 * 1000; %% [gDM / m^3] Root density  Jackson et al., 1997
-    R_C = 0.488; %% [gC/gDM] Ratio Carbon-Dry Matter in root   Jackson et al.,  1997
+    root_den = 250 * 1000; % [gDM / m^3] Root density  Jackson et al., 1997
+    R_C = 0.488; % [gC/gDM] Ratio Carbon-Dry Matter in root   Jackson et al.,  1997
     nn = numel(Rl);
+
+    % This is used to simulate the root growth
     if (~isnan(Ac)) || (Ac > 0)
         Rl = Rl .* DeltZ0;
         Delta_Rltot = BR / R_C / root_den / (pi * (rroot^2)); %% %% root length index [m root / m^2 PFT]
-        % for i=1:nn
-        %   if Rl(i)>50000
-        %       frac(i)=0;
-        %    end
-        % end
         if ~isnan(frac)
             frac = frac / sum(sum(frac));
         else
@@ -48,11 +39,6 @@
         Delta_Rl = frac .* bbx * Delta_Rltot;
         Rl = Rl + Delta_Rl;
         Rl = Rl ./ DeltZ0;
-    else
     end
-    % for i=1:nn
-    %    if Rl(i)>50000
-    %       Rl(i)=50000;
-    %    end
-    % end
 end
+

src/ebal.m

@@ -461,11 +461,11 @@
     Rntot_PAR       = LAI * (Fc * Rnh_PAR  + equations.meanleaf(canopy, Rnu_PAR, 'angles_and_layers', Ps)); % net PAR leaves
     aPAR_Cab_eta        = LAI * (Fc * (profiles.etah .* Rnh_PAR) + equations.meanleaf(canopy, profiles.etau .* Rnu_PAR, 'angles_and_layers', Ps));
     % ... green ePAR * relative fluorescence emission efficiency
-    %%%%%%%%%%%%%%%%%%% [Delta_Rltot] = Root_properties(Actot,rroot);
-    %%%%%%%%%%%%%%%%%%% Delta_Rl = fc*Delta_Rltot;
-    %%%%%%%%%%%%%%%%%%% Rl = Rl + Delta_Rl;
-    %%%%%%%%%%%%%%%%%%% Rltot = sum(sum(Rl));
-    %%%%%%%%%%%%%%%%%%% fc = Rl./Rltot;
+    % [Delta_Rltot] = Root_properties(Rl, Ac, rroot, frac, bbx, KT, DeltZ, sfactor, LAI_msr);
+    % Delta_Rl = fc*Delta_Rltot;
+    % Rl = Rl + Delta_Rl;
+    % Rltot = sum(sum(Rl));
+    % fc = Rl./Rltot;
     % sum of soil fluxes and average temperature
     %   (note that averaging temperature is physically not correct...)
     Rnstot          = Fs * Rns;           %                   Net radiation soil