Solving mach number for 5000 different iterations
22 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
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
Antworten (1)
Walter Roberson
am 28 Mär. 2020
%% 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);
3 Kommentare
Walter Roberson
am 29 Mär. 2020
Bearbeitet: Walter Roberson
am 29 Mär. 2020
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 ?
Siehe auch
Kategorien
Mehr zu Applications finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!