Beam Stress from Strains#

../../../../_images/beam_stress_from_strains.png

Modelling an L-frame under a distributed load and evaluating the normal stress from strains:

  • Create an L-frame from two beams with a rectangular cross-section

  • Pin both far ends, leaving the corner between them rigid

  • Apply a uniform distributed load to the horizontal beam

  • Retrieve Young’s modulus from the material properties

  • Evaluate the maximum normal stress at the outer fibres from member strains

Keywords:
L-frame rectangular cross-section distributed load nodal support material properties member strains normal stress
from math import inf
from dlubal.api import rfem, common

# -------------------------------------------------------
# This example demonstrates how to model two perpendicular
# beams with a rectangular cross-section. The top
# (horizontal) beam is loaded by a distributed load.
# While neglecting self-weight, the maximum normal stress
# on the horizontal beam is evaluated from member strains.
# -------------------------------------------------------

# Editable parameters (SI units)
MODEL_NAME = "beam_stress_from_strains"
MATERIAL = "S235"
CROSS_SECTION_WIDTH = 0.025        # b [m]
CROSS_SECTION_HEIGHT = 0.05        # h [m]
VERTICAL_BEAM_LENGTH = 1.0         # [m]
HORIZONTAL_BEAM_LENGTH = 1.0       # [m]
DISTRIBUTED_LOAD = 10000.0         # q [N/m] = 10 kN/m


def define_structure() -> list:
    """Define and return a list of structural objects."""

    return [
        rfem.structure_core.Material(
            no=1,
            name=MATERIAL,
        ),

        rfem.structure_core.CrossSection(
            no=1,
            type=rfem.structure_core.CrossSection.TYPE_PARAMETRIC_MASSIVE_I,
            material=1,
            b=CROSS_SECTION_WIDTH,
            h=CROSS_SECTION_HEIGHT,
            shear_stiffness_deactivated=True,
        ),

        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=0.0, coordinate_2=0.0, coordinate_3=-VERTICAL_BEAM_LENGTH),
        rfem.structure_core.Node(no=3, coordinate_1=HORIZONTAL_BEAM_LENGTH, coordinate_2=0.0, coordinate_3=-VERTICAL_BEAM_LENGTH),

        rfem.structure_core.Line(no=1, definition_nodes=[1, 2]),
        rfem.structure_core.Line(no=2, definition_nodes=[2, 3]),

        # The two members meet rigidly at node 2. That corner is what makes the
        # frame statically indeterminate: it carries part of the span moment
        # into the vertical member instead of leaving the beam simply supported.
        rfem.structure_core.Member(no=1, line=1, cross_section_start=1),
        rfem.structure_core.Member(no=2, line=2, cross_section_start=1),

        # Pinned at both far ends: translations restrained, free to rotate about
        # the global Y axis, which is the axis this frame bends about.
        rfem.types_for_nodes.NodalSupport(
            no=1,
            nodes=[1, 3],
            spring=common.Vector3d(x=inf, y=inf, z=inf),
            rotational_restraint=common.Vector3d(x=inf, y=0, z=inf),
        )
    ]


def define_loading() -> list:
    """Define and return a list of loading objects."""

    return [
        rfem.loading.StaticAnalysisSettings(no=1),

        rfem.loading.LoadCase(
            no=1,
            static_analysis_settings=1,
            self_weight_active=False,
        ),

        rfem.loads.MemberLoad(
            no=1,
            members=[2],
            load_case=1,
            load_type=rfem.loads.MemberLoad.LOAD_TYPE_FORCE,
            load_distribution=rfem.loads.MemberLoad.LOAD_DISTRIBUTION_UNIFORM,
            magnitude=DISTRIBUTED_LOAD,
        ),
    ]


with rfem.Application() as rfem_app:

    app_info = rfem_app.get_application_info()
    print(f"\nApplication Info:\n{app_info}")

    # Modelling
    rfem_app.close_all_models(save_changes=False)
    rfem_app.create_model(name=MODEL_NAME)

    rfem_app.delete_all_objects()
    rfem_app.create_object_list(define_structure() + define_loading())

    # Calculation
    calculation_info = rfem_app.calculate_all(skip_warnings=True)
    print(f"\nCalculation Info:\n{calculation_info}")

    # Results
    member_strains_df = rfem_app.get_results(
        results_type=rfem.results.STATIC_ANALYSIS_MEMBERS_STRAINS,
        filters=[
            rfem.results.ResultsFilter(column_id="member_no", filter_expression="2"),
            rfem.results.ResultsFilter(column_id="loading", filter_expression="LC1"),
        ],
    ).data

    # Material Properties
    material = rfem_app.get_object(
        rfem.structure_core.Material(no=1)
    )
    material_values_tree = material.material_values.rows[0].material_values_tree
    young_modulus = common.tree_table.get_values_by_key(
        tree=material_values_tree,
        key='e', # young_modulus
        occurrence=1
    )
    if young_modulus is None:
        raise RuntimeError("Young's modulus could not be read from the material values.")
    print(f"young_modulus: {young_modulus}")

    # Stress evaluation (horizontal beam).
    # The load acts in the global X-Z plane, so the beam bends about its local y
    # axis and the normal stress varies over the section HEIGHT. Only the axial
    # strain and that curvature contribute:
    #
    #     sigma(z) = E * (epsilon_x -+ kappa_y * h/2)   at the outer fibres
    #
    # The sign convention does not matter here, because the maximum is taken
    # over the absolute values of both fibres.
    required_columns = ["epsilon_x", "kappa_y"]
    missing = [c for c in required_columns if c not in member_strains_df.columns]
    if missing:
        raise RuntimeError(
            f"Expected columns {missing} not found. "
            f"Available columns: {list(member_strains_df.columns)}"
        )

    half_height = CROSS_SECTION_HEIGHT / 2
    sigma_top = young_modulus * (member_strains_df["epsilon_x"] - member_strains_df["kappa_y"] * half_height)
    sigma_bottom = young_modulus * (member_strains_df["epsilon_x"] + member_strains_df["kappa_y"] * half_height)

    sigma_max = max(sigma_top.abs().max(), sigma_bottom.abs().max())
    print(f"sigma_max: {sigma_max} Pa ({sigma_max / 1e6:.2f} MPa)")