Repository: fast-aircraft-design/FAST-OAD_CS25 (language: python)

Issue:
Add options to default load factor calculation
The current load factor calculation in [cs25.py](https://github.com/fast-aircraft-design/FAST-OAD_CS25/blob/main/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py) is the maximum between a fixed 2.5g manoeuvre load and the calculated gust load. The sizing loads then take into account the fuel load alleviation.

We'd like to be able to let the user modify the manoeuvre load and manipulate the gust load through a tuning coefficient so as to simulate Gust Load Alleviation and Manoeuvre load Alleviation functions. Additionally, for dry wing, we'd like the option to deactivate the fuel load alleviation.

This can be done by adding an option and two user inputs.

Error in gust load calculation
In class [cs_25.py](https://github.com/fast-aircraft-design/FAST-OAD_CS25/blob/main/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py#L128) line [128 ](https://github.com/fast-aircraft-design/FAST-OAD_CS25/blob/f1b0156a53b57ac9b99efd7ca0bcb441d7b5a445/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py#L128) there is an error in the formula for calculating the gust load:
`
n_gust = 1 + (sea_level_density / 2 / 9.81) * k_g * u_gust * (
            vc_eas * cl_alpha / mass / wing_area
        )
`

**Expected behavior**
The last two terms should be in parenthesis. The wing loading is at the denominator as in the book AP avion p85 (french edition). The fix is:

`
n_gust = 1 + (sea_level_density / 2 / 9.81) * k_g * u_gust * (
            vc_eas * cl_alpha / (mass / wing_area)
        )
`

The error was giving loads very close to 1.0 for an A320 like while it should be around 2.5. As the effective minimum load cannot be lower than 2.5, this error has effect only for low wing loading.


Error in gust load calculation
In class [cs_25.py](https://github.com/fast-aircraft-design/FAST-OAD_CS25/blob/main/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py#L128) line [128 ](https://github.com/fast-aircraft-design/FAST-OAD_CS25/blob/f1b0156a53b57ac9b99efd7ca0bcb441d7b5a445/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py#L128) there is an error in the formula for calculating the gust load:
`
n_gust = 1 + (sea_level_density / 2 / 9.81) * k_g * u_gust * (
            vc_eas * cl_alpha / mass / wing_area
        )
`

**Expected behavior**
The last two terms should be in parenthesis. The wing loading is at the denominator as in the book AP avion p85 (french edition). The fix is:

`
n_gust = 1 + (sea_level_density / 2 / 9.81) * k_g * u_gust * (
            vc_eas * cl_alpha / (mass / wing_area)
        )
`

The error was giving loads very close to 1.0 for an A320 like while it should be around 2.5. As the effective minimum load cannot be lower than 2.5, this error has effect only for low wing loading.




Discussion / hints:

Fix patch:
diff --git a/src/fastoad_cs25/models/variable_descriptions.txt b/src/fastoad_cs25/models/variable_descriptions.txt
index 14b7039..a0e0d38 100644
--- a/src/fastoad_cs25/models/variable_descriptions.txt
+++ b/src/fastoad_cs25/models/variable_descriptions.txt
@@ -181,6 +181,8 @@ data:load_case:lc1:altitude || altitude for sizing load case 1 (gust with minimu
 data:load_case:lc2:U_gust || gust vertical speed for sizing load case 2 (gust with maximum aircraft mass)
 data:load_case:lc2:Vc_EAS || equivalent air speed for sizing load case 2 (gust with maximum aircraft mass)
 data:load_case:lc2:altitude || altitude for sizing load case 2 (gust with maximum aircraft mass)
+data:load_case:manoeuvre_load_factor || limit positive manoeuvring load factor
+data:load_case:gust_intensity || a factor for gust load alleviation, 1.0: full gust load is applied, 0.0: gust load is fully alleviated
 data:propulsion:MTO_thrust || maximum thrust of one engine at sea level
 data:propulsion:SFC || Specific Fuel Consumption (can be a vector)
 data:propulsion:altitude || altitude for propulsion calculation (can be a vector)
diff --git a/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py b/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py
index c9a15ac..cc2c489 100644
--- a/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py
+++ b/src/fastoad_cs25/models/weight/mass_breakdown/cs25.py
@@ -15,15 +15,15 @@
 #  along with this program.  If not, see <https://www.gnu.org/licenses/>.
 
 import numpy as np
+import openmdao.api as om
 from fastoad.module_management.service_registry import RegisterSubmodel
-from openmdao.core.explicitcomponent import ExplicitComponent
 from stdatm import Atmosphere
 
 from .constants import SERVICE_GUST_LOADS
 
 
 @RegisterSubmodel(SERVICE_GUST_LOADS, "fastoad.submodel.gust_loads")
-class Loads(ExplicitComponent):
+class Loads(om.ExplicitComponent):
     """
     Computes gust load cases
 
@@ -34,6 +34,15 @@ class Loads(ExplicitComponent):
 
     """
 
+    def initialize(self):
+        self.options.declare(
+            "fuel_load_alleviation",
+            types=bool,
+            default=True,
+            desc="If False this simulates a dry wing,"
+            "i.e. the sizing load 2 does not take into account the fuel weight.",
+        )
+
     def setup(self):
         self.add_input("data:geometry:wing:area", val=np.nan, units="m**2")
         self.add_input("data:geometry:wing:span", val=np.nan, units="m")
@@ -47,6 +56,8 @@ def setup(self):
         self.add_input("data:load_case:lc2:U_gust", val=np.nan, units="m/s")
         self.add_input("data:load_case:lc2:altitude", val=np.nan, units="ft")
         self.add_input("data:load_case:lc2:Vc_EAS", val=np.nan, units="m/s")
+        self.add_input("data:load_case:manoeuvre_load_factor", val=2.5)
+        self.add_input("data:load_case:gust_intensity", val=1.0)
 
         self.add_output("data:mission:sizing:cs25:load_factor_1")
         self.add_output("data:mission:sizing:cs25:load_factor_2")
@@ -72,6 +83,8 @@ def compute(self, inputs, outputs, discrete_inputs=None, discrete_outputs=None):
         u_gust2 = inputs["data:load_case:lc2:U_gust"]
         alt_2 = inputs["data:load_case:lc2:altitude"]
         vc_eas2 = inputs["data:load_case:lc2:Vc_EAS"]
+        n_manoeuvre = inputs["data:load_case:manoeuvre_load_factor"]
+        gust_intensity = inputs["data:load_case:gust_intensity"]
 
         # calculation of mean geometric chord
         chord_geom = wing_area / span
@@ -88,7 +101,7 @@ def compute(self, inputs, outputs, discrete_inputs=None, discrete_outputs=None):
             cl_alpha,
             u_gust1,
         )
-        n1 = 1.5 * max(2.5, n_gust_1)
+        n1 = 1.5 * max(n_manoeuvre, n_gust_1 * gust_intensity)
         n1m1 = n1 * m1
 
         # load case #2
@@ -102,9 +115,13 @@ def compute(self, inputs, outputs, discrete_inputs=None, discrete_outputs=None):
             cl_alpha,
             u_gust2,
         )
-        n2 = 1.5 * max(2.5, n_gust_2)
-        mcv = min(0.8 * mfw, mtow - mzfw)
-        n2m2 = n2 * (mtow - 0.55 * mcv)
+        n2 = 1.5 * max(n_manoeuvre, n_gust_2 * gust_intensity)
+
+        if not self.options["fuel_load_alleviation"]:
+            n2m2 = n2 * mtow
+        else:
+            mcv = min(0.8 * mfw, mtow - mzfw)
+            n2m2 = n2 * (mtow - 0.55 * mcv)
 
         outputs["data:mission:sizing:cs25:load_factor_1"] = n1
         outputs["data:mission:sizing:cs25:load_factor_2"] = n2
@@ -129,7 +146,7 @@ def __n_gust(mass, wing_area, rho, sea_level_density, chord_geom, vc_eas, cl_alp
         mu_g = 2 * mass / rho / wing_area / chord_geom / cl_alpha
         k_g = 0.88 * mu_g / (5.3 + mu_g)  # attenuation factor
         n_gust = 1 + (sea_level_density / 2 / 9.81) * k_g * u_gust * (
-            vc_eas * cl_alpha / mass / wing_area
+            vc_eas * cl_alpha / (mass / wing_area)
         )
 
         return n_gust
diff --git a/src/fastoad_cs25/notebooks/01_tutorial/data/CeRAS01_baseline.xml b/src/fastoad_cs25/notebooks/01_tutorial/data/CeRAS01_baseline.xml
index 9bddfd6..33ef5d7 100644
--- a/src/fastoad_cs25/notebooks/01_tutorial/data/CeRAS01_baseline.xml
+++ b/src/fastoad_cs25/notebooks/01_tutorial/data/CeRAS01_baseline.xml
@@ -334,6 +334,7 @@
             </taxi>
         </propulsion>
         <load_case>
+            <gust_intensity>0.0</gust_intensity>
             <lc1>
                 <U_gust units="m/s">
                     15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>
diff --git a/src/fastoad_cs25/notebooks/02_CeRAS_case_study/data/CeRAS_reference_data.xml b/src/fastoad_cs25/notebooks/02_CeRAS_case_study/data/CeRAS_reference_data.xml
index 1639f9e..ef59349 100644
--- a/src/fastoad_cs25/notebooks/02_CeRAS_case_study/data/CeRAS_reference_data.xml
+++ b/src/fastoad_cs25/notebooks/02_CeRAS_case_study/data/CeRAS_reference_data.xml
@@ -331,6 +331,7 @@
             </taxi>
         </propulsion>
         <load_case>
+            <gust_intensity>0.0</gust_intensity>
             <lc1>
                 <U_gust units="m/s" is_input="True">
                     15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>


Test patch:
diff --git a/src/fastoad_cs25/models/weight/mass_breakdown/tests/data/mass_breakdown_inputs.xml b/src/fastoad_cs25/models/weight/mass_breakdown/tests/data/mass_breakdown_inputs.xml
index e7b8681..367e7eb 100644
--- a/src/fastoad_cs25/models/weight/mass_breakdown/tests/data/mass_breakdown_inputs.xml
+++ b/src/fastoad_cs25/models/weight/mass_breakdown/tests/data/mass_breakdown_inputs.xml
@@ -91,6 +91,7 @@
       <MTO_thrust units="N">117880.0<!--maximum thrust of one engine at sea level--></MTO_thrust>
     </propulsion>
     <load_case>
+      <gust_intensity>0.0</gust_intensity>
       <lc1>
         <U_gust units="m/s">15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>
         <Vc_EAS units="kn">375.0<!--equivalent air speed for sizing load case 1 (gust with minimum aircraft mass)--></Vc_EAS>
diff --git a/src/fastoad_cs25/models/weight/mass_breakdown/tests/test_mass_breakdown.py b/src/fastoad_cs25/models/weight/mass_breakdown/tests/test_mass_breakdown.py
index 4a11e66..273794e 100644
--- a/src/fastoad_cs25/models/weight/mass_breakdown/tests/test_mass_breakdown.py
+++ b/src/fastoad_cs25/models/weight/mass_breakdown/tests/test_mass_breakdown.py
@@ -96,6 +96,7 @@ def test_compute_loads():
         "data:load_case:lc2:Vc_EAS",
     ]
     ivc = get_indep_var_comp(input_list)
+    ivc.add_output("data:load_case:gust_intensity", val=0.5)
     problem = run_system(Loads(), ivc)
 
     n1m1 = problem["data:mission:sizing:cs25:sizing_load_1"]
@@ -107,6 +108,37 @@ def test_compute_loads():
     assert n1 == pytest.approx(3.75, abs=0.01)
     assert n2 == pytest.approx(3.75, abs=0.01)
 
+    # Now without fuel alleviation
+    problem = run_system(Loads(fuel_load_alleviation=False), ivc)
+
+    n1m1 = problem["data:mission:sizing:cs25:sizing_load_1"]
+    n2m2 = problem["data:mission:sizing:cs25:sizing_load_2"]
+
+    assert n1m1 == pytest.approx(240968, abs=10)
+    assert n2m2 == pytest.approx(284242, abs=10)
+
+    # Now with different MLA
+    ivc.add_output("data:load_case:manoeuvre_load_factor", val=2.25)
+    problem = run_system(Loads(), ivc)
+
+    n1 = problem["data:mission:sizing:cs25:load_factor_1"]
+    n2 = problem["data:mission:sizing:cs25:load_factor_2"]
+    assert n1 == pytest.approx(3.375, abs=0.01)
+    assert n2 == pytest.approx(3.375, abs=0.01)
+
+    # Now without gust load alleviation
+    ivc = get_indep_var_comp(input_list)
+    problem = run_system(Loads(), ivc)
+
+    n1 = problem["data:mission:sizing:cs25:load_factor_1"]
+    n2 = problem["data:mission:sizing:cs25:load_factor_2"]
+    n1m1 = problem["data:mission:sizing:cs25:sizing_load_1"]
+    n2m2 = problem["data:mission:sizing:cs25:sizing_load_2"]
+    assert n1 == pytest.approx(4.198, abs=0.01)
+    assert n2 == pytest.approx(3.81, abs=0.01)
+    assert n1m1 == pytest.approx(269789, abs=10)
+    assert n2m2 == pytest.approx(258500, abs=10)
+
 
 def test_compute_wing_weight():
     """Tests wing weight computation from sample XML data"""
diff --git a/tests/integration_tests/oad_process/data/CeRAS01_legacy.xml b/tests/integration_tests/oad_process/data/CeRAS01_legacy.xml
index 6b0cdb7..81b5b16 100644
--- a/tests/integration_tests/oad_process/data/CeRAS01_legacy.xml
+++ b/tests/integration_tests/oad_process/data/CeRAS01_legacy.xml
@@ -226,6 +226,7 @@
       </rubber_engine>
     </propulsion>
     <load_case>
+      <gust_intensity>0.5</gust_intensity>
       <lc1>
         <U_gust units="m/s">15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>
         <Vc_EAS units="m/s">375.0<!--equivalent air speed for sizing load case 1 (gust with minimum aircraft mass)--></Vc_EAS>
diff --git a/tests/integration_tests/oad_process/data/CeRAS01_legacy_breguet_result.xml b/tests/integration_tests/oad_process/data/CeRAS01_legacy_breguet_result.xml
index a679b84..674140a 100644
--- a/tests/integration_tests/oad_process/data/CeRAS01_legacy_breguet_result.xml
+++ b/tests/integration_tests/oad_process/data/CeRAS01_legacy_breguet_result.xml
@@ -222,6 +222,7 @@
       </rubber_engine>
     </propulsion>
     <load_case>
+      <gust_intensity>0.5</gust_intensity>
       <lc1>
         <U_gust units="m/s">15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>
         <Vc_EAS units="m/s">375.0<!--equivalent air speed for sizing load case 1 (gust with minimum aircraft mass)--></Vc_EAS>
diff --git a/tests/integration_tests/oad_process/data/CeRAS01_legacy_mission_result.xml b/tests/integration_tests/oad_process/data/CeRAS01_legacy_mission_result.xml
index d1dfdb1..a5559eb 100644
--- a/tests/integration_tests/oad_process/data/CeRAS01_legacy_mission_result.xml
+++ b/tests/integration_tests/oad_process/data/CeRAS01_legacy_mission_result.xml
@@ -226,6 +226,7 @@
       </rubber_engine>
     </propulsion>
     <load_case>
+      <gust_intensity>0.5</gust_intensity>
       <lc1>
         <U_gust units="m/s">15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>
         <Vc_EAS units="m/s">375.0<!--equivalent air speed for sizing load case 1 (gust with minimum aircraft mass)--></Vc_EAS>
diff --git a/tests/integration_tests/oad_process/data/no_kink_breguet_result.xml b/tests/integration_tests/oad_process/data/no_kink_breguet_result.xml
index 53110b2..f5b5db8 100644
--- a/tests/integration_tests/oad_process/data/no_kink_breguet_result.xml
+++ b/tests/integration_tests/oad_process/data/no_kink_breguet_result.xml
@@ -318,6 +318,7 @@
             </rubber_engine>
         </propulsion>
         <load_case>
+            <gust_intensity>0.0</gust_intensity>
             <lc1>
                 <U_gust units="m/s" is_input="True">
                     15.25<!--gust vertical speed for sizing load case 1 (gust with minimum aircraft mass)--></U_gust>