diff --git a/README.md b/README.md
index cfaa5572c..b117f0fb2 100644
--- a/README.md
+++ b/README.md
@@ -116,6 +116,7 @@ warning: : The module does not work completely as expected. It is not a bug, but
- [Optimisation](ceasiompy/Optimisation/README.md) :x:
- [SMTrain](ceasiompy/SMTrain/README.md) :x:
- [SMUse](ceasiompy/SMUse/README.md) :x:
+- [ThermoData](ceasiompy/ThermoData/README.md) :heavy_check_mark:
diff --git a/ceasiompy/CPACS2GMSH/README.md b/ceasiompy/CPACS2GMSH/README.md
index 37b8fc312..f80e80ab2 100644
--- a/ceasiompy/CPACS2GMSH/README.md
+++ b/ceasiompy/CPACS2GMSH/README.md
@@ -40,7 +40,7 @@ Multiple options are available with `CPACS2GMSH`, you can see these options if y
General options:
-* `Display mesh with GMSHl : False`
+* `Display mesh with GMSH : False`
Open the gmsh GUI after the generation of the surface mesh (2D) and the domain mesh (3D). This option is usefully to control the quality of the automated generated mesh or make extra operation with gmsh GUI.
Domain:
diff --git a/ceasiompy/SU2Run/__specs__.py b/ceasiompy/SU2Run/__specs__.py
index ecd05f417..243ad926b 100644
--- a/ceasiompy/SU2Run/__specs__.py
+++ b/ceasiompy/SU2Run/__specs__.py
@@ -28,6 +28,7 @@
SU2_TARGET_CL_XPATH,
SU2_UPDATE_WETTED_AREA_XPATH,
SU2MESH_XPATH,
+ SU2_CONFIG_RANS_XPATH,
)
from ceasiompy.utils.moduleinterfaces import CPACSInOut
@@ -203,6 +204,18 @@
gui_group="CPU",
)
+cpacs_inout.add_input(
+ var_name="RANS calculation",
+ var_type=list,
+ default_value=["Euler", "RANS"],
+ unit="1",
+ descr="Running an Euler or a RANS calculation",
+ xpath=SU2_CONFIG_RANS_XPATH,
+ gui=True,
+ gui_name="Euler or RANS simulation",
+ gui_group="SU2 Parameters",
+)
+
cpacs_inout.add_input(
var_name="max_iter",
var_type=int,
diff --git a/ceasiompy/SU2Run/files/config_template_rans.cfg b/ceasiompy/SU2Run/files/config_template_rans.cfg
new file mode 100644
index 000000000..1af0ceee0
--- /dev/null
+++ b/ceasiompy/SU2Run/files/config_template_rans.cfg
@@ -0,0 +1,198 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% %
+% SU2 configuration file %
+% Case description: VKI turbine %
+% Author: Francisco Palacios, Thomas D. Economon %
+% Institution: Stanford University %
+% Date: Feb 18th, 2013 %
+% File Version 8.0.0 "Harrier" %
+% %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------%
+%
+SOLVER= RANS
+%
+KIND_TURB_MODEL= SA
+%
+MATH_PROBLEM= DIRECT
+%
+RESTART_SOL= NO
+%
+SYSTEM_MEASUREMENTS= SI
+
+% -------------------- COMPRESSIBLE FREE-STREAM DEFINITION --------------------%
+%
+MACH_NUMBER= 0.3
+%
+AOA= 0.0
+%
+FREESTREAM_PRESSURE= 1.013E5
+%
+FREESTREAM_OPTION= DENSITY_FS
+%
+FREESTREAM_DENSITY = 1.255
+%
+FREESTREAM_TEMPERATURE =288.15
+%
+REYNOLDS_NUMBER= 4E6
+%
+REYNOLDS_LENGTH= 1.0
+
+% ---------------------- REFERENCE VALUE DEFINITION ---------------------------%
+%
+REF_ORIGIN_MOMENT_X = 0.00
+%
+REF_ORIGIN_MOMENT_Y = 0.00
+%
+REF_ORIGIN_MOMENT_Z = 0.00
+%
+REF_LENGTH= 45.0
+%
+REF_AREA= 1
+% Compressible flow non-dimensionalization (DIMENSIONAL, FREESTREAM_PRESS_EQ_ONE,
+% FREESTREAM_VEL_EQ_MACH, FREESTREAM_VEL_EQ_ONE)
+REF_DIMENSIONALIZATION= FREESTREAM_PRESS_EQ_ONE
+%
+% ----------------------- BOUNDARY CONDITION DEFINITION -----------------------%
+%
+MARKER_HEATFLUX = (S4CreatedbyGmsh, 0.0)
+%
+MARKER_FAR= (S3Farfield)
+%
+%MARKER_SYM =()
+%
+% ------------------------ SURFACES IDENTIFICATION ----------------------------%
+%
+% Marker(s) of the surface to be plotted or designed
+MARKER_PLOTTING= (S4CreatedbyGmsh)
+%
+% Marker(s) of the surface where the functional (Cd, Cl, etc.) will be evaluated
+MARKER_MONITORING= (S4CreatedbyGmsh)
+
+% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------%
+%
+NUM_METHOD_GRAD= GREEN_GAUSS
+%
+CFL_NUMBER= 1
+%
+CFL_ADAPT= NO
+%
+CFL_ADAPT_PARAM= ( 0.5, 1, 1.0, 100.0 )
+%
+ITER = 5000
+
+% Objective function in gradient evaluation (DRAG, LIFT, SIDEFORCE, MOMENT_X,
+% MOMENT_Y, MOMENT_Z, EFFICIENCY, BUFFET,
+% EQUIVALENT_AREA, NEARFIELD_PRESSURE,
+% FORCE_X, FORCE_Y, FORCE_Z, THRUST,
+% TORQUE, TOTAL_HEATFLUX,
+% MAXIMUM_HEATFLUX, INVERSE_DESIGN_PRESSURE,
+% INVERSE_DESIGN_HEATFLUX, SURFACE_TOTAL_PRESSURE,
+% SURFACE_MASSFLOW, SURFACE_STATIC_PRESSURE, SURFACE_MACH)
+% For a weighted sum of objectives: separate by commas, add OBJECTIVE_WEIGHT and MARKER_MONITORING in matching order.
+OBJECTIVE_FUNCTION= DRAG
+
+% ------------------------ LINEAR SOLVER DEFINITION ---------------------------%
+% Linear solver or smoother for implicit formulations:
+% BCGSTAB, FGMRES, RESTARTED_FGMRES, CONJUGATE_GRADIENT (self-adjoint problems only), SMOOTHER.
+LINEAR_SOLVER= FGMRES
+%
+% Preconditioner of the Krylov linear solver or type of smoother (ILU, LU_SGS, LINELET, JACOBI)
+LINEAR_SOLVER_PREC= ILU
+%
+% Minimum error of the linear solver for implicit formulations
+LINEAR_SOLVER_ERROR= 1E-12
+%
+% Max number of iterations of the linear solver for the implicit formulation
+LINEAR_SOLVER_ITER= 3
+%
+% Number of elements to apply the criteria
+CONV_CAUCHY_ELEMS= 1000
+%
+% Epsilon to control the series convergence
+CONV_CAUCHY_EPS= 1E-10
+% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------%
+% Convective numerical method (JST, LAX-FRIEDRICH, CUSP, ROE, AUSM, AUSMPLUSUP,
+% AUSMPLUSUP2, HLLC, TURKEL_PREC, MSW, FDS, SLAU, SLAU2)
+CONV_NUM_METHOD_FLOW= JST
+%
+% Time discretization (RUNGE-KUTTA_EXPLICIT, EULER_IMPLICIT, EULER_EXPLICIT)
+TIME_DISCRE_FLOW= EULER_IMPLICIT
+%
+% -------------------- TURBULENT NUMERICAL METHOD DEFINITION ------------------%
+% Convective numerical method (SCALAR_UPWIND)
+CONV_NUM_METHOD_TURB= SCALAR_UPWIND
+%
+% Monotonic Upwind Scheme for Conservation Laws (TVD) in the turbulence equations.
+% Required for 2nd order upwind schemes (NO, YES)
+MUSCL_TURB= NO
+%
+% Slope limiter (NONE, VENKATAKRISHNAN, VENKATAKRISHNAN_WANG,
+% BARTH_JESPERSEN, VAN_ALBADA_EDGE)
+SLOPE_LIMITER_TURB= VENKATAKRISHNAN
+%
+% Time discretization (EULER_IMPLICIT)
+TIME_DISCRE_TURB= EULER_IMPLICIT
+% -------------------------- MULTIGRID PARAMETERS -----------------------------%
+%
+% Multi-Grid Levels (0 = no multi-grid)
+MGLEVEL = 3
+%
+% Multi-grid cycle (V_CYCLE, W_CYCLE, FULLMG_CYCLE)
+MGCYCLE = W_CYCLE
+%
+% Multi-Grid PreSmoothing Level
+MG_PRE_SMOOTH = ( 1.0, 2.0, 3.0, 3.0 )
+%
+% Multi-Grid PostSmoothing Level
+MG_POST_SMOOTH = ( 0.0, 0.0, 0.0, 0.0 )
+%
+% Jacobi implicit smoothing of the correction
+MG_CORRECTION_SMOOTH = ( 0.0, 0.0, 0.0, 0.0 )
+%
+% Damping factor for the residual restriction
+MG_DAMP_RESTRICTION = 0.9
+%
+% Damping factor for the correction prolongation
+MG_DAMP_PROLONGATION = 0.9
+% --------------------------- CONVERGENCE PARAMETERS --------------------------%
+%
+% Min value of the residual (log10 of the residual)
+CONV_RESIDUAL_MINVAL= -12
+%
+% Start convergence criteria at iteration number
+CONV_STARTITER= 10
+%
+
+% ------------------------- INPUT/OUTPUT INFORMATION --------------------------%
+%
+MESH_FILENAME= labair_penta.su2
+%
+MESH_FORMAT= SU2
+%
+MESH_OUT_FILENAME= mesh_out.su2
+%
+SOLUTION_FILENAME= solution_flow.dat
+%
+HISTORY_OUTPUT = (INNER_ITER, RMS_RES, AERO_COEFF)
+%
+TABULAR_FORMAT= CSV
+%
+CONV_FILENAME= history
+%
+RESTART_FILENAME= restart_flow.dat
+%
+RESTART_ADJ_FILENAME= restart_adj.dat
+%
+VOLUME_FILENAME= flow
+%
+SURFACE_FILENAME= surface_flow
+%
+OUTPUT_WRT_FREQ= 100
+%
+BREAKDOWN_FILENAME = forces_breakdown_def.dat
+%
+WRT_FORCES_BREAKDOWN = YES
+%
+SCREEN_OUTPUT= (INNER_ITER, RMS_DENSITY, RMS_ENERGY, LIFT, DRAG)
\ No newline at end of file
diff --git a/ceasiompy/SU2Run/func/su2config.py b/ceasiompy/SU2Run/func/su2config.py
index 7b77a9bc6..8e08de5ed 100644
--- a/ceasiompy/SU2Run/func/su2config.py
+++ b/ceasiompy/SU2Run/func/su2config.py
@@ -64,6 +64,8 @@
SU2_ROTATION_RATE_XPATH,
SU2_TARGET_CL_XPATH,
SU2MESH_XPATH,
+ ENGINE_TYPE_XPATH,
+ ENGINE_BC,
)
from ceasiompy.utils.configfiles import ConfigFile
from cpacspy.cpacsfunctions import (
@@ -294,6 +296,34 @@ def add_actuator_disk(cfg, cpacs, case_dir_path, actuator_disk_file, mesh_marker
f.close()
+# adding the results of ThermoData to the config file of su2
+
+
+def add_thermodata(cfg, cpacs, alt, case_nb, alt_list):
+
+ if cpacs.tixi.checkElement(ENGINE_TYPE_XPATH):
+ log.info("adding engine BC to the SU2 config file")
+ engine_type = get_value(cpacs.tixi, ENGINE_TYPE_XPATH)
+ log.info(f"engine type {engine_type}")
+ alt = alt_list[case_nb]
+ Atm = Atmosphere(alt)
+ tot_temp_in = Atm.temperature[0]
+ tot_pressure_in = Atm.pressure[0]
+ if len(alt_list) > 1:
+ tot_temp_out = get_value(cpacs.tixi, ENGINE_BC + "/temperatureOutlet").split(";")
+ tot_pressure_out = get_value(cpacs.tixi, ENGINE_BC + "/pressureOutlet").split(";")
+ tot_temp_out = tot_temp_out[case_nb]
+ tot_pressure_out = tot_pressure_out[case_nb]
+ else:
+ tot_temp_out = get_value(cpacs.tixi, ENGINE_BC + "/temperatureOutlet")
+ tot_pressure_out = get_value(cpacs.tixi, ENGINE_BC + "/pressureOutlet")
+ cfg["INLET_TYPE"] = "TOTAL_CONDITIONS"
+ cfg["MARKER_INLET"] = (
+ f"(INLET_ENGINE, {tot_temp_in}, {tot_pressure_in}, {1},{0},{0}, "
+ f"OUTLET_ENGINE,{tot_temp_out},{tot_pressure_out}, {1},{0},{0})"
+ )
+
+
def generate_su2_cfd_config(cpacs_path, cpacs_out_path, wkdir):
"""Function to create SU2 config file.
@@ -341,7 +371,7 @@ def generate_su2_cfd_config(cpacs_path, cpacs_out_path, wkdir):
if aeromap_list:
aeromap_default = aeromap_list[0]
- log.info(f'The aeromap is {aeromap_default}')
+ log.info(f"The aeromap is {aeromap_default}")
aeromap_uid = get_value_or_default(cpacs.tixi, SU2_AEROMAP_UID_XPATH, aeromap_default)
@@ -448,6 +478,7 @@ def generate_su2_cfd_config(cpacs_path, cpacs_out_path, wkdir):
cfg["HISTORY_OUTPUT"] = "(INNER_ITER, RMS_RES, AERO_COEFF)"
# Parameters which will vary for the different cases (alt,mach,aoa,aos)
+
for case_nb in range(len(alt_list)):
cfg["MESH_FILENAME"] = str(su2_mesh)
@@ -470,6 +501,8 @@ def generate_su2_cfd_config(cpacs_path, cpacs_out_path, wkdir):
f"_aoa{round(aoa, 1)}_aos{round(aos, 1)}"
)
+ add_thermodata(cfg, cpacs, alt, case_nb, alt_list)
+
case_dir_path = Path(wkdir, case_dir_name)
if not case_dir_path.exists():
case_dir_path.mkdir()
diff --git a/ceasiompy/SU2Run/func/su2config_rans.py b/ceasiompy/SU2Run/func/su2config_rans.py
new file mode 100644
index 000000000..77b6bbeb6
--- /dev/null
+++ b/ceasiompy/SU2Run/func/su2config_rans.py
@@ -0,0 +1,582 @@
+"""
+CEASIOMpy: Conceptual Aircraft Design Software
+
+Developed by CFS ENGINEERING, 1015 Lausanne, Switzerland
+
+Function generate or modify SU2 configuration files
+
+Python version: >=3.8
+
+| Author: Aidan Jungo
+| Creation: 2020-02-24
+
+TODO:
+
+ * add and test control surface functions
+
+"""
+
+# =================================================================================================
+# IMPORTS
+# =================================================================================================
+
+from pathlib import Path
+from shutil import copyfile
+
+from ambiance import Atmosphere
+from ceasiompy.SU2Run.func.su2actuatordiskfile import (
+ get_advanced_ratio,
+ get_radial_stations,
+ save_plots,
+ thrust_calculator,
+ write_actuator_disk_data,
+ write_header,
+)
+from ceasiompy.CPACS2GMSH.func.mesh_sizing import wings_size
+from ceasiompy.SU2Run.func.su2utils import get_mesh_markers, get_su2_config_template_rans
+from ceasiompy.utils.ceasiomlogger import get_logger
+from ceasiompy.utils.commonnames import (
+ ACTUATOR_DISK_FILE_NAME,
+ ACTUATOR_DISK_INLET_SUFFIX,
+ ACTUATOR_DISK_OUTLET_SUFFIX,
+ CONFIG_CFD_NAME,
+ SU2_FORCES_BREAKDOWN_NAME,
+)
+from ceasiompy.utils.commonxpath import (
+ GMSH_SYMMETRY_XPATH,
+ PROP_XPATH,
+ RANGE_XPATH,
+ SU2_ACTUATOR_DISK_XPATH,
+ SU2_AEROMAP_UID_XPATH,
+ SU2_BC_FARFIELD_XPATH,
+ SU2_BC_WALL_XPATH,
+ SU2_CFL_NB_XPATH,
+ SU2_CFL_ADAPT_XPATH,
+ SU2_CFL_ADAPT_PARAM_DOWN_XPATH,
+ SU2_CFL_ADAPT_PARAM_UP_XPATH,
+ SU2_CFL_MAX_XPATH,
+ SU2_CFL_MIN_XPATH,
+ SU2_CONTROL_SURF_XPATH,
+ SU2_DAMPING_DER_XPATH,
+ SU2_DEF_MESH_XPATH,
+ SU2_FIXED_CL_XPATH,
+ SU2_MAX_ITER_XPATH,
+ SU2_MG_LEVEL_XPATH,
+ SU2_ROTATION_RATE_XPATH,
+ SU2_TARGET_CL_XPATH,
+ SU2MESH_XPATH,
+ ENGINE_TYPE_XPATH,
+ ENGINE_BC,
+)
+from ceasiompy.utils.configfiles import ConfigFile
+from cpacspy.cpacsfunctions import (
+ create_branch,
+ get_string_vector,
+ get_value,
+ get_value_or_default,
+)
+from cpacspy.cpacspy import CPACS
+
+log = get_logger()
+
+MODULE_DIR = Path(__file__).parent
+
+
+# =================================================================================================
+# CLASSES
+# =================================================================================================
+
+
+# =================================================================================================
+# FUNCTIONS
+# =================================================================================================
+
+
+def add_damping_derivatives(cfg, wkdir, case_dir_name, rotation_rate):
+ """Add damping derivatives parameter to the config file and save them to their respective
+ directory.
+
+ Args:
+ cfg (ConfigFile): ConfigFile object.
+ wkdir (Path): Path to the working directory
+ case_dir_name (str): Name of the case directory
+ rotation_rate (float): Rotation rate that will be impose to calculate damping derivatives
+ """
+
+ cfg["GRID_MOVEMENT"] = "ROTATING_FRAME"
+
+ cfg["ROTATION_RATE"] = f"{rotation_rate} 0.0 0.0"
+ case_dir = Path(wkdir, f"{case_dir_name}_dp")
+ case_dir.mkdir()
+ cfg.write_file(Path(case_dir, CONFIG_CFD_NAME), overwrite=True)
+
+ cfg["ROTATION_RATE"] = f"0.0 {rotation_rate} 0.0"
+ case_dir = Path(wkdir, f"{case_dir_name}_dq")
+ case_dir.mkdir()
+ cfg.write_file(Path(case_dir, CONFIG_CFD_NAME), overwrite=True)
+
+ cfg["ROTATION_RATE"] = f"0.0 0.0 {rotation_rate}"
+ case_dir = Path(wkdir, f"{case_dir_name}_dr")
+ case_dir.mkdir()
+ cfg.write_file(Path(case_dir, CONFIG_CFD_NAME), overwrite=True)
+
+ log.info("Damping derivatives cases directories has been created.")
+
+
+def add_actuator_disk(cfg, cpacs, case_dir_path, actuator_disk_file, mesh_markers, alt, mach):
+ """Add actuator disk parameter to the config file.
+
+ Args:
+ cfg (ConfigFile): ConfigFile object.
+ cpacs (CPACS): CPACS object from cpacspy library
+ case_dir_path (Path): Path object to the current case directory
+ actuator_disk_file (Path): Path to the actuator disk file
+ mesh_markers (dict): Dictionary containing all the mesh markers
+
+ Returns:
+ cfg (ConfigFile): ConfigFile object.
+ """
+
+ ad_inlet_marker = sorted(mesh_markers["actuator_disk_inlet"])
+ ad_outlet_marker = sorted(mesh_markers["actuator_disk_outlet"])
+
+ if "None" in ad_inlet_marker or "None" in ad_outlet_marker:
+ return
+
+ if len(ad_inlet_marker) != len(ad_outlet_marker):
+ raise ValueError(
+ "The number of inlet and outlet markers of the actuator disk must be the same."
+ )
+
+ # Get rotorcraft configuration (propeller)
+ try:
+ rotorcraft_config = cpacs.rotorcraft.configuration
+ except AttributeError:
+ raise ValueError(
+ "The actuator disk is defined but no rotorcraft configuration is defined in "
+ "the CPACS file."
+ )
+
+ rotor_uid_pos = {}
+ for i in range(1, rotorcraft_config.get_rotor_count() + 1):
+ rotor = rotorcraft_config.get_rotor(i)
+
+ rotor_uid = rotor.get_uid()
+ pos_x = rotor.get_translation().x
+ pos_y = rotor.get_translation().y
+ pos_z = rotor.get_translation().z
+ radius = rotor.get_radius()
+ hub_radius = 0.0 # TODO: get correctly from CPACS
+
+ rotor_xpath = cpacs.tixi.uIDGetXPath(rotor_uid)
+
+ number_of_blades_xpath = (
+ rotor_xpath + "/rotorHub/rotorBladeAttachments/rotorBladeAttachment/numberOfBlades"
+ )
+ number_of_blades = get_value_or_default(cpacs.tixi, number_of_blades_xpath, 3)
+
+ # TODO: this is the nominal speed, how to get a speed which correspond to each flight cond.
+ rotational_velocity_xpath = rotor_xpath + "/nominalRotationsPerMinute"
+ rotational_velocity = (
+ get_value_or_default(cpacs.tixi, rotational_velocity_xpath, 3000) / 60.0
+ )
+
+ rotor_uid_pos[rotor_uid] = (
+ pos_x,
+ pos_y,
+ pos_z,
+ radius,
+ hub_radius,
+ number_of_blades,
+ rotational_velocity,
+ )
+
+ cfg["ACTDISK_DOUBLE_SURFACE"] = "YES"
+ cfg["ACTDISK_TYPE"] = "VARIABLE_LOAD"
+ cfg["ACTDISK_FILENAME"] = ACTUATOR_DISK_FILE_NAME
+ cfg["MGLEVEL"] = 0 # Calculation diverges if multigrid is used with a disk actuator
+
+ actdisk_markers = []
+
+ f = open(actuator_disk_file, "w")
+ f = write_header(f)
+
+ for maker_inlet, marker_outlet in zip(ad_inlet_marker, ad_outlet_marker):
+ inlet_uid = maker_inlet.split(ACTUATOR_DISK_INLET_SUFFIX)[0]
+ outlet_uid = marker_outlet.split(ACTUATOR_DISK_OUTLET_SUFFIX)[0]
+
+ if inlet_uid != outlet_uid:
+ raise ValueError("The inlet and outlet markers of the actuator disk must be the same.")
+
+ if "_mirrored" in maker_inlet:
+ uid = inlet_uid.split("_mirrored")[0]
+ sym = -1
+ else:
+ uid = inlet_uid
+ sym = 1
+
+ center = []
+ center.append(round(rotor_uid_pos[uid][0], 5))
+ center.append(round(sym * rotor_uid_pos[uid][1], 5))
+ center.append(round(rotor_uid_pos[uid][2], 5))
+
+ axis = (1.0, 0.0, 0.0) # TODO: get the axis by applying the rotation matrix
+ radius = round(rotor_uid_pos[uid][3], 5)
+ hub_radius = round(rotor_uid_pos[uid][4], 5)
+ number_of_blades = round(rotor_uid_pos[uid][5], 5)
+ rotational_velocity = round(rotor_uid_pos[uid][6], 5)
+
+ actdisk_markers.append(maker_inlet)
+ actdisk_markers.append(marker_outlet)
+ actdisk_markers.append(str(center[0]))
+ actdisk_markers.append(str(center[1]))
+ actdisk_markers.append(str(center[2]))
+ actdisk_markers.append(str(center[0]))
+ actdisk_markers.append(str(center[1]))
+ actdisk_markers.append(str(center[2]))
+
+ Atm = Atmosphere(alt)
+ free_stream_velocity = mach * Atm.speed_of_sound[0]
+
+ radial_stations = get_radial_stations(radius, hub_radius)
+ advanced_ratio = get_advanced_ratio(free_stream_velocity, rotational_velocity, radius)
+
+ prandtl_correction_xpath = PROP_XPATH + "/propeller/blade/loss"
+ prandtl_correction = get_value_or_default(cpacs.tixi, prandtl_correction_xpath, True)
+
+ thrust_xpath = PROP_XPATH + "/propeller/thrust"
+ thrust = get_value_or_default(cpacs.tixi, thrust_xpath, 3000)
+ total_thrust_coefficient = float(
+ thrust / (Atm.density * rotational_velocity**2 * (radius * 2) ** 4)
+ )
+
+ (
+ radial_thrust_coefs,
+ radial_power_coefs,
+ non_dimensional_radius,
+ optimal_axial_interference_factor,
+ optimal_rotational_interference_factor,
+ prandtl_correction_values,
+ ) = thrust_calculator(
+ radial_stations,
+ total_thrust_coefficient,
+ radius,
+ free_stream_velocity,
+ prandtl_correction,
+ number_of_blades,
+ rotational_velocity,
+ )
+
+ save_plots(
+ radial_stations,
+ radial_thrust_coefs,
+ radial_power_coefs,
+ non_dimensional_radius,
+ optimal_axial_interference_factor,
+ optimal_rotational_interference_factor,
+ prandtl_correction_values,
+ case_dir_path,
+ inlet_uid,
+ )
+
+ f = write_actuator_disk_data(
+ file=f,
+ inlet_marker=maker_inlet,
+ outlet_marker=marker_outlet,
+ center=center,
+ axis=axis,
+ radius=radius,
+ advanced_ratio=advanced_ratio,
+ radial_stations=radial_stations,
+ radial_thrust_coefs=radial_thrust_coefs,
+ radial_power_coefs=radial_power_coefs,
+ )
+
+ cfg["MARKER_ACTDISK"] = " (" + ", ".join(actdisk_markers) + " )"
+
+ f.close()
+
+
+def add_thermodata(cfg, cpacs, alt, case_nb, alt_list):
+ if cpacs.tixi.checkElement(ENGINE_TYPE_XPATH):
+ log.info("adding engine BC to the SU2 config file")
+ engine_type = get_value(cpacs.tixi, ENGINE_TYPE_XPATH)
+ log.info(f"engine type {engine_type}")
+ alt = alt_list[case_nb]
+ Atm = Atmosphere(alt)
+ tot_temp_in = Atm.temperature[0]
+ tot_pressure_in = Atm.pressure[0]
+ if len(alt_list) > 1:
+ tot_temp_out = get_value(cpacs.tixi, ENGINE_BC + "/temperatureOutlet").split(";")
+ tot_pressure_out = get_value(cpacs.tixi, ENGINE_BC + "/pressureOutlet").split(";")
+ tot_temp_out = tot_temp_out[case_nb]
+ tot_pressure_out = tot_pressure_out[case_nb]
+ else:
+ tot_temp_out = get_value(cpacs.tixi, ENGINE_BC + "/temperatureOutlet")
+ tot_pressure_out = get_value(cpacs.tixi, ENGINE_BC + "/pressureOutlet")
+ cfg["INLET_TYPE"] = "TOTAL_CONDITIONS"
+ cfg["MARKER_INLET"] = (
+ f"(INLET_ENGINE, {tot_temp_in}, {tot_pressure_in}, {1},{0},{0}, "
+ f"OUTLET_ENGINE,{tot_temp_out},{tot_pressure_out}, {1},{0},{0})"
+ )
+
+
+def add_reynods_number(alt, mach, cfg, cpacs_path):
+
+ Atm = Atmosphere(alt)
+
+ # Get speed from Mach Number
+ speed = mach * Atm.speed_of_sound[0]
+
+ ref_chord = wings_size(cpacs_path)[0] / 0.15
+ print(ref_chord)
+
+ # Reynolds number based on the mean chord
+ reynolds = ref_chord * speed / Atm.kinematic_viscosity[0]
+ log.info(f"Reynolds number= {int(reynolds)}")
+ cfg["REYNOLDS_NUMBER"] = int(reynolds)
+
+
+def generate_su2_cfd_config_rans(cpacs_path, cpacs_out_path, wkdir):
+ """Function to create SU2 config file.
+
+ Function 'generate_su2_cfd_config' reads data in the CPACS file and generate configuration
+ files for one or multiple flight conditions (alt,mach,aoa,aos)
+
+ Source:
+ * SU2 config template: https://github.com/su2code/SU2/blob/master/config_template.cfg
+
+ Args:
+ cpacs_path (Path): Path to CPACS file
+ cpacs_out_path (Path):Path to CPACS output file
+ wkdir (Path): Path to the working directory
+
+ """
+
+ cpacs = CPACS(cpacs_path)
+
+ # creare delle nuove xpath per la mesh su2
+
+ su2_mesh = Path(get_value(cpacs.tixi, SU2MESH_XPATH))
+ if not su2_mesh.is_file():
+ raise FileNotFoundError(f"SU2 mesh file {su2_mesh} not found")
+
+ mesh_markers = get_mesh_markers(su2_mesh)
+
+ create_branch(cpacs.tixi, SU2_BC_WALL_XPATH)
+ bc_wall_str = ";".join(mesh_markers["wall"])
+ cpacs.tixi.updateTextElement(SU2_BC_WALL_XPATH, bc_wall_str)
+
+ create_branch(cpacs.tixi, SU2_BC_FARFIELD_XPATH)
+ bc_farfiled_str = ";".join(mesh_markers["engine_intake"] + mesh_markers["engine_exhaust"])
+ cpacs.tixi.updateTextElement(SU2_BC_FARFIELD_XPATH, bc_farfiled_str)
+
+ create_branch(cpacs.tixi, SU2_ACTUATOR_DISK_XPATH)
+ bc_actuator_disk_str = ";".join(
+ mesh_markers["actuator_disk_inlet"] + mesh_markers["actuator_disk_outlet"]
+ )
+ cpacs.tixi.updateTextElement(SU2_ACTUATOR_DISK_XPATH, bc_actuator_disk_str)
+
+ fixed_cl = get_value_or_default(cpacs.tixi, SU2_FIXED_CL_XPATH, "NO")
+ target_cl = get_value_or_default(cpacs.tixi, SU2_TARGET_CL_XPATH, 1.0)
+
+ if fixed_cl == "NO":
+ # Get the first aeroMap as default one or create automatically one
+ aeromap_list = cpacs.get_aeromap_uid_list()
+
+ if aeromap_list:
+ aeromap_default = aeromap_list[0]
+ log.info(f"The aeromap is {aeromap_default}")
+
+ aeromap_uid = get_value_or_default(cpacs.tixi, SU2_AEROMAP_UID_XPATH, aeromap_default)
+
+ activate_aeromap = cpacs.get_aeromap_by_uid(aeromap_uid)
+ alt_list = activate_aeromap.get("altitude").tolist()
+ mach_list = activate_aeromap.get("machNumber").tolist()
+ aoa_list = activate_aeromap.get("angleOfAttack").tolist()
+ aos_list = activate_aeromap.get("angleOfSideslip").tolist()
+
+ else:
+ default_aeromap = cpacs.create_aeromap("DefaultAeromap")
+ default_aeromap.description = "AeroMap created automatically"
+
+ mach = get_value_or_default(cpacs.tixi, RANGE_XPATH + "/cruiseMach", 0.3)
+ alt = get_value_or_default(cpacs.tixi, RANGE_XPATH + "/cruiseAltitude", 10000)
+
+ default_aeromap.add_row(alt=alt, mach=mach, aos=0.0, aoa=0.0)
+ default_aeromap.save()
+
+ alt_list = [alt]
+ mach_list = [mach]
+ aoa_list = [0.0]
+ aos_list = [0.0]
+
+ aeromap_uid = get_value_or_default(cpacs.tixi, SU2_AEROMAP_UID_XPATH, "DefaultAeromap")
+ log.info(f"{aeromap_uid} has been created")
+
+ else: # if fixed_cl == 'YES':
+ log.info("Configuration file for fixed CL calculation will be created.")
+
+ fixed_cl_aeromap = cpacs.create_aeromap("aeroMap_fixedCL_SU2")
+ fixed_cl_aeromap.description = f"AeroMap created for SU2 fixed CL value of {target_cl}"
+
+ mach = get_value_or_default(cpacs.tixi, RANGE_XPATH + "/cruiseMach", 0.78)
+ alt = get_value_or_default(cpacs.tixi, RANGE_XPATH + "/cruiseAltitude", 12000)
+
+ fixed_cl_aeromap.add_row(alt=alt, mach=mach, aos=0.0, aoa=0.0)
+ fixed_cl_aeromap.save()
+
+ alt_list = [alt]
+ mach_list = [mach]
+ aoa_list = [0.0]
+ aos_list = [0.0]
+
+ cfg = ConfigFile(get_su2_config_template_rans())
+
+ # Check if symmetry plane is defined (Default: False)
+ sym_factor = 1.0
+ if get_value_or_default(cpacs.tixi, GMSH_SYMMETRY_XPATH, False):
+ log.info("Symmetry plane is defined. The reference area will be divided by 2.")
+ sym_factor = 2.0
+
+ # General parameters
+ cfg["RESTART_SOL"] = "NO"
+ cfg["REF_LENGTH"] = cpacs.aircraft.ref_length
+ cfg["REF_AREA"] = cpacs.aircraft.ref_area / sym_factor
+ cfg["REF_ORIGIN_MOMENT_X"] = cpacs.aircraft.ref_point_x
+ cfg["REF_ORIGIN_MOMENT_Y"] = cpacs.aircraft.ref_point_y
+ cfg["REF_ORIGIN_MOMENT_Z"] = cpacs.aircraft.ref_point_z
+
+ # Settings
+
+ cfl_down = get_value_or_default(cpacs.tixi, SU2_CFL_ADAPT_PARAM_DOWN_XPATH, 0.5)
+ cfl_up = get_value_or_default(cpacs.tixi, SU2_CFL_ADAPT_PARAM_UP_XPATH, 1.5)
+ cfl_min = get_value_or_default(cpacs.tixi, SU2_CFL_MIN_XPATH, 0.5)
+ cfl_max = get_value_or_default(cpacs.tixi, SU2_CFL_MAX_XPATH, 100)
+
+ if get_value_or_default(cpacs.tixi, SU2_CFL_ADAPT_XPATH, True):
+ cfg["CFL_ADAPT"] = "YES"
+
+ else:
+ cfg["CFL_ADAPT"] = "NO"
+
+ cfg["INNER_ITER"] = int(get_value_or_default(cpacs.tixi, SU2_MAX_ITER_XPATH, 200))
+ cfg["CFL_NUMBER"] = get_value_or_default(cpacs.tixi, SU2_CFL_NB_XPATH, 1.0)
+ cfg["CFL_ADAPT_PARAM"] = f"( {cfl_down}, {cfl_up}, {cfl_min}, {cfl_max} )"
+ cfg["MGLEVEL"] = int(get_value_or_default(cpacs.tixi, SU2_MG_LEVEL_XPATH, 3))
+
+ # Fixed CL mode (AOA will not be taken into account)
+ cfg["FIXED_CL_MODE"] = fixed_cl
+ cfg["TARGET_CL"] = target_cl
+ cfg["DCL_DALPHA"] = "0.1"
+ cfg["UPDATE_AOA_ITER_LIMIT"] = "50"
+ cfg["ITER_DCL_DALPHA"] = "80"
+
+ # Mesh Marker
+ bc_wall_str = f"( {','.join(mesh_markers['wall'])} )"
+
+ cfg["MARKER_EULER"] = bc_wall_str
+ farfield_bc = (
+ mesh_markers["farfield"] + mesh_markers["engine_intake"] + mesh_markers["engine_exhaust"]
+ )
+ cfg["MARKER_FAR"] = f"( {','.join(farfield_bc)} )"
+ cfg["MARKER_SYM"] = f"( {','.join(mesh_markers['symmetry'])} )"
+ cfg["MARKER_PLOTTING"] = bc_wall_str
+ cfg["MARKER_MONITORING"] = bc_wall_str
+ cfg["MARKER_MOVING"] = "( NONE )" # TODO: when do we need to define MARKER_MOVING?
+ cfg["DV_MARKER"] = bc_wall_str
+
+ # Output
+ cfg["WRT_FORCES_BREAKDOWN"] = "YES"
+ cfg["BREAKDOWN_FILENAME"] = SU2_FORCES_BREAKDOWN_NAME
+ cfg["OUTPUT_FILES"] = "(RESTART, PARAVIEW, SURFACE_PARAVIEW)"
+ cfg["HISTORY_OUTPUT"] = "(INNER_ITER, RMS_RES, AERO_COEFF)"
+
+ # Parameters which will vary for the different cases (alt,mach,aoa,aos)
+
+ for case_nb in range(len(alt_list)):
+ cfg["MESH_FILENAME"] = str(su2_mesh)
+
+ alt = alt_list[case_nb]
+ mach = mach_list[case_nb]
+ aoa = aoa_list[case_nb]
+ aos = aos_list[case_nb]
+
+ Atm = Atmosphere(alt)
+
+ cfg["MACH_NUMBER"] = mach
+ cfg["AOA"] = aoa
+ cfg["SIDESLIP_ANGLE"] = aos
+ cfg["FREESTREAM_PRESSURE"] = Atm.pressure[0]
+ cfg["FREESTREAM_TEMPERATURE"] = Atm.temperature[0]
+ cfg["ROTATION_RATE"] = "0.0 0.0 0.0"
+
+ case_dir_name = (
+ f"Case{str(case_nb).zfill(2)}_alt{alt}_mach{round(mach, 2)}"
+ f"_aoa{round(aoa, 1)}_aos{round(aos, 1)}"
+ )
+
+ add_thermodata(cfg, cpacs, alt, case_nb, alt_list)
+
+ add_reynods_number(alt, mach, cfg, cpacs_path)
+
+ case_dir_path = Path(wkdir, case_dir_name)
+ if not case_dir_path.exists():
+ case_dir_path.mkdir()
+
+ if get_value_or_default(cpacs.tixi, SU2_ACTUATOR_DISK_XPATH, False):
+ actuator_disk_file = Path(wkdir, ACTUATOR_DISK_FILE_NAME)
+ add_actuator_disk(
+ cfg, cpacs, case_dir_path, actuator_disk_file, mesh_markers, alt, mach
+ )
+
+ if actuator_disk_file.exists():
+ case_actuator_disk_file = Path(case_dir_path, ACTUATOR_DISK_FILE_NAME)
+ copyfile(actuator_disk_file, case_actuator_disk_file)
+
+ bc_wall_str = (
+ "("
+ + ",".join(
+ mesh_markers["wall"]
+ + mesh_markers["actuator_disk_inlet"]
+ + mesh_markers["actuator_disk_outlet"]
+ )
+ + ")"
+ )
+
+ cfg["MARKER_PLOTTING"] = bc_wall_str
+ cfg["MARKER_MONITORING"] = bc_wall_str
+
+ config_output_path = Path(case_dir_path, CONFIG_CFD_NAME)
+ cfg.write_file(config_output_path, overwrite=True)
+
+ if get_value_or_default(cpacs.tixi, SU2_DAMPING_DER_XPATH, False):
+ rotation_rate = get_value_or_default(cpacs.tixi, SU2_ROTATION_RATE_XPATH, 1.0)
+ add_damping_derivatives(cfg, wkdir, case_dir_name, rotation_rate)
+
+ # Control surfaces deflections (TODO: create a subfunctions)
+ if get_value_or_default(cpacs.tixi, SU2_CONTROL_SURF_XPATH, False):
+ # Get deformed mesh list
+ if cpacs.tixi.checkElement(SU2_DEF_MESH_XPATH):
+ su2_def_mesh_list = get_string_vector(cpacs.tixi, SU2_DEF_MESH_XPATH)
+ else:
+ log.warning("No SU2 deformed mesh has been found!")
+ su2_def_mesh_list = []
+
+ for su2_def_mesh in su2_def_mesh_list:
+ mesh_path = Path(wkdir, "MESH", su2_def_mesh)
+ config_dir_path = Path(wkdir, case_dir_name + "_" + su2_def_mesh.split(".")[0])
+ config_dir_path.mkdir()
+ cfg["MESH_FILENAME"] = mesh_path
+
+ cfg.write_file(Path(config_dir_path, CONFIG_CFD_NAME), overwrite=True)
+
+ cpacs.save_cpacs(cpacs_out_path, overwrite=True)
+
+
+# =================================================================================================
+# MAIN
+# =================================================================================================
+
+if __name__ == "__main__":
+ log.info("Nothing to execute!")
diff --git a/ceasiompy/SU2Run/func/su2utils.py b/ceasiompy/SU2Run/func/su2utils.py
index ab08a43ef..24db788bf 100644
--- a/ceasiompy/SU2Run/func/su2utils.py
+++ b/ceasiompy/SU2Run/func/su2utils.py
@@ -170,6 +170,14 @@ def get_su2_config_template():
return su2_config_template_path
+def get_su2_config_template_rans():
+
+ su2_dir = get_module_path("SU2Run")
+ su2_config_template_path_rans = Path(su2_dir, "files", "config_template_rans.cfg")
+
+ return su2_config_template_path_rans
+
+
def get_su2_aerocoefs(force_file):
"""Get aerodynamic coefficients and velocity from SU2 forces file (forces_breakdown.dat)
diff --git a/ceasiompy/SU2Run/su2run.py b/ceasiompy/SU2Run/su2run.py
index 4b103b22c..ea8499714 100644
--- a/ceasiompy/SU2Run/su2run.py
+++ b/ceasiompy/SU2Run/su2run.py
@@ -26,6 +26,7 @@
from pathlib import Path
from ceasiompy.SU2Run.func.su2config import generate_su2_cfd_config
+from ceasiompy.SU2Run.func.su2config_rans import generate_su2_cfd_config_rans
from ceasiompy.SU2Run.func.su2results import get_su2_results
from ceasiompy.utils.ceasiomlogger import get_logger
from ceasiompy.utils.ceasiompyutils import (
@@ -34,7 +35,7 @@
run_software,
)
from ceasiompy.utils.commonnames import CONFIG_CFD_NAME, SU2_FORCES_BREAKDOWN_NAME
-from ceasiompy.utils.commonxpath import SU2_NB_CPU_XPATH
+from ceasiompy.utils.commonxpath import SU2_NB_CPU_XPATH, SU2_CONFIG_RANS_XPATH
from ceasiompy.utils.moduleinterfaces import get_toolinput_file_path, get_tooloutput_file_path
from cpacspy.cpacsfunctions import get_value_or_default, open_tixi
@@ -72,7 +73,6 @@ def run_SU2_multi(wkdir, nb_proc=1):
raise OSError(f"No Case directory has been found in the working directory: {wkdir}")
for config_dir in sorted(case_dir_list):
-
config_cfd = [c for c in config_dir.iterdir() if c.name == CONFIG_CFD_NAME]
if not config_cfd:
@@ -103,7 +103,6 @@ def run_SU2_multi(wkdir, nb_proc=1):
def main(cpacs_path, cpacs_out_path):
-
log.info("----- Start of " + MODULE_NAME + " -----")
tixi = open_tixi(cpacs_path)
@@ -114,7 +113,15 @@ def main(cpacs_path, cpacs_out_path):
# Temporary CPACS to be stored after "generate_su2_cfd_config"
cpacs_tmp_cfg = Path(cpacs_out_path.parent, "ConfigTMP.xml")
- generate_su2_cfd_config(cpacs_path, cpacs_tmp_cfg, results_dir)
+ config_file_type = get_value_or_default(tixi, SU2_CONFIG_RANS_XPATH, "Euler")
+
+ if config_file_type == "RANS":
+ log.info("RANS simulation")
+ generate_su2_cfd_config_rans(cpacs_path, cpacs_tmp_cfg, results_dir)
+ else:
+ log.info("Euler simulation")
+ generate_su2_cfd_config(cpacs_path, cpacs_tmp_cfg, results_dir)
+
run_SU2_multi(results_dir, nb_proc)
get_su2_results(cpacs_tmp_cfg, cpacs_out_path, results_dir)
@@ -122,7 +129,6 @@ def main(cpacs_path, cpacs_out_path):
if __name__ == "__main__":
-
cpacs_path = get_toolinput_file_path(MODULE_NAME)
cpacs_out_path = get_tooloutput_file_path(MODULE_NAME)
diff --git a/ceasiompy/ThermoData/README.md b/ceasiompy/ThermoData/README.md
new file mode 100644
index 000000000..a4bd6082f
--- /dev/null
+++ b/ceasiompy/ThermoData/README.md
@@ -0,0 +1,50 @@
+
+
+# ThermoData
+
+**Categories:** Engine Boundary conditions
+
+**State**: :heavy_check_mark:
+
+
+
+`ThermoData` is a module to provide the outlet conditions of a given engine. It can calculate different operating conditions and save the results in a text file. This module is derived starting from the OpenSource code [pyCycle](https://github.com/OpenMDAO/pycycle) developed by Eric S. Hendricks and Justin S. Gray. It can perform calculations on both turbojet and turbofan engines, with the possibility to customize the parameters. The main parameters that can be modified are: the rotational speed of the shaft/s, the temperature at the inlet of the turbine/s, the efficiency of compressor/s and turbine/s (HP and LP) and the compressor/s pressure ratio. The results are automatically written inside the configuration file of the module `SU2Run` to be able to perform the CFD calculations if needed.
+
+## Inputs
+`ThermoData` can be run on is own by giving an Aeromap that contains Altitude and Mach number with the addition of the Net force that needs to be chosen. Otherwise it can take as an input the values from `CPACS2GMSH` module.
+
+## Analyses
+`ThermoData` compute the values obtained at the engine outlet giving a "EngineBC.dat" file as an output. If the workflow continues with the `SU2Run` run module the results are added to the config file of SU2 to perform the simulation.
+
+## Outputs
+`ThermoData` output is the "EngineBC.dat" file with stored inside for the turbojet engine:
+* T_tot_out= Nozzle total temperature outlet[K],
+* T_stat_out= Nozzle static temperature outlet[K],
+* P_tot_out= Nozzle total pressure outlet [Pa],
+* P_stat_out= Nozzle static pressure outlet [Pa].
+* V_stat_out= Nozzle static velocity outlet[m/s],
+* MN_out= Nozzle Mach number outlet [adim],
+* massflow_stat_out= Nozzle massflow outlet [Kg/s]\\
+
+For the turbofan engine are added also the values at the exit of the bypass nozzle.
+
+
+
+## Installation or requirements
+`ThermoData` needs the installation of the openMDAO and pycycle suite that are included in the python environment of CEASIOMpy.
+
+## Limitations
+
+To be able to change the engine parameters other than those given as input (altitude, mach, net force) the turbojet and turbofan functions in the module must be modified by coding.
+
+## More information
+
+* [pyCycle Github repository](https://github.com/OpenMDAO/pycycle)
+
+* [OpenMDAO documentation ](https://openmdao.org/newdocs/versions/latest/main.html)
+
+* [SU2 website](https://su2code.github.io/)
+
+* [SU2 Github repository](https://github.com/su2code/SU2)
+
+
diff --git a/ceasiompy/ThermoData/__specs__.py b/ceasiompy/ThermoData/__specs__.py
new file mode 100644
index 000000000..7fe6bb9b5
--- /dev/null
+++ b/ceasiompy/ThermoData/__specs__.py
@@ -0,0 +1,67 @@
+from pathlib import Path
+from ceasiompy.utils.moduleinterfaces import CPACSInOut
+from ceasiompy.utils.commonxpath import (
+ SU2_FIXED_CL_XPATH,
+ SU2_TARGET_CL_XPATH,
+ ENGINE_TYPE_XPATH,
+ RANGE_XPATH,
+)
+
+# ===== Module Status =====
+# True if the module is active
+# False if the module is disabled (not working or not ready)
+module_status = True
+
+# ===== Results directory path =====
+
+RESULTS_DIR = Path("Results", "Thermodata")
+
+# ===== CPACS inputs and outputs =====
+
+cpacs_inout = CPACSInOut()
+
+# ===== Input =====
+
+
+cpacs_inout.add_input(
+ var_name="net_force",
+ var_type=float,
+ default_value=2000,
+ unit="N",
+ descr="Engine net force",
+ xpath=RANGE_XPATH + "/NetForce",
+ gui=True,
+ gui_name="NetForce",
+ gui_group="Cruise",
+)
+
+cpacs_inout.add_input(
+ var_name="engine type",
+ var_type=list,
+ default_value=[0, 1],
+ unit=None,
+ descr="0: TBJ, 1: TBF ",
+ xpath=ENGINE_TYPE_XPATH,
+ gui=True,
+ gui_name="0 for Turbojet 1 for Turbofan",
+ gui_group="User inputs",
+)
+
+
+# ===== Output =====
+
+cpacs_inout.add_output(
+ var_name="target_cl",
+ default_value=None,
+ unit="1",
+ descr="Value of CL to achieve to have a level flight with the given conditions",
+ xpath=SU2_TARGET_CL_XPATH,
+)
+
+cpacs_inout.add_output(
+ var_name="fixed_cl",
+ default_value=None,
+ unit="-",
+ descr="FIXED_CL_MODE parameter for SU2",
+ xpath=SU2_FIXED_CL_XPATH,
+)
diff --git a/ceasiompy/ThermoData/func/turbofan_func.py b/ceasiompy/ThermoData/func/turbofan_func.py
new file mode 100644
index 000000000..27481771f
--- /dev/null
+++ b/ceasiompy/ThermoData/func/turbofan_func.py
@@ -0,0 +1,424 @@
+"""
+CEASIOMpy: Conceptual Aircraft Design Software
+
+Developed by CFS ENGINEERING, 1015 Lausanne, Switzerland
+
+Function to run the PyCycle code for the turbofan engine
+
+Python version: >=3.8
+
+| Author: Francesco Marcucci
+| Creation: 2023-12-12
+
+"""
+import sys
+
+import openmdao.api as om
+
+import pycycle.api as pyc
+
+from scipy.constants import convert_temperature
+
+from ceasiompy.utils.ceasiomlogger import get_logger
+
+log = get_logger()
+
+# =================================================================================================
+# FUNCTIONS
+# =================================================================================================
+
+
+def turbofan_analysis(alt, MN, Fn):
+ class HBTF(pyc.Cycle):
+ def setup(self):
+
+ # Setup the problem by including all the relevant components here
+
+ # Create any relevant short hands here:
+
+ design = self.options["design"]
+
+ USE_TABULAR = False
+ if USE_TABULAR:
+ self.options["thermo_method"] = "TABULAR"
+ self.options["thermo_data"] = pyc.AIR_JETA_TAB_SPEC
+ FUEL_TYPE = "FAR"
+ else:
+ self.options["thermo_method"] = "CEA"
+ self.options["thermo_data"] = pyc.species_data.janaf
+ FUEL_TYPE = "Jet-A(g)"
+
+ # Add subsystems to build the engine deck:
+ self.add_subsystem("fc", pyc.FlightConditions())
+ self.add_subsystem("inlet", pyc.Inlet())
+
+ self.add_subsystem(
+ "fan",
+ pyc.Compressor(map_data=pyc.FanMap, bleed_names=[], map_extrap=True),
+ promotes_inputs=[("Nmech", "LP_Nmech")],
+ )
+ self.add_subsystem("splitter", pyc.Splitter())
+ self.add_subsystem("duct4", pyc.Duct())
+ self.add_subsystem(
+ "lpc",
+ pyc.Compressor(map_data=pyc.LPCMap, map_extrap=True),
+ promotes_inputs=[("Nmech", "LP_Nmech")],
+ )
+ self.add_subsystem("duct6", pyc.Duct())
+ self.add_subsystem(
+ "hpc",
+ pyc.Compressor(
+ map_data=pyc.HPCMap,
+ bleed_names=["cool1", "cool2", "cust"],
+ map_extrap=True,
+ ),
+ promotes_inputs=[("Nmech", "HP_Nmech")],
+ )
+ self.add_subsystem("bld3", pyc.BleedOut(bleed_names=["cool3", "cool4"]))
+ self.add_subsystem("burner", pyc.Combustor(fuel_type=FUEL_TYPE))
+ self.add_subsystem(
+ "hpt",
+ pyc.Turbine(map_data=pyc.HPTMap, bleed_names=["cool3", "cool4"], map_extrap=True),
+ promotes_inputs=[("Nmech", "HP_Nmech")],
+ )
+ self.add_subsystem("duct11", pyc.Duct())
+ self.add_subsystem(
+ "lpt",
+ pyc.Turbine(map_data=pyc.LPTMap, bleed_names=["cool1", "cool2"], map_extrap=True),
+ promotes_inputs=[("Nmech", "LP_Nmech")],
+ )
+ self.add_subsystem("duct13", pyc.Duct())
+ self.add_subsystem("core_nozz", pyc.Nozzle(nozzType="CV", lossCoef="Cv"))
+
+ self.add_subsystem("byp_bld", pyc.BleedOut(bleed_names=["bypBld"]))
+ self.add_subsystem("duct15", pyc.Duct())
+ self.add_subsystem("byp_nozz", pyc.Nozzle(nozzType="CV", lossCoef="Cv"))
+
+ # Create shaft instances. Note that LP shaft has 3 ports! => no gearbox
+ self.add_subsystem(
+ "lp_shaft",
+ pyc.Shaft(num_ports=3),
+ promotes_inputs=[("Nmech", "LP_Nmech")],
+ )
+ self.add_subsystem(
+ "hp_shaft",
+ pyc.Shaft(num_ports=2),
+ promotes_inputs=[("Nmech", "HP_Nmech")],
+ )
+ self.add_subsystem("perf", pyc.Performance(num_nozzles=2, num_burners=1))
+
+ # Now use the explicit connect method to make connections -- connect(, )
+
+ # Connect the inputs to perf group
+ self.connect("inlet.Fl_O:tot:P", "perf.Pt2")
+ self.connect("hpc.Fl_O:tot:P", "perf.Pt3")
+ self.connect("burner.Wfuel", "perf.Wfuel_0")
+ self.connect("inlet.F_ram", "perf.ram_drag")
+ self.connect("core_nozz.Fg", "perf.Fg_0")
+ self.connect("byp_nozz.Fg", "perf.Fg_1")
+
+ # LP-shaft connections
+ self.connect("fan.trq", "lp_shaft.trq_0")
+ self.connect("lpc.trq", "lp_shaft.trq_1")
+ self.connect("lpt.trq", "lp_shaft.trq_2")
+ # HP-shaft connections
+ self.connect("hpc.trq", "hp_shaft.trq_0")
+ self.connect("hpt.trq", "hp_shaft.trq_1")
+ # Ideally expanding flow by conneting flight condition static
+ # pressure to nozzle exhaust pressure
+ self.connect("fc.Fl_O:stat:P", "core_nozz.Ps_exhaust")
+ self.connect("fc.Fl_O:stat:P", "byp_nozz.Ps_exhaust")
+
+ balance = self.add_subsystem("balance", om.BalanceComp())
+ if design:
+ balance.add_balance("W", units="lbm/s", eq_units="lbf")
+ # Here balance.W is implicit state variable that is the OUTPUT of balance object
+ self.connect(
+ "balance.W", "fc.W"
+ ) # Connect the output of balance to the relevant input
+ self.connect(
+ "perf.Fn", "balance.lhs:W"
+ ) # This statement makes perf.Fn the LHS of the balance eqn.
+ self.promotes("balance", inputs=[("rhs:W", "Fn_DES")])
+
+ balance.add_balance("FAR", eq_units="degR", lower=1e-4, val=0.017)
+ self.connect("balance.FAR", "burner.Fl_I:FAR")
+ self.connect("burner.Fl_O:tot:T", "balance.lhs:FAR")
+ self.promotes("balance", inputs=[("rhs:FAR", "T4_MAX")])
+
+ # Note that for the following two balances
+ # the mult val is set to -1 so that the NET torque is zero
+ balance.add_balance(
+ "lpt_PR",
+ val=1.5,
+ lower=1.001,
+ upper=8,
+ eq_units="hp",
+ use_mult=True,
+ mult_val=-1,
+ )
+ self.connect("balance.lpt_PR", "lpt.PR")
+ self.connect("lp_shaft.pwr_in_real", "balance.lhs:lpt_PR")
+ self.connect("lp_shaft.pwr_out_real", "balance.rhs:lpt_PR")
+
+ balance.add_balance(
+ "hpt_PR",
+ val=1.5,
+ lower=1.001,
+ upper=8,
+ eq_units="hp",
+ use_mult=True,
+ mult_val=-1,
+ )
+ self.connect("balance.hpt_PR", "hpt.PR")
+ self.connect("hp_shaft.pwr_in_real", "balance.lhs:hpt_PR")
+ self.connect("hp_shaft.pwr_out_real", "balance.rhs:hpt_PR")
+
+ # Set up all the flow connections:
+ self.pyc_connect_flow("fc.Fl_O", "inlet.Fl_I")
+ self.pyc_connect_flow("inlet.Fl_O", "fan.Fl_I")
+ self.pyc_connect_flow("fan.Fl_O", "splitter.Fl_I")
+ self.pyc_connect_flow("splitter.Fl_O1", "duct4.Fl_I")
+ self.pyc_connect_flow("duct4.Fl_O", "lpc.Fl_I")
+ self.pyc_connect_flow("lpc.Fl_O", "duct6.Fl_I")
+ self.pyc_connect_flow("duct6.Fl_O", "hpc.Fl_I")
+ self.pyc_connect_flow("hpc.Fl_O", "bld3.Fl_I")
+ self.pyc_connect_flow("bld3.Fl_O", "burner.Fl_I")
+ self.pyc_connect_flow("burner.Fl_O", "hpt.Fl_I")
+ self.pyc_connect_flow("hpt.Fl_O", "duct11.Fl_I")
+ self.pyc_connect_flow("duct11.Fl_O", "lpt.Fl_I")
+ self.pyc_connect_flow("lpt.Fl_O", "duct13.Fl_I")
+ self.pyc_connect_flow("duct13.Fl_O", "core_nozz.Fl_I")
+ self.pyc_connect_flow("splitter.Fl_O2", "byp_bld.Fl_I")
+ self.pyc_connect_flow("byp_bld.Fl_O", "duct15.Fl_I")
+ self.pyc_connect_flow("duct15.Fl_O", "byp_nozz.Fl_I")
+
+ # Bleed flows:
+ self.pyc_connect_flow("hpc.cool1", "lpt.cool1", connect_stat=False)
+ self.pyc_connect_flow("hpc.cool2", "lpt.cool2", connect_stat=False)
+ self.pyc_connect_flow("bld3.cool3", "hpt.cool3", connect_stat=False)
+ self.pyc_connect_flow("bld3.cool4", "hpt.cool4", connect_stat=False)
+
+ # Specify solver settings:
+ newton = self.nonlinear_solver = om.NewtonSolver()
+ newton.options["atol"] = 1e-8
+
+ # set this very small, so it never activates and we rely on atol
+ newton.options["rtol"] = 1e-99
+ newton.options["iprint"] = 2
+ newton.options["maxiter"] = 50
+ newton.options["solve_subsystems"] = True
+ newton.options["max_sub_solves"] = 1000
+ newton.options["reraise_child_analysiserror"] = False
+ # ls = newton.linesearch = BoundsEnforceLS()
+ ls = newton.linesearch = om.ArmijoGoldsteinLS()
+ ls.options["maxiter"] = 3
+ ls.options["rho"] = 0.75
+ # ls.options['print_bound_enforce'] = True
+
+ self.linear_solver = om.DirectSolver()
+
+ super().setup()
+
+ class MPhbtf(pyc.MPCycle):
+ def setup(self):
+
+ self.pyc_add_pnt(
+ "DESIGN", HBTF(thermo_method="CEA")
+ ) # Create an instance of the High Bypass ratio Turbofan
+
+ super().setup()
+
+ def viewer(prob, pt, file=sys.stdout):
+ """
+ print a report of all the relevant cycle properties
+ """
+
+ if pt == "DESIGN":
+ MN = prob["DESIGN.fc.Fl_O:stat:MN"]
+ # LPT_PR = prob["DESIGN.balance.lpt_PR"]
+ # HPT_PR = prob["DESIGN.balance.hpt_PR"]
+ # FAR = prob["DESIGN.balance.FAR"]
+ else:
+ MN = prob[pt + ".fc.Fl_O:stat:MN"]
+ # LPT_PR = prob[pt + ".lpt.PR"]
+ # HPT_PR = prob[pt + ".hpt.PR"]
+ # FAR = prob[pt + ".balance.FAR"]
+
+ summary_data = (
+ MN,
+ prob[pt + ".fc.alt"],
+ prob[pt + ".inlet.Fl_O:stat:W"],
+ prob[pt + ".perf.Fn"],
+ prob[pt + ".perf.Fg"],
+ prob[pt + ".inlet.F_ram"],
+ prob[pt + ".perf.OPR"],
+ prob[pt + ".perf.TSFC"],
+ prob[pt + ".splitter.BPR"],
+ )
+
+ print(file=file, flush=True)
+ print(file=file, flush=True)
+ print(file=file, flush=True)
+ print(
+ "----------------------------------------------------------------------------",
+ file=file,
+ flush=True,
+ )
+ print(" POINT:", pt, file=file, flush=True)
+ print(
+ "----------------------------------------------------------------------------",
+ file=file,
+ flush=True,
+ )
+ print(" PERFORMANCE CHARACTERISTICS", file=file, flush=True)
+ print(
+ " Mach Alt W Fn Fg Fram OPR TSFC BPR ",
+ file=file,
+ flush=True,
+ )
+ print(
+ " %7.5f %7.1f %7.3f %7.1f %7.1f %7.1f %7.3f %7.5f %7.3f" % summary_data,
+ file=file,
+ flush=True,
+ )
+
+ fs_names = [
+ "fc.Fl_O",
+ "core_nozz.Fl_O",
+ "byp_nozz.Fl_O",
+ ]
+ fs_full_names = [f"{pt}.{fs}" for fs in fs_names]
+ pyc.print_flow_station(prob, fs_full_names, file=file)
+
+ prob = om.Problem()
+
+ prob.model = MPhbtf()
+
+ prob.setup()
+
+ prob.set_val("DESIGN.fan.PR", 1.685)
+ prob.set_val("DESIGN.fan.eff", 0.8948)
+
+ prob.set_val("DESIGN.lpc.PR", 1.935)
+ prob.set_val("DESIGN.lpc.eff", 0.9243)
+
+ prob.set_val("DESIGN.hpc.PR", 9.369)
+ prob.set_val("DESIGN.hpc.eff", 0.8707)
+
+ prob.set_val("DESIGN.hpt.eff", 0.8888)
+ prob.set_val("DESIGN.lpt.eff", 0.8996)
+
+ prob.set_val("DESIGN.fc.alt", alt * 3.2808399, units="ft")
+ prob.set_val("DESIGN.fc.MN", MN)
+
+ prob.set_val("DESIGN.T4_MAX", 2857, units="degR")
+ prob.set_val("DESIGN.Fn_DES", Fn * 0.2248089431, units="lbf") # 1 N = 0.2248089431 lbf
+
+ # Set initial guesses for balances
+ prob["DESIGN.balance.FAR"] = 0.1
+ prob["DESIGN.balance.W"] = 10.0
+ prob["DESIGN.balance.lpt_PR"] = 2
+ prob["DESIGN.balance.hpt_PR"] = 2.0
+ prob["DESIGN.fc.balance.Pt"] = 2
+ prob["DESIGN.fc.balance.Tt"] = 500.0
+
+ prob.set_solver_print(level=-1)
+ prob.set_solver_print(level=2, depth=1)
+
+ viewer_file = open("hbtf_des_view.out", "w")
+
+ # for PC in [1, 0.9]:
+ # viewer(prob, "DESIGN", file=viewer_file)
+ # prob.run_model()
+ # print()
+
+ viewer(prob, "DESIGN", file=viewer_file)
+ prob.run_model()
+
+ # Obtaining variables names
+
+ # variable_names = open("variables_hptf.txt", "w")
+ # print(
+ # prob.model.list_outputs(val=True, units=True, implicit=False),
+ # file=variable_names,
+ # )
+
+ # BYPASS VARIABLES
+ T_tot_out_byp = convert_temperature(
+ (prob.get_val("DESIGN.byp_nozz.throat_total.flow.Fl_O:tot:T")),
+ "Rankine",
+ "Celsius",
+ )
+ V_stat_out_byp = prob.get_val("DESIGN.byp_nozz.mux.Fl_O:stat:V") * 0.3048
+ MN_out_byp = prob.get_val("DESIGN.byp_nozz.mux.Fl_O:stat:MN")
+ P_tot_out_byp = prob.get_val("DESIGN.byp_nozz.throat_total.flow.Fl_O:tot:P") * 6894.7573 # Pa
+ massflow_stat_out_byp = prob.get_val("DESIGN.byp_nozz.mux.Fl_O:stat:W") * 0.45359237 # kg/s
+ T_stat_out_byp = convert_temperature(
+ (prob.get_val("DESIGN.byp_nozz.mux.Fl_O:stat:T")), "Rankine", "Celsius"
+ ) # celsius
+
+ # CORE VARIABLES
+ T_tot_out_core = convert_temperature(
+ (prob.get_val("DESIGN.core_nozz.throat_total.flow.Fl_O:tot:T")),
+ "Rankine",
+ "Kelvin",
+ )
+ V_stat_out_core = prob.get_val("DESIGN.core_nozz.mux.Fl_O:stat:V") * 0.3048
+ MN_out_core = prob.get_val("DESIGN.core_nozz.mux.Fl_O:stat:MN")
+ P_tot_out_core = (
+ prob.get_val("DESIGN.core_nozz.throat_total.flow.Fl_O:tot:P") * 6894.7573
+ ) # Pa
+ massflow_stat_out_core = prob.get_val("DESIGN.core_nozz.mux.Fl_O:stat:W") * 0.45359237 # kg/s
+ T_stat_out_core = convert_temperature(
+ (prob.get_val("DESIGN.core_nozz.mux.Fl_O:stat:T")), "Rankine", "Kelvin"
+ ) # celsius
+
+ return (
+ T_tot_out_byp,
+ V_stat_out_byp,
+ MN_out_byp,
+ P_tot_out_byp,
+ massflow_stat_out_byp,
+ T_stat_out_byp,
+ T_tot_out_core,
+ V_stat_out_core,
+ MN_out_core,
+ P_tot_out_core,
+ massflow_stat_out_core,
+ T_stat_out_core,
+ )
+
+
+def write_hbtf_file(
+ file,
+ T_tot_out_byp,
+ V_stat_out_byp,
+ MN_out_byp,
+ P_tot_out_byp,
+ massflow_stat_out_byp,
+ T_stat_out_byp,
+ T_tot_out_core,
+ V_stat_out_core,
+ MN_out_core,
+ P_tot_out_core,
+ massflow_stat_out_core,
+ T_stat_out_core,
+):
+
+ file.write(f"T_tot_out_core = {T_tot_out_core} [K]\n")
+ file.write(f"V_stat_out_core = {V_stat_out_core} [m/s]\n")
+ file.write(f"MN_out_core = {MN_out_core} [adim]\n")
+ file.write(f"P_tot_out_core = {P_tot_out_core} [Pa]\n")
+ file.write(f"massflow_out_core = {massflow_stat_out_core} [kg/s]\n")
+ file.write(f"T_stat_out_core = {T_stat_out_core} [K]\n")
+ file.write(f"T_tot_out_byp = {T_tot_out_byp} [K]\n")
+ file.write(f"V_stat_out_byp = {V_stat_out_byp} [m/s]\n")
+ file.write(f"MN_out_byp = {MN_out_byp} [adim]\n")
+ file.write(f"P_tot_out_byp = {P_tot_out_byp} [Pa]\n")
+ file.write(f"massflow_stat_out_byp = {massflow_stat_out_byp} [kg/s]\n")
+ file.write(f"T_stat_out_byp = {T_stat_out_byp} [K]\n")
+
+ log.info("hbtf.dat file generated!")
+
+ return file
diff --git a/ceasiompy/ThermoData/func/turbojet_func.py b/ceasiompy/ThermoData/func/turbojet_func.py
new file mode 100644
index 000000000..1e02b3550
--- /dev/null
+++ b/ceasiompy/ThermoData/func/turbojet_func.py
@@ -0,0 +1,222 @@
+"""
+CEASIOMpy: Conceptual Aircraft Design Software
+
+Developed by CFS ENGINEERING, 1015 Lausanne, Switzerland
+
+Function to run the PyCycle code for the turbojet engine
+
+Python version: >=3.8
+
+| Author: Francesco Marcucci
+| Creation: 2023-12-12
+
+"""
+
+import openmdao.api as om
+
+import pycycle.api as pyc
+
+from scipy.constants import convert_temperature
+
+from ceasiompy.utils.ceasiomlogger import get_logger
+
+log = get_logger()
+
+# =================================================================================================
+# FUNCTIONS
+# =================================================================================================
+
+
+def turbojet_analysis(alt, MN, Fn):
+ class Turbojet(pyc.Cycle):
+ def setup(self):
+
+ USE_TABULAR = True
+
+ if USE_TABULAR:
+ self.options["thermo_method"] = "TABULAR"
+ self.options["thermo_data"] = pyc.AIR_JETA_TAB_SPEC
+ FUEL_TYPE = "FAR"
+
+ design = self.options["design"]
+
+ # Add engine elements
+ self.add_subsystem("fc", pyc.FlightConditions())
+ self.add_subsystem("inlet", pyc.Inlet())
+ self.add_subsystem(
+ "comp",
+ pyc.Compressor(map_data=pyc.AXI5, map_extrap=True),
+ promotes_inputs=["Nmech"],
+ )
+ self.add_subsystem("burner", pyc.Combustor(fuel_type=FUEL_TYPE))
+ self.add_subsystem(
+ "turb", pyc.Turbine(map_data=pyc.LPT2269), promotes_inputs=["Nmech"]
+ )
+ self.add_subsystem("nozz", pyc.Nozzle(nozzType="CD", lossCoef="Cv"))
+ self.add_subsystem("shaft", pyc.Shaft(num_ports=2), promotes_inputs=["Nmech"])
+ self.add_subsystem("perf", pyc.Performance(num_nozzles=1, num_burners=1))
+
+ # Connect flow stations
+ self.pyc_connect_flow("fc.Fl_O", "inlet.Fl_I", connect_w=False)
+ self.pyc_connect_flow("inlet.Fl_O", "comp.Fl_I")
+ self.pyc_connect_flow("comp.Fl_O", "burner.Fl_I")
+ self.pyc_connect_flow("burner.Fl_O", "turb.Fl_I")
+ self.pyc_connect_flow("turb.Fl_O", "nozz.Fl_I")
+
+ # Make other non-flow connections
+ # Connect turbomachinery elements to shaft
+ self.connect("comp.trq", "shaft.trq_0")
+ self.connect("turb.trq", "shaft.trq_1")
+
+ # Connnect nozzle exhaust to freestream static conditions
+ self.connect("fc.Fl_O:stat:P", "nozz.Ps_exhaust")
+
+ # Connect outputs to perfomance element
+ self.connect("inlet.Fl_O:tot:P", "perf.Pt2")
+ self.connect("comp.Fl_O:tot:P", "perf.Pt3")
+ self.connect("burner.Wfuel", "perf.Wfuel_0")
+ self.connect("inlet.F_ram", "perf.ram_drag")
+ self.connect("nozz.Fg", "perf.Fg_0")
+
+ # Add balances for design and off-design
+ balance = self.add_subsystem("balance", om.BalanceComp())
+ if design:
+
+ balance.add_balance("W", units="lbm/s", eq_units="lbf", rhs_name="Fn_target")
+ self.connect("balance.W", "inlet.Fl_I:stat:W")
+ self.connect("perf.Fn", "balance.lhs:W")
+
+ balance.add_balance(
+ "FAR", eq_units="degR", lower=1e-4, val=0.017, rhs_name="T4_target"
+ )
+ self.connect("balance.FAR", "burner.Fl_I:FAR")
+ self.connect("burner.Fl_O:tot:T", "balance.lhs:FAR")
+
+ balance.add_balance(
+ "turb_PR", val=1.5, lower=1.001, upper=8, eq_units="hp", rhs_val=0.0
+ )
+ self.connect("balance.turb_PR", "turb.PR")
+ self.connect("shaft.pwr_net", "balance.lhs:turb_PR")
+
+ newton = self.nonlinear_solver = om.NewtonSolver()
+ newton.options["atol"] = 1e-6
+ newton.options["rtol"] = 1e-6
+ newton.options["iprint"] = 2
+ newton.options["maxiter"] = 15
+ newton.options["solve_subsystems"] = True
+ newton.options["max_sub_solves"] = 100
+ newton.options["reraise_child_analysiserror"] = False
+
+ self.linear_solver = om.DirectSolver()
+
+ super().setup()
+
+ class MPTurbojet(pyc.MPCycle):
+ def setup(self):
+ self.pyc_add_pnt("DESIGN", Turbojet())
+
+ self.set_input_defaults("DESIGN.Nmech", 8070.0, units="rpm")
+ self.set_input_defaults("DESIGN.inlet.MN", 0.60)
+ self.set_input_defaults("DESIGN.comp.MN", 0.020) # .2
+ self.set_input_defaults("DESIGN.burner.MN", 0.020) # .2
+ self.set_input_defaults("DESIGN.turb.MN", 0.4)
+
+ self.pyc_add_cycle_param("burner.dPqP", 0.03)
+ self.pyc_add_cycle_param("nozz.Cv", 0.99)
+
+ super().setup()
+
+ # defining the optimization problem
+
+ prob = om.Problem()
+
+ prob.model = MPTurbojet()
+ prob.setup(check=False)
+
+ # define input values
+ prob.set_val("DESIGN.fc.alt", alt * 3.2808399, units="ft")
+ prob.set_val("DESIGN.fc.MN", MN)
+ prob.set_val(
+ "DESIGN.balance.Fn_target", Fn * 0.2248089431, units="lbf"
+ ) # 1 N = 0.2248089431 lbf
+ prob.set_val("DESIGN.balance.T4_target", 2370.0, units="degR")
+ prob.set_val("DESIGN.comp.PR", 13.5)
+ prob.set_val("DESIGN.comp.eff", 0.83)
+ prob.set_val("DESIGN.turb.eff", 0.86)
+
+ if Fn > 2700:
+ prob["DESIGN.balance.FAR"] = 0.0175506829934
+ prob["DESIGN.balance.W"] = 75.453135137
+ prob["DESIGN.balance.turb_PR"] = 4.46138725662
+ prob["DESIGN.fc.balance.Pt"] = 14.6955113159
+ prob["DESIGN.fc.balance.Tt"] = 518.665288153
+
+ elif 2700 <= Fn < 850:
+ prob["DESIGN.balance.FAR"] = 0.0175506829934
+ prob["DESIGN.balance.W"] = 50.453135137
+ prob["DESIGN.balance.turb_PR"] = 4.46138725662
+ prob["DESIGN.fc.balance.Pt"] = 14.6955113159
+ prob["DESIGN.fc.balance.Tt"] = 518.665288153
+
+ else:
+ prob["DESIGN.balance.FAR"] = 0.0175506829934
+ prob["DESIGN.balance.W"] = 10.453135137
+ prob["DESIGN.balance.turb_PR"] = 4.46138725662
+ prob["DESIGN.fc.balance.Pt"] = 14.6955113159
+ prob["DESIGN.fc.balance.Tt"] = 518.665288153
+
+ prob.set_solver_print(level=-1)
+ # prob.set_solver_print(level=2, depth=1)
+
+ prob.run_model()
+
+ print()
+
+ # command to visualize all the output of the system
+ # a = prob.model.list_outputs(val=True, units=True)
+
+ T_tot_out = convert_temperature(
+ (prob.get_val("DESIGN.nozz.throat_total.flow.Fl_O:tot:T")), "Rankine", "Kelvin"
+ )
+ V_stat_out = (prob.get_val("DESIGN.nozz.mux.Fl_O:stat:V")) * 0.3048
+ MN_out = prob.get_val("DESIGN.nozz.mux.Fl_O:stat:MN")
+ T_stat_out = convert_temperature(
+ prob.get_val("DESIGN.nozz.mux.Fl_O:stat:T"), "Rankine", "Kelvin"
+ ) # celsius
+ massflow_stat_out = prob.get_val("DESIGN.nozz.mux.Fl_O:stat:W") * 0.45359237 # kg/s
+ P_tot_out = prob.get_val("DESIGN.nozz.throat_total.flow.Fl_O:tot:P") * 6894.7573 # Pa
+ P_stat_out = prob.get_val("DESIGN.nozz.mux.Fl_O:stat:P") * 6894.7573 # Pa
+
+ return (
+ T_tot_out,
+ V_stat_out,
+ MN_out,
+ P_tot_out,
+ massflow_stat_out,
+ T_stat_out,
+ P_stat_out,
+ )
+
+
+def write_turbojet_file(
+ file,
+ T_tot_out,
+ V_stat_out,
+ MN_out,
+ P_tot_out,
+ massflow_stat_out,
+ T_stat_out,
+ P_stat_out,
+):
+
+ file.write(f"T_tot_out = {T_tot_out} [K]\n")
+ file.write(f"V_stat_out = {V_stat_out} [m/s]\n")
+ file.write(f"MN_out = {MN_out} [adim]\n")
+ file.write(f"P_tot_out = {P_tot_out} [Pa]\n")
+ file.write(f"massflow_out = {massflow_stat_out} [kg/s]\n")
+ file.write(f"T_stat_out = {T_stat_out} [K]\n")
+ file.write(f"P_stat_out = {P_stat_out} [Pa]\n")
+
+ log.info("turbojet.dat file generated!")
+
+ return file
diff --git a/ceasiompy/ThermoData/tests/test_thermodata.py b/ceasiompy/ThermoData/tests/test_thermodata.py
new file mode 100644
index 000000000..38ba491a7
--- /dev/null
+++ b/ceasiompy/ThermoData/tests/test_thermodata.py
@@ -0,0 +1,177 @@
+"""
+CEASIOMpy: Conceptual Aircraft Design Software
+
+Developed by CFS ENGINEERING, 1015 Lausanne, Switzerland
+
+Test functions of 'ceasiompy/ThermoData/func/turbojet_func.py'
+
+Python version: >=3.8
+
+
+| Author : Francesco Marcucci
+| Creation: 2024-02-09
+
+"""
+# =================================================================================================
+# IMPORTS
+# =================================================================================================
+
+from pathlib import Path
+
+import sys
+
+import numpy as np
+
+from ceasiompy.ThermoData.func.turbojet_func import (
+ turbojet_analysis,
+ write_turbojet_file,
+)
+
+from ceasiompy.ThermoData.func.turbofan_func import (
+ turbofan_analysis,
+ write_hbtf_file,
+)
+
+sys.path.append("/home/cfse/Stage_Francesco/Thermodata")
+
+# =================================================================================================
+# CLASSES
+# =================================================================================================
+
+
+# =================================================================================================
+# FUNCTIONS
+# =================================================================================================
+
+
+def test_turbojet_func():
+ """Test function 'turbojet_analysis'"""
+ alt = 1000
+ MN = 0.3
+ Fn = 2000
+ new_sol = turbojet_analysis(alt, MN, Fn)
+ correct_sol = np.array(
+ [
+ 1006.6296749,
+ 798.79704119,
+ 1.50478324,
+ 326886.19429314,
+ 2.89168807,
+ 727.76800186,
+ 89874.50518856,
+ ]
+ ).reshape((7, 1))
+ np.testing.assert_almost_equal(new_sol, correct_sol, 3)
+
+
+def test_write_turbojet_file(tmp_path):
+ """Test function 'write_turbojet_file'"""
+ T_tot_out = 300
+ V_stat_out = 200
+ MN_out = 0.5
+ P_tot_out = 100000
+ massflow_stat_out = 50
+ T_stat_out = 400
+ P_stat_out = 200000
+
+ test_thermodata_path = Path(tmp_path, "EngineBC.dat")
+
+ with open(test_thermodata_path, "w") as file:
+ write_turbojet_file(
+ file,
+ T_tot_out,
+ V_stat_out,
+ MN_out,
+ P_tot_out,
+ massflow_stat_out,
+ T_stat_out,
+ P_stat_out,
+ )
+
+ with open(test_thermodata_path, "r") as file:
+ content = [line.strip() for line in file.readlines()]
+ content.append("")
+
+ # print("content=", content)
+ # expected_content = test_thermodata_path.read_text().split("\n")
+ # print("expected_content=", expected_content)
+
+ assert test_thermodata_path.read_text().split("\n") == content
+
+
+def test_turbofan_func():
+ """Test function 'turbofan_analysis'"""
+ alt = 1000
+ MN = 0.3
+ Fn = 2000
+ new_sol_tuple = turbofan_analysis(alt, MN, Fn)
+ new_sol = np.concatenate(new_sol_tuple)
+ correct_sol = np.array(
+ [
+ 65.04136793,
+ 323.23292298,
+ 0.95295201,
+ 161198.12289882,
+ 1.78761453,
+ 13.08931717,
+ 1270.61685498,
+ 724.51083329,
+ 1.0,
+ 964349.32127573,
+ 1.42155524,
+ 1139.22870449,
+ ]
+ )
+ np.testing.assert_almost_equal(new_sol, correct_sol, 3)
+
+
+def test_write_hbtf_file(tmp_path):
+ """Test function 'write_hbtf_file'"""
+ T_tot_out_core = 300
+ V_stat_out_core = 200
+ MN_out_core = 0.5
+ P_tot_out_core = 100000
+ massflow_out_core = 50
+ T_stat_out_core = 400
+ T_tot_out_byp = 200000
+ V_stat_out_byp = 300
+ MN_out_byp = 1.2
+ P_tot_out_byp = 100000
+ massflow_stat_out_byp = 1.6
+ T_stat_out_byp = 13
+
+ test_thermodata_path = Path(tmp_path, "EngineBC.dat")
+
+ with open(test_thermodata_path, "w") as file:
+ write_hbtf_file(
+ file,
+ T_tot_out_core,
+ V_stat_out_core,
+ MN_out_core,
+ P_tot_out_core,
+ massflow_out_core,
+ T_stat_out_core,
+ T_tot_out_byp,
+ V_stat_out_byp,
+ MN_out_byp,
+ P_tot_out_byp,
+ massflow_stat_out_byp,
+ T_stat_out_byp,
+ )
+
+ with open(test_thermodata_path, "r") as file:
+ content = [line.strip() for line in file.readlines()]
+ content.append("")
+
+ assert test_thermodata_path.read_text().split("\n") == content
+
+
+# =================================================================================================
+# MAIN
+# =================================================================================================
+
+if __name__ == "__main__":
+
+ print("Test ThermoData")
+ print("To run test use the following command:")
+ print(">> pytest -v")
diff --git a/ceasiompy/ThermoData/thermodata.py b/ceasiompy/ThermoData/thermodata.py
new file mode 100644
index 000000000..9b31e3b5c
--- /dev/null
+++ b/ceasiompy/ThermoData/thermodata.py
@@ -0,0 +1,217 @@
+"""
+CEASIOMpy: Conceptual Aircraft Design Software
+
+Developed by CFS ENGINEERING, 1015 Lausanne, Switzerland
+
+Calculate outlet conditions fot turbojet and turbofan engines by using the open source code PyCycle
+and saving those conditions in a text file
+
+| Author: Francesco Marcucci
+| Creation: 2023-12-12
+"""
+
+# =================================================================================================
+# IMPORTS
+# =================================================================================================
+
+import shutil
+from pathlib import Path
+
+from ceasiompy.ThermoData.func.turbofan_func import (
+ turbofan_analysis,
+ write_hbtf_file,
+)
+
+from ceasiompy.ThermoData.func.turbojet_func import (
+ turbojet_analysis,
+ write_turbojet_file,
+)
+
+from ceasiompy.utils.ceasiomlogger import get_logger
+
+from ceasiompy.utils.moduleinterfaces import (
+ get_toolinput_file_path,
+ get_tooloutput_file_path,
+)
+from ceasiompy.utils.ceasiompyutils import get_results_directory
+from ceasiompy.utils.commonxpath import (
+ ENGINE_TYPE_XPATH,
+ ENGINE_BC,
+ RANGE_XPATH,
+ SU2_AEROMAP_UID_XPATH,
+)
+from ceasiompy.utils.commonnames import (
+ ENGINE_BOUNDARY_CONDITIONS,
+)
+from cpacspy.cpacsfunctions import (
+ get_value_or_default,
+ add_float_vector,
+ create_branch,
+)
+from cpacspy.cpacspy import CPACS
+
+
+log = get_logger()
+
+MODULE_DIR = Path(__file__).parent
+MODULE_NAME = MODULE_DIR.name
+
+
+# =================================================================================================
+# CLASSES
+# =================================================================================================
+
+
+# =================================================================================================
+# FUNCTIONS
+# =================================================================================================
+
+
+def thermo_data_run(cpacs_path, cpacs_out_path, wkdir):
+ """Running the PyCycle code by choosing between turbojet or turbofan engine
+
+ Args
+ cpacs_path (Path): Path to CPACS file
+ cpacs_out_path (Path):Path to CPACS output file
+ wkdir (str): Path to the working directory
+ """
+
+ if not wkdir.exists():
+ raise OSError(f"The working directory : {wkdir} does not exit!")
+
+ cpacs = CPACS(cpacs_path)
+ tixi = cpacs.tixi
+
+ Fn = get_value_or_default(tixi, RANGE_XPATH + "/NetForce", 2000)
+
+ default_aeromap = cpacs.create_aeromap("DefaultAeromap")
+ default_aeromap.description = "AeroMap created automatically"
+ mach = get_value_or_default(cpacs.tixi, RANGE_XPATH + "/cruiseMach", 0.3)
+ alt = get_value_or_default(cpacs.tixi, RANGE_XPATH + "/cruiseAltitude", 10000)
+ default_aeromap.add_row(alt=alt, mach=mach, aos=0.0, aoa=0.0)
+ default_aeromap.save()
+ alt_list = [alt]
+ mach_list = [mach]
+ aeromap_uid = get_value_or_default(cpacs.tixi, SU2_AEROMAP_UID_XPATH, "DefaultAeromap")
+ log.info(f"{aeromap_uid} has been created")
+
+ aeromap_list = cpacs.get_aeromap_uid_list()
+
+ if aeromap_list:
+ aeromap_default = aeromap_list[0]
+ log.info(f"The aeromap is {aeromap_default}")
+ aeromap_uid = get_value_or_default(cpacs.tixi, SU2_AEROMAP_UID_XPATH, aeromap_default)
+ activate_aeromap = cpacs.get_aeromap_by_uid(aeromap_uid)
+ alt_list = activate_aeromap.get("altitude").tolist()
+ mach_list = activate_aeromap.get("machNumber").tolist()
+ T_tot_out_array = []
+ P_tot_out_array = []
+
+ for case_nb, alt in enumerate(alt_list):
+ alt = alt_list[case_nb]
+ MN = mach_list[case_nb]
+ case_dir_name = f"Case{str(case_nb).zfill(2)}_alt{alt}_mach{round(MN, 2)}"
+ case_dir_path = Path(wkdir, case_dir_name)
+
+ if not case_dir_path.exists():
+ case_dir_path.mkdir()
+
+ EngineBC = Path(case_dir_path, ENGINE_BOUNDARY_CONDITIONS)
+
+ f = open(EngineBC, "w")
+
+ engine_type = get_value_or_default(tixi, ENGINE_TYPE_XPATH, 0)
+ create_branch(cpacs.tixi, ENGINE_BC)
+
+ if engine_type == 0:
+ (
+ T_tot_out,
+ V_stat_out,
+ MN_out,
+ P_tot_out,
+ massflow_stat_out,
+ T_stat_out,
+ P_stat_out,
+ ) = turbojet_analysis(alt, MN, Fn)
+
+ T_tot_out_array.append(T_tot_out)
+ P_tot_out_array.append(P_tot_out)
+
+ f = write_turbojet_file(
+ file=f,
+ T_tot_out=T_tot_out,
+ V_stat_out=V_stat_out,
+ MN_out=MN_out,
+ P_tot_out=P_tot_out,
+ massflow_stat_out=massflow_stat_out,
+ T_stat_out=T_stat_out,
+ P_stat_out=P_stat_out,
+ )
+
+ else:
+
+ (
+ T_tot_out_byp,
+ V_stat_out_byp,
+ MN_out_byp,
+ P_tot_out_byp,
+ massflow_stat_out_byp,
+ T_stat_out_byp,
+ T_tot_out_core,
+ V_stat_out_core,
+ MN_out_core,
+ P_tot_out_core,
+ massflow_stat_out_core,
+ T_stat_out_core,
+ ) = turbofan_analysis(alt, MN, Fn)
+
+ T_tot_out_array.append(T_tot_out_core)
+ P_tot_out_array.append(P_tot_out_core)
+
+ f = write_hbtf_file(
+ file=f,
+ T_tot_out_byp=T_tot_out_byp,
+ V_stat_out_byp=V_stat_out_byp,
+ MN_out_byp=MN_out_byp,
+ P_tot_out_byp=P_tot_out_byp,
+ massflow_stat_out_byp=massflow_stat_out_byp,
+ T_stat_out_byp=T_stat_out_byp,
+ T_tot_out_core=T_tot_out_core,
+ V_stat_out_core=V_stat_out_core,
+ MN_out_core=MN_out_core,
+ P_tot_out_core=P_tot_out_core,
+ massflow_stat_out_core=massflow_stat_out_core,
+ T_stat_out_core=T_stat_out_core,
+ )
+ add_float_vector(tixi, ENGINE_BC + "/temperatureOutlet", T_tot_out_array)
+ add_float_vector(tixi, ENGINE_BC + "/pressureOutlet", P_tot_out_array)
+
+ cpacs.save_cpacs(cpacs_out_path, overwrite=True)
+
+
+# =================================================================================================
+# MAIN
+# =================================================================================================
+
+
+def main(cpacs_path, cpacs_out_path):
+
+ log.info("----- Start of " + MODULE_NAME + " -----")
+
+ results_dir = get_results_directory("ThermoData")
+
+ thermo_data_run(cpacs_path, cpacs_out_path, results_dir)
+
+ folder_name = "reports"
+
+ shutil.rmtree(folder_name)
+
+ log.info("----- End of " + MODULE_NAME + " -----")
+
+
+if __name__ == "__main__":
+
+ cpacs_path = get_toolinput_file_path(MODULE_NAME)
+ cpacs_out_path = get_tooloutput_file_path(MODULE_NAME)
+
+ main(cpacs_path, cpacs_out_path)
diff --git a/ceasiompy/utils/commonnames.py b/ceasiompy/utils/commonnames.py
index 2fcdabd97..0de79f1c1 100644
--- a/ceasiompy/utils/commonnames.py
+++ b/ceasiompy/utils/commonnames.py
@@ -43,3 +43,6 @@
# WEIGHT & BALANCE
MTOM_FIGURE_NAME = "MTOM_Prediction.png"
+
+# PYCYCLE
+ENGINE_BOUNDARY_CONDITIONS = "EngineBC.dat"
diff --git a/ceasiompy/utils/commonxpath.py b/ceasiompy/utils/commonxpath.py
index 0b9e28a57..a297e734e 100644
--- a/ceasiompy/utils/commonxpath.py
+++ b/ceasiompy/utils/commonxpath.py
@@ -138,6 +138,7 @@
SU2_DEF_MESH_XPATH = SU2_XPATH + "/availableDeformedMesh"
SU2_ACTUATOR_DISK_XPATH = SU2_XPATH + "/options/includeActuatorDisk"
+SU2_CONFIG_RANS_XPATH = SU2_XPATH + "/options/config_type"
# RANGE
RANGE_LD_RATIO_XPATH = CEASIOMPY_XPATH + "/ranges/lDRatio"
@@ -180,3 +181,7 @@
CHECK_LONGITUDINAL_STABILITY_XPATH = STABILITY_XPATH + "/stabilityToCheck/longitudinal"
CHECK_DIRECTIONAL_STABILITY_XPATH = STABILITY_XPATH + "/stabilityToCheck/directional"
CHECK_LATERAL_STABILITY_XPATH = STABILITY_XPATH + "/stabilityToCheck/lateral"
+
+# PYCYCLE
+ENGINE_TYPE_XPATH = CEASIOMPY_XPATH + "/ThermoData"
+ENGINE_BC = CEASIOMPY_XPATH + "/BC"
diff --git a/environment.yml b/environment.yml
index b191c3297..d798c246f 100644
--- a/environment.yml
+++ b/environment.yml
@@ -47,3 +47,4 @@ dependencies:
- smt==1.2.0
- streamlit==1.14.0
- streamlit-autorefresh
+ - om-pycycle
diff --git "a/src/streamlit/pages/01_ \342\236\241_Workflow.py" "b/src/streamlit/pages/01_ \342\236\241_Workflow.py"
index 27ac8306e..ac987512f 100644
--- "a/src/streamlit/pages/01_ \342\236\241_Workflow.py"
+++ "b/src/streamlit/pages/01_ \342\236\241_Workflow.py"
@@ -66,6 +66,7 @@ def section_predefined_workflow():
predefine_workflows = [
["PyTornado", "WeightConventional"],
["CPACS2GMSH", "SU2Run", "SkinFriction"],
+ ["CPACS2GMSH", "ThermoData", "SU2Run"],
["CPACS2SUMO", "SUMOAutoMesh", "SU2Run", "ExportCSV"],
]