Wind Simulation from Wizard#

../../../../../_images/wind_simulation_user_defined.png

Letting the ‘Wind Simulation’ load wizard generate the wind load cases from a range of wind directions:

Activate the wind simulation add-on and orient the global Z-axis upwards
Model a closed box building of four walls and a flat roof on a fixed base
Define a horizontal terrain plane as the floor of the CFD wind tunnel
Create the shrink wrapping, surface roughness, wind profile and analysis settings
Create the wind simulation load wizard over a range of wind directions
Calculate an unrelated load case, which makes the wizard generate the wind load cases
Read the generated load cases from the wizard and calculate them
Read the wind force resultants and the deformations per wind direction

Keywords:

wind simulation load wizard CFD wind direction generated load case wind tunnel

from math import inf, radians

from dlubal.api import rfem, common

# -------------------------------------------------------
# This example runs the wind simulation add-on on a simple closed box building,
# but the wind load cases are not written by hand - they are generated by the
# 'Wind Simulation' load wizard from a range of wind directions.
#
# The wizard generates during the preparation of a calculation. There is no
# separate generation call on the API, so the workflow has two calculation
# steps:
#
#   1. calculate an unrelated, already existing load case. The preparation
#      processes every active wizard in the model regardless of the selection,
#      so the wind load cases are created without running the CFD solver.
#   2. read the generated load cases from the wizard and calculate them.
#
# 'WindSimulation.generate_into_load_cases' is the read-only result of step 1 -
# one row per wind direction, holding the direction, the generated load case
# and the wind profile. Writing into it has no effect.
#
# Every wind direction is a separate CFD run. With the settings below one
# direction takes roughly 3-4 minutes on a regular workstation, so the whole
# script runs for several minutes.
#
# See 'wind_simulation_user_defined.py' for the same model with the wind load
# cases defined explicitly.
# -------------------------------------------------------

# Building dimensions [m]
LENGTH = 10.0     # along global X
WIDTH = 8.0       # along global Y
HEIGHT = 6.0      # along global Z
THICKNESS = 0.2   # wall and roof thickness

# Range of wind directions around the global Z-axis, clockwise [deg]
WIND_DIRECTION_START = 0.0
WIND_DIRECTION_END = 90.0
WIND_DIRECTION_STEP = 90.0


def define_structure() -> list:
    """Define a closed box building - four walls and a flat roof."""

    return [

        # Material
        rfem.structure_core.Material(
            no=1,
            name='C30/37 | EN 1992-1-1:2004/A1:2014',
        ),

        # Nodes - base
        rfem.structure_core.Node(no=1, coordinate_1=0.0, coordinate_2=0.0, coordinate_3=0.0),
        rfem.structure_core.Node(no=2, coordinate_1=LENGTH, coordinate_2=0.0, coordinate_3=0.0),
        rfem.structure_core.Node(no=3, coordinate_1=LENGTH, coordinate_2=WIDTH, coordinate_3=0.0),
        rfem.structure_core.Node(no=4, coordinate_1=0.0, coordinate_2=WIDTH, coordinate_3=0.0),

        # Nodes - top
        rfem.structure_core.Node(no=5, coordinate_1=0.0, coordinate_2=0.0, coordinate_3=HEIGHT),
        rfem.structure_core.Node(no=6, coordinate_1=LENGTH, coordinate_2=0.0, coordinate_3=HEIGHT),
        rfem.structure_core.Node(no=7, coordinate_1=LENGTH, coordinate_2=WIDTH, coordinate_3=HEIGHT),
        rfem.structure_core.Node(no=8, coordinate_1=0.0, coordinate_2=WIDTH, coordinate_3=HEIGHT),

        # Lines - base
        rfem.structure_core.Line(no=1, definition_nodes=[1, 2]),
        rfem.structure_core.Line(no=2, definition_nodes=[2, 3]),
        rfem.structure_core.Line(no=3, definition_nodes=[3, 4]),
        rfem.structure_core.Line(no=4, definition_nodes=[4, 1]),

        # Lines - vertical
        rfem.structure_core.Line(no=5, definition_nodes=[1, 5]),
        rfem.structure_core.Line(no=6, definition_nodes=[2, 6]),
        rfem.structure_core.Line(no=7, definition_nodes=[3, 7]),
        rfem.structure_core.Line(no=8, definition_nodes=[4, 8]),

        # Lines - top
        rfem.structure_core.Line(no=9, definition_nodes=[5, 6]),
        rfem.structure_core.Line(no=10, definition_nodes=[6, 7]),
        rfem.structure_core.Line(no=11, definition_nodes=[7, 8]),
        rfem.structure_core.Line(no=12, definition_nodes=[8, 5]),

        # Surfaces - four walls and the roof
        rfem.structure_core.Surface(no=1, boundary_lines=[1, 6, 9, 5]),
        rfem.structure_core.Surface(no=2, boundary_lines=[2, 7, 10, 6]),
        rfem.structure_core.Surface(no=3, boundary_lines=[3, 8, 11, 7]),
        rfem.structure_core.Surface(no=4, boundary_lines=[4, 5, 12, 8]),
        rfem.structure_core.Surface(no=5, boundary_lines=[9, 10, 11, 12]),

        # Thickness
        rfem.structure_core.Thickness(
            no=1,
            material=1,
            uniform_thickness=THICKNESS,
            assigned_to_surfaces=[1, 2, 3, 4, 5],
        ),

        # Fixed line support along the base
        rfem.types_for_lines.LineSupport(
            no=1,
            lines=[1, 2, 3, 4],
            spring=common.Vector3d(x=inf, y=inf, z=inf),
        ),

        # Static analysis settings of the load cases
        rfem.loading.StaticAnalysisSettings(
            no=1,
            analysis_type=rfem.loading.StaticAnalysisSettings.ANALYSIS_TYPE_GEOMETRICALLY_LINEAR,
        ),
    ]


def define_terrain() -> rfem.base_data_objects.Terrain:
    """Define a horizontal terrain plane - the floor of the CFD wind tunnel."""

    return rfem.base_data_objects.Terrain(
        no=1,
        type=rfem.base_data_objects.Terrain.TYPE_HORIZONTAL_PLANE,
        center_of_terrain_z=0.0,
        bounding_box_offset_x=30.0,
        bounding_box_offset_y=30.0,
    )


def define_wind_simulation() -> list:
    """Define the wind simulation objects.

    The CFD model is built from a shrink-wrapped envelope of the structure, so
    a shrink wrapping is required for the main structure and one for the
    surrounding objects. Both are referenced by the mesh settings further down.
    The wind profile describes the atmospheric boundary layer at the tunnel
    inlet, and the load wizard turns the range of wind directions into one wind
    simulation load case per direction.

    Load case no. 1 is not a wind case - it only gives the wizard something to
    calculate, because the generation runs in the preparation of a calculation.
    """

    return [

        # Shrink wrapping - main structure
        rfem.types_for_wind_simulation.ShrinkWrapping(
            no=1,
            simplification_defined_by=rfem.types_for_wind_simulation.ShrinkWrapping.SIMPLIFICATION_DEFINED_BY_LEVEL_OF_DETAILS,
            level_of_detail=2,
            small_openings_closure_type=rfem.types_for_wind_simulation.ShrinkWrapping.SMALL_OPENINGS_CLOSURE_TYPE_PERCENT_OF_MODEL_DIAMETER,
            closure_relative_to_model_parameter=0.05,
            deactivate_shrink_wrapping=False,
        ),

        # Shrink wrapping - surrounding objects
        rfem.types_for_wind_simulation.ShrinkWrapping(
            no=2,
            simplification_defined_by=rfem.types_for_wind_simulation.ShrinkWrapping.SIMPLIFICATION_DEFINED_BY_LEVEL_OF_DETAILS,
            level_of_detail=0,
            small_openings_closure_type=rfem.types_for_wind_simulation.ShrinkWrapping.SMALL_OPENINGS_CLOSURE_TYPE_PERCENT_OF_MODEL_DIAMETER,
            closure_relative_to_model_parameter=0.2,
            deactivate_shrink_wrapping=False,
        ),

        # Surface roughness of the model in the wind tunnel
        rfem.types_for_wind_simulation.RoughnessAndPermeability(
            no=1,
            type_of_surface=rfem.types_for_wind_simulation.RoughnessAndPermeability.TYPE_OF_SURFACE_SMOOTH,
        ),

        # Wind profile at the tunnel inlet
        rfem.load_wizards.WindProfile(
            no=1,
            type=rfem.load_wizards.WindProfile.TYPE_USER_DEFINED,
            user_defined_input_type=rfem.load_wizards.WindProfile.USER_DEFINED_INPUT_TYPE_CONSTANT,
            user_defined_inlet_variables=rfem.load_wizards.WindProfile.USER_DEFINED_INLET_VARIABLES_I_TUL,
            user_defined_wind_profile=rfem.load_wizards.WindProfile.UserDefinedWindProfileTable(
                rows=[
                    # height [m], velocity [m/s], turbulence intensity [-]
                    rfem.load_wizards.WindProfile.UserDefinedWindProfileRow(
                        no=1, height=0.0, velocity=15.0, turbulence_intensity=0.20),
                    rfem.load_wizards.WindProfile.UserDefinedWindProfileRow(
                        no=2, height=5.0, velocity=20.0, turbulence_intensity=0.18),
                    rfem.load_wizards.WindProfile.UserDefinedWindProfileRow(
                        no=3, height=10.0, velocity=23.0, turbulence_intensity=0.16),
                    rfem.load_wizards.WindProfile.UserDefinedWindProfileRow(
                        no=4, height=20.0, velocity=26.0, turbulence_intensity=0.14),
                    rfem.load_wizards.WindProfile.UserDefinedWindProfileRow(
                        no=5, height=40.0, velocity=29.0, turbulence_intensity=0.12),
                ]
            ),
        ),

        # Steady-flow CFD simulation with a reduced iteration count
        rfem.loading.WindSimulationAnalysisSettings(
            no=1,
            simulation_type=rfem.loading.WindSimulationAnalysisSettings.SIMULATION_TYPE_STEADY_FLOW,
            turbulence_model_type=rfem.loading.WindSimulationAnalysisSettings.TURBULENCE_MODEL_TYPE_EPSILON,
            finite_volume_mesh_density=0.2,
            minimum_number_of_iterations=60,
            maximum_number_of_iterations=300,
            residual_type=rfem.loading.WindSimulationAnalysisSettings.RESIDUAL_TYPE_PRESSURE,
            residual_pressure=0.01,
            consider_turbulence=True,
            snap_to_model_edges=True,
            use_potential_flow_solver_for_initial_condition=True,
            member_load_distribution=rfem.loading.WindSimulationAnalysisSettings.MEMBER_LOAD_DISTRIBUTION_CONCENTRATED,
        ),

        # Load case that triggers the generation of the wind load cases
        rfem.loading.LoadCase(
            no=1,
            name='Self-weight',
            analysis_type=rfem.loading.LoadCase.ANALYSIS_TYPE_STATIC_ANALYSIS,
            action_category=rfem.loading.LoadCase.ACTION_CATEGORY_PERMANENT_G,
            static_analysis_settings=1,
            self_weight_active=True,
            self_weight_factor_z=-1.0,
            to_solve=True,
        ),

        # Wind simulation load wizard over the range of wind directions. It has
        # to be active, otherwise the preparation of the calculation skips it
        # and no load cases are generated.
        rfem.load_wizards.WindSimulation(
            no=1,
            active=True,
            wind_definition_type=rfem.load_wizards.WindSimulation.WIND_DEFINITION_TYPE_UNIFORM_WIND_PROFILE,
            wind_profile=1,
            wind_simulation_analysis_settings=1,
            wind_direction_type=rfem.load_wizards.WindSimulation.WIND_DIRECTION_TYPE_UNIFORM,
            uniform_wind_direction_range_start=radians(WIND_DIRECTION_START),
            uniform_wind_direction_range_end=radians(WIND_DIRECTION_END),
            uniform_wind_direction_step=radians(WIND_DIRECTION_STEP),
        ),
    ]


with rfem.Application() as rfem_app:

    # Create an empty model
    rfem_app.close_all_models(save_changes=False)
    rfem_app.create_model(name='wind_simulation_from_wizard')

    # Activate the wind simulation add-on and orient the global Z-axis upwards
    base_data: rfem.BaseData = rfem_app.get_base_data()
    base_data.addons.wind_simulation_active = True
    base_data.general_settings.global_axes_orientation = \
        rfem.BaseData.GeneralSettings.GLOBAL_AXES_ORIENTATION_ZUP
    rfem_app.set_base_data(base_data=base_data)

    # Clean up the objects pre-created by the wind simulation add-on
    rfem_app.delete_all_objects()

    # Create the structure and the wind simulation objects
    rfem_app.create_object_list(
        objs=define_structure() + define_wind_simulation()
    )

    # Terrain no. 1 always exists in the model and is not removed by
    # 'delete_all_objects', so it is updated instead of created
    rfem_app.update_object(
        obj=define_terrain()
    )

    # Mesh settings of the wind simulation - shrink wrapping, terrain, solver run
    mesh_settings = rfem_app.get_mesh_settings()
    mesh_settings.wind_simulation.consider_terrain_enabled = True
    mesh_settings.wind_simulation.shrink_wrapping_main_structure = 1
    mesh_settings.wind_simulation.shrink_wrapping_surrounding_objects = 2
    mesh_settings.wind_simulation.consider_surface_thickness_above_enabled = True
    mesh_settings.wind_simulation.consider_surface_thickness_above_value = 0.005
    mesh_settings.wind_simulation.run_rwind_in_background_enabled = True
    rfem_app.set_mesh_settings(mesh_settings=mesh_settings)

    # Right after creation the wizard has generated nothing yet
    wizard = rfem_app.get_object(rfem.load_wizards.WindSimulation(no=1))
    print(f'Wizard: {wizard.name}')
    print(f'Generated load cases before the calculation: '
          f'{len(wizard.generate_into_load_cases.rows)}')

    # Step 1 | Generate the wind load cases. Calculating the self-weight load
    # case no. 1 is enough: the preparation of the calculation processes every
    # active wizard in the model, and since load case no. 1 is not a wind
    # simulation case the CFD solver is not started.
    rfem_app.calculate_specific(
        loadings=[rfem.ObjectId(no=1, object_type=rfem.ObjectType.OBJECT_TYPE_LOAD_CASE)],
        skip_warnings=True,
    )

    # Read the generated wind load cases from the wizard
    wizard = rfem_app.get_object(rfem.load_wizards.WindSimulation(no=1))
    wind_load_cases = [
        rfem.ObjectId(no=row.load_case, object_type=rfem.ObjectType.OBJECT_TYPE_LOAD_CASE)
        for row in wizard.generate_into_load_cases.rows
    ]

    print(f'\nGenerated load cases after the calculation: '
          f'{len(wizard.generate_into_load_cases.rows)}')
    for row in wizard.generate_into_load_cases.rows:
        load_case = rfem_app.get_object(rfem.loading.LoadCase(no=row.load_case))
        print(f"  LC{row.load_case} | {row.direction * 180.0 / 3.141592653589793:5.1f} deg | "
              f"wind profile {row.wind_profile} | '{load_case.name}' | "
              f"generated: {load_case.is_generated} | wind tunnel "
              f"{load_case.wind_simulation_wind_tunnel_depth:.1f} x "
              f"{load_case.wind_simulation_wind_tunnel_width:.1f} x "
              f"{load_case.wind_simulation_wind_tunnel_height:.1f} m")

    # Step 2 | Calculate the generated wind load cases
    print('\nCalculating - this takes a few minutes per wind direction ...')
    print(rfem_app.calculate_specific(
        loadings=wind_load_cases,
        skip_warnings=True,
    ).message)

    # ----------------------------------------------------------------------
    # Results
    # ----------------------------------------------------------------------

    # Retrieve results per wind direction. The result tables are wider than the
    # terminal, so they are printed with 'to_string()' to keep them untruncated.
    for loading in wind_load_cases:

        load_case = rfem_app.get_object(rfem.loading.LoadCase(no=loading.no))
        print(f'\n=== LC{loading.no} | {load_case.name} ===')

        # The wind simulation results of a load case can be missing even though
        # the calculation reported success, so this is checked explicitly.
        if not rfem_app.has_results(loading=loading).value:
            print('No results - the wind simulation of this load case did not finish.')
            continue

        # Solver messages, e.g. on the convergence of the residual pressure
        errors_and_warnings = rfem_app.get_result_table(
            table=rfem.results.ResultTable.ERRORS_AND_WARNINGS_TABLE,
            loading=loading,
        )
        print(f'\nErrors and Warnings:\n{errors_and_warnings.data.to_string(index=False)}')

        # Summary of the static analysis - the wind force resultants handed over
        # by the CFD solver, checked against the support forces. The 'value'
        # column holds SI base units (N, Nm, m).
        summary = rfem_app.get_result_table(
            table=rfem.results.ResultTable.STATIC_ANALYSIS_SUMMARY_TABLE,
            loading=loading,
        )
        totals = summary.data[summary.data['description'].str.contains(
            'sum of loads|sum of support forces', na=False)]
        print(f'\nWind Force Resultants and Support Forces [N]:\n'
              f'{totals[["description", "value"]].to_string(index=False)}')

    # Deformations of the walls and the roof in the mesh nodes [m, rad].
    # 'get_results' returns every calculated loading, so the self-weight load
    # case is filtered out and only the wind directions are kept.
    deformations = rfem_app.get_results(
        results_type=rfem.results.ResultsType.STATIC_ANALYSIS_SURFACES_GLOBAL_DEFORMATIONS_MESH_NODES
    )
    wind_deformations = deformations.data[deformations.data['loading'].isin(
        [f'LC{loading.no}' for loading in wind_load_cases])]
    print(f'\nSurfaces Global Deformations:\n{wind_deformations.head(10).to_string(index=False)}')

    # Governing deformation per wind direction
    governing = wind_deformations.loc[
        wind_deformations.groupby('loading')['u_abs'].idxmax(),
        ['loading', 'surface_no', 'mesh_node_no', 'u_abs', 'u_x', 'u_y', 'u_z'],
    ]
    print(f'\nMaximum Absolute Deformation per Load Case [m]:\n{governing.to_string(index=False)}')

using System.Globalization;
using Google.Protobuf;
using Common = Dlubal.Api.Common;
using Rfem = Dlubal.Api.Rfem;

// -------------------------------------------------------
// This example runs the wind simulation add-on on a simple closed box building,
// but the wind load cases are not written by hand - they are generated by the
// 'Wind Simulation' load wizard from a range of wind directions.
//
// The wizard generates during the preparation of a calculation. There is no
// separate generation call on the API, so the workflow has two calculation
// steps:
//
//   1. calculate an unrelated, already existing load case. The preparation
//      processes every active wizard in the model regardless of the selection,
//      so the wind load cases are created without running the CFD solver.
//   2. read the generated load cases from the wizard and calculate them.
//
// 'WindSimulation.GenerateIntoLoadCases' is the read-only result of step 1 -
// one row per wind direction, holding the direction, the generated load case
// and the wind profile. Writing into it has no effect.
//
// Every wind direction is a separate CFD run. With the settings below one
// direction takes roughly 3-4 minutes on a regular workstation, so the whole
// script runs for several minutes.
//
// See 'WindSimulationUserDefined.cs' for the same model with the wind load
// cases defined explicitly.
// -------------------------------------------------------

// Building dimensions [m]
const double Length = 10.0;     // along global X
const double Width = 8.0;       // along global Y
const double Height = 6.0;      // along global Z
const double Thickness = 0.2;   // wall and roof thickness

// Range of wind directions around the global Z-axis, clockwise [deg]
const double WindDirectionStart = 0.0;
const double WindDirectionEnd = 90.0;
const double WindDirectionStep = 90.0;

static double Radians(double degrees) => degrees * Math.PI / 180.0;

static double Degrees(double radians) => radians * 180.0 / Math.PI;

// Define a closed box building - four walls and a flat roof.
static List<IMessage> DefineStructure()
{
    return new List<IMessage>
    {
        // Material
        new Rfem.StructureCore.Material { No = 1, Name = "C30/37 | EN 1992-1-1:2004/A1:2014" },

        // Nodes - base
        new Rfem.StructureCore.Node { No = 1, Coordinate1 = 0.0, Coordinate2 = 0.0, Coordinate3 = 0.0 },
        new Rfem.StructureCore.Node { No = 2, Coordinate1 = Length, Coordinate2 = 0.0, Coordinate3 = 0.0 },
        new Rfem.StructureCore.Node { No = 3, Coordinate1 = Length, Coordinate2 = Width, Coordinate3 = 0.0 },
        new Rfem.StructureCore.Node { No = 4, Coordinate1 = 0.0, Coordinate2 = Width, Coordinate3 = 0.0 },

        // Nodes - top
        new Rfem.StructureCore.Node { No = 5, Coordinate1 = 0.0, Coordinate2 = 0.0, Coordinate3 = Height },
        new Rfem.StructureCore.Node { No = 6, Coordinate1 = Length, Coordinate2 = 0.0, Coordinate3 = Height },
        new Rfem.StructureCore.Node { No = 7, Coordinate1 = Length, Coordinate2 = Width, Coordinate3 = Height },
        new Rfem.StructureCore.Node { No = 8, Coordinate1 = 0.0, Coordinate2 = Width, Coordinate3 = Height },

        // Lines - base
        new Rfem.StructureCore.Line { No = 1, DefinitionNodes = { 1, 2 } },
        new Rfem.StructureCore.Line { No = 2, DefinitionNodes = { 2, 3 } },
        new Rfem.StructureCore.Line { No = 3, DefinitionNodes = { 3, 4 } },
        new Rfem.StructureCore.Line { No = 4, DefinitionNodes = { 4, 1 } },

        // Lines - vertical
        new Rfem.StructureCore.Line { No = 5, DefinitionNodes = { 1, 5 } },
        new Rfem.StructureCore.Line { No = 6, DefinitionNodes = { 2, 6 } },
        new Rfem.StructureCore.Line { No = 7, DefinitionNodes = { 3, 7 } },
        new Rfem.StructureCore.Line { No = 8, DefinitionNodes = { 4, 8 } },

        // Lines - top
        new Rfem.StructureCore.Line { No = 9, DefinitionNodes = { 5, 6 } },
        new Rfem.StructureCore.Line { No = 10, DefinitionNodes = { 6, 7 } },
        new Rfem.StructureCore.Line { No = 11, DefinitionNodes = { 7, 8 } },
        new Rfem.StructureCore.Line { No = 12, DefinitionNodes = { 8, 5 } },

        // Surfaces - four walls and the roof
        new Rfem.StructureCore.Surface { No = 1, BoundaryLines = { 1, 6, 9, 5 } },
        new Rfem.StructureCore.Surface { No = 2, BoundaryLines = { 2, 7, 10, 6 } },
        new Rfem.StructureCore.Surface { No = 3, BoundaryLines = { 3, 8, 11, 7 } },
        new Rfem.StructureCore.Surface { No = 4, BoundaryLines = { 4, 5, 12, 8 } },
        new Rfem.StructureCore.Surface { No = 5, BoundaryLines = { 9, 10, 11, 12 } },

        // Thickness
        new Rfem.StructureCore.Thickness
        {
            No = 1,
            Material = 1,
            UniformThickness = Thickness,
            AssignedToSurfaces = { 1, 2, 3, 4, 5 },
        },

        // Fixed line support along the base
        new Rfem.TypesForLines.LineSupport
        {
            No = 1,
            Lines = { 1, 2, 3, 4 },
            Spring = new Common.Vector3d
            {
                X = double.PositiveInfinity,
                Y = double.PositiveInfinity,
                Z = double.PositiveInfinity,
            },
        },

        // Static analysis settings of the load cases
        new Rfem.Loading.StaticAnalysisSettings
        {
            No = 1,
            AnalysisType = Rfem.Loading.StaticAnalysisSettings.Types.AnalysisType.GeometricallyLinear,
        },
    };
}

// Define a horizontal terrain plane - the floor of the CFD wind tunnel.
static Rfem.BaseDataObjects.Terrain DefineTerrain()
{
    return new Rfem.BaseDataObjects.Terrain
    {
        No = 1,
        Type = Rfem.BaseDataObjects.Terrain.Types.Type.HorizontalPlane,
        CenterOfTerrainZ = 0.0,
        BoundingBoxOffsetX = 30.0,
        BoundingBoxOffsetY = 30.0,
    };
}

// Define the wind simulation objects.
//
// The CFD model is built from a shrink-wrapped envelope of the structure, so
// a shrink wrapping is required for the main structure and one for the
// surrounding objects. Both are referenced by the mesh settings further down.
// The wind profile describes the atmospheric boundary layer at the tunnel
// inlet, and the load wizard turns the range of wind directions into one wind
// simulation load case per direction.
//
// Load case no. 1 is not a wind case - it only gives the wizard something to
// calculate, because the generation runs in the preparation of a calculation.
static List<IMessage> DefineWindSimulation()
{
    return new List<IMessage>
    {
        // Shrink wrapping - main structure
        new Rfem.TypesForWindSimulation.ShrinkWrapping
        {
            No = 1,
            SimplificationDefinedBy = Rfem.TypesForWindSimulation.ShrinkWrapping.Types.SimplificationDefinedBy.LevelOfDetails,
            LevelOfDetail = 2,
            SmallOpeningsClosureType = Rfem.TypesForWindSimulation.ShrinkWrapping.Types.SmallOpeningsClosureType.PercentOfModelDiameter,
            ClosureRelativeToModelParameter = 0.05,
            DeactivateShrinkWrapping = false,
        },

        // Shrink wrapping - surrounding objects
        new Rfem.TypesForWindSimulation.ShrinkWrapping
        {
            No = 2,
            SimplificationDefinedBy = Rfem.TypesForWindSimulation.ShrinkWrapping.Types.SimplificationDefinedBy.LevelOfDetails,
            LevelOfDetail = 0,
            SmallOpeningsClosureType = Rfem.TypesForWindSimulation.ShrinkWrapping.Types.SmallOpeningsClosureType.PercentOfModelDiameter,
            ClosureRelativeToModelParameter = 0.2,
            DeactivateShrinkWrapping = false,
        },

        // Surface roughness of the model in the wind tunnel
        new Rfem.TypesForWindSimulation.RoughnessAndPermeability
        {
            No = 1,
            TypeOfSurface = Rfem.TypesForWindSimulation.RoughnessAndPermeability.Types.TypeOfSurface.Smooth,
        },

        // Wind profile at the tunnel inlet
        new Rfem.LoadWizards.WindProfile
        {
            No = 1,
            Type = Rfem.LoadWizards.WindProfile.Types.Type.UserDefined,
            UserDefinedInputType = Rfem.LoadWizards.WindProfile.Types.UserDefinedInputType.Constant,
            UserDefinedInletVariables = Rfem.LoadWizards.WindProfile.Types.UserDefinedInletVariables.ITul,
            UserDefinedWindProfile = new Rfem.LoadWizards.WindProfile.Types.UserDefinedWindProfileTable
            {
                Rows =
                {
                    // height [m], velocity [m/s], turbulence intensity [-]
                    new Rfem.LoadWizards.WindProfile.Types.UserDefinedWindProfileRow
                        { No = 1, Height = 0.0, Velocity = 15.0, TurbulenceIntensity = 0.20 },
                    new Rfem.LoadWizards.WindProfile.Types.UserDefinedWindProfileRow
                        { No = 2, Height = 5.0, Velocity = 20.0, TurbulenceIntensity = 0.18 },
                    new Rfem.LoadWizards.WindProfile.Types.UserDefinedWindProfileRow
                        { No = 3, Height = 10.0, Velocity = 23.0, TurbulenceIntensity = 0.16 },
                    new Rfem.LoadWizards.WindProfile.Types.UserDefinedWindProfileRow
                        { No = 4, Height = 20.0, Velocity = 26.0, TurbulenceIntensity = 0.14 },
                    new Rfem.LoadWizards.WindProfile.Types.UserDefinedWindProfileRow
                        { No = 5, Height = 40.0, Velocity = 29.0, TurbulenceIntensity = 0.12 },
                },
            },
        },

        // Steady-flow CFD simulation with a reduced iteration count
        new Rfem.Loading.WindSimulationAnalysisSettings
        {
            No = 1,
            SimulationType = Rfem.Loading.WindSimulationAnalysisSettings.Types.SimulationType.SteadyFlow,
            TurbulenceModelType = Rfem.Loading.WindSimulationAnalysisSettings.Types.TurbulenceModelType.Epsilon,
            FiniteVolumeMeshDensity = 0.2,
            MinimumNumberOfIterations = 60,
            MaximumNumberOfIterations = 300,
            ResidualType = Rfem.Loading.WindSimulationAnalysisSettings.Types.ResidualType.Pressure,
            ResidualPressure = 0.01,
            ConsiderTurbulence = true,
            SnapToModelEdges = true,
            UsePotentialFlowSolverForInitialCondition = true,
            MemberLoadDistribution = Rfem.Loading.WindSimulationAnalysisSettings.Types.MemberLoadDistribution.Concentrated,
        },

        // Load case that triggers the generation of the wind load cases
        new Rfem.Loading.LoadCase
        {
            No = 1,
            Name = "Self-weight",
            AnalysisType = Rfem.Loading.LoadCase.Types.AnalysisType.StaticAnalysis,
            ActionCategory = Rfem.Loading.LoadCase.Types.ActionCategory.PermanentG,
            StaticAnalysisSettings = 1,
            SelfWeightActive = true,
            SelfWeightFactorZ = -1.0,
            ToSolve = true,
        },

        // Wind simulation load wizard over the range of wind directions. It has
        // to be active, otherwise the preparation of the calculation skips it
        // and no load cases are generated.
        new Rfem.LoadWizards.WindSimulation
        {
            No = 1,
            Active = true,
            WindDefinitionType = Rfem.LoadWizards.WindSimulation.Types.WindDefinitionType.UniformWindProfile,
            WindProfile = 1,
            WindSimulationAnalysisSettings = 1,
            WindDirectionType = Rfem.LoadWizards.WindSimulation.Types.WindDirectionType.Uniform,
            UniformWindDirectionRangeStart = Radians(WindDirectionStart),
            UniformWindDirectionRangeEnd = Radians(WindDirectionEnd),
            UniformWindDirectionStep = Radians(WindDirectionStep),
        },
    };
}

static double ParseDoubleInvariant(object? value)
{
    if (value is null) return double.NaN;

    if (value is double d) return d;

    return double.TryParse(
        Convert.ToString(value, CultureInfo.InvariantCulture),
        NumberStyles.Any,
        CultureInfo.InvariantCulture,
        out var parsed) ? parsed : double.NaN;
}

// -------------------------------------------------------
// MAIN SCRIPT
// -------------------------------------------------------

ApplicationRfem? rfemApp = null;

try
{
    rfemApp = new ApplicationRfem();

    // Create an empty model
    rfemApp.close_all_models(saveChanges: false);
    rfemApp.create_model(name: "wind_simulation_from_wizard");

    // Activate the wind simulation add-on and orient the global Z-axis upwards
    var baseData = rfemApp.get_base_data();
    baseData.Addons.WindSimulationActive = true;
    baseData.GeneralSettings.GlobalAxesOrientation =
        Rfem.BaseData.Types.GeneralSettings.Types.GlobalAxesOrientation.Zup;
    rfemApp.set_base_data(baseData: baseData);

    // Clean up the objects pre-created by the wind simulation add-on
    rfemApp.delete_all_objects();

    // Create the structure and the wind simulation objects
    var objects = DefineStructure();
    objects.AddRange(DefineWindSimulation());
    rfemApp.create_object_list(objects);

    // Terrain no. 1 always exists in the model and is not removed by
    // 'delete_all_objects', so it is updated instead of created
    rfemApp.update_object(DefineTerrain());

    // Mesh settings of the wind simulation - shrink wrapping, terrain, solver run
    var meshSettings = rfemApp.get_mesh_settings();
    meshSettings.WindSimulation.ConsiderTerrainEnabled = true;
    meshSettings.WindSimulation.ShrinkWrappingMainStructure = 1;
    meshSettings.WindSimulation.ShrinkWrappingSurroundingObjects = 2;
    meshSettings.WindSimulation.ConsiderSurfaceThicknessAboveEnabled = true;
    meshSettings.WindSimulation.ConsiderSurfaceThicknessAboveValue = 0.005;
    meshSettings.WindSimulation.RunRwindInBackgroundEnabled = true;
    rfemApp.set_mesh_settings(meshSettings: meshSettings);

    // Right after creation the wizard has generated nothing yet
    var wizard = rfemApp.get_object<Rfem.LoadWizards.WindSimulation>(
        new Rfem.LoadWizards.WindSimulation { No = 1 });
    Console.WriteLine($"Wizard: {wizard.Name}");
    Console.WriteLine("Generated load cases before the calculation: " +
                      $"{wizard.GenerateIntoLoadCases.Rows.Count}");

    // Step 1 | Generate the wind load cases. Calculating the self-weight load
    // case no. 1 is enough: the preparation of the calculation processes every
    // active wizard in the model, and since load case no. 1 is not a wind
    // simulation case the CFD solver is not started.
    rfemApp.calculate_specific(
        loadings: new List<Rfem.ObjectId>
        {
            new Rfem.ObjectId { No = 1, ObjectType = Rfem.ObjectType.LoadCase },
        },
        skipWarnings: true);

    // Read the generated wind load cases from the wizard
    wizard = rfemApp.get_object<Rfem.LoadWizards.WindSimulation>(
        new Rfem.LoadWizards.WindSimulation { No = 1 });
    var windLoadCases = wizard.GenerateIntoLoadCases.Rows
        .Select(row => new Rfem.ObjectId { No = row.LoadCase, ObjectType = Rfem.ObjectType.LoadCase })
        .ToList();

    Console.WriteLine("\nGenerated load cases after the calculation: " +
                      $"{wizard.GenerateIntoLoadCases.Rows.Count}");
    foreach (var row in wizard.GenerateIntoLoadCases.Rows)
    {
        var generatedCase = rfemApp.get_object<Rfem.Loading.LoadCase>(
            new Rfem.Loading.LoadCase { No = row.LoadCase });
        Console.WriteLine(
            $"  LC{row.LoadCase} | {Degrees(row.Direction),5:F1} deg | " +
            $"wind profile {row.WindProfile} | '{generatedCase.Name}' | " +
            $"generated: {generatedCase.IsGenerated} | wind tunnel " +
            $"{generatedCase.WindSimulationWindTunnelDepth:F1} x " +
            $"{generatedCase.WindSimulationWindTunnelWidth:F1} x " +
            $"{generatedCase.WindSimulationWindTunnelHeight:F1} m");
    }

    // Step 2 | Calculate the generated wind load cases
    Console.WriteLine("\nCalculating - this takes a few minutes per wind direction ...");
    Console.WriteLine(rfemApp.calculate_specific(
        loadings: windLoadCases,
        skipWarnings: true).Message);

    // ----------------------------------------------------------------------
    // Results
    // ----------------------------------------------------------------------

    // Retrieve results per wind direction.
    foreach (var loading in windLoadCases)
    {
        var loadCase = rfemApp.get_object<Rfem.Loading.LoadCase>(
            new Rfem.Loading.LoadCase { No = loading.No });
        Console.WriteLine($"\n=== LC{loading.No} | {loadCase.Name} ===");

        // The wind simulation results of a load case can be missing even though
        // the calculation reported success, so this is checked explicitly.
        if (!rfemApp.has_results(loading: loading).Value)
        {
            Console.WriteLine("No results - the wind simulation of this load case did not finish.");
            continue;
        }

        // Solver messages, e.g. on the convergence of the residual pressure
        var errorsAndWarnings = rfemApp.get_result_table(
            table: Rfem.Results.ResultTable.ErrorsAndWarningsTable,
            loading: loading);
        Console.WriteLine("\nErrors and Warnings:");
        errorsAndWarnings.Print();

        // Summary of the static analysis - the wind force resultants handed over
        // by the CFD solver, checked against the support forces. The 'value'
        // column holds SI base units (N, Nm, m).
        var summary = rfemApp.get_result_table(
            table: Rfem.Results.ResultTable.StaticAnalysisSummaryTable,
            loading: loading);
        Console.WriteLine("\nWind Force Resultants and Support Forces [N]:");
        foreach (var row in summary.Data.Rows)
        {
            var description = Convert.ToString(row["description"], CultureInfo.InvariantCulture) ?? "";
            if (description.Contains("sum of loads") || description.Contains("sum of support forces"))
            {
                Console.WriteLine($"  {description}: {row["value"]}");
            }
        }
    }

    // Deformations of the walls and the roof in the mesh nodes [m, rad].
    // 'get_results' returns every calculated loading, so the self-weight load
    // case is filtered out and only the wind directions are kept.
    var windLoadCaseNames = windLoadCases.Select(loading => $"LC{loading.No}").ToHashSet();
    var deformations = rfemApp.get_results(
        resultsType: Rfem.Results.ResultsType.StaticAnalysisSurfacesGlobalDeformationsMeshNodes);
    var windDeformations = deformations.Data.Rows
        .Select(row => new
        {
            Loading = Convert.ToString(row["loading"], CultureInfo.InvariantCulture) ?? "",
            SurfaceNo = row["surface_no"],
            MeshNodeNo = row["mesh_node_no"],
            UAbs = ParseDoubleInvariant(row["u_abs"]),
        })
        .Where(r => windLoadCaseNames.Contains(r.Loading))
        .ToList();

    Console.WriteLine("\nSurfaces Global Deformations:");
    foreach (var row in windDeformations.Take(10))
    {
        Console.WriteLine(
            $"  {row.Loading} | surface {row.SurfaceNo} | mesh node {row.MeshNodeNo} | " +
            $"u_abs = {row.UAbs:E3}");
    }

    // Governing deformation per wind direction
    Console.WriteLine("\nMaximum Absolute Deformation per Load Case [m]:");
    var governing = windDeformations
        .Where(r => !double.IsNaN(r.UAbs))
        .GroupBy(r => r.Loading)
        .Select(group => group.OrderByDescending(r => r.UAbs).First())
        .OrderBy(r => r.Loading);

    foreach (var row in governing)
    {
        Console.WriteLine(
            $"  {row.Loading} | surface {row.SurfaceNo} | mesh node {row.MeshNodeNo} | " +
            $"u_abs = {row.UAbs:E3}");
    }
}
catch (Exception ex)
{
    Console.WriteLine($"Error: {ex.Message}");
}
finally
{
    if (rfemApp != null) rfemApp.close_connection();
}