diff --git a/aviary/subsystems/mass/flops_based/test/test_wing_detailed.py b/aviary/subsystems/mass/flops_based/test/test_wing_detailed.py index f7981c1a3..dfd750fe4 100644 --- a/aviary/subsystems/mass/flops_based/test/test_wing_detailed.py +++ b/aviary/subsystems/mass/flops_based/test/test_wing_detailed.py @@ -6,6 +6,7 @@ from openmdao.utils.testing_utils import use_tempdirs from parameterized import parameterized + from aviary.subsystems.mass.flops_based.wing_detailed import ( BWBDetailedWingBendingFact, DetailedWingBendingFact, @@ -21,6 +22,7 @@ print_case, ) from aviary.variable_info.functions import setup_model_options +from aviary.variable_info.options import get_option_defaults from aviary.variable_info.variables import Aircraft, Mission, Settings omit_cases = ['LargeSingleAisle2FLOPS'] @@ -34,7 +36,7 @@ class DetailedWingBendingTest(unittest.TestCase): def setUp(self): self.prob = om.Problem() - # Skip model that doesn't use detailed wing and BWB cases. + # Skip models that don't use detailed wing and skip BWB cases. @parameterized.expand(get_flops_case_names(omit=omit_cases), name_func=print_case) def test_case(self, case_name): prob = self.prob @@ -420,6 +422,101 @@ def test_extreme_engine_loc(self): pod_inertia_factor = prob.get_val(Aircraft.Wing.ENG_POD_INERTIA_FACTOR) assert_near_equal(pod_inertia_factor, 1.0, tolerance=1e-10) + def test_intensity_factor(self): + """ + data taken from high_wing_single_aisle.csv + Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL = 3, 1, 1.5, 2.5 + """ + options = get_option_defaults() + options.set_val(Aircraft.Engine.NUM_ENGINES, val=[2], units='unitless') + options.set_val(Aircraft.Engine.NUM_WING_ENGINES, val=[2], units='unitless') + options.set_val(Aircraft.Propulsion.TOTAL_NUM_WING_ENGINES, val=2, units='unitless') + options.set_val( + Aircraft.Wing.INPUT_STATION_DISTRIBUTION, val=[0.0, 0.5, 1.0], units='unitless' + ) + options.set_val(Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL, val=3, units='unitless') + options.set_val(Aircraft.Wing.NUM_INTEGRATION_STATIONS, val=46, units='unitless') + + prob = self.prob + self.prob.model.add_subsystem( + 'wing', + DetailedWingBendingFact(), + promotes_inputs=['*'], + promotes_outputs=['*'], + ) + + prob.model.set_input_defaults( + Aircraft.Wing.LOAD_PATH_SWEEP_DISTRIBUTION, val=[23, 22.0], units='deg' + ) + prob.model.set_input_defaults( + Aircraft.Wing.THICKNESS_TO_CHORD_DISTRIBUTION, val=[0.1, 0.1, 0.1], units='unitless' + ) + prob.model.set_input_defaults( + Aircraft.Wing.CHORD_PER_SEMISPAN_DISTRIBUTION, val=[0.13, 0.115, 0.06], units='unitless' + ) + prob.model.set_input_defaults(Aircraft.Design.GROSS_MASS, val=150000.0, units='lbm') + prob.model.set_input_defaults(Aircraft.Wing.ASPECT_RATIO, val=9.40817168, units='unitless') + prob.model.set_input_defaults( + Aircraft.Wing.ASPECT_RATIO_REFERENCE, val=9.40817168, units='unitless' + ) # it was 0.01 in high_wing_single_aisle.csv + prob.model.set_input_defaults( + Aircraft.Wing.STRUT_BRACING_FACTOR, val=1.01, units='unitless' + ) + prob.model.set_input_defaults( + Aircraft.Wing.AEROELASTIC_TAILORING_FACTOR, val=0.01, units='unitless' + ) + prob.model.set_input_defaults(Aircraft.Wing.THICKNESS_TO_CHORD, val=0.1, units='unitless') + prob.model.set_input_defaults( + Aircraft.Wing.THICKNESS_TO_CHORD_REFERENCE, val=0.01, units='unitless' + ) + + setup_model_options(prob, options) + prob.setup(check=False, force_alloc_complex=True) + + prob.set_val(Aircraft.Engine.POD_MASS, np.array([5587.40886136]), units='lbm') + prob.set_val(Aircraft.Engine.WING_LOCATIONS, [0.22]) + + prob.run_model() + + bending_mat_factor = prob.get_val(Aircraft.Wing.BENDING_MATERIAL_FACTOR) + assert_near_equal(bending_mat_factor, 2.9367664396, tolerance=1e-9) + + pod_inertia_factor = prob.get_val(Aircraft.Wing.ENG_POD_INERTIA_FACTOR) + assert_near_equal(pod_inertia_factor, 0.9772998541, tolerance=1e-9) + + options.set_val(Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL, val=1, units='unitless') + setup_model_options(prob, options) + prob.setup(check=False, force_alloc_complex=True) + prob.run_model() + + bending_mat_factor = prob.get_val(Aircraft.Wing.BENDING_MATERIAL_FACTOR) + assert_near_equal(bending_mat_factor, 1.51247369, tolerance=1e-9) + + pod_inertia_factor = prob.get_val(Aircraft.Wing.ENG_POD_INERTIA_FACTOR) + assert_near_equal(pod_inertia_factor, 1.0, tolerance=1e-9) + + options.set_val(Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL, val=1.5, units='unitless') + setup_model_options(prob, options) + prob.setup(check=False, force_alloc_complex=True) + prob.run_model() + + bending_mat_factor = prob.get_val(Aircraft.Wing.BENDING_MATERIAL_FACTOR) + assert_near_equal(bending_mat_factor, 1.94203494, tolerance=1e-9) + + pod_inertia_factor = prob.get_val(Aircraft.Wing.ENG_POD_INERTIA_FACTOR) + assert_near_equal(pod_inertia_factor, 1.0, tolerance=1e-9) + + options.set_val(Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL, val=2.5, units='unitless') + setup_model_options(prob, options) + prob.setup(check=False, force_alloc_complex=True) + prob.run_model() + + bending_mat_factor = prob.get_val(Aircraft.Wing.BENDING_MATERIAL_FACTOR) + assert_near_equal(bending_mat_factor, 2.61652282, tolerance=1e-9) + + pod_inertia_factor = prob.get_val(Aircraft.Wing.ENG_POD_INERTIA_FACTOR) + assert_near_equal(pod_inertia_factor, 1.0, tolerance=1e-9) + def test_IO(self): assert_match_varnames(self.prob.model) @@ -607,7 +704,8 @@ def test_case1(self): pod_inertia_expected = 1.0 assert_near_equal(BENDING_MATERIAL_FACTOR, BENDING_MATERIAL_FACTOR_expected, tolerance=1e-9) assert_near_equal(prob.get_val('calculated_wing_area'), 5399.4057051, tolerance=1e-9) - # current BWB data set does not check the following + # For BWB with non wing engines, inertia_factor_prod = 1 always. + # Need a different dataset with wing engines to check this one assert_near_equal(pod_inertia, pod_inertia_expected, tolerance=1e-9) def test_case2(self): @@ -693,7 +791,6 @@ def test_case2(self): pod_inertia_expected = 1.0 assert_near_equal(BENDING_MATERIAL_FACTOR, BENDING_MATERIAL_FACTOR_expected, tolerance=1e-9) assert_near_equal(prob.get_val('calculated_wing_area'), 4151.88659141, tolerance=1e-9) - # current BWB data set does not check the following assert_near_equal(pod_inertia, pod_inertia_expected, tolerance=1e-9) @@ -701,4 +798,4 @@ def test_case2(self): unittest.main() # test = DetailedWingBendingTest() # test.setUp() - # test.test_extreme_engine_loc() + # test.test_intensity_factor() diff --git a/aviary/subsystems/mass/flops_based/wing_detailed.py b/aviary/subsystems/mass/flops_based/wing_detailed.py index e9e6f1c9b..6611d0c0c 100644 --- a/aviary/subsystems/mass/flops_based/wing_detailed.py +++ b/aviary/subsystems/mass/flops_based/wing_detailed.py @@ -5,7 +5,45 @@ from aviary.variable_info.enums import Verbosity from aviary.variable_info.functions import add_aviary_input, add_aviary_option, add_aviary_output -from aviary.variable_info.variables import Aircraft, Mission, Settings +from aviary.variable_info.variables import Aircraft, Settings + + +def load_intensity_by_factor(load_dist_factor, intn_stations): + """ + Calculate load intensity at wing integration stations for the given load_dist_factor. + + Parameters: + load_dist_factor (float): 0 or 1 <= load_dist_factor <= 3. + 1.0 : triangular distribution + 2.0 : elliptical distribution (default) + 3.0 : rectangular distribution + 1.0-2.0 : blend of triangular and elliptical + 2.0-3.0 : blend of elliptical and rectangular + intn_stations (array): integration stations + + Note: In FLOPS, load_dist_factor can be 0 in which case load intensities are computed + based on input pressure distribution. + """ + if load_dist_factor == 1: + load_intensity = 1.0 - intn_stations + elif load_dist_factor == 2: + load_intensity = np.sqrt(1.0 - intn_stations**2) + elif load_dist_factor == 3: + load_intensity = np.ones(len(intn_stations)) + elif 1 < load_dist_factor < 2: + load_intensity = (load_dist_factor - 1.0) * np.sqrt(1.0 - intn_stations**2) + ( + 2.0 - load_dist_factor + ) * (1.0 - intn_stations) + elif 2 < load_dist_factor < 3: + load_intensity = (3.0 - load_dist_factor) * np.sqrt(1.0 - intn_stations**2) + ( + load_dist_factor - 2.0 + ) * np.ones(len(intn_stations)) + else: + raise om.AnalysisError( + f'{load_dist_factor} is not a valid value for ' + f'{Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL}, it must be between 1 and 3 (inclusive).' + ) + return load_intensity class DetailedWingBendingFact(om.ExplicitComponent): @@ -85,15 +123,6 @@ def compute(self, inputs, outputs): num_wing_engines = self.options[Aircraft.Engine.NUM_WING_ENGINES] num_engine_type = len(num_wing_engines) - # TODO: Support all options for this parameter. - # 0.0 : input distribution - # 1.0 : triangular distribution - # 2.0 : elliptical distribution (default) - # 3.0 : rectangular distribution - # 1.0-2.0 : blend of triangular and elliptical - # 2.0-3.0 : blend of elliptical and rectangular - load_distribution_factor = self.options[Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL] - load_path_sweep = inputs[Aircraft.Wing.LOAD_PATH_SWEEP_DISTRIBUTION] thickness_to_chord = inputs[Aircraft.Wing.THICKNESS_TO_CHORD_DISTRIBUTION] chord = inputs[Aircraft.Wing.CHORD_PER_SEMISPAN_DISTRIBUTION] @@ -138,18 +167,8 @@ def compute(self, inputs, outputs): (dy[1:] + 2.0 * integration_stations[1:-1]) * dy[1:] * sweep_int_stations[1:-1] ) - # TODO: add all load_distribution_factor options - if load_distribution_factor == 1: - load_intensity = 1.0 - integration_stations - elif load_distribution_factor == 2: - load_intensity = np.sqrt(1.0 - integration_stations**2) - elif load_distribution_factor == 3: - load_intensity = np.ones(num_integration_stations + 1) - else: - raise om.AnalysisError( - f'{load_distribution_factor} is not a valid value for ' - f'{Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL}, it must be "1", "2", or "3".' - ) + load_distribution_factor = self.options[Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL] + load_intensity = load_intensity_by_factor(load_distribution_factor, integration_stations) chord_interp = InterpND( method='slinear', points=(inp_stations), x_interp=integration_stations @@ -341,8 +360,6 @@ def setup_partials(self): def compute(self, inputs, outputs): verbosity = self.options[Settings.VERBOSITY] num_integration_stations = self.options[Aircraft.Wing.NUM_INTEGRATION_STATIONS] - num_wing_engines = self.options[Aircraft.Engine.NUM_WING_ENGINES] - num_engine_type = len(num_wing_engines) width = inputs[Aircraft.Fuselage.MAX_WIDTH][0] wingspan = inputs[Aircraft.Wing.SPAN][0] rate_span = (wingspan - width) / wingspan @@ -370,15 +387,6 @@ def compute(self, inputs, outputs): # For BWB, always start from inp_stations_mod[1], not inp_stations_mod[0] inp_stations_mod = inp_stations_mod[1:] - # TODO: Support all options for this parameter. - # 0.0 : input distribution - # 1.0 : triangular distribution - # 2.0 : elliptical distribution (default) - # 3.0 : rectangular distribution - # 1.0-2.0 : blend of triangular and elliptical - # 2.0-3.0 : blend of elliptical and rectangular - load_distribution_factor = self.options[Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL] - load_path_sweep = inputs['BWB_LOAD_PATH_SWEEP_DISTRIBUTION'] load_path_sweep_mod = np.array(load_path_sweep[1:]) @@ -439,18 +447,8 @@ def compute(self, inputs, outputs): (dy[0:] + 2.0 * integration_stations[0:-1]) * dy[0:] * sweep_int_stations[0:-1] ) - # TODO: add all load_distribution_factor options - if load_distribution_factor == 1: - load_intensity = 1.0 - integration_stations - elif load_distribution_factor == 2: - load_intensity = np.sqrt(1.0 - integration_stations**2) - elif load_distribution_factor == 3: - load_intensity = np.ones(num_integration_stations + 1) - else: - raise om.AnalysisError( - f'{load_distribution_factor} is not a valid value for ' - f'{Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL}, it must be "1", "2", or "3".' - ) + load_distribution_factor = self.options[Aircraft.Wing.LOAD_DISTRIBUTION_CONTROL] + load_intensity = load_intensity_by_factor(load_distribution_factor, integration_stations) chord_interp = InterpND( method='slinear', points=(inp_stations_mod), x_interp=integration_stations @@ -523,13 +521,14 @@ def compute(self, inputs, outputs): outputs[Aircraft.Wing.BENDING_MATERIAL_FACTOR] = bt + num_wing_engines = self.options[Aircraft.Engine.NUM_WING_ENGINES] + num_engine_type = len(num_wing_engines) engine_locations = inputs[Aircraft.Engine.WING_LOCATIONS] gross_mass = inputs[Aircraft.Design.GROSS_MASS] # NOTE pod mass assumed the same for wing/non-wing mounted engines, only using # wing mounted pods here pod_mass = inputs[Aircraft.Engine.POD_MASS] if np.sum(num_wing_engines) > 0: - # TODO: the rest is not checked. inertia_factor = np.zeros(num_engine_type, dtype=chord.dtype) eel = np.zeros(len(dy) + 1, dtype=chord.dtype) diff --git a/aviary/variable_info/variable_meta_data.py b/aviary/variable_info/variable_meta_data.py index 4ef9a792b..fcd24b9b8 100644 --- a/aviary/variable_info/variable_meta_data.py +++ b/aviary/variable_info/variable_meta_data.py @@ -5690,7 +5690,7 @@ 'FLOPS': 'WTIN.PDIST', # ['&DEFINE.WTIN.PDIST', 'WDEF.PDIST'], }, units='unitless', - desc='controls spatial distribution of integration stations for detailed wing', + desc='controls spatial distribution of integration stations for detailed wing, in [1, 3]', default_value=2.0, option=True, )