vpasolve function incorrect solution

Hello

I am trying to solve this set of equations, however am getting incorrect solution:

syms Q1 Q2 f1 f2 Re1 Re2
eq1= Re1==(31.77E6)*Q1;
eq2= Re2==(15.88E6)*Q2;
eq3= 1/(f1^(1/2))==-1.8*log10(6.9/(Re1)+((1.125E-3)/3.7)^(1.11));
eq4= 1/(f2^(1/2))==-1.8*log10(6.9/(Re2)+((5.625E-4)/3.7)^(1.11));
eq5= 0.572==8*Q1*(1+sqrt(32*f1/f2))+(29.05E6)*f1*(Q1^2);
eq6= Q2==sqrt(32*f1/f2)*Q1;
S=vpasolve(eq1,eq2,eq3,eq4,eq5,eq6);

Matlab answer for f1 for example is 415671.45, however the correct answer should be approx 0.022. What am i doing wrong?

[madhan ravi]

You can specify the range into the solver.

[Samer Ali]

How can i specify the range for just one variable, lets say i know the solution of f1 to be between 0 and 1.

[Walter Roberson]

There seem to be at least 3 valid solutions. The solution that MATLAB produces is one of them. Another has f1 near 0.02653046284, and a third has f1 near 384805.6353. The f2 values are fairly similar (but not exactly the same) for all three of those.

[Walter Roberson]

S = vpasolve([eq1,eq2,eq3,eq4,eq5,eq6],[f1 f2 Q1 Q2 Re1 Re2], [0 1; [-inf(5,1), inf(5,1)]] );

[Samer Ali]

Hello Walter

executing your command:

S = vpasolve([eq1,eq2,eq3,eq4,eq5,eq6],[f1 f2 Q1 Q2 Re1 Re2], [0 1; [-inf(5,1), inf(5,1)]] );

Matlab gives me: Empty sym: 0-by-1?

The solution associated with f1 = 0.02653 should be the correct one but how i would get Matlab to search for this solution?

[Alex Sha]

1:
q1: 0.000823993387369594
q2: 0.00526536192506034
f1: 0.0265304628607812
f2: 0.0207914989607889
re1: 26178.269916732
re2: 83613.9473699583

feval:
-3.63797880709171E-12
1.16415321826935E-10
-1.77635683940025E-15
7.99360577730113E-15
2.22044604925031E-16
1.11022302462516E-16

2:
q1: 0.000859390439731193
q2: 4.34840412474342E-7
f1: 0.0263399277193314
f2: 3292196.85150343
re1: 27302.83427026
re2: 6.90526575009254

feval:
-3.63797880709171E-12
-1.06581410364015E-14
3.5527136788005E-15
4.41487124636097E-16
0
-1.70167860312205E-17

Answer 1

[Walter Roberson]

S1 = solve([eq1,eq2,eq5,eq6],[Q1,Q2,Re1,Re2]);
E2 = subs([eq3,eq4], S1);
sol1 = vpasolve(E2(1,:), [0 1; 0 1]);

The other branch, by solving E2(2,:), eventually gives up.

[Samer Ali]

Thank you Walter for the help. I solved it using Iterative Gauss Seidel method, which i guess was more reliable.