How to plot values of a matrix?

I have attached my code below. I want to plot the values of 'steadystate' during the while loop as 4 different lines on the same plot. How do I do that?

while abs(steadystate-concoutGuess)>0.0000001
      bigmatrix=alt_bigmatrix;
      concoutGuess=steadystate;
      concout=steadystate;
      fA= ((flowinA*concinA-flowout*concout(1))/Vol) - k1*concout(1)*concout(2);
      fB= ((flowinB*concinB-flowout*concout(2))/Vol) - k1*concout(1)*concout(2) - k2*concout(2)*concout(3);
      fC= ((-flowout*concout(3))/Vol) + k1*concout(1)*concout(2) - k2*concout(2)*concout(3);
      fD= ((-flowout*concout(4))/Vol) + k1*concout(2)*concout(3);
      molflow = [ fA; fB; fC; fD ];
      bigmatrix= subs(bigmatrix);
      differential=inv(bigmatrix);
      differential=vpa(differential,8);
      molflow=vpa(molflow,8);
      steadystate = concoutGuess - (differential*molflow);
      steadystate = vpa(steadystate,8)
end

Answer 1

[Image Analyst]

Have you tried the plot() function?

What did you intialize steadystate and concoutGuess with, so we can try your code?

To post formatted code, read this

[Remston Martis]

Oh, I understand what you meant by your question. I initialized concoutGuess as [0.1;0.1;0.1;0.1] steadystate was then calculated using concoutGuess and the first iteration gave steadystate = [0.1504734;-0.049281077; 0.14977213;0.34987724]

[Image Analyst]

This does not work:

concoutGuess = [0.1;0.1;0.1;0.1]
while abs(steadystate-concoutGuess)>0.0000001
	bigmatrix=alt_bigmatrix;
	concoutGuess=steadystate;
	concout=steadystate;
	fA= ((flowinA*concinA-flowout*concout(1))/Vol) - k1*concout(1)*concout(2);
	fB= ((flowinB*concinB-flowout*concout(2))/Vol) - k1*concout(1)*concout(2) - k2*concout(2)*concout(3);
	fC= ((-flowout*concout(3))/Vol) + k1*concout(1)*concout(2) - k2*concout(2)*concout(3);
	fD= ((-flowout*concout(4))/Vol) + k1*concout(2)*concout(3);
	molflow = [ fA; fB; fC; fD ];
	bigmatrix= subs(bigmatrix);
	differential=inv(bigmatrix);
	differential=vpa(differential,8);
	molflow=vpa(molflow,8);
	steadystate = concoutGuess - (differential*molflow);
	steadystate = vpa(steadystate,8)
end

because you didn't give the intial values for steadystate. Again, what did you initialize steadystate to?

[Remston Martis]

concoutGuess is a user input but steadystate is not a user input.

This is my code:

syms concoutA concoutB concoutC concoutD fA fB fC fD %Declaring variables in symbolic to perform operations i.e. jacobian and differentiation
syms flowinA flowinB flowinC flowinD Vol concinA concinB concinC concinD 
syms flowout k1 k2
fA= ((flowinA*concinA-flowout*concoutA)/Vol) - k1*concoutA*concoutB; %Equations derived from steady state mass balance of a CSTR
fB= ((flowinB*concinB-flowout*concoutB)/Vol) - k1*concoutA*concoutB - k2*concoutB*concoutC; %All the variables are symbolic for now
fC= ((-flowout*concoutC)/Vol) + k1*concoutA*concoutB - k2*concoutB*concoutC;
fD= ((-flowout*concoutD)/Vol) + k1*concoutB*concoutC;
molflow = [ fA; fB; fC; fD ]; %molar flow rate as a vector
concout= [concoutA,concoutB, concoutC, concoutD];
bigmatrix = jacobian(molflow,concout); %using jacobian inbuilt function which allows molflow to be differentiated by concout w.r.t. each of the chemical species
alt_bigmatrix=bigmatrix; %storing bigmatrix as symbolic to simplify Newton Rhapson procedures later on
concinA = input('Enter the concentration of Chemical A flowing into the reactor [kmol/m^3] : '); %User inputs/initial conditions necessary to solve the problem
concinB = input('Enter the concentration of Chemical B flowing into the reactor [kmol/m^3] : ');
k1=input('Enter desired rate of reaction 1 A+B->C [m^3/kmol*sec] : ');
k2=input('Enter desired rate of reaction 2 B+C->D [m^3/kmol*sec] : ');
flowinA=input('Enter the flowrate of Chemical A into the reactor [m^3/sec] : ');
flowinB=input('Enter the flowrate of Chemical B into the reactor [m^3/sec] : ');
Vol=input('Enter desired volume of reactor vessel [m^3] : ');
%IMPORTANT!!:For input of 'concoutGuess', kindly make sure to enter vector in correct
%format: "[A;B;C;D]" as a column vector and not a row vector else
%calculations will fail
concoutGuess= input('Enter initial guess for final concentrations of chemical species A, B, C, D respectively in the format "[A; B; C; D]" kmol/m^3 : ');
concoutA=concoutGuess(1);
concoutB=concoutGuess(2);
concoutC=concoutGuess(3);
concoutD=concoutGuess(4);
flowout=flowinA+flowinB; %total volumetric flow out of CSTR
fA= ((flowinA*concinA-flowout*concoutA)/Vol) - k1*concoutA*concoutB;
fB= ((flowinB*concinB-flowout*concoutB)/Vol) - k1*concoutA*concoutB - k2*concoutB*concoutC;
fC= ((-flowout*concoutC)/Vol) + k1*concoutA*concoutB - k2*concoutB*concoutC;
fD= ((-flowout*concoutD)/Vol) + k1*concoutB*concoutC;
molflow = [ fA; fB; fC; fD ];
bigmatrix= subs(bigmatrix); %substitute and perform arithmetic operations on symbolic matrix with user inputs
differential=inv(bigmatrix); %inverse of bigmatrix required for Newton Rhapson
%At this point the matrix is still in symbolic form but sustituted with the user
%inputs and inverted
differential=vpa(differential,8); %Rounding to specified significant figures
molflow=vpa(molflow,8);%Rounding to specified significant figures
steadystate = concoutGuess - (differential*molflow); %Newton Rhapson Equation to caltulate steady state concentrations
steadystate = vpa(steadystate,8)
while abs(steadystate-concoutGuess)>0.0000001 %Convergence error criteria for the Newton Rhapson loop
bigmatrix=alt_bigmatrix; %Using the stored symbolic matrix for the loop calculations
concoutGuess=steadystate; %Updating values
concout=steadystate;
fA= ((flowinA*concinA-flowout*concout(1))/Vol) - k1*concout(1)*concout(2); %Recalcualting values similar to previous calculations
fB= ((flowinB*concinB-flowout*concout(2))/Vol) - k1*concout(1)*concout(2) - k2*concout(2)*concout(3);
fC= ((-flowout*concout(3))/Vol) + k1*concout(1)*concout(2) - k2*concout(2)*concout(3);
fD= ((-flowout*concout(4))/Vol) + k1*concout(2)*concout(3);
molflow = [ fA; fB; fC; fD ];
bigmatrix= subs(bigmatrix); %Sustituting updated values and calculating in symbolic form
differential=inv(bigmatrix);
differential=vpa(differential,8);%Rounding to specified significant figures
molflow=vpa(molflow,8);
steadystate = concoutGuess - (differential*molflow); %Newton Rhapson formula for steady state
steadystate = vpa(steadystate,8)
end

Hope this helps.

[Image Analyst]

Sorry, I don't have the Symbolic Math toolbox. I've added it to the Products list on the upper right so others will know.