[Draft] FLOPS based BWB mission bench test

aviary/mission/height_energy_problem_configurator.py

@@ -394,7 +394,11 @@ def add_post_mission_systems(self, aviary_group):
             )
 
         if aviary_group.post_mission_info['include_landing']:
-            self._add_landing_systems(aviary_group)
+            if 'aircraft:wing:area' in aviary_group.aviary_inputs:
+                self._add_landing_systems(aviary_group)
+            else:
+                print('Aircraft.Wing.AREA is not given. Set include_landing = False')
+                aviary_group.post_mission_info['include_landing'] = False
 
         aviary_group.add_subsystem(
             'range_constraint',

aviary/models/aircraft/blended_wing_body/bwb_detailed_FLOPS.csv

@@ -152,7 +152,6 @@ aircraft:wing:var_sweep_mass_penalty,0.0,unitless
 aircraft:wing:wetted_area_scaler,1.0,unitless
 mission:constraints:max_mach,0.85,unitless
 mission:design:gross_mass,874099,lbm
-mission:design:lift_coefficient,-1.0,unitless
 mission:design:range,7750,NM
 mission:design:thrust_takeoff_per_eng,0.25,lbf
 mission:landing:initial_velocity,140,ft/s
@@ -164,6 +163,7 @@ mission:summary:fuel_flow_scaler,1.0,unitless
 settings:aerodynamics_method,FLOPS,unitless
 settings:equations_of_motion,height_energy,unitless
 settings:mass_method,FLOPS,unitless
+settings:verbosity,1,unitless
 
 # Unconverted Values
 AERIN.FLLDG,11000

aviary/models/aircraft/blended_wing_body/bwb_simple_FLOPS.csv

@@ -151,7 +151,6 @@ aircraft:wing:var_sweep_mass_penalty,0.0,unitless
 aircraft:wing:wetted_area_scaler,1.0,unitless
 mission:constraints:max_mach,0.85,unitless
 mission:design:gross_mass,874099,lbm
-mission:design:lift_coefficient,-1.0,unitless
 mission:design:range,7750,NM
 mission:design:thrust_takeoff_per_eng,0.25,lbf
 mission:landing:initial_velocity,140,ft/s
@@ -163,6 +162,7 @@ mission:summary:fuel_flow_scaler,1.0,unitless
 settings:aerodynamics_method,FLOPS,unitless
 settings:equations_of_motion,height_energy,unitless
 settings:mass_method,FLOPS,unitless
+settings:verbosity,1,unitless
 
 # Unconverted Values
 AERIN.FLLDG,11000

aviary/subsystems/aerodynamics/aerodynamics_builder.py

@@ -494,7 +494,15 @@ def get_parameters(self, aviary_inputs=None, **kwargs):
                         params[Aircraft.Design.LIFT_DEPENDENT_DRAG_POLAR] = opts
 
             if method == 'computed':
-                for var in COMPUTED_CORE_INPUTS:
+                try:
+                    design_type = aviary_inputs.get_val(Aircraft.Design.TYPE)
+                except:
+                    design_type = AircraftTypes.TRANSPORT
+                if design_type is AircraftTypes.BLENDED_WING_BODY:
+                    core_inputs_computed = COMPUTED_CORE_INPUTS_BWB
+                else:
+                    core_inputs_computed = COMPUTED_CORE_INPUTS
+                for var in core_inputs_computed:
                     meta = _MetaData[var]
 
                     val = meta['default_value']
@@ -622,7 +630,10 @@ def get_parameters(self, aviary_inputs=None, **kwargs):
                     'tabular_cruise, low_speed, tabular_low_speed)'
                 )
 
-            design_type = aviary_inputs.get_val(Aircraft.Design.TYPE)
+            try:
+                design_type = aviary_inputs.get_val(Aircraft.Design.TYPE)
+            except:
+                design_type = AircraftTypes.TRANSPORT
 
             if design_type is AircraftTypes.BLENDED_WING_BODY:
                 all_vars.add(Aircraft.Fuselage.LIFT_CURVE_SLOPE_MACH0)
@@ -729,6 +740,37 @@ def report(self, prob, reports_folder, **kwargs):
     Mission.Design.MACH,
 ]
 
+COMPUTED_CORE_INPUTS_BWB = [
+    Aircraft.Design.BASE_AREA,
+    Aircraft.Design.LIFT_DEPENDENT_DRAG_COEFF_FACTOR,
+    Aircraft.Design.SUBSONIC_DRAG_COEFF_FACTOR,
+    Aircraft.Design.SUPERSONIC_DRAG_COEFF_FACTOR,
+    Aircraft.Design.ZERO_LIFT_DRAG_COEFF_FACTOR,
+    Aircraft.Fuselage.CHARACTERISTIC_LENGTH,
+    Aircraft.Fuselage.CROSS_SECTION,
+    Aircraft.Fuselage.DIAMETER_TO_WING_SPAN,
+    Aircraft.Fuselage.FINENESS,
+    Aircraft.Fuselage.LAMINAR_FLOW_LOWER,
+    Aircraft.Fuselage.LAMINAR_FLOW_UPPER,
+    Aircraft.Fuselage.LENGTH_TO_DIAMETER,
+    Aircraft.Fuselage.WETTED_AREA,
+    Aircraft.Wing.AREA,
+    Aircraft.Wing.ASPECT_RATIO,
+    Aircraft.Wing.CHARACTERISTIC_LENGTH,
+    Aircraft.Wing.FINENESS,
+    Aircraft.Wing.LAMINAR_FLOW_LOWER,
+    Aircraft.Wing.LAMINAR_FLOW_UPPER,
+    Aircraft.Wing.MAX_CAMBER_AT_70_SEMISPAN,
+    Aircraft.Wing.SPAN_EFFICIENCY_FACTOR,
+    Aircraft.Wing.SWEEP,
+    Aircraft.Wing.TAPER_RATIO,
+    Aircraft.Wing.THICKNESS_TO_CHORD,
+    Aircraft.Wing.WETTED_AREA,
+    # Mission.Summary.GROSS_MASS,
+    Mission.Design.LIFT_COEFFICIENT,
+    Mission.Design.MACH,
+]
+
 TABULAR_CORE_INPUTS = [
     Aircraft.Wing.AREA,
     Aircraft.Design.SUBSONIC_DRAG_COEFF_FACTOR,

aviary/subsystems/aerodynamics/flops_based/skin_friction.py

@@ -1,6 +1,7 @@
 import numpy as np
 import openmdao.api as om
 
+from aviary.variable_info.enums import AircraftTypes
 from aviary.variable_info.functions import add_aviary_input, add_aviary_option
 from aviary.variable_info.variables import Aircraft, Dynamic
 
@@ -36,17 +37,23 @@ def initialize(self):
             desc='The number of points at which the cross product is computed.',
         )
 
+        add_aviary_option(self, Aircraft.Design.TYPE)
         add_aviary_option(self, Aircraft.Engine.NUM_ENGINES)
         add_aviary_option(self, Aircraft.Fuselage.NUM_FUSELAGES)
         add_aviary_option(self, Aircraft.VerticalTail.NUM_TAILS)
 
     def setup(self):
         nn = self.options['num_nodes']
+        design_type = self.options[Aircraft.Design.TYPE]
         num_engines = self.options[Aircraft.Engine.NUM_ENGINES]
         num_fuselages = self.options[Aircraft.Fuselage.NUM_FUSELAGES]
         num_tails = self.options[Aircraft.VerticalTail.NUM_TAILS]
 
-        self.nc = nc = 2 + num_tails + num_fuselages + int(sum(num_engines))
+        if design_type is AircraftTypes.BLENDED_WING_BODY:
+            # No horizontal tail for BWB
+            self.nc = nc = 1 + num_tails + num_fuselages + int(sum(num_engines))
+        else:
+            self.nc = nc = 2 + num_tails + num_fuselages + int(sum(num_engines))
 
         # Simulation inputs
         add_aviary_input(self, Dynamic.Atmosphere.TEMPERATURE, shape=nn, units='degR')

aviary/subsystems/aerodynamics/flops_based/skin_friction_drag.py

@@ -1,6 +1,7 @@
 import numpy as np
 import openmdao.api as om
 
+from aviary.variable_info.enums import AircraftTypes
 from aviary.variable_info.functions import add_aviary_input, add_aviary_option, get_units
 from aviary.variable_info.variables import Aircraft
 
@@ -33,6 +34,7 @@ def initialize(self):
             desc='The number of points at which the cross product is computed.',
         )
 
+        add_aviary_option(self, Aircraft.Design.TYPE)
         add_aviary_option(self, Aircraft.Engine.NUM_ENGINES)
         add_aviary_option(self, Aircraft.Fuselage.NUM_FUSELAGES)
         add_aviary_option(self, Aircraft.VerticalTail.NUM_TAILS)
@@ -47,12 +49,17 @@ def initialize(self):
 
     def setup(self):
         nn = self.options['num_nodes']
+        design_type = self.options[Aircraft.Design.TYPE]
 
         nvtail = self.options[Aircraft.VerticalTail.NUM_TAILS]
         nfuse = self.options[Aircraft.Fuselage.NUM_FUSELAGES]
         num_engines = self.options[Aircraft.Engine.NUM_ENGINES]
 
-        self.nc = nc = 2 + nvtail + nfuse + int(sum(num_engines))
+        if design_type is AircraftTypes.BLENDED_WING_BODY:
+            # No horizontal tail for BWB
+            self.nc = nc = 1 + nvtail + nfuse + int(sum(num_engines))
+        else:
+            self.nc = nc = 2 + nvtail + nfuse + int(sum(num_engines))
 
         # Computed by other components in drag group.
         self.add_input('skin_friction_coeff', np.zeros((nn, nc)), units='unitless')

aviary/subsystems/geometry/flops_based/canard.py

@@ -12,6 +12,7 @@
 class Canard(om.ExplicitComponent):
     """Calculate the wetted area of canard."""
 
+    # TODO: what is it for?
     def initialize(self):
         self.options.declare(
             'aviary_options',

aviary/subsystems/geometry/flops_based/fuselage.py

@@ -167,7 +167,7 @@ def compute(self, inputs, outputs):
             if verbosity > Verbosity.BRIEF:
                 raise UserWarning(
                     'Passenger compartment lenght is longer than recommended maximum'
-                    ' length (of 190 ft). Suggest using detailed layout algorithm.'
+                    ' length. Suggest using detailed layout algorithm.'
                 )
 
         outputs[Aircraft.Fuselage.PASSENGER_COMPARTMENT_LENGTH] = pax_compart_length
@@ -480,8 +480,8 @@ def compute(self, inputs, outputs):
 
         # Enforce maximum number of bays
         num_bays_max = self.options[Aircraft.BWB.MAX_NUM_BAYS]
-        num_bays = int(0.5 + max_width.real / bay_width_max.real)
-        if num_bays > num_bays_max and num_bays_max > 0:
+        num_bays = int(0.5 + max_width / bay_width_max)
+        if num_bays.real > num_bays_max and num_bays_max > 0:
             num_bays = num_bays_max
         outputs[Aircraft.BWB.NUM_BAYS] = smooth_int_tanh(num_bays, mu=20.0)
 
@@ -657,27 +657,30 @@ def compute(self, inputs, outputs):
             pax_compart_length = root_chord + tan_sweep * max_width / 2.0
 
             # Enforce maximum number of bays
-            z = 0.5 + max_width / bay_width_max
-            z = z[0]
-            num_bays = int(z.real)
-            if num_bays > num_bays_max and num_bays_max > 0:
+            num_bays_tmp = 0.5 + max_width / bay_width_max
+            if num_bays_tmp[0].real > num_bays_max and num_bays_max > 0:
                 num_bays = num_bays_max
+            else:
+                num_bays = int(num_bays_tmp[0].real)
 
             # Enforce maximum bay width
             bay_width = max_width / num_bays
             if bay_width > bay_width_max:
                 bay_width = bay_width_max
-                num_bays = int(0.999 + max_width / bay_width)
-                if num_bays > num_bays_max and num_bays_max > 0:
+                num_bays_tmp = 0.999 + max_width / bay_width
+                if num_bays_tmp.real > num_bays_max and num_bays_max > 0:
                     num_bays = num_bays_max
                     max_width = bay_width_max * bay_width
                     pax_compart_length = area_cabin / max_width + tan_sweep * max_width / 4.0
                     root_chord = pax_compart_length - tan_sweep * max_width / 2.0
+                else:
+                    num_bays = smooth_int_tanh(num_bays_tmp, mu=40.0)
+
             # If number of bays has changed, recalculate cabin area
 
         length = pax_compart_length / rear_spar_percent_chord
         max_height = height_to_width * length
-        outputs[Aircraft.BWB.NUM_BAYS] = smooth_int_tanh(num_bays, mu=20.0)
+        outputs[Aircraft.BWB.NUM_BAYS] = num_bays
 
         outputs[Aircraft.Fuselage.LENGTH] = length
         outputs[Aircraft.Fuselage.PASSENGER_COMPARTMENT_LENGTH] = pax_compart_length

aviary/subsystems/test/test_flops_based_premission.py

@@ -503,7 +503,7 @@ def test_case_geom(self):
         assert_near_equal(prob[Aircraft.Wing.ROOT_CHORD], 63.96, tol)
         assert_near_equal(prob[Aircraft.Fuselage.CABIN_AREA], 5173.187202504683, tol)
         assert_near_equal(prob[Aircraft.Fuselage.MAX_HEIGHT], 15.125, tol)
-        assert_near_equal(prob[Aircraft.BWB.NUM_BAYS], 5.0, tol)
+        assert_near_equal(prob[Aircraft.BWB.NUM_BAYS], 5.0, 1e-4)
         # BWBFuselagePrelim
         assert_near_equal(prob[Aircraft.Fuselage.REF_DIAMETER], 39.8525, tol)
         assert_near_equal(prob[Aircraft.Fuselage.PLANFORM_AREA], 7390.267432149546, tol)

aviary/utils/math.py

@@ -215,7 +215,7 @@ def smooth_int_tanh(x, mu=10.0):
     """
     Smooth approximation of int(x) using tanh.
     """
-    f = np.floor(x)
+    f = np.floor(x.real) + x.imag * 1j
     frac = x - f
     t = np.tanh(mu * (frac - 0.5))
     s = 0.5 * (t + 1)
@@ -228,7 +228,7 @@ def d_smooth_int_tanh(x, mu=10.0):
     Smooth approximation of int(x) using tanh.
     Returns (y, dy_dx).
     """
-    f = np.floor(x)
+    f = np.floor(x) + x.imag * 1j
     frac = x - f
     t = np.tanh(mu * (frac - 0.5))
     dy_dx = 0.5 * mu * (1 - t**2)