What is the best way to change the matrix(dMdt in my code) dimension to plot it?

function mu

% 6 ODEs:

% dMsh/dt = kG*Mshc*Mshn/Msh - sigmaSS*klitt*Msh/(1 + KMlitt/Msh);

% dmrt/dt = kG*Mrtc*Mrtn)/Mrt - sigmaSS*klitt*Mrt/(1+KMlitt/Mrt);

% dMshC/dt = kc*Msh/((1 + sigmaSS*Msh/KM)*(1 + sigmaPI*Mshc/(Msh*Jc)))...

% -fc*kG*Mshc*Mshn/Msh - (Mshc/Msh - Mrtc/Mrt)/(pc/Msh.^q + pc/Mrt.^q);

% dmrtC/dt = (Mshc/Msh - Mrtc/Mrt)/((pc/Msh.^q + pc/Mrt.^q) - fc*kG*Mrtc*Mrtn/Mrt;

% dMshN/dt = (Mrtn/Mrt - Mshn/Msh)/(pn/Mrt.^q + pn/Msh.^q) - fc*kG*Mshc*Mshn/Msh;

% dmrtN/dt = kn*Mrt/((1+sigmaSS*Mrt/KM) * (1+sigmaPI*Mrtn/(Mrt*Jc)))...

%- fn*kG*Mrtc*Mrtn/Mrt - (Mrtn/Mrt - Mshn/Msh)/(pn/Mrt.^q + pn/Msh.^q);

%Define Constants

kG = 200; % Growth Rate

klitt = 0.05; % Litter Rate

KMlitt = 0.5; % Litter Parameter

Jc = 0.1; % Inhibition of C assimilation

Jn = 0.01; % Inhibition of N uptake

fc = 0.5; % Fractions of C in structural

fn = 0.025; %% fN 0.025 is better % Fractions of N in structural

pc = 1; % Transport resistance coefficients

pn = 1; % Transport resistance coefficients

q = 1; % Transport resistance scaling parameter

kc = 0.1; % C assimilation parameter

kn = 0.02; % N uptake parameter

KM = 1; % Parameter giving asymtotic values of photosynthesis and N uptake

%Switches

sigmaPI_vec = [0,1]; % Product Inhibition 1(on) 0(off)

sigmaSS = 0; % Steady state growth

for i = 1:2

sigmaPI = sigmaPI_vec(i);

%Initial values defined in problem

Msh0 = 0.8;

Mrt0 = 0.9;

Mshc0 = 0;

Mrtc0 = 0.03;

Mshn0 = 0.04;

Mrtn0 = 0.001;

x0 = [Msh0 , Mrt0 , Mshc0 , Mrtc0 , Mshn0 , Mrtn0 ];

%% First solve from 0 to 400

tspan1 = [0, 400];

x01 = x0;

[t1,x1] = ode15s(@fun_call,tspan1,x01);

mass_old = x1(:,1)+x1(:,2);

%% Second solve from 400 to 500 (reduce shoot by 75%)

x02 = x1(end,:);

% Apply 75% reduction

x02(1) = 0.25*x02(1); % Msh

x02(3) = 0.25*x02(3); % MshC

x02(5) = 0.25*x02(5); % MshN

tspan2 = [400, 500];

[t2 , x2] = ode15s(@fun_call,tspan2,x02);

mass_new = x2(:,1) + x2(:,2);

%% Combine results for plotting

t_all = [t1;t2];

x_all = [x1;x2];

M = [mass_old; mass_new]; % For figure 2_Part A Specific Growth Rate

dMdt = diff(M)./diff(t_all);

%dMdt(end+1) = NaN; % Add length

mu = dMdt./M; % Specific Growth Rate of Plant

figure (2)

hold on

plot(t_all, mu, 'LineWidth', 2)

%axis([380 500 0.00 0.14])

grid on

legend('No Product Inhibition','Product Inhibition');

xlabel('Time')

ylabel('Specific Growth Rate, mu')

title('A, Specific Growth Rate')

end

function u = fun_call(~,x)

u = zeros(6,1);

Q1 = x(3)/x(1) - x(4)/x(2);

Q2 = pc/x(1)^q + pc/x(2)^q;

R1 = x(6)/x(2) - x(5)/x(1);

R2 = pn/x(2)^q + pn/x(1)^q;

S = sigmaSS*klitt;

u(1) = kG*x(3)*x(5)/x(1) - S*x(1)/(1 + KMlitt/x(1));

u(2) = kG*x(4)*x(6)/x(2) - S*x(2)/(1+KMlitt/x(2));

u(3) = kc*x(1)/((1 + sigmaSS*x(1)/KM)*(1 + sigmaPI*x(3)/(x(1)*Jc)))...

-fc*kG*x(3)*x(5)/x(1) - Q1/Q2;

u(4) = Q1/Q2 - fc*kG*x(4)*x(6)/x(2);

u(5) = R1/R2 - fn*kG*x(3)*x(5)/x(1);

u(6) = kn*x(2)/((1 + sigmaSS*x(2)/KM) * (1 + sigmaPI*x(6)/(x(2)*Jn)))...

- fn*kG*x(4)*x(6)/x(2) - R1/R2;

end

end

[AT]

I have attached the correct plot part. Mine is not like this.

Answer 1

[Gowtham]

Hello

I understand you want to change the dimension of the matrix dMdt in order to plot it.

I attempted to run the provided code at my end and found out the dimensions of dMdt and M are different (using 'size' function). This size mismatch between dMdt and M is resulting in an error during the step where mu is defined.

To address this issue, a possible solution is to adjust the size of the M array so that each element in dMdt aligns with a corresponding element in M.

Assuming that the M array represents the initial mass for calculating the change dMdt, one approach is to remove the last element from the M array. This adjustment ensures that the remaining elements in M align with the corresponding elements in dMdt in the specified order. Consequently, when plotting the data, it is also necessary to trim the t_all array in a similar manner to maintain consistency.

Below is the sample code snippet that could be updated in the code to address the issue:

M = M(1:end-1);
mu = dMdt./M; % Specific Growth Rate of Plant

figure (2)
hold on
plot(t_all(1:end-1), mu, 'LineWidth', 2)   

Kindly refer to the following MATLAB documentation for further information on how to use 'size' and 'array indexing'.

• https://www.mathworks.com/help/matlab/ref/size.html
• https://www.mathworks.com/help/matlab/math/array-indexing.html

Hope this helps.

Best Regards,

Gowtham