Where is this error coming from?

I have relatively limited MATLAB experience so please, bare with me for this question. When I try and run the following code:

if true
function [y1, y2, y3, x1, x2, x3] = flash1(T, P, z1, z2, z3)
T = 250;
P=1520;
z1=0.5;
z2=0.3;
z3=0.2;

%Inlet Compositions - 1) Pentane, 2) Hexane, 3) Cyclohexane
%Determining Pisat
%Antoine Coefficients
A1=6.87362;
B1=1075.78;
C1=233.205;
A2=6.87024;
B2=1168.72;
C2=224.21;
A3=6.8413;
B3=1201.531;
C3=222.647;
%Calculating saturation pressure of each component, Temp is in Celsius
P1s=A1-B1/(T+C1);
P2s=A2-B2/(T+C2);
P3s=A3-B3/(T+C3);
% Calculating the equilibrium coefficient for each component
K1 = P1s/P;
K2 = P2s/P;
K3 = P3s/P;

% Rachford Rice Equation
S = solve('z1*(K1-1)/(q*(K1-1)+1)+z2*(K2-1)/(q*(K2-1)+1)+z3*(K3-1)/(q*(K3-1)+1)','q');
x1=z1/(S*(K1-1)+1); <-------This Line
x2=z2/(S*(K2-1)+1);
x3=z3/(S*(K3-1)+1);
y1=K1*x1;
y2=K2*x2;
y3=K3*x3;
end

I get the following errors at the indicated line.

??? Error using ==> maple at 129 Error, (in linalg:-linsolve) 2nd argument fails to evaluate to a vector or a matrix

Error in ==> sym.mldivide at 30 X = maple('linsolve',char(A),char(B),'''_rank''');

Error in ==> sym.mrdivide at 29 X = (A.'\B.').';

Error in ==> flash1 at 35 x1=z1/(S*(K1-1)+1);

Thanks for the help

Answer 1

[Matt Fig]

That:

if true

is not part of the code, right? If it is, then you should be not even able to run the code....

Did you look at the output S? S has two solutions; it is a 2-by-1 one symbolic variable. To turn it into a 2-by-one double value, do:

S = subs(S);

Now use the element-by-element division:

x1=z1./(S*(K1-1)+1);
x2=z2./(S*(K2-1)+1);
x3=z3./(S*(K3-1)+1);
y1=K1*x1;
y2=K2*x2;
y3=K3*x3;

[KeJoSa]

Hi Matt, As a follow up to this question is there a way to only return one of the two solutions to this equation? There are obviously two roots but, only one of them makes any physical sense (the one that is less than or equal to 1 but, greater than 0).

Thanks

[Matt Fig]

Yes, do it with logical indexing.

A = [-9 0 9 12 56];
A = A(A>0 & A<10)

Answer 2

[Walter Roberson]

You must be using an older version; newer versions would use MuPAD instead of Maple.

ANyhow, you need to put

syms q

before your solve() call, and remove the quotes in your expression in that call.

S = solve(z1*(K1-1)/(q*(K1-1)+1)+z2*(K2-1)/(q*(K2-1)+1)+z3*(K3-1)/(q*(K3-1)+1),q)

[Matt Fig]

Actually, this is an extension of a previous question. KeJoSa indicated that S returns o.k. in that question. Proceeding with the remainder of the code as I did above works in 2011b.

[KeJoSa]

Indeed, this is an extension of that question. Thanks again Matt!