From fbf8ebf7dd98bf1cd9266144b7a7b4206f4dd06b Mon Sep 17 00:00:00 2001
From: Thomas Lambert <dev@tlambert.be>
Date: Tue, 12 Mar 2024 14:45:38 +0100
Subject: [PATCH] feat(Polar): add interpolated polars
This allows the user to create polars for interpolated airfoils.
When the shape of the local airfoil is unknown but is between sections
of known airfoils, this allows to create interpolated airfoils polars
based on the known polars of the bounding airfoils.
---
+af_tools/@Polar/Polar.m | 2 +
+af_tools/@Polar/interppolar.m | 145 +++++++++++++++++++++++++++++++++
2 files changed, 147 insertions(+)
create mode 100644 +af_tools/@Polar/interppolar.m
diff --git a/+af_tools/@Polar/Polar.m b/+af_tools/@Polar/Polar.m
index 196137a..69fbc62 100644
--- a/+af_tools/@Polar/Polar.m
+++ b/+af_tools/@Polar/Polar.m
@@ -180,6 +180,8 @@ classdef Polar < handle
self = loadpolar(self, polarFile) % Load polar from file
self = polypolar(polyCl, polyCd) % Create polar using polynomial coefficients
+ % Create a polar by interpolating between two different polar objects
+ self = interppolar(Polar1, Polar2, position)
end
end
diff --git a/+af_tools/@Polar/interppolar.m b/+af_tools/@Polar/interppolar.m
new file mode 100644
index 0000000..054ff05
--- /dev/null
+++ b/+af_tools/@Polar/interppolar.m
@@ -0,0 +1,145 @@
+function self = interppolar(Polar1, Polar2, weightPolar2)
+ % INTERPPOLAR Create a new polar object by interpolation of two other polars.
+ % This static method can be used in place of a normal constructor.
+ % This method can be used when the user wants to interpolate the polars bewteen two known
+ % airfoils. Consider that the airfoil at x=0 is a NACA 0012, and the one at x=1 is a CLARY-Y.
+ % Rather than assuming that all segments in [0,1] are NACA 0012 and all segments after 1 are
+ % CLARK-Y, this function allows to create polars for the [0,1] interval that are a mix between
+ % the two bounds. If we are at 0.1, the newly created Polar object, will then be 90% of NACA
+ % 0012 and 10% of CLARK-Y.
+ % Note:
+ % The two input Polars should have the same range of AOA, Reynolds, Mach, and nCrit.
+ % -----
+ %
+ % Usage:
+ % MiedPolar = Polar.INTERPPOLAR(Polar1, Polar2, weightPolar2) create an new Polar object that
+ % is a linear interpolation between Polar1 and Polar2, with a given weight for Polar2.
+ %
+ % Inputs:
+ % Polar1 : Polar object for one part of the interpolation.
+ % Polar2 : Polar object for the second part of the interpolation. Muse have the same Re, Mach
+ % and nCrit values as Polar1.
+ % weightPolar2: Scalar in [0-1] to skew the interpolation properly between Polar1 and Polar2.
+ %
+ % Output:
+ % self : Polar object that results from the interpolation.
+ %
+ % Example:
+ % MixedPolar = af_tools.interppolar(PolarNACA0012, PolarClarkY, 0.6); Creates an interpolated
+ % Polar object MixedPolar that results from the interpolation of PolarNACA0012
+ % and PolarClarkY, asusming that the MixedPolar is made of 60% ClarkY and 40%
+ % NACA0012.
+ %
+ % See also: af_tools.Polar.
+ %
+ % -----
+ % <a href="https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox">Documentation (online)</a>
+
+ % ----------------------------------------------------------------------------------------------
+ % (c) Copyright 2022-2024 University of Liege
+ % Author: Thomas Lambert <t.lambert@uliege.be>
+ % ULiege - Aeroelasticity and Experimental Aerodynamics
+ % MIT License
+ % Repo: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox
+ % Issues: https://gitlab.uliege.be/am-dept/matlab_airfoil_toolbox/-/issues
+ % ----------------------------------------------------------------------------------------------
+
+ % --- Defaults and constants
+ DEBUG = true;
+
+ % --- Input checks
+ % Check if both Polars are Polars objects and weightPolar2 is in [0,1]
+ validateattributes(Polar1, {'af_tools.Polar'}, {'scalar'}, mfilename(), 'Polar1', 1);
+ validateattributes(Polar2, {'af_tools.Polar'}, {'scalar'}, mfilename(), 'Polar2', 2);
+ validateattributes(weightPolar2, {'double'}, {'scalar', 'nonnegative', '<=', 1}, ...
+ mfilename(), 'weightPolar2', 3);
+
+ % Check that both Polars have the same AOA, Reynolds, Mach and nCrit
+ equalfieldscheck(Polar1, Polar2, 'aoa');
+ equalfieldscheck(Polar1, Polar2, 'reynolds');
+ equalfieldscheck(Polar1, Polar2, 'mach');
+ equalfieldscheck(Polar1, Polar2, 'nCrit');
+
+ % --- Interpolate the polars
+
+ self = af_tools.Polar;
+
+ wPol2 = weightPolar2;
+ wPol1 = 1 - wPol2;
+ nPolars = numel(Polar1.reynolds);
+
+ % Set base polar information
+ afName = {sprintf('MIX-%d%%%s-%d%%%s', round(wPol1 * 100), Polar1.airfoil{1}, ...
+ round(wPol2 * 100), Polar2.airfoil{2})};
+ self.airfoil = repmat(afName, 1, nPolars);
+
+ self.reynolds = Polar1.reynolds;
+ self.mach = Polar1.mach;
+ self.nCrit = Polar1.nCrit;
+
+ % Interpolate the angle of attack, cl, cd and cm using the base polars
+ if Polar1.aoa ~= Polar2.aoa
+ % Reinterpolate the original polars so they share the same range of AOA.
+ aoaMin = min(min(Polar1.aoa(:, 1)), min(Polar2.aoa(:, 1)));
+ aoaMax = max(max(Polar1.aoa(:, 1)), max(Polar2.aoa(:, 1)));
+ aoaStep = min(diff(Polar1.aoa(1:2, 1)), diff(Polar2.aoa(1:2, 1)));
+ aoaVect = aoaMin:aoaStep:aoaMax;
+
+ aoa1 = repmapt(aoaVect', 1, nPolars);
+ aoa2 = repmapt(aoaVect', 1, nPolars);
+ cl1 = interp1(Polar1.aoa(:, 1), Polar1.cl, aoaVect');
+ cl2 = interp1(Polar2.aoa(:, 1), Polar2.cl, aoaVect');
+ cd1 = interp1(Polar1.aoa(:, 1), Polar1.cd, aoaVect');
+ cd2 = interp1(Polar2.aoa(:, 1), Polar2.cd, aoaVect');
+ cm1 = interp1(Polar1.aoa(:, 1), Polar1.cm, aoaVect');
+ cm2 = interp1(Polar2.aoa(:, 1), Polar2.cm, aoaVect');
+ else
+ % Use direclty the polar data
+ aoa1 = Polar1.aoa;
+ aoa2 = Polar2.aoa;
+ cl1 = Polar1.cl;
+ cl2 = Polar2.cl;
+ cd1 = Polar1.cd;
+ cd2 = Polar2.cd;
+ cm1 = Polar1.cm;
+ cm2 = Polar2.cm;
+ end
+ self.aoa = lininterparray(wPol1, aoa1, wPol2, aoa2);
+ self.cl = lininterparray(wPol1, cl1, wPol2, cl2);
+ self.cd = lininterparray(wPol1, cd1, wPol2, cd2);
+ self.cm = lininterparray(wPol1, cm1, wPol2, cm2);
+
+ % --- DEBUGGING
+ if DEBUG
+ % Plot the original polars for each Reynolds as well as the interpolated one
+ for j = 1:nPolars
+ figure('Name', 'Combined Polars');
+ hold on;
+ plot(Polar1.aoa(:, j), Polar1.cl(:, j));
+ plot(Polar2.aoa(:, j), Polar2.cl(:, j));
+ plot(self.aoa(:, j), self.cl(:, j));
+ legend(Polar1.airfoil{1}, Polar2.airfoil{1}, self.airfoil{1}, 'Location', 'NorthWest');
+ grid on;
+ end
+ end
+
+end
+
+function interpArray = lininterparray(weight1, array1, weight2, array2)
+ % LININTERPARRAY Linear interpolation between two arrays
+
+ interpArray = weight1 * array1 + weight2 * array2;
+
+end
+
+function equalfieldscheck(Polar1, Polar2, field)
+ % FIELDCHECK Check if some fields are equal in both Polar objects
+
+ if Polar1.(field) ~= Polar2.(field)
+ error('Polar:interppolar:different%s', ...
+ ['Impossible to interpolate between the two Polar objects as they have '...
+ 'different %s.\n'...
+ 'Please use two Polar objects with identical AOA, Reynolds, Mach and nCrit'], ...
+ field, field);
+ end
+end
--
GitLab