I am having a problem, I have a function popmesh1 that calculates a matrix P( , ). I then use sens\_analysis to store a element of that matrix and continue. How do I output P( , ) so that it is kept track of? I keep getting that the matrix is size (0,0). Also, how do I pass the matrix into another function? Sorry to post whole code, I want to concrete and clear, I'm pretty new to MATLAB

```
    function pdemeshpop(varargin)
global par;
global re;
global P;
global par_n;
global a1;
clc
clear

%INIITIALIZE MESH:- Can change time and age for refining of mesh
time=linspace(0,800,4000);
age=linspace(0,time(end),4000);

dt=time(2)-time(1);
dtao=dt;

P=zeros(length(time),9); % State matrix over time.

P1=zeros(length(time),length(age)); % Mesh for population of P1mod.

prodrev=zeros(length(time),length(age));
p1tot=zeros(length(time));
p2tot=zeros(length(time));
f=zeros(length(time));
A_1=zeros(length(time),1);

%Parameters
G=log(2)/30; %This growth rate had been set to nthroot(2,20), but I think it should be log(2)/20 for a doubling time of 20 mins. Units 1/min
R=.75; %Reduction in growth rate due to viral production, range from 0.5-0.75
global A_s; %Number of virus produced each minute from one cell? Units 1/min
A_s = 35; %Source for this?
global re;
re = varargin(1); %Reduction in efficiency of virus production in P1mod
c=1.5e9; %Concentration of cells in saturated culture. Units 1/cm^3 Source: http://bionumbers.hms.harvard.edu/bionumber.aspx?&id=100984&ver=2
K=3e-11; %Adsorption rate. Units cm^3/min. Source: Tzagoloff, H., and D. Pratt. 1964. The initial steps in infection with coliphage M13. Virology 24:372?380.
i = 1; %Is flipping of switch induced or not induced? if i==1 then switch is induced.

if i==1
    S_i=0.5; %S_i is probability that a switch will flip during a timestep in a p1ori cell. Units pure number (a probablility). Ranges from 0 to 1.
elseif i==0
    S_i= 0.005;
end

%IC and BC implementation for the 9 dependent variables <<<10?
P0=zeros(9,1);
P0(1)=100; %Initial concentration of senders. Units: cells/ml.
P0(2)=10000; %Initial concentration of primary receivers. Units: cells/ml.
P0(3)=10000;  %Initial concentration of secondary receivers. Units: cells/ml.

%The loop below covers the initial conditions and BC of t=0,all ages
for i=1:9
    P(1,i)=P0(i); 
end

%Iterative solution

for m=1:length(time)-1 % m is timestep

        %Simplifications
        p1tot(m)=P(m,2)+P(m,4)+P(m,6);
        p2tot(m)=P(m,3)+P(m,5)+P(m,7);
        f(m)=1-(P(m,1)+p1tot(m)+p2tot(m))/c; 

        %Senders
        P(m+1,1)=dt*(P(m,1)*G*R*f(m))+P(m,1);

        %Primary Receivers
        P(m+1,2)=dt*((P(m,2)*G*f(m))-K*(P(m,8)+P(m,9))*P(m,2))+P(m,2); 

        %Secondary Receivers
        P(m+1,3)=dt*((P(m,3)*G*f(m))-K*(P(m,8)+P(m,9))*P(m,3))+P(m,3);

        %Primary Original
        P(m+1,4)=dt*((P(m,4)*G*f(m))+K*P(m,8)*P(m,2)-S_i*P(m,4))+P(m,4);

        %Secondary Original
        P(m+1,5)=dt*((P(m,5)*G*f(m))+K*P(m,8)*P(m,3))+P(m,5);

    for n=1:m
        t=(m-1)*dt; %Why not t=m*dt?
        tao=(n-1)*(dtao);%Determines current age basket
        prodrev(m,n)=rate(t-tao); %Calculates corresponding rate of production of phage (reversed)   

        %Primary Modified
        if n==1
        P1(m+1,n)=dt*(K*P(m,2)*P(m,9)+S_i*P(m,4)); %Left hand side boundary (New cells at age zero)  
        else
        P1(m+1,n)=dt*(-((P1(m,n)-P1(m,n-1))/dtao)+P1(m,n)*G*R*f(m))+P1(m,n);    
        end

    end

        P(m+1,6)=sum(P1(m+1,:)); %phi1mod

        %Secondary Modified
        P(m+1,7)=dt*((P(m,7)*G*f(m))+K*P(m,9)*P(m,3))+P(m,7);

        %Original
        P(m+1,8)=dt*(A_s*P(m,1))+P(m,8);

        if m<2
            A_1(m)=0;
        else
            convolution(m,:)=prodrev(m,:).*P1(m,:); 

           %A_1(m)=dtao*trapz(conv(prod1(m,:), P1(m,:))); 
           A_1(m)=dtao*trapz([convolution(m,:) 0]);
           %A_1 obtained by convolving the discrete vectors of P1 and prod1
           %then finding the area under the curve

        end    

        %Modified
        P(m+1,9)=dt*A_1(m)+P(m,9);

end

P1;

end

function prod=rate(tao)
%Function generates production rate values of the infected cells based on their age 
global A_s;
global re;
a=re*A_s; %Max production rate

ageofcell=tao;
if ageofcell<=10
        prod=0;
    elseif ageofcell<=50
        prod=(a/40)*(ageofcell-10);
    else
        prod=a;

end
end

```

and my other code that calls the above function, and that I want to pass P(time(length),7) to:

```
function sens_analysis

global par;
global re;
global par_n;
global P;
time=linspace(0,800,4000);

pdemeshpop_final_re_sens;
re_0 = re;
par = re;
s_nom_ss = a1;
delta = 0.05;
par_n = par*(1+delta);
pdemeshpop_final_re_sens_par(par_n);  % similar to pdemeshpop_final_re_sens
s_pert_ss = P(length(time),7);
abs_sens = (s_pert_ss - s_nom_ss)/(delta*re_0);
rel_sens = abs_sens*(re_value/s_nom_ss);
end

```

Again, sorry to post whole code, felt it was a necessary evil. The global variables might also be unnecessary. Could be something obvious. I might need to store P first somehow. Can someone please explain this carefully? Thank you!