% Example 18-Open: Concentration profile upstream and in reaction zone
%

Pe=1;
Pe1=5;
Da=1;
lam=linspace(-1,1,100);
set (gcf,'Position',[0,0.1,500,300])
for i=1:100
q=(1+4*Da/Pe)^0.5;
A(i)=(1+q)*exp(q*Pe*(1-lam(i))/2);
B(i)=(1-q)*exp(-q*Pe*(1-lam(i))/2);
A0=(1+q)*exp(q*Pe*(1-0)/2);
B0=(1-q)*exp(-q*Pe*(1-0)/2);
C=(1+q)^2*exp(q*Pe/2);
D=(1-q)^2*exp(-q*Pe/2);
Sai(i)=(2*exp(Pe*lam(i)/2).*(A(i)-B(i)))/(C-D);
Sai0=(2*exp(Pe*0/2)*(A0-B0))/(C-D);
Sai1(i)=1-(1-Sai0)*exp((lam(i)*Pe1));
if (lam(i) >=0) 
    Sa(i) = Sai(i);
else
    Sa(i) = Sai1(i);
end 
end
plot(lam(1:100),Sa(1:100),'LineWidth',2)
xlabel('Dimensionless distance, \lambda =z/L')
ylabel('Dimensionless concentration, \psi =C_A/C_{A0}')
axis([-1 1 0 1])
grid on

%