Skip to content
Snippets Groups Projects
Verified Commit 4faf74fa authored by Thomas Lambert's avatar Thomas Lambert :helicopter:
Browse files

refactor(coeffs): prevent confict with cd command

parent 34dc2698
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 205 additions and 200 deletions
+af_tools/+utils/parsepolarinputs.m
+ 23
− 23
View file @ 4faf74fa
function [alpha, cl, cd, cm, idxOpts] = parsepolarinputs(optList, varargin)
function [alpha, c_l, c_d, c_m, idxOpts] = parsepolarinputs(optList, varargin)
% PARSEPOLARINPUTS Parses the input and checks their validity.
% Returns the arrays for the main polar arrays (alpha, cl, cd, cm), as well as the index where
% the remaining options start.
......@@ -38,7 +38,7 @@ function [alpha, cl, cd, cm, idxOpts] = parsepolarinputs(optList, varargin)
% Prompt user for polar
[filename, path] = uigetfile('*.mat', 'Select a Polar file');
tmp = load([path, filename]);
[alpha, cl, cd, cm] = extractpolar(tmp.Polar);
[alpha, c_l, c_d, c_m] = extractpolar(tmp.Polar);
clear tmp;
else
if isstruct(varargin{1})
......@@ -46,29 +46,29 @@ function [alpha, cl, cd, cm, idxOpts] = parsepolarinputs(optList, varargin)
% Extract polar values
Polar = varargin{1};
validateattributes(Polar, {'struct'}, {'nonempty'}, mfilename(), 'Polar', 1);
[alpha, cl, cd, cm] = extractpolar(varargin{1});
[alpha, c_l, c_d, c_m] = extractpolar(varargin{1});
else
[alpha, cl, cd, cm, idxOpts] = parsearrays(optList, varargin);
[alpha, c_l, c_d, c_m, idxOpts] = parsearrays(optList, varargin);
end
end
% Transform row vectors into column vectors
alpha = vecttocol(alpha);
cl = vecttocol(cl);
cd = vecttocol(cd);
cm = vecttocol(cm);
c_l = vecttocol(c_l);
c_d = vecttocol(c_d);
c_m = vecttocol(c_m);
% Validate input arrays or inputs coming from Polar structure
validateattributes(alpha, {'numeric'}, {'real', '2d', 'nonempty', 'increasing'}, ...
mfilename(), 'alpha', 1);
validatearray(cl, size(alpha, 1), 2);
validatearray(cd, size(alpha, 1), 3);
validatearray(cm, size(alpha, 1), 4);
validatearray(c_l, size(alpha, 1), 2);
validatearray(c_d, size(alpha, 1), 3);
validatearray(c_m, size(alpha, 1), 4);
% Standardize alpha size so we have one column per Polar as well
if size(alpha, 2) == 1 && ~isempty(cl)
alpha = repmat(alpha, [1, size(cl, 2)]);
if size(alpha, 2) == 1 && ~isempty(c_l)
alpha = repmat(alpha, [1, size(c_l, 2)]);
end
end
......@@ -81,13 +81,13 @@ function bool = isoption(optList, var)
end
% --------------------------------------------------------------------------------------------------
function [alpha, cl, cd, cm] = extractpolar(Polar)
function [alpha, c_l, c_d, c_m] = extractpolar(Polar)
% EXTRACTPOLAR Extract alpha, cl, cd and cm fields Polar structure.
alpha = Polar.alpha;
cl = Polar.cl;
cd = Polar.cd;
cm = Polar.cm;
c_l = Polar.cl;
c_d = Polar.cd;
c_m = Polar.cm;
end
% --------------------------------------------------------------------------------------------------
......@@ -108,22 +108,22 @@ function [var, idxOpts, foundOpts] = assignvar(optList, data, idx)
end
% --------------------------------------------------------------------------------------------------
function [alpha, cl, cd, cm, idxOpts] = parsearrays(optList, data)
function [alpha, c_l, c_d, c_m, idxOpts] = parsearrays(optList, data)
% PARSEARRAYS Parse user-provided arrays.
alpha = [];
cl = [];
cd = [];
cm = [];
c_l = [];
c_d = [];
c_m = [];
if length(data) >= 1
[alpha, idxOpts, foundOpts] = assignvar(optList, data, 1);
if length(data) >= 2 && ~foundOpts
[cl, idxOpts, foundOpts] = assignvar(optList, data, 2);
[c_l, idxOpts, foundOpts] = assignvar(optList, data, 2);
if length(data) >= 3 && ~foundOpts
[cd, idxOpts, foundOpts] = assignvar(optList, data, 3);
[c_d, idxOpts, foundOpts] = assignvar(optList, data, 3);
if length(data) >= 4 && ~foundOpts
[cm, idxOpts, ~] = assignvar(optList, data, 4);
[c_m, idxOpts, ~] = assignvar(optList, data, 4);
end
end
end
......
+af_tools/+utils/trimtorange.m
+ 10
− 10
View file @ 4faf74fa
function [alphaTrim, clTrim, cdTrim] = trimtorange(range, alpha, cl, cd)
function [alpha_trim, c_l_trim, c_d_trim] = trimtorange(range, alpha, c_l, c_d)
% TRIMTORANGE Trims polars to a given range of angles of attack.
% -----
% (c) Copyright 2022 University of Liege
......@@ -12,23 +12,23 @@ function [alphaTrim, clTrim, cdTrim] = trimtorange(range, alpha, cl, cd)
narginchk(2, 4);
for i = 1:size(alpha, 2)
if nargin >= 3
clTrim(:, i) = cl(alpha(:, i) >= range(1) & alpha(:, i) <= range(end), i);
if nargin == 4 && ~isempty(cd)
cdTrim(:, i) = cd(alpha(:, i) >= range(1) & alpha(:, i) <= range(end), i);
c_l_trim(:, i) = c_l(alpha(:, i) >= range(1) & alpha(:, i) <= range(end), i);
if nargin == 4 && ~isempty(c_d)
c_d_trim(:, i) = c_d(alpha(:, i) >= range(1) & alpha(:, i) <= range(end), i);
else
cdTrim = [];
c_d_trim = [];
end
else
nargoutchk(1, 1);
end
alphaTrim(:, i) = alpha(alpha(:, i) >= range(1) & alpha(:, i) <= range(end), i);
alpha_trim(:, i) = alpha(alpha(:, i) >= range(1) & alpha(:, i) <= range(end), i);
end
% If impossible to trim, just return the input
if isempty(alphaTrim)
alphaTrim = alpha;
clTrim = cl;
cdTrim = cd;
if isempty(alpha_trim)
alpha_trim = alpha;
c_l_trim = c_l;
c_d_trim = c_d;
end
end
+af_tools/extendpolar.m
+ 46
− 45
View file @ 4faf74fa
function [alphaExt, clExt, cdExt] = extendpolar(varargin)
function [alpha_ext, c_l_ext, c_d_ext] = extendpolar(varargin)
% EXTENDPOLAR Extends polars to the full range of AOA ([-180, 180] deg).
% This function implements different empirical methods to extend the polars outside of their
% usual range.
......@@ -82,8 +82,8 @@ function [alphaExt, clExt, cdExt] = extendpolar(varargin)
OPTION_LIST = {'method', 'limit'}; % No options for this function
% Extract and validate the inputs
[alpha, cl, cd, ~, idxOpts] = parsepolarinputs(OPTION_LIST, varargin{:});
if isempty(cd)
[alpha, c_l, c_d, ~, idxOpts] = parsepolarinputs(OPTION_LIST, varargin{:});
if isempty(c_d)
error('MATLAB:extendpolar:noCD', ...
'To be extended, the polars must specify alpha, cl and cd.');
end
......@@ -96,7 +96,7 @@ function [alphaExt, clExt, cdExt] = extendpolar(varargin)
end
try
findstall(alpha, cl, cd);
findstall(alpha, c_l, c_d);
catch
error('MATLAB:extendpolar:noStall', ...
'To be extended, a polar must at least contain the stall point.');
......@@ -107,14 +107,14 @@ function [alphaExt, clExt, cdExt] = extendpolar(varargin)
if alpha(1, 1) == -pi || alpha(end, 1) == pi
% If input already ranges from [-pi,pi], just send it back
alphaExt = alpha;
clExt = cl;
cdExt = cd;
alpha_ext = alpha;
c_l_ext = c_l;
c_d_ext = c_d;
else
% Extend using the selected method
switch method
case 'Viterna'
[alphaExt, clExt, cdExt] = viterna(limit, alpha, cl, cd);
[alpha_ext, c_l_ext, c_d_ext] = viterna(limit, alpha, c_l, c_d);
end
end
......@@ -147,7 +147,7 @@ function [method, limit] = parseoptioninputs(varargin)
end
% --------------------------------------------------------------------------------------------------
function [alphaExt, clExt, cdExt] = viterna(limit, alpha, cl, cd)
function [alpha_ext, c_l_ext, c_d_ext] = viterna(limit, alpha, c_l, c_d)
% VITERNA Implementation of the original Viterna model
% Source:
% Viterna & Corrigan. 1982. "Fixed Pitch Rotor Performance of Large Horizontal Axis Wind
......@@ -176,14 +176,14 @@ function [alphaExt, clExt, cdExt] = viterna(limit, alpha, cl, cd)
SCALING_CL = 0.7; % Suggested in AeroDyn's manual
% Loop for every polar data (every column of cl)
for i = 1:size(cl, 2)
for i = 1:size(c_l, 2)
% Trim the polars to nonNaN values only
% (some Cl, Cd may be NaN if they come from xf2mat with trimAoas = false)
id = find(~isnan(cl(:, i)) == 1);
id = find(~isnan(c_l(:, i)) == 1);
trimmedAlpha = alpha(id, i);
trimmedCl = cl(id, i);
trimmedCd = cd(id, i);
trimmedCl = c_l(id, i);
trimmedCd = c_d(id, i);
% Get limit point
if strcmp(limit, 'stall')
......@@ -203,73 +203,74 @@ function [alphaExt, clExt, cdExt] = viterna(limit, alpha, cl, cd)
sas = sin(alpha_l);
cas = cos(alpha_l);
cd_max = 1.11 + 0.018 * DEFAULT_AR;
VitCoeffs.a1 = cd_max / 2;
VitCoeffs.b1 = cd_max;
VitCoeffs.a2 = (cl_l - cd_max .* sas .* cas) .* sas ./ cas.^2;
VitCoeffs.b2 = (cd_l - cd_max .* sas.^2) ./ cas;
c_d_max = 1.11 + 0.018 * DEFAULT_AR;
VitCoeffs.a1 = c_d_max / 2;
VitCoeffs.b1 = c_d_max;
VitCoeffs.a2 = (cl_l - c_d_max .* sas .* cas) .* sas ./ cas.^2;
VitCoeffs.b2 = (cd_l - c_d_max .* sas.^2) ./ cas;
% Minimum step between 2 AOA
alpha_step = (min(diff(trimmedAlpha)));
% 1. Extend Cl/Cd to ]limit, 90]
alpha1 = (alpha_l + alpha_step:alpha_step:pi / 2)';
[cl1, cd1] = viternaeq(VitCoeffs, alpha1, 1);
alpha_1 = (alpha_l + alpha_step:alpha_step:pi / 2)';
[cl1, cd1] = viternaeq(VitCoeffs, alpha_1, 1);
% 2. Extend to ]90, 180-limit]
alpha2 = (pi / 2 + alpha_step:alpha_step:pi - alpha_l)';
[cl2, cd2] = viternaeq(VitCoeffs, pi - alpha2, -SCALING_CL);
alpha_2 = (pi / 2 + alpha_step:alpha_step:pi - alpha_l)';
[c_l_2, c_d_2] = viternaeq(VitCoeffs, pi - alpha_2, -SCALING_CL);
% 3. Extend to ]180-limit,180]
alpha3 = (pi - alpha_l + alpha_step:alpha_step:pi)';
[~, cd3] = viternaeq(VitCoeffs, pi - alpha3, 1);
alpha_3 = (pi - alpha_l + alpha_step:alpha_step:pi)';
[~, c_d_3] = viternaeq(VitCoeffs, pi - alpha_3, 1);
% Correct lift with linear variation
cl3 = (alpha3 - pi) / alpha_l * cl_l * SCALING_CL;
c_l_3 = (alpha_3 - pi) / alpha_l * cl_l * SCALING_CL;
% 4. Extend to [-limit, -alpha(1)[ (by linear interpolation)
if trimmedAlpha(1) <= -alpha_l
alpha4 = [];
cl4 = [];
cd4 = [];
alpha_4 = [];
c_l_4 = [];
c_d_4 = [];
alpha_lneg = trimmedAlpha(1);
else
alpha4 = (-alpha_l:alpha_step:trimmedAlpha(1) - alpha_step)';
cl4 = interp1([-alpha_l, trimmedAlpha(1)], [-cl_l * SCALING_CL, trimmedCl(1)], alpha4);
cd4 = interp1([-alpha_l, trimmedAlpha(1)], [cd_l, trimmedCd(1)], alpha4);
alpha_4 = (-alpha_l:alpha_step:trimmedAlpha(1) - alpha_step)';
c_l_4 = interp1([-alpha_l, trimmedAlpha(1)], [-cl_l * SCALING_CL, trimmedCl(1)], ...
alpha_4);
c_d_4 = interp1([-alpha_l, trimmedAlpha(1)], [cd_l, trimmedCd(1)], alpha_4);
alpha_lneg = -alpha_l;
end
% 5. Extend to [-90,-limit[
alpha5 = (-pi / 2:alpha_step:alpha_lneg - alpha_step)';
[cl5, cd5] = viternaeq(VitCoeffs, alpha5, SCALING_CL);
alpha_5 = (-pi / 2:alpha_step:alpha_lneg - alpha_step)';
[c_l_5, c_d_5] = viternaeq(VitCoeffs, alpha_5, SCALING_CL);
% 6. Remaining part is simply a relfexion
alpha6 = (-pi:alpha_step:-pi / 2 - alpha_step)';
cl6 = -flipud([cl2; cl3]);
cd6 = flipud([cd2; cd3]);
alpha_6 = (-pi:alpha_step:-pi / 2 - alpha_step)';
c_l_6 = -flipud([c_l_2; c_l_3]);
c_d_6 = flipud([c_d_2; c_d_3]);
% Concatenate everything to get the full extended polars
alphaExt(:, i) = [alpha6; alpha5; alpha4; trimmedAlpha; alpha1; alpha2; alpha3];
clExt(:, i) = [cl6; cl5; cl4; trimmedCl; cl1; cl2; cl3];
cdExt(:, i) = [cd6; cd5; cd4; trimmedCd; cd1; cd2; cd3];
alpha_ext(:, i) = [alpha_6; alpha_5; alpha_4; trimmedAlpha; alpha_1; alpha_2; alpha_3];
c_l_ext(:, i) = [c_l_6; c_l_5; c_l_4; trimmedCl; cl1; c_l_2; c_l_3];
c_d_ext(:, i) = [c_d_6; c_d_5; c_d_4; trimmedCd; cd1; c_d_2; c_d_3];
% Trim negatives cd to cdmin
cdmin = min(trimmedCd);
cdExt(cdExt <= cdmin, i) = cdmin;
c_d_min = min(trimmedCd);
c_d_ext(c_d_ext <= c_d_min, i) = c_d_min;
end
end
% --------------------------------------------------------------------------------------------------
function [cl, cd] = viternaeq(VitCoeffs, alpha, clScaling)
function [c_l, c_d] = viternaeq(VitCoeffs, alpha, clScaling)
% VITERNAEQ Viterna's equation for the extension of the polars
% Source:
% Viterna & Corrigan. 1982. "Fixed Pitch Rotor Performance of Large Horizontal Axis Wind
% Turbines".
cl = VitCoeffs.a1 * sin(2 * alpha) + VitCoeffs.a2 * cos(alpha).^2 ./ sin(alpha);
cl = cl * clScaling;
c_l = VitCoeffs.a1 * sin(2 * alpha) + VitCoeffs.a2 * cos(alpha).^2 ./ sin(alpha);
c_l = c_l * clScaling;
cd = VitCoeffs.b1 * sin(alpha).^2 + VitCoeffs.b2 * cos(alpha);
c_d = VitCoeffs.b1 * sin(alpha).^2 + VitCoeffs.b2 * cos(alpha);
end
+af_tools/findcllinearrange.m
+ 20
− 20
View file @ 4faf74fa
function [alphaRange, clRange, clSlope] = findcllinearrange(varargin)
function [alphaRange, c_lRange, c_lSlope] = findcllinearrange(varargin)
% FINDCLLINEARRANGE Finds the linear range of the lift coefficient its slope.
% This function approximates the lift coefficient linear range by looking at its first and
% second derivatives (using hard-coded tolerance values).
......@@ -67,21 +67,21 @@ function [alphaRange, clRange, clSlope] = findcllinearrange(varargin)
OPTION_LIST = {}; % No options for this function
% Extract and validate the inputs
[alpha, cl] = parsepolarinputs(OPTION_LIST, varargin{:});
[alpha, c_l] = parsepolarinputs(OPTION_LIST, varargin{:});
% Trim extended polars to original range
% (avoid issues if extended polars are used as input)
[alpha, cl] = trimtorange(ORIGINAL_RANGE, alpha, cl);
[alpha, c_l] = trimtorange(ORIGINAL_RANGE, alpha, c_l);
% -------------------------------------------
% Find linear range for each polar
[alphaRange, clRange, clSlope] = linearrangefinder(alpha, cl);
[alphaRange, c_lRange, c_lSlope] = linearrangefinder(alpha, c_l);
end
% --------------------------------------------------------------------------------------------------
function [alphaRange, clRange, clSlope] = linearrangefinder(alpha, cl)
function [alphaRange, c_lRange, c_lSlope] = linearrangefinder(alpha, c_l)
% FINDLINEARRANGEVECT Finds the linearrange and returns the corresponding range and slope.
% We consider that the linear range meets the two following criterions:
% - First derivative close to 2*pi (TOL_DCL defines how close)
......@@ -93,17 +93,17 @@ function [alphaRange, clRange, clSlope] = linearrangefinder(alpha, cl)
TOL_DDCL = 0.01; % Tolerance on the curvature
% Initialization
alphaRange = zeros(2, size(cl, 2));
clRange = zeros(2, size(cl, 2));
clSlope = zeros(1, size(cl, 2));
alphaRange = zeros(2, size(c_l, 2));
c_lRange = zeros(2, size(c_l, 2));
c_lSlope = zeros(1, size(c_l, 2));
% Loop for every polar data (every column of cl)
for i = 1:size(cl, 2)
[tmpRange, rangeFound] = findrange(i, alpha, cl, TOL_DCL, TOL_DDCL, IDEAL_LIFT_SLOPE);
for i = 1:size(c_l, 2)
[tmpRange, rangeFound] = findrange(i, alpha, c_l, TOL_DCL, TOL_DDCL, IDEAL_LIFT_SLOPE);
% Check issues in case linear range was not found
if ~rangeFound
ploterror(alpha, cl, i, tmpRange);
ploterror(alpha, c_l, i, tmpRange);
error('MATLAB:findcllinearrange:noRangeFound', ...
['findcllinearrange was not able to determine the linear range of the given '...
'polar.\nMaybe the polar does not have enough points, or it is malformed.\n '...
......@@ -111,9 +111,9 @@ function [alphaRange, clRange, clSlope] = linearrangefinder(alpha, cl)
end
clRange(:, i) = cl(tmpRange, i);
c_lRange(:, i) = c_l(tmpRange, i);
alphaRange(:, i) = alpha(tmpRange, i);
clSlope(i) = diff(clRange(:, i)) / diff(alphaRange(:, i));
c_lSlope(i) = diff(c_lRange(:, i)) / diff(alphaRange(:, i));
end
......@@ -141,17 +141,17 @@ function bool = iscurvaturezero(ddCl, tol)
end
% --------------------------------------------------------------------------------------------------
function ploterror(alpha, cl, i, tmpRange)
function ploterror(alpha, c_l, i, tmpRange)
% PLOTERROR Plots the polar in case of error
figure('Name', 'linearrange:Debug');
hold on;
plot(alpha(:, i), cl(:, i), '--');
plot(alpha(tmpRange(1):tmpRange(2), i), cl(tmpRange(1):tmpRange(2), i), ...
plot(alpha(:, i), c_l(:, i), '--');
plot(alpha(tmpRange(1):tmpRange(2), i), c_l(tmpRange(1):tmpRange(2), i), ...
'color', 'r', 'lineWidth', 1.2);
plot(alpha(tmpRange(1), i), cl(tmpRange(1), i), ...
plot(alpha(tmpRange(1), i), c_l(tmpRange(1), i), ...
's', 'MarkerFaceColor', 'b', 'Markersize', 12);
plot(alpha(tmpRange(2), i), cl(tmpRange(2), i), ...
plot(alpha(tmpRange(2), i), c_l(tmpRange(2), i), ...
's', 'MarkerFaceColor', 'b', 'Markersize', 12);
grid on;
legend('Polar', 'Attempt at linear range', 'Missing point(s)');
......@@ -159,14 +159,14 @@ function ploterror(alpha, cl, i, tmpRange)
end
% --------------------------------------------------------------------------------------------------
function [tmpRange, rangeFound] = findrange(i, alpha, cl, TOL_DCL, TOL_DDCL, IDEAL_LIFT_SLOPE)
function [tmpRange, rangeFound] = findrange(i, alpha, c_l, TOL_DCL, TOL_DDCL, IDEAL_LIFT_SLOPE)
% FINDRANGE Finds the range iteratively.
rangeFound = false;
% First and second derivatives
dalpha = diff(alpha(:, i));
dcl = diff(cl(:, i)) ./ dalpha(1);
dcl = diff(c_l(:, i)) ./ dalpha(1);
ddcl = diff(dcl) ./ dalpha(1);
% Adapt tolerance to angle of attack
......
+af_tools/findstall.m
+ 16
− 16
View file @ 4faf74fa
function [alpha_s, cl_s, cd_s] = findstall(varargin)
function [alpha_s, c_l_s, c_d_s] = findstall(varargin)
% FINDSTALL Finds the stall point and returns the corresponding AOA, CL and CD.
% This function uses Matlab's findpeak in order to determine the position of the stall point.
%
......@@ -69,36 +69,36 @@ function [alpha_s, cl_s, cd_s] = findstall(varargin)
OPTION_LIST = {}; % No options for this function
% Extract and validate the inputs
[alpha, cl, cd] = parsepolarinputs(OPTION_LIST, varargin{:});
[alpha, c_l, c_d] = parsepolarinputs(OPTION_LIST, varargin{:});
% Trim extended polars to original range
% (avoid issues if extended polars are used as input)
[alpha, cl, cd] = trimtorange(ORIGINAL_RANGE, alpha, cl, cd);
[alpha, c_l, c_d] = trimtorange(ORIGINAL_RANGE, alpha, c_l, c_d);
% -------------------------------------------
% Find stall for each polar
[alpha_s, cl_s, cd_s] = stallfinder(alpha, cl, cd);
[alpha_s, c_l_s, c_d_s] = stallfinder(alpha, c_l, c_d);
end
% --------------------------------------------------------------------------------------------------
function [alpha_s, cl_s, cd_s] = stallfinder(alpha, cl, cd)
function [alpha_s, c_l_s, c_d_s] = stallfinder(alpha, c_l, c_d)
% STALLFINDER Finds the stall point and returns the corresponding AOA, CL and CD
PLOT_POLAR_IF_ERROR = true; % Plot the polar if error so user can see why it fails
% Initialization
alpha_s = zeros(1, size(cl, 2));
cl_s = zeros(1, size(cl, 2));
if isempty(cd)
cd_s = cd;
alpha_s = zeros(1, size(c_l, 2));
c_l_s = zeros(1, size(c_l, 2));
if isempty(c_d)
c_d_s = c_d;
else
cd_s = zeros(1, size(cl, 2));
c_d_s = zeros(1, size(c_l, 2));
end
% Loop for every polar data (every column of cl)
for i = 1:size(cl, 2)
for i = 1:size(c_l, 2)
% Find stall for each polar using findpeaks
maxFound = false;
......@@ -106,7 +106,7 @@ function [alpha_s, cl_s, cd_s] = stallfinder(alpha, cl, cd)
% Iterate in order to find only one peak (largest one should indicate the stall)
while ~maxFound && minWidth <= length(alpha(:, i))
[~, locs] = findpeaks(cl(:, i), 'MinPeakWidth', minWidth);
[~, locs] = findpeaks(c_l(:, i), 'MinPeakWidth', minWidth);
if numel(locs) == 1
maxFound = true;
else
......@@ -118,7 +118,7 @@ function [alpha_s, cl_s, cd_s] = stallfinder(alpha, cl, cd)
% If findpeaks does not work
if PLOT_POLAR_IF_ERROR
figure('Name', 'Debug: Polar with no stall');
plot(rad2deg(alpha(:, i)), cl(:, i));
plot(rad2deg(alpha(:, i)), c_l(:, i));
title('Polar with no clear stall peak');
xlabel('AOA [deg]');
ylabel('C_l [-]');
......@@ -131,9 +131,9 @@ function [alpha_s, cl_s, cd_s] = stallfinder(alpha, cl, cd)
end
alpha_s(i) = alpha(locs, i);
cl_s(i) = cl(locs, i);
if ~isempty(cd)
cd_s(i) = cd(locs, i);
c_l_s(i) = c_l(locs, i);
if ~isempty(c_d)
c_d_s(i) = c_d(locs, i);
end
end
......
+af_tools/findzerolift.m
+ 28
− 28
View file @ 4faf74fa
function [alpha_zeroL, cd_zeroL] = findzerolift(varargin)
function [alpha_0, c_d_0] = findzerolift(varargin)
% FINDZEROLIFT Finds the zero-lift angle and the associated cd.
% This function finds the zero-lift angle by interpolating around the last negative and first
% positive values for the lift coefficient.
......@@ -20,19 +20,19 @@ function [alpha_zeroL, cd_zeroL] = findzerolift(varargin)
% -----
%
% Usage:
% [ALPHA_ZEROL, CD_ZEROL] = FINDZEROLIFT prompts the user to input a Polar structure and
% determines the zero-lift angle of attack (ALPHA_ZEROL) and the drag coefficient at this
% angle (CD_ZEROL).
% [ALPHA_0, CD_0] = FINDZEROLIFT prompts the user to input a Polar structure and
% determines the zero-lift angle of attack (ALPHA_0) and the drag coefficient at this
% angle (CD_0).
%
% [ALPHA_ZEROL, CD_ZEROL] = FINDZEROLIFT(POLAR) uses the POLAR strcture given as input instead
% [ALPHA_0, CD_0] = FINDZEROLIFT(POLAR) uses the POLAR strcture given as input instead
% of prompting the user to select the file.
%
% [ALPHA_ZEROL, CD_ZEROL] = FINDZEROLIFT(ALPHA, CL) uses directly the values given in in ALPHA
% and CL to determine the zero-lift angle. As no CD is provided, CD_ZEROL will be an empty
% [ALPHA_0, CD_0] = FINDZEROLIFT(ALPHA, CL) uses directly the values given in in ALPHA
% and CL to determine the zero-lift angle. As no CD is provided, CD_0 will be an empty
% vector.
%
% [ALPHA_ZEROL, CD_ZEROL] = FINDZEROLIFT(ALPHA, CL, CD) uses directly the values given in in
% ALPHA and CL to determine the zero-lift angle. As CD is provided, CD_ZEROL will be returned
% [ALPHA_0, CD_0] = FINDZEROLIFT(ALPHA, CL, CD) uses directly the values given in in
% ALPHA and CL to determine the zero-lift angle. As CD is provided, CD_0 will be returned
% as well.
%
% Inputs:
......@@ -42,8 +42,8 @@ function [alpha_zeroL, cd_zeroL] = findzerolift(varargin)
% cd : Array with drag coefficients [aoas x nbPolar]
%
% Output:
% alpha_zeroL : Zero-lift angle of attack, IN RADIANS
% cd_zeroL : Drag at zero-lift angle of attack
% alpha_0 : Zero-lift angle of attack, IN RADIANS
% cd_0 : Drag at zero-lift angle of attack
%
% Example:
% FINDZEROLIFT
......@@ -72,21 +72,21 @@ function [alpha_zeroL, cd_zeroL] = findzerolift(varargin)
OPTION_LIST = {}; % No options for this function
% Extract and validate the inputs
[alpha, cl, cd] = parsepolarinputs(OPTION_LIST, varargin{:});
[alpha, c_l, c_d] = parsepolarinputs(OPTION_LIST, varargin{:});
% Trim extended polars to original range
% (avoid issues if extended polars are used as input)
[alpha, cl, cd] = trimtorange(ORIGINAL_RANGE, alpha, cl, cd);
[alpha, c_l, c_d] = trimtorange(ORIGINAL_RANGE, alpha, c_l, c_d);
% -------------------------------------------
% Find zero lift for each polar
[alpha_zeroL, cd_zeroL] = zerofinder(alpha, cl, cd);
[alpha_0, c_d_0] = zerofinder(alpha, c_l, c_d);
end
% --------------------------------------------------------------------------------------------------
function [alpha_zeroL, cd_zeroL] = zerofinder(alpha, cl, cd)
function [alpha_0, c_d_0] = zerofinder(alpha, c_l, c_d)
% FINDZEROLIFTVECT Finds the zero-lift point and returns the corresponding AOA and CD
% For vectors!
......@@ -97,23 +97,23 @@ function [alpha_zeroL, cd_zeroL] = zerofinder(alpha, cl, cd)
PLOT_POLAR_IF_ERROR = true; % Plot the polar if error so user can see why it fails
% Initialization
alpha_zeroL = zeros(1, size(cl, 2));
if isempty(cd)
cd_zeroL = cd;
alpha_0 = zeros(1, size(c_l, 2));
if isempty(c_d)
c_d_0 = c_d;
else
cd_zeroL = zeros(1, size(cl, 2));
c_d_0 = zeros(1, size(c_l, 2));
end
% Loop for every polar data (every column of cl)
for i = 1:size(cl, 2)
for i = 1:size(c_l, 2)
idx = find(cl(2:end, i) > 0, 1); % Find first cl > 0
idx = find(c_l(2:end, i) > 0, 1); % Find first cl > 0
idx = idx + 1; % Needs a +1 because we looked for cl(2:end)
if isempty(idx)
if PLOT_POLAR_IF_ERROR
figure('Name', 'Debug: Polar with no zero lift');
plot(rad2deg(alpha(:, i)), cl(:, i));
plot(rad2deg(alpha(:, i)), c_l(:, i));
title('Polar with no clear zero lift angle');
xlabel('AOA [deg]');
ylabel('C_l [-]');
......@@ -123,10 +123,10 @@ function [alpha_zeroL, cd_zeroL] = zerofinder(alpha, cl, cd)
'positive Cl']);
end
if sign(cl(idx, i)) == sign(cl(idx - 1, i))
if sign(c_l(idx, i)) == sign(c_l(idx - 1, i))
if PLOT_POLAR_IF_ERROR
figure('Name', 'Debug: Polar with no zero lift');
plot(rad2deg(alpha(:, i)), cl(:, i));
plot(rad2deg(alpha(:, i)), c_l(:, i));
title('Polar with no clear zero lift angle');
xlabel('AOA [deg]');
ylabel('C_l [-]');
......@@ -137,11 +137,11 @@ function [alpha_zeroL, cd_zeroL] = zerofinder(alpha, cl, cd)
end
dummy(:, 1) = alpha(idx - 1:idx, i); % AOA of the two points around cl = 0
dummy(:, 2) = cl(idx - 1:idx, i); % Cl of the two points around cl = 0
alpha_zeroL(i) = interp1(dummy(:, 2), dummy(:, 1), 0);
dummy(:, 2) = c_l(idx - 1:idx, i); % Cl of the two points around cl = 0
alpha_0(i) = interp1(dummy(:, 2), dummy(:, 1), 0);
if ~isempty(cd)
cd_zeroL(i) = interp1(alpha(:, i), cd(:, i), alpha_zeroL(:, i));
if ~isempty(c_d)
c_d_0(i) = interp1(alpha(:, i), c_d(:, i), alpha_0(:, i));
end
end
......
+af_tools/plotpolars.m
+ 6
− 6
View file @ 4faf74fa
......@@ -74,22 +74,22 @@ function plotpolars(varargin)
OPTION_LIST = {}; % No options for this function
% Extract and validate the inputs
[alpha, cl, cd, cm] = parsepolarinputs(OPTION_LIST, varargin{:});
[alpha, c_l, c_d, c_m] = parsepolarinputs(OPTION_LIST, varargin{:});
% -------------------------------------------
% Plots
set(0, 'defaultTextInterpreter', 'latex'); % Default interpreter for all plots
plotsinglepolar(rad2deg(alpha), cl, '$C_l - \alpha$', '$\alpha$ [deg]', '$C_l$ [-]');
plotsinglepolar(rad2deg(alpha), c_l, '$C_l - \alpha$', '$\alpha$ [deg]', '$C_l$ [-]');
if nargin > 2 || nargin <= 1
plotsinglepolar(rad2deg(alpha), cd, '$C_d - \alpha$', '$\alpha$ [deg]', '$C_d$ [-]');
plotsinglepolar(rad2deg(alpha), cl ./ cd, '$C_l/C_d - \alpha$', '$\alpha$ [deg]', ...
plotsinglepolar(rad2deg(alpha), c_d, '$C_d - \alpha$', '$\alpha$ [deg]', '$C_d$ [-]');
plotsinglepolar(rad2deg(alpha), c_l ./ c_d, '$C_l/C_d - \alpha$', '$\alpha$ [deg]', ...
'$C_l/C_d$ [-]');
plotsinglepolar(cl, cd, '$C_d - C_l$', '$C_l$ [-]', '$C_d$ [-]');
plotsinglepolar(c_l, c_d, '$C_d - C_l$', '$C_l$ [-]', '$C_d$ [-]');
end
if nargin > 3 || nargin <= 1
plotsinglepolar(rad2deg(alpha), cm, '$C_m - \alpha$', '$\alpha$ [deg]', '$C_m$ [-]');
plotsinglepolar(rad2deg(alpha), c_m, '$C_m - \alpha$', '$\alpha$ [deg]', '$C_m$ [-]');
end
end
......
+af_tools/xf2mat.m
+ 7
− 7
View file @ 4faf74fa
......@@ -243,20 +243,20 @@ function aoaRange = findmaxrange(aoaRange, resArray)
end
% --------------------------------------------------------------------------------------------------
function [alpha, cl, cd, cm] = standardizepolar(Input, aoaRange)
function [alpha, c_l, c_d, c_m] = standardizepolar(Input, aoaRange)
% STANDARDIZEPOLAR Standardize polar output.
% Re-interpolate everything over the same range of AOAs
alpha = zeros(length(aoaRange), length(Input));
cl = zeros(length(aoaRange), length(Input));
cd = zeros(length(aoaRange), length(Input));
cm = zeros(length(aoaRange), length(Input));
c_l = zeros(length(aoaRange), length(Input));
c_d = zeros(length(aoaRange), length(Input));
c_m = zeros(length(aoaRange), length(Input));
for i = 1:length(Input)
alpha(:, i) = deg2rad(aoaRange');
cl(:, i) = interp1(Input(i).resultArray(:, 1), Input(i).resultArray(:, 2), aoaRange);
cd(:, i) = interp1(Input(i).resultArray(:, 1), Input(i).resultArray(:, 3), aoaRange);
cm(:, i) = interp1(Input(i).resultArray(:, 1), Input(i).resultArray(:, 5), aoaRange);
c_l(:, i) = interp1(Input(i).resultArray(:, 1), Input(i).resultArray(:, 2), aoaRange);
c_d(:, i) = interp1(Input(i).resultArray(:, 1), Input(i).resultArray(:, 3), aoaRange);
c_m(:, i) = interp1(Input(i).resultArray(:, 1), Input(i).resultArray(:, 5), aoaRange);
end
end
......
CHANGELOG.md
+ 2
− 0
View file @ 4faf74fa
......@@ -17,6 +17,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
- Copyright notice: change main copyright holder to the University
- Doc: Revamp `CONTRIBUTING.md` to put emphasis on MISS_HIT and conventions
- Refactor: Rework functions to lower their complexity using mh_metrics
- Refactor: Rename lift, drag and moment coefficient variables to avoid conflict
with `cd` command.
### Deprecated
......
CONTRIBUTING.md
+ 4
− 2
View file @ 4faf74fa
......@@ -154,11 +154,13 @@ still apply:
The following naming convention is adopted:
- File names must be the same as the main function name.
- Functions should be _alllowercase_: `calculatestuff()`, `getallvariables()`
- Functions must be _alllowercase_: `calculatestuff()`, `getallvariables()`
- Constants should be all caps, with underscore between words: `MAX_ITERATIONS`
- In general, variables should use _camelCase_: `myVariable`, `liftCoefficient`
- In general, variables should use _camelCase_: `myVariable`
- Matlab structures should be _PascalCase_: `MyStructure`
- Indices are allowed if they make sense: `alpha_0`, `cDrag_0`, `alpha_s`
- Aerodynamic coefficients must have an underscore: `c_d` (prevents conflict
with `cd` command)
## Merge Request Process
......
tests/test_findcllinearrange.m
+ 5
− 5
View file @ 4faf74fa
......@@ -121,9 +121,9 @@ function test_noLinearRange(testCase)
% Ensure that the function throws an exception if no linear range was found
alpha = deg2rad(0:10)';
cl = rand(size(alpha));
c_l = rand(size(alpha));
verifyError(testCase, @() af_tools.findcllinearrange(alpha, cl), 'MATLAB:findcllinearrange:noRangeFound');
verifyError(testCase, @() af_tools.findcllinearrange(alpha, c_l), 'MATLAB:findcllinearrange:noRangeFound');
end
......@@ -133,11 +133,11 @@ function test_correctResults(testCase)
% Ensure that the function throws an exception if no linear range was found
alpha = deg2rad(-10:10);
cl = 2 * pi * alpha;
c_l = 2 * pi * alpha;
[alphaRange, clRange, clSlope] = af_tools.findcllinearrange(alpha, cl);
[alphaRange, clRange, clSlope] = af_tools.findcllinearrange(alpha, c_l);
verifyEqual(testCase, alphaRange, [alpha(1), alpha(end)]');
verifyEqual(testCase, clRange, [cl(1), cl(end)]');
verifyEqual(testCase, clRange, [c_l(1), c_l(end)]');
verifyEqual(testCase, clSlope, 2 * pi, 'AbsTol', 1e-9);
end
tests/test_findstall.m
+ 3
− 3
View file @ 4faf74fa
......@@ -140,10 +140,10 @@ function test_correctResults(testCase)
% Check if function works as intended on basic cases
alpha = deg2rad(0:120) / 4;
cl = sin(4 * alpha);
cd = cos(4 * alpha);
c_l = sin(4 * alpha);
c_d = cos(4 * alpha);
[alpha_s, cl_s, cd_s] = af_tools.findstall(alpha, cl, cd);
[alpha_s, cl_s, cd_s] = af_tools.findstall(alpha, c_l, c_d);
verifyEqual(testCase, alpha_s, (pi / 2) / 4, 'AbsTol', 1e-9);
verifyEqual(testCase, cl_s, 1, 'AbsTol', 1e-9);
......
tests/test_findzerolift.m
+ 5
− 5
View file @ 4faf74fa
......@@ -148,11 +148,11 @@ function test_correctResults(testCase)
alphad = (-10:10)';
alpha = deg2rad(alphad);
cl = alphad - 5;
cd = alphad + 5;
c_l = alphad - 5;
c_d = alphad + 5;
[alpha_zeroL, cd_zeroL] = af_tools.findzerolift(alpha, cl, cd);
[alpha_0, cd_0] = af_tools.findzerolift(alpha, c_l, c_d);
verifyEqual(testCase, alpha_zeroL, deg2rad(5));
verifyEqual(testCase, cd_zeroL, 10);
verifyEqual(testCase, alpha_0, deg2rad(5));
verifyEqual(testCase, cd_0, 10);
end
tests/test_parsepolarinputs.m
+ 11
− 11
View file @ 4faf74fa
......@@ -274,12 +274,12 @@ function test_emptyInputsInPolar(testCase)
AOA = 1:10;
EmptyPolar = struct('alpha', AOA, 'cl', [], 'cd', [], 'cm', []);
[alpha, cl, cd, cm, idxOpts] = af_tools.utils.parsepolarinputs(EMPTY_OPTS, EmptyPolar);
[alpha, c_l, c_d, c_m, idxOpts] = af_tools.utils.parsepolarinputs(EMPTY_OPTS, EmptyPolar);
verifyEqual(testCase, alpha, AOA');
verifyEqual(testCase, cl, []);
verifyEqual(testCase, cd, []);
verifyEqual(testCase, cm, []);
verifyEqual(testCase, c_l, []);
verifyEqual(testCase, c_d, []);
verifyEqual(testCase, c_m, []);
verifyEqual(testCase, idxOpts, 2);
end
......@@ -289,16 +289,16 @@ function test_emptyInputsArrays(testCase)
EMPTY_OPTS = {};
AOA = 1:10;
cl = [];
cd = [];
cm = [];
c_l = [];
c_d = [];
c_m = [];
[alpha, cl, cd, cm, idxOpts] = af_tools.utils.parsepolarinputs(EMPTY_OPTS, AOA, cl, cd, cm);
[alpha, c_l, c_d, c_m, idxOpts] = af_tools.utils.parsepolarinputs(EMPTY_OPTS, AOA, c_l, c_d, c_m);
verifyEqual(testCase, alpha, AOA');
verifyEqual(testCase, cl, []);
verifyEqual(testCase, cd, []);
verifyEqual(testCase, cm, []);
verifyEqual(testCase, c_l, []);
verifyEqual(testCase, c_d, []);
verifyEqual(testCase, c_m, []);
verifyEqual(testCase, idxOpts, 5);
end
......
tests/test_trimtorange.m
+ 19
− 19
View file @ 4faf74fa
......@@ -70,11 +70,11 @@ function test_rangeInAOA(testCase)
range = [-1, 1]';
alpha = (-10:10)';
cl = alpha + 10;
cd = alpha - 10;
c_l = alpha + 10;
c_d = alpha - 10;
expAlphaTrim = (-1:1)';
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, cl, cd);
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, c_l, c_d);
verifyEqual(testCase, alphaTrim, expAlphaTrim);
verifyEqual(testCase, clTrim, expAlphaTrim + 10);
......@@ -87,14 +87,14 @@ function test_rangeOutAOA(testCase)
range = [-20, -15]';
alpha = (-10:10)';
cl = alpha + 10;
cd = alpha - 10;
c_l = alpha + 10;
c_d = alpha - 10;
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, cl, cd);
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, c_l, c_d);
verifyEqual(testCase, alphaTrim, alpha);
verifyEqual(testCase, clTrim, cl);
verifyEqual(testCase, cdTrim, cd);
verifyEqual(testCase, clTrim, c_l);
verifyEqual(testCase, cdTrim, c_d);
end
......@@ -103,11 +103,11 @@ function test_rangePartOutAOA(testCase)
range = (-20:0)';
alpha = (-10:10)';
cl = alpha + 10;
cd = alpha - 10;
c_l = alpha + 10;
c_d = alpha - 10;
expAlphaTrim = (-10:0)';
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, cl, cd);
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, c_l, c_d);
verifyEqual(testCase, alphaTrim, expAlphaTrim);
verifyEqual(testCase, clTrim, expAlphaTrim + 10);
......@@ -120,14 +120,14 @@ function test_rangeEqualAOA(testCase)
alpha = (-10:10)';
range = (-10:10)';
cl = (20:40)';
cd = (-20:0)';
c_l = (20:40)';
c_d = (-20:0)';
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, cl, cd);
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, c_l, c_d);
verifyEqual(testCase, alphaTrim, alpha);
verifyEqual(testCase, clTrim, cl);
verifyEqual(testCase, cdTrim, cd);
verifyEqual(testCase, clTrim, c_l);
verifyEqual(testCase, cdTrim, c_d);
end
......@@ -136,11 +136,11 @@ function test_scalarRange(testCase)
range = 0;
alpha = (-10:10)';
cl = alpha + 10;
cd = alpha - 10;
c_l = alpha + 10;
c_d = alpha - 10;
expAlphaTrim = 0;
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, cl, cd);
[alphaTrim, clTrim, cdTrim] = af_tools.utils.trimtorange(range, alpha, c_l, c_d);
verifyEqual(testCase, alphaTrim, expAlphaTrim);
verifyEqual(testCase, clTrim, expAlphaTrim + 10);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment