Area Mach Number Relation

10 次查看(过去 30 天)
Steven Castrillon
Steven Castrillon 2019-9-30
回答: Chris 2019-10-1
The following code shows me how to get converged Mach number solutions for both subsonic and supersonic given Area ratio (ARatio), however, how can i input a range of ARatio instead of just one value so that the solutions are an array of converged mach numbers?
clear;
clc;
%% INPUTS
% Define some paramters
g = 1.4;
gm1 = g-1;
gp1 = g+1;
% Define anonymous function with two inputs (M and ARatio)
% - Will be used in the methods below
% - Pass M and ARatio as arguments to AM_EQN to get function value
% - funVal = AM_EQN(M,ARatio)
AM_EQN = @(M,ARatio) sqrt((1/M^2)*(((2+gm1*M^2)/gp1)^...
(gp1/gm1)))-ARatio;
% Solve for Msub and Msup using this area ratio (A/A*)
ARatio = 1.5;
% Error tolerance
errTol = 1e-4;
% Flags for printing iterations to screen
verboseBisection = 0;
verboseIncremental = 0;
%% SUBSONIC INCREMENTAL SEARCH
% Initial values
dM = 0.1; % Initial M step size
M = 1e-6; % Initial M value
iConvSub = 0; % Initial converge index
if (verboseIncremental == 1)
fprintf('Incremental Search Method: Subsonic\n');
fprintf('-----------------------------------\n');
end
% Iterate to solve for root
iterMax = 100; % Maximum iterations
stepMax = 100; % Maximum step iterations
for i = 1:1:iterMax
for j = 1:1:stepMax
% Evaluate function at j and j+1
fj = AM_EQN(M,ARatio);
fjp1 = AM_EQN(M+dM,ARatio);
% Print iterations to command window
if (verboseIncremental == 1)
fprintf('fj | fjp1: %3.4f\t%3.4f\n',fj,fjp1);
end
% Update M depending on sign change or not
% - If no sign change, then we are not bounding root yet
% - If sign change, then we are bounding the root, update dM
if (fj*fjp1 > 0)
M = M + dM; % Update M
elseif (fj*fjp1 < 0)
dM = dM*0.1; % Refine the M increment
break; % Break out of j loop
end
end % END: j Loop
% Check for convergence
if (abs(fj-fjp1) <= errTol) % If converged
iConvSub = i; % Set converged index
break; % Exit loop
end
end % END: i Loop
% Set subsonic Mach number to final M from iterations
Msub = M;
%% SUPERSONIC INCREMENTAL SEARCH
% Initial values
dM = 1; % Initial M step size
M = 1+1e-6; % Initial M value
iConvSup = 0; % Initial converge index
if (verboseIncremental == 1)
fprintf('\nIncremental Search Method: Supersonic\n');
fprintf('-------------------------------------\n');
end
% Iterate to solve for root
iterMax = 100; % Maximum iterations
stepMax = 100; % Maximum step iterations
for i = 1:1:iterMax
for j = 1:1:stepMax
% Evaluate function at j and j+1
fj = AM_EQN(M,ARatio);
fjp1 = AM_EQN(M+dM,ARatio);
% Print iterations to command window
if (verboseIncremental == 1)
fprintf('fj | fjp1: %3.4f\t%3.4f\n',fj,fjp1);
end
% Update M depending on sign change or not
% - If no sign change, then we are not bounding root yet
% - If sign change, then we are bounding the root, update dM
if (fj*fjp1 > 0)
M = M + dM; % Update M
elseif (fj*fjp1 < 0)
dM = dM*0.1; % Refine the M increment
break; % Break out of j loop
end
end % END: j Loop
% Check for convergence
if (abs(fj-fjp1) <= errTol) % If converged
iConvSup = i; % Set converged index
break; % Exit loop
end
end % END: i Loop
% Set supersonic Mach number to final M from iterations
Msup = M;
% Print solutions to command window
fprintf('==== INCREMENTAL SEARCH SOLVER ====\n');
fprintf('Msub: %3.4f after %i iterations\n',Msub,iConvSub);
fprintf('Msup: %3.4f after %i iterations\n',Msup,iConvSup);
fprintf('===================================\n\n');

回答(2 个)

Chris
Chris 2019-9-30
You can format posts to look llike ML text and make them easier to read. Without reading any of your code it sounds like you can write a wraper function if you can not eaisly vectorize the code.
  1 个评论
Steven Castrillon
Steven Castrillon 2019-9-30
编辑:Steven Castrillon 2019-9-30
Hopefully this is better formatting. Here my input ARatio is: ARatio = 1.5
However if i input ARatio as: ARatio = (0.1:0.1:10) i get the following error
Error in AREAMACH2 (line 62)
if (fj*fjp1 > 0)
If this can be solved by using a wraper function could you please help?
clear;
clc;
%% INPUTS
% Define some paramters
g = 1.4;
gm1 = g-1;
gp1 = g+1;
% Define anonymous function with two inputs (M and ARatio)
% - Will be used in the methods below
% - Pass M and ARatio as arguments to AM_EQN to get function value
% - funVal = AM_EQN(M,ARatio)
AM_EQN = @(M,ARatio) sqrt((1/M^2)*(((2+gm1*M^2)/gp1)^...
(gp1/gm1)))-ARatio;
% Solve for Msub and Msup using this area ratio (A/A*)
ARatio = 1.5;
% Error tolerance
errTol = 1e-4;
% Flags for printing iterations to screen
verboseBisection = 0;
verboseIncremental = 0;
%% SUBSONIC INCREMENTAL SEARCH
% Initial values
dM = 0.1; % Initial M step size
M = 1e-6; % Initial M value
iConvSub = 0; % Initial converge index
if (verboseIncremental == 1)
fprintf('Incremental Search Method: Subsonic\n');
fprintf('-----------------------------------\n');
end
% Iterate to solve for root
iterMax = 100; % Maximum iterations
stepMax = 100; % Maximum step iterations
for i = 1:1:iterMax
for j = 1:1:stepMax
% Evaluate function at j and j+1
fj = AM_EQN(M,ARatio);
fjp1 = AM_EQN(M+dM,ARatio);
% Print iterations to command window
if (verboseIncremental == 1)
fprintf('fj | fjp1: %3.4f\t%3.4f\n',fj,fjp1);
end
% Update M depending on sign change or not
% - If no sign change, then we are not bounding root yet
% - If sign change, then we are bounding the root, update dM
if (fj*fjp1 > 0)
M = M + dM; % Update M
elseif (fj*fjp1 < 0)
dM = dM*0.1; % Refine the M increment
break; % Break out of j loop
end
end % END: j Loop
% Check for convergence
if (abs(fj-fjp1) <= errTol) % If converged
iConvSub = i; % Set converged index
break; % Exit loop
end
end % END: i Loop
% Set subsonic Mach number to final M from iterations
Msub = M;
%% SUPERSONIC INCREMENTAL SEARCH
% Initial values
dM = 1; % Initial M step size
M = 1+1e-6; % Initial M value
iConvSup = 0; % Initial converge index
if (verboseIncremental == 1)
fprintf('\nIncremental Search Method: Supersonic\n');
fprintf('-------------------------------------\n');
end
% Iterate to solve for root
iterMax = 100; % Maximum iterations
stepMax = 100; % Maximum step iterations
for i = 1:1:iterMax
for j = 1:1:stepMax
% Evaluate function at j and j+1
fj = AM_EQN(M,ARatio);
fjp1 = AM_EQN(M+dM,ARatio);
% Print iterations to command window
if (verboseIncremental == 1)
fprintf('fj | fjp1: %3.4f\t%3.4f\n',fj,fjp1);
end
% Update M depending on sign change or not
% - If no sign change, then we are not bounding root yet
% - If sign change, then we are bounding the root, update dM
if (fj*fjp1 > 0)
M = M + dM; % Update M
elseif (fj*fjp1 < 0)
dM = dM*0.1; % Refine the M increment
break; % Break out of j loop
end
end % END: j Loop
% Check for convergence
if (abs(fj-fjp1) <= errTol) % If converged
iConvSup = i; % Set converged index
break; % Exit loop
end
end % END: i Loop
% Set supersonic Mach number to final M from iterations
Msup = M;
% Print solutions to command window
fprintf('==== INCREMENTAL SEARCH SOLVER ====\n');
fprintf('Msub: %3.4f after %i iterations\n',Msub,iConvSub);
fprintf('Msup: %3.4f after %i iterations\n',Msup,iConvSup);
fprintf('===================================\n\n');

请先登录,再进行评论。


Chris
Chris 2019-10-1
Put all that code info a function with inputs/outputs; call it for each ARatio you want. https://en.wikipedia.org/wiki/Wrapper_function
Control statements with vectors are commonly inapropreate.

类别

Help CenterFile Exchange 中查找有关 Dates and Time 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by