Solving mach number for 5000 different iterations

조회 수: 42 (최근 30일)
Zachary Kubala
Zachary Kubala 2020년 3월 28일
편집: Walter Roberson 2020년 3월 29일
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
  댓글 수: 1
Zachary Kubala
Zachary Kubala 2020년 3월 28일
It only returns a 2x1 symbolic, I would like a V to be 1x5000 for each 6,8,10,12, and 14

댓글을 달려면 로그인하십시오.

답변 (1개)

Walter Roberson
Walter Roberson 2020년 3월 28일
%% 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
Zachary Kubala
Zachary Kubala 2020년 3월 29일
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
Walter Roberson 2020년 3월 29일
편집: Walter Roberson 2020년 3월 29일
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 ?

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Applications에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by