MATLAB Answers

0

Area Mach Number Relation

Asked by Steven Castrillon on 30 Sep 2019
Latest activity Answered by Chris
on 1 Oct 2019
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');

  0 Comments

Sign in to comment.

2 Answers

Chris
Answer by Chris
on 30 Sep 2019

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 Comment

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');

Sign in to comment.


Chris
Answer by Chris
on 1 Oct 2019

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.

  0 Comments

Sign in to comment.