% solve PFR for a given air temperature at v = 0
% for simplified model with constant heat capacities
function taf = solvePFRSimplified(Ta0)
global R
R = 8.31446; % kJ/kgmol/K  also  m3*kPa/kgmol/K
Tref = 298.15;

% Basic data
MWA = 58.08;
MWAir = 28.96;

% Reactor
NoTubes = 1000;
TotalVolume = 2; % m3
VolPerTube = TotalVolume/NoTubes; % m3
TubeID = 26.7e-3; % m
TubeXC = pi*TubeID^2/4; % m2
TubeLength = VolPerTube/TubeXC; % m

% Heat Transfer
U = 400; % kJ/m2/hr/K
A = 4/TubeID; % m2/m3

% Acetone Feed
ReactorMassFeed = 7850; % kg/hr
FAM = ReactorMassFeed/NoTubes; % kg/hr per tube
FAF = FAM/MWA; % kgmol/hr
TF = 1035; % K
P = 162; % kPa
CAF = P/R/TF; % kgmol/m3

% Air Feed
AirFeed = 11088*8; % kg/hr
FaM = AirFeed/NoTubes; % kg/hr per tube
Fa = FaM/MWAir; % kgmol/hr

% Solution span
vspan = linspace(0,VolPerTube,200);

% Initial conditions

y0 = [ FAF ; TF ; Ta0  ];

% anonymous function
PFRanon = @(v,y) PFRAcetoneSimplified(v,y,FAF,CAF,TF,Fa,U,A);

% solve ODEs
[ v, yout ] = ode45(PFRanon,vspan,y0);

% unpack output
Taout = yout(:,3);
n = length(Taout);
taf = Taout(n);
