Solving mach number for 5000 different iterations

I am trying to solve an equation for Mach number for the Rayleigh pitot tube equation over 5000 iterations for each time I acquired data.

%% Equations
%p02/p1 ratio
p02_p1_6 = p02_6./p1_6;
p02_p1_8 = p02_8./p1_8;
p02_p1_10 = p02_10./p1_10;
p02_p1_12 = p02_12./p1_12;
p02_p1_14 = p02_14./p1_14;


%% Find M with rayleigh's pitot tube formula
syms M6 M8 M10 M12 M14
g = 1.4;

for i = 1:5000
eqn1 = (((((g+1)^2)*M6^2)/(4*g*(M6^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M6^2)/(g+1))-p02_p1_6(i) == 0;
V6 = vpasolve(eqn1,M6);
eqn2 = (((((g+1)^2)*M8^2)/(4*g*(M8^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M8^2)/(g+1))-p02_p1_8(i) == 0;
V8 = vpasolve(eqn2,M8);
eqn3 = (((((g+1)^2)*M10^2)/(4*g*(M10^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M10^2)/(g+1))-p02_p1_10(i) == 0;
V10 = vpasolve(eqn3,M10);
eqn4 = (((((g+1)^2)*M12^2)/(4*g*(M12^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M12^2)/(g+1))-p02_p1_12(i) == 0;
V12 = vpasolve(eqn4,M12);
eqn5 = (((((g+1)^2)*M14^2)/(4*g*(M14^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M14^2)/(g+1))-p02_p1_14(i) == 0;
V14 = vpasolve(eqn5,M14);

end

It takes about 15 minutes to get an output but every time it returns symbolics

[Zachary Kubala]

It only returns a 2x1 symbolic, I would like a V to be 1x5000 for each 6,8,10,12, and 14

Answer 1

[Walter Roberson]

%% Equations
%p02/p1 ratio
p02_p1_6 = p02_6./p1_6;
p02_p1_8 = p02_8./p1_8;
p02_p1_10 = p02_10./p1_10;
p02_p1_12 = p02_12./p1_12;
p02_p1_14 = p02_14./p1_14;


%% Find M with rayleigh's pitot tube formula
syms M6 M8 M10 M12 M14
g = 1.4;

eqn1 = (((((g+1)^2)*M6^2)/(4*g*(M6^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M6^2)/(g+1));
eqn2 = (((((g+1)^2)*M8^2)/(4*g*(M8^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M8^2)/(g+1));
eqn3 = (((((g+1)^2)*M10^2)/(4*g*(M10^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M10^2)/(g+1));
eqn4 = (((((g+1)^2)*M12^2)/(4*g*(M12^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M12^2)/(g+1));
eqn5 = (((((g+1)^2)*M14^2)/(4*g*(M14^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M14^2)/(g+1));

%hidden loops
V6 = arrayfun(@(v6) vpasolve(eqn1 == v6, M6), p02_p1_6);
V8 = arrayfun(@(v8) vpasolve(eqn2 == v8, M8), p02_p1_8);
V10 = arrayfun(@(v10) vpasolve(eqn3 == v10, M10), p02_p1_10);
V12 = arrayfun(@(v12) vpasolve(eqn4 == v12, M12), p02_p1_12);
V14 = arrayfun(@(v14) vpasolve(eqn5 == v14, M14), p02_p1_14);

[Zachary Kubala]

Thank you!

[Zachary Kubala]

Walter, for "arrayfun(@(v6) ... what kind of function did you use? I recieved the error:

Error using arrayfun

Non-scalar in Uniform output, at index 1, output 1.

Set 'UniformOutput' to false.

[Walter Roberson]

V6 = arrayfun(@(v6) vpasolve(eqn1 == v6, M6), p02_p1_6, 'uniform', 0);
V8 = arrayfun(@(v8) vpasolve(eqn2 == v8, M8), p02_p1_8, 'uniform', 0);
V10 = arrayfun(@(v10) vpasolve(eqn3 == v10, M10), p02_p1_10, 'uniform', 0);
V12 = arrayfun(@(v12) vpasolve(eqn4 == v12, M12), p02_p1_12, 'uniform', 0);
V14 = arrayfun(@(v14) vpasolve(eqn5 == v14, M14), p02_p1_14, 'uniform', 0);

The message is telling you that vpasolve() did not return exactly one output for some value of p02_p1_6.

Indeed, looking at the structure of the formulas, you have a quadratic and there will be exactly two solutions, which can be easily found through a closed-form formula.

All of your formula are the same, just with different variable names. There is no reason to use separate equations.

V6 = arrayfun(@(v6) vpasolve(eqn1 == v6, M6), p02_p1_6, 'uniform', 0);
V8 = arrayfun(@(v8) vpasolve(eqn1 == v8, M6), p02_p1_8, 'uniform', 0);
V10 = arrayfun(@(v10) vpasolve(eqn1 == v10, M6), p02_p1_10, 'uniform', 0);
V12 = arrayfun(@(v12) vpasolve(eqn1 == v12, M6), p02_p1_12, 'uniform', 0);
V14 = arrayfun(@(v14) vpasolve(eqn1 == v14, M6), p02_p1_14, 'uniform', 0);

%% Find M with rayleigh's pitot tube formula
eqn1 = (((((g+1)^2)*M6^2)/(4*g*(M6^2)-2*(g-1)))^(g/g-1)*...
    (1-g+2*g*M6^2)/(g+1))-p02_p1_6(i) == 0;

Look more carefully. That is not rayleigh's pitot tube formula, and the mistake in the formula is why you are getting a quadratic instead of a single value.

Hint: what is 5/5-1 ?