"""Section class for conducting cross-section analyses.
Also contains the MaterialGroup class for grouping elements of different material
properties.
"""
from __future__ import annotations
import warnings
from dataclasses import dataclass, field
from typing import TYPE_CHECKING, Any, Callable, cast
import matplotlib.patches as mpatches
import matplotlib.tri as tri
import numpy as np
import numpy.typing as npt
from matplotlib.colors import ListedColormap
from rich.console import Console
from rich.live import Live
from rich.panel import Panel
from rich.progress import Progress, TaskID
from rich.table import Table
from scipy.sparse import coo_matrix, csc_matrix, linalg
from shapely import Point, Polygon
from shapely.strtree import STRtree
import sectionproperties.analysis.fea as fea
import sectionproperties.analysis.plastic_section as sp_plastic
import sectionproperties.analysis.solver as solver
import sectionproperties.post.post as post
import sectionproperties.post.stress_post as sp_stress_post
import sectionproperties.pre.geometry as sp_geom
import sectionproperties.pre.pre as pre
if TYPE_CHECKING:
import matplotlib.axes
[docs]
class Section:
"""Class for structural cross-sections.
Stores the finite element geometry, mesh and material information and provides
methods to compute the cross-section properties.
Attributes:
geometry (Geometry | CompoundGeometry): Cross-section geometry
materials (list[Material]): List of materials in the geometry
material_groups (list[MaterialGroup]): List of
class:`~sectionproperties.pre.pre.MaterialGroup` objects, which contain the
finite elements and stress results (if applicable) for each defined material
mesh (dict[str, Any]): Finite element mesh generated by ``triangle``
num_nodes (int): Number of nodes in the finite element mesh
elements (list[Tri6]): List of finite element objects describing the
cross-section mesh
mesh_search_tree (STRtree): STRtree to enable efficient lookup of elements in
the mesh
section_props (SectionProperties): Class to store calculated section properties
"""
[docs]
def __init__(
self,
geometry: sp_geom.Geometry | sp_geom.CompoundGeometry,
time_info: bool = False,
) -> None:
"""Inits the Section class.
The constructor extracts information from the provided mesh object and creates
and stores the corresponding finite element objects.
Args:
geometry: Cross-section geometry object used to generate the mesh
time_info: If set to True, a detailed description of the computation and the
time cost is printed to the terminal for every computation performed
Raises:
ValueError: If geometry does not contain a mesh
AssertionError: If the number of materials does not equal the number of
regions
"""
if not hasattr(geometry, "mesh") or not geometry.mesh:
raise ValueError(
"Selected Geometry or CompoundGeometry "
"object does not contain a mesh.\n"
"Try running {geometry}.create_mesh() before adding to "
"a Section object for analysis."
)
self.geometry = geometry
self.time_info = time_info
self.mesh = geometry.mesh
self.materials: list[pre.Material] = []
if isinstance(self.geometry, sp_geom.CompoundGeometry):
self.materials = [geom.material for geom in self.geometry.geoms]
else:
self.materials = [self.geometry.material]
# extract mesh data
self._mesh_nodes = np.array(self.mesh["vertices"], dtype=float)
self._mesh_elements = np.array(self.mesh["triangles"], dtype=int)
self._mesh_attributes = np.array(
self.mesh["triangle_attributes"].T[0], dtype=int
)
# swap mid-node order to retain node ordering consistency
self._mesh_elements[:, [3, 4, 5]] = self._mesh_elements[:, [5, 3, 4]]
# save total number of nodes in mesh
self.num_nodes = len(self._mesh_nodes)
# initialise material_groups variable
self.material_groups: list[MaterialGroup] = []
# check that the right number of material properties are specified
msg = f"Number of materials ({len(self.materials)}), "
msg += f"should match the number of regions ({max(self._mesh_attributes) + 1})."
assert len(self.materials) == max(self._mesh_attributes) + 1, msg
# add a MaterialGroup object to the material_groups list for each uniquely
# encountered material
for idx, material in enumerate(self.materials):
# add the first material to the list
if idx == 0:
self.material_groups.append(
MaterialGroup(num_nodes=self.num_nodes, material=material)
)
else:
# if the material hasn't been encountered
if material not in self.materials[:idx]:
self.material_groups.append(
MaterialGroup(num_nodes=self.num_nodes, material=material)
)
# initialise list holding all element objects
self.elements: list[fea.Tri6] = []
# build the mesh one element at a time
for i, node_ids in enumerate(self._mesh_elements):
# create a list containing the vertex and mid-node coordinates
coords = self._mesh_nodes[node_ids, :].transpose()
# if materials are specified, get the material
if self.materials:
# get attribute index of current element
att_el = self._mesh_attributes[i]
# fetch the material
material = self.materials[att_el]
# if there are no materials specified, use a default material
else: # Should not happen but included as failsafe
material = pre.DEFAULT_MATERIAL
# add tri6 elements to the mesh
new_element = fea.Tri6(i, coords, node_ids, material)
self.elements.append(new_element)
# add element to relevant MaterialGroup
for group in self.material_groups:
if material is group.material:
group.add_element(element=new_element)
break
# create the search tree
if self.geometry.mesh is not None:
p_mesh = [
Polygon(self.geometry.mesh["vertices"][tri][0:3])
for tri in self.geometry.mesh["triangles"]
]
self.mesh_search_tree = STRtree(p_mesh)
# initialise class storing section properties
self.section_props = post.SectionProperties()
[docs]
def calculate_geometric_properties(self) -> None:
r"""Calculates geometric (area) properties.
Calculates the geometric properties of the cross-section and stores them in the
:class:`~sectionproperties.post.post.SectionProperties` object contained in the
``section_props`` variable.
Note:
If materials are specified for the cross-section, the moments of area and
section moduli are elastic modulus weighted.
.. admonition:: Geometric Properties
Below is a list of the section properties calculated by the
``calculate_geometric_properties()`` method:
- Cross-sectional area
- Cross-sectional perimeter
- Cross-sectional mass
- Area weighted material properties (composite only) :math:`E_{eff}`,
:math:`G_{eff}`, :math:`{\nu}_{eff}`
- Modulus weighted area (axial rigidity)
- First moments of area
- Second moments of area about the global axis
- Second moments of area about the centroidal axis
- Elastic centroid
- Centroidal section moduli
- Radii of gyration
- Principal axis properties
"""
def calculate_geom(progress: Progress | None = None) -> None:
# initialise properties
self.section_props.area = 0.0
self.section_props.perimeter = 0.0
self.section_props.mass = 0.0
self.section_props.ea = 0.0
self.section_props.ga = 0.0
self.section_props.qx = 0.0
self.section_props.qy = 0.0
self.section_props.ixx_g = 0.0
self.section_props.iyy_g = 0.0
self.section_props.ixy_g = 0.0
# calculate perimeter
self.section_props.perimeter = self.geometry.calculate_perimeter()
if progress:
task = progress.add_task(
description="[red]Calculating geometric properties",
total=len(self.elements),
)
else:
task = None
# calculate global geometric properties
for el in self.elements:
(
area,
qx,
qy,
ixx_g,
iyy_g,
ixy_g,
e,
g,
rho,
) = el.geometric_properties()
self.section_props.area += area
self.section_props.mass += area * rho
self.section_props.ea += area * e
self.section_props.ga += area * g
self.section_props.qx += qx * e
self.section_props.qy += qy * e
self.section_props.ixx_g += ixx_g * e
self.section_props.iyy_g += iyy_g * e
self.section_props.ixy_g += ixy_g * e
if progress and task is not None:
progress.update(task_id=task, advance=1)
self.section_props.nu_eff = (
self.section_props.ea / (2 * self.section_props.ga) - 1
)
self.section_props.e_eff = self.section_props.ea / self.section_props.area
self.section_props.g_eff = self.section_props.ga / self.section_props.area
self.section_props.calculate_elastic_centroid()
self.section_props.calculate_centroidal_properties(
node_list=self.mesh["vertices"]
)
if progress and task is not None:
msg = "[bold green]:white_check_mark: Geometric analysis complete"
progress.update(task_id=task, description=msg)
# conduct geometric analysis
if self.time_info:
progress = solver.create_progress() # create progress bar
progress_table = Table.grid() # create table to print progress in
panel = Panel.fit(
renderable=progress,
title="[bold]Geometric Analysis",
border_style="red",
padding=(1, 1),
)
progress_table.add_row(panel)
with Live(progress_table, refresh_per_second=10) as live:
calculate_geom(progress=progress) # calculate geometric properties
panel.border_style = "green"
live.refresh()
print() # create a blank line
else:
calculate_geom() # calculate geometric properties
[docs]
def calculate_warping_properties(
self,
solver_type: str = "direct",
) -> None:
"""Calculates warping properties.
Calculates all the warping properties of the cross-section and stores them in
the :class:`~sectionproperties.post.post.SectionProperties` object contained in
the ``section_props`` variable.
Args:
solver_type: Solver used for solving systems of linear equations, either
using the ``"direct"`` method or ``"cgs"`` iterative method
Raises:
RuntimeError: If the geometric properties have not been calculated prior to
calling this method
Note:
If materials are specified, the values calculated for the torsion constant,
warping constant and shear areas are elastic modulus weighted.
Warning:
The geometric properties must be calculated prior to the calculation of the
warping properties.
.. admonition:: Warping Properties
Below is a list of the section properties calculated by the
``calculate_warping_properties()`` method:
- Torsion constant
- Shear centre
- Shear area
- Warping constant
- Monosymmetry constant
"""
# check that a geometric analysis has been performed
if None in [
self.section_props.area,
self.section_props.ixx_c,
self.section_props.cx,
]:
err = "Calculate geometric properties before performing a warping analysis."
raise RuntimeError(err)
# check if any geometry are disjoint
if isinstance(self.geometry, sp_geom.CompoundGeometry):
polygons = [sec_geom.geom for sec_geom in self.geometry.geoms]
disjoint_regions = sp_geom.check_geometry_disjoint(lop=polygons)
if disjoint_regions:
msg = "\nThe section geometry contains disjoint regions which is "
msg += "invalid for warping analysis.\n Please revise your geometry to "
msg += "ensure there is connectivity between all regions.\n Please see "
msg += "https://sectionproperties.rtfd.io/en/stable/user_guide/"
msg += "analysis.html#warping-analysis for more information."
warnings.warn(msg)
# create a new Section with the origin shifted to the centroid for calculation
# of the warping properties such that the Lagrangian multiplier approach can be
# utilised
warping_section = Section(geometry=self.geometry)
# shift the coordinates of each element
# N.B. the mesh class attribute remains unshifted!
cx, cy = self.get_c()
for el in warping_section.elements:
el.coords[0, :] -= cx
el.coords[1, :] -= cy
# entire warping analysis is contained in following definition
def warping_analysis(progress: Progress | None = None) -> None:
# assemble stiffness matrix and load vector for warping function
if progress:
msg = f"[red]Assembling {self.num_nodes}x{self.num_nodes} stiffness "
msg += "matrix"
task = progress.add_task(
description=msg,
total=len(warping_section.elements),
)
k_lg, f_torsion = warping_section.assemble_torsion(
progress=progress, task=task
)
msg = f"[green]:white_check_mark: {self.num_nodes}x{self.num_nodes} "
msg += "stiffness matrix assembled"
progress.update(task_id=task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
k_lg, f_torsion = warping_section.assemble_torsion()
# ILU decomposition of stiffness matrices
def ilu_decomp(
progress: Progress | None = None,
task: TaskID | None = None,
) -> linalg.LinearOperator:
# ILU decomposition on Lagrangian stiffness matrix
k_lg_precond = linalg.LinearOperator(
(self.num_nodes + 1, self.num_nodes + 1), linalg.spilu(k_lg).solve
)
if progress and task:
progress.update(task_id=task, advance=1)
return k_lg_precond
# if the cgs method is used, perform ILU decomposition
if solver_type == "cgs":
if progress:
task = progress.add_task(
description="[red]Performing ILU decomposition",
total=2,
)
k_lg_precond = ilu_decomp(progress=progress, task=task)
msg = "[green]:white_check_mark: ILU decomposition complete"
progress.update(task_id=task, description=msg)
else:
k_lg_precond = ilu_decomp()
# solve for warping function
def solve_warping() -> npt.NDArray[np.float64]:
if solver_type == "cgs" and k_lg_precond:
omega = solver.solve_cgs_lagrange(
k_lg=k_lg, f=f_torsion, m=k_lg_precond
)
elif solver_type == "direct":
omega = solver.solve_direct_lagrange(k_lg=k_lg, f=f_torsion)
else:
raise ValueError("solver_type must be 'cgs' or 'direct'.")
return omega
if progress:
msg = f"[red]Solving warping function ({solver_type})"
task = progress.add_task(description=msg, total=1)
omega = solve_warping()
progress.update(task, advance=1)
msg = "[green]:white_check_mark: Warping function solved "
msg += f"({solver_type})"
progress.update(task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
omega = solve_warping()
# save the warping function
self.section_props.omega = omega
# determine the torsion constant
if self.is_composite():
ixx_c, iyy_c, ixy_c = self.get_eic()
else:
ixx_c, iyy_c, ixy_c = self.get_ic()
self.section_props.j = ixx_c + iyy_c - omega.dot(f_torsion)
# assemble shear function load vectors
if self.is_composite():
nu_eff = self.get_nu_eff()
else:
nu_eff = 0.0
def assemble_shear_load(
progress: Progress | None = None,
task: TaskID | None = None,
) -> tuple[npt.NDArray[np.float64], npt.NDArray[np.float64]]:
f_psi = np.zeros(self.num_nodes)
f_phi = np.zeros(self.num_nodes)
for el in warping_section.elements:
f_psi_el, f_phi_el = el.shear_load_vectors(
ixx=ixx_c,
iyy=iyy_c,
ixy=ixy_c,
nu=nu_eff,
)
f_psi[el.node_ids] += f_psi_el
f_phi[el.node_ids] += f_phi_el
if progress and task:
progress.update(task_id=task, advance=1)
return f_psi, f_phi
if progress:
task = progress.add_task(
description="[red]Assembling shear function vectors",
total=len(warping_section.elements),
)
f_psi, f_phi = assemble_shear_load(progress=progress, task=task)
msg = "[green]:white_check_mark: Shear function vectors assembled"
progress.update(task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
f_psi, f_phi = assemble_shear_load()
# solve for shear functions psi and phi
def solve_shear_functions(
progress: Progress | None = None,
task: TaskID | None = None,
) -> tuple[npt.NDArray[np.float64], npt.NDArray[np.float64]]:
if solver_type == "cgs" and k_lg_precond:
psi_shear = solver.solve_cgs_lagrange(
k_lg=k_lg, f=f_psi, m=k_lg_precond
)
if progress and task:
progress.update(task_id=task, advance=1)
phi_shear = solver.solve_cgs_lagrange(
k_lg=k_lg, f=f_phi, m=k_lg_precond
)
if progress and task:
progress.update(task_id=task, advance=1)
else:
psi_shear = solver.solve_direct_lagrange(k_lg=k_lg, f=f_psi)
if progress and task:
progress.update(task_id=task, advance=1)
phi_shear = solver.solve_direct_lagrange(k_lg=k_lg, f=f_phi)
if progress and task:
progress.update(task_id=task, advance=1)
return psi_shear, phi_shear
if progress:
task = progress.add_task(
description=f"[red]Solving shear functions ({solver_type})",
total=2,
)
psi_shear, phi_shear = solve_shear_functions(
progress=progress, task=task
)
msg = "[green]:white_check_mark: Shear functions solved "
msg += f"({solver_type})"
progress.update(task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
psi_shear, phi_shear = solve_shear_functions()
# save the shear functions
self.section_props.psi_shear = psi_shear
self.section_props.phi_shear = phi_shear
# assemble shear centre and warping moment integrals
def assemble_sc_warping_integrals(
progress: Progress | None = None,
task: TaskID | None = None,
) -> tuple[float, float, float, float, float, float]:
sc_xint = 0.0
sc_yint = 0.0
q_omega = 0.0
i_omega = 0.0
i_xomega = 0.0
i_yomega = 0.0
for el in warping_section.elements:
(
sc_xint_el,
sc_yint_el,
q_omega_el,
i_omega_el,
i_xomega_el,
i_yomega_el,
) = el.shear_warping_integrals(
ixx=ixx_c,
iyy=iyy_c,
ixy=ixy_c,
omega=omega[el.node_ids],
)
sc_xint += sc_xint_el
sc_yint += sc_yint_el
q_omega += q_omega_el
i_omega += i_omega_el
i_xomega += i_xomega_el
i_yomega += i_yomega_el
if progress and task:
progress.update(task, advance=1)
return sc_xint, sc_yint, q_omega, i_omega, i_xomega, i_yomega
if progress:
task = progress.add_task(
description="[red]Assembling shear and warping integrals",
total=len(warping_section.elements),
)
(
sc_xint,
sc_yint,
q_omega,
i_omega,
i_xomega,
i_yomega,
) = assemble_sc_warping_integrals(progress=progress, task=task)
msg = "[green]:white_check_mark: Shear and warping integrals assembled"
progress.update(task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
(
sc_xint,
sc_yint,
q_omega,
i_omega,
i_xomega,
i_yomega,
) = assemble_sc_warping_integrals()
# calculate shear centres (elasticity approach)
phi = self.get_phi()
delta_s = 2 * (1 + nu_eff) * (ixx_c * iyy_c - ixy_c**2)
x_se = (1 / delta_s) * ((nu_eff / 2 * sc_xint) - f_torsion.dot(phi_shear))
y_se = (1 / delta_s) * ((nu_eff / 2 * sc_yint) + f_torsion.dot(psi_shear))
x11_se, y22_se = fea.principal_coordinate(phi=phi, x=x_se, y=y_se)
# calculate shear centres (Trefftz's approach)
x_st = (ixy_c * i_xomega - iyy_c * i_yomega) / (ixx_c * iyy_c - ixy_c**2)
y_st = (ixx_c * i_xomega - ixy_c * i_yomega) / (ixx_c * iyy_c - ixy_c**2)
# save shear centres
self.section_props.delta_s = delta_s
self.section_props.x_se = x_se
self.section_props.y_se = y_se
self.section_props.x11_se = x11_se
self.section_props.y22_se = y22_se
self.section_props.x_st = x_st
self.section_props.y_st = y_st
# calculate warping constant
if self.is_composite():
ea = self.get_ea()
else:
ea = self.get_area()
self.section_props.gamma = (
i_omega - q_omega**2 / ea - y_se * i_xomega + x_se * i_yomega
)
def assemble_shear_deformation(
progress: Progress | None = None,
task: TaskID | None = None,
) -> tuple[float, float, float]:
# assemble shear deformation coefficients
kappa_x = 0.0
kappa_y = 0.0
kappa_xy = 0.0
for el in warping_section.elements:
(kappa_x_el, kappa_y_el, kappa_xy_el) = el.shear_coefficients(
ixx=ixx_c,
iyy=iyy_c,
ixy=ixy_c,
psi_shear=psi_shear[el.node_ids],
phi_shear=phi_shear[el.node_ids],
nu=nu_eff,
)
kappa_x += kappa_x_el
kappa_y += kappa_y_el
kappa_xy += kappa_xy_el
if progress and task:
progress.update(task_id=task, advance=1)
return kappa_x, kappa_y, kappa_xy
if progress:
task = progress.add_task(
description="[red]Assembling shear deformation coefficients",
total=len(warping_section.elements),
)
kappa_x, kappa_y, kappa_xy = assemble_shear_deformation(
progress=progress, task=task
)
msg = "[green]:white_check_mark: Shear deformation coefficients "
msg += "assembled"
progress.update(task_id=task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
kappa_x, kappa_y, kappa_xy = assemble_shear_deformation()
# calculate shear areas wrt global axis
self.section_props.a_sx = delta_s**2 / kappa_x
self.section_props.a_sy = delta_s**2 / kappa_y
self.section_props.a_sxy = delta_s**2 / kappa_xy
# calculate shear areas wrt principal bending axis:
area = self.get_area()
alpha_xx = kappa_x * area / delta_s**2
alpha_yy = kappa_y * area / delta_s**2
alpha_xy = kappa_xy * area / delta_s**2
# rotate the tensor by the principal axis angle
phi_rad = phi * np.pi / 180
r = np.array(
[
[np.cos(phi_rad), np.sin(phi_rad)],
[-np.sin(phi_rad), np.cos(phi_rad)],
]
)
rotated_alpha = r.dot(
np.array([[alpha_xx, alpha_xy], [alpha_xy, alpha_yy]])
).dot(np.transpose(r))
# recalculate the shear area based on the rotated alpha value
self.section_props.a_s11 = area / rotated_alpha[0, 0]
self.section_props.a_s22 = area / rotated_alpha[1, 1]
# calculate the monosymmetry consants
def calculate_monosymmetry_integrals(
progress: Progress | None = None,
task: TaskID | None = None,
) -> tuple[float, float, float, float]:
int_x = 0.0
int_y = 0.0
int_11 = 0.0
int_22 = 0.0
for el in warping_section.elements:
(
int_x_el,
int_y_el,
int_11_el,
int_22_el,
) = el.monosymmetry_integrals(phi=phi)
int_x += int_x_el
int_y += int_y_el
int_11 += int_11_el
int_22 += int_22_el
if progress and task:
progress.update(task_id=task, advance=1)
return (int_x, int_y, int_11, int_22)
if progress:
task = progress.add_task(
description="[red]Assembling monosymmetry integrals",
total=len(self.elements),
)
int_x, int_y, int_11, int_22 = calculate_monosymmetry_integrals(
progress=progress, task=task
)
msg = "[green]:white_check_mark: Monosymmetry integrals assembled"
progress.update(task_id=task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
int_x, int_y, int_11, int_22 = calculate_monosymmetry_integrals()
# calculate the monosymmetry constants
if self.is_composite():
i11_c, i22_c = self.get_eip()
else:
i11_c, i22_c = self.get_ip()
self.section_props.beta_x_plus = -int_x / ixx_c + 2 * y_se
self.section_props.beta_x_minus = int_x / ixx_c - 2 * y_se
self.section_props.beta_y_plus = -int_y / iyy_c + 2 * x_se
self.section_props.beta_y_minus = int_y / iyy_c - 2 * x_se
self.section_props.beta_11_plus = -int_11 / i11_c + 2 * y22_se
self.section_props.beta_11_minus = int_11 / i11_c - 2 * y22_se
self.section_props.beta_22_plus = -int_22 / i22_c + 2 * x11_se
self.section_props.beta_22_minus = int_22 / i22_c - 2 * x11_se
# conduct warping analysis
if self.time_info:
# create warping progress
progress = solver.create_progress()
total_task = progress.add_task(
description="[bold red]Conducting warping analysis...", total=7
)
progress_table = Table.grid()
panel = Panel.fit(
renderable=progress,
title="[bold]Warping Analysis",
border_style="red",
padding=(1, 1),
)
progress_table.add_row(panel)
with Live(progress_table, refresh_per_second=10) as live:
warping_analysis(progress=progress)
msg = "[bold green]:white_check_mark: Warping analysis completed"
progress.update(task_id=total_task, description=msg)
panel.border_style = "green"
live.refresh()
print()
else:
warping_analysis()
[docs]
def calculate_frame_properties(
self,
solver_type: str = "direct",
) -> tuple[float, float, float, float, float, float]:
"""Calculates section properties to be used for frame analysis.
Calculates and returns the properties required for a frame analysis. The
properties are also stored in the
:class:`~sectionproperties.post.post.SectionProperties` object contained in the
``section_props`` variable.
This method is more efficient than running a geometric and warping analysis as
unnecessary values are not calculated.
Args:
solver_type: Solver used for solving systems of linear equations, either
using the ``"direct"`` method or ``"cgs"`` iterative method
Returns:
Cross-section properties to be used for a frame analysis (``area``, ``ixx``,
``iyy``, ``ixy``, ``j``, ``phi``)
Note:
If materials are specified, the values calculated for the moments of area
and the torsion constant are elastic modulus weighted.
.. admonition:: Warping Properties
Below is a list of the section properties calculated by the
``calculate_frame_properties()`` method:
- Cross-sectional area
- Second moments of area about the centroidal axis
- Torsion constant
- Principal axis angle
"""
def calculate_geom(progress: Progress | None = None) -> None:
# initialise geometric properties
self.section_props.area = 0.0
self.section_props.ea = 0.0
self.section_props.qx = 0.0
self.section_props.qy = 0.0
self.section_props.ixx_g = 0.0
self.section_props.iyy_g = 0.0
self.section_props.ixy_g = 0.0
self.section_props.ixx_c = 0.0
self.section_props.iyy_c = 0.0
self.section_props.ixy_c = 0.0
self.section_props.j = 0.0
self.section_props.phi = 0.0
if progress:
task = progress.add_task(
description="[red]Calculating geometric properties",
total=len(self.elements),
)
else:
task = None
# calculate global geometric properties
for el in self.elements:
area, qx, qy, ixx_g, iyy_g, ixy_g, e, _, _ = el.geometric_properties()
self.section_props.area += area
self.section_props.ea += area * e
self.section_props.qx += qx * e
self.section_props.qy += qy * e
self.section_props.ixx_g += ixx_g * e
self.section_props.iyy_g += iyy_g * e
self.section_props.ixy_g += ixy_g * e
if progress and task is not None:
progress.update(task_id=task, advance=1)
# calculate elastic centroid location
self.section_props.calculate_elastic_centroid()
# calculate second moments of area about the centroidal xy axis
if self.is_composite():
ea = self.get_ea()
qx, qy = self.get_eq()
ixx_g, iyy_g, ixy_g = self.get_eig()
else:
ea = self.get_area()
qx, qy = self.get_q()
ixx_g, iyy_g, ixy_g = self.get_ig()
self.section_props.ixx_c = ixx_g - qx**2 / ea
self.section_props.iyy_c = iyy_g - qy**2 / ea
self.section_props.ixy_c = ixy_g - qx * qy / ea
# calculate the principal axis angle
if self.is_composite():
ixx_c, iyy_c, ixy_c = self.get_eic()
else:
ixx_c, iyy_c, ixy_c = self.get_ic()
delta = (((ixx_c - iyy_c) / 2) ** 2 + ixy_c**2) ** 0.5
i11_c = (ixx_c + iyy_c) / 2 + delta
# calculate initial principal axis angle
if abs(ixx_c - i11_c) < 1e-12 * i11_c:
self.section_props.phi = 0.0
else:
self.section_props.phi = np.arctan2(ixx_c - i11_c, ixy_c) * 180 / np.pi
if progress and task is not None:
msg = "[bold green]:white_check_mark: Geometric analysis complete"
progress.update(task_id=task, description=msg)
# conduct geometric analysis
if self.time_info:
progress = solver.create_progress() # create progress bar
progress_table = Table.grid() # create table to print progress in
panel = Panel.fit(
renderable=progress,
title="[bold]Geometric (Frame) Analysis",
border_style="red",
padding=(1, 1),
)
progress_table.add_row(panel)
with Live(progress_table, refresh_per_second=10) as live:
calculate_geom(progress=progress) # calculate geometric properties
panel.border_style = "green"
live.refresh()
print() # create a blank line
else:
calculate_geom() # calculate geometric properties
# entire warping analysis is contained in following definition
def warping_analysis(progress: Progress | None = None) -> None:
# create a new Section with the origin shifted to the centroid for
# calculation of the warping properties
warping_section = Section(geometry=self.geometry)
# shift the coordinates of each element N.B. the mesh class attribute
# remains unshifted
cx, cy = self.get_c()
for el in warping_section.elements:
el.coords[0, :] -= cx
el.coords[1, :] -= cy
# assemble stiffness matrix and load vector for warping function
if progress:
msg = f"[red]Assembling {self.num_nodes}x{self.num_nodes} stiffness "
msg += "matrix"
task = progress.add_task(
description=msg,
total=len(warping_section.elements),
)
k_lg, f_torsion = warping_section.assemble_torsion(
progress=progress, task=task
)
msg = f"[green]:white_check_mark: {self.num_nodes}x{self.num_nodes} "
msg += "stiffness matrix assembled"
progress.update(task_id=task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
k_lg, f_torsion = warping_section.assemble_torsion()
# ILU decomposition of stiffness matrices
def ilu_decomp(
progress: Progress | None = None,
task: TaskID | None = None,
) -> linalg.LinearOperator:
# ILU decomposition on Lagrangian stiffness matrix
k_lg_precond = linalg.LinearOperator(
(self.num_nodes + 1, self.num_nodes + 1), linalg.spilu(k_lg).solve
)
if progress and task:
progress.update(task_id=task, advance=1)
return k_lg_precond
# if the cgs method is used, perform ILU decomposition
if solver_type == "cgs":
if progress:
task = progress.add_task(
description="[red]Performing ILU decomposition",
total=2,
)
k_lg_precond = ilu_decomp(progress=progress, task=task)
msg = "[green]:white_check_mark: ILU decomposition complete"
progress.update(task_id=task, description=msg)
else:
k_lg_precond = ilu_decomp()
# solve for warping function
def solve_warping() -> npt.NDArray[np.float64]:
if solver_type == "cgs" and k_lg_precond:
omega = solver.solve_cgs_lagrange(
k_lg=k_lg, f=f_torsion, m=k_lg_precond
)
elif solver_type == "direct":
omega = solver.solve_direct_lagrange(k_lg=k_lg, f=f_torsion)
else:
raise ValueError("solver_type must be 'cgs' or 'direct'.")
return omega
if progress:
msg = f"[red]Solving warping function ({solver_type})"
task = progress.add_task(description=msg, total=1)
omega = solve_warping()
progress.update(task, advance=1)
msg = "[green]:white_check_mark: Warping function solved "
msg += f"({solver_type})"
progress.update(task, description=msg)
progress.update(cast(TaskID, 0), advance=1)
else:
omega = solve_warping()
# save the warping function
self.section_props.omega = omega
# determine the torsion constant
if self.is_composite():
ixx_c, iyy_c, _ = self.get_eic()
else:
ixx_c, iyy_c, _ = self.get_ic()
self.section_props.j = ixx_c + iyy_c - omega.dot(f_torsion)
# conduct warping analysis
if self.time_info:
# create warping progress
progress = solver.create_progress()
total_task = progress.add_task(
description="[bold red]Conducting warping analysis...", total=2
)
progress_table = Table.grid()
panel = Panel.fit(
renderable=progress,
title="[bold]Warping (Frame) Analysis",
border_style="red",
padding=(1, 1),
)
progress_table.add_row(panel)
with Live(progress_table, refresh_per_second=10) as live:
warping_analysis(progress=progress)
msg = "[bold green]:white_check_mark: Warping analysis completed"
progress.update(task_id=total_task, description=msg)
panel.border_style = "green"
live.refresh()
print()
else:
warping_analysis()
# check all properties were calculated
assert self.section_props.area is not None
assert self.section_props.ixx_c is not None
assert self.section_props.iyy_c is not None
assert self.section_props.ixy_c is not None
assert self.section_props.j is not None
assert self.section_props.phi is not None
return (
self.section_props.area,
self.section_props.ixx_c,
self.section_props.iyy_c,
self.section_props.ixy_c,
self.section_props.j,
self.section_props.phi,
)
[docs]
def calculate_plastic_properties(
self,
verbose: bool = False,
) -> None:
"""Calculates plastic properties.
Calculates the plastic properties of the cross-section and stores them in the
:class:`~sectionproperties.post.post.SectionProperties` object contained in
the ``section_props`` variable.
Args:
verbose: If set to True, the number of iterations required for each plastic
axis is printed to the terminal.
Raises:
RuntimeError: If the geometric properties have not been calculated prior to
calling this method
Note:
If materials are specified, the values calculated for the plastic section
moduli are displayed as plastic moments (i.e :math:`M_p = f_y S`) and the
shape factors are not calculated.
Warning:
The geometric properties must be calculated prior to the calculation of the
plastic properties.
.. admonition:: Plastic Properties
Below is a list of the section properties calculated by the
``calculate_plastic_properties()`` method:
- Plastic centroids (centroidal and principal axes)
- Plastic section moduli (centroidal and principal axes)
- Shape factors, non-composite only (centroidal and principal axe)
"""
# check that a geometric analysis has been performed
if self.section_props.cx is None:
err = "Calculate geometric properties before performing a plastic analysis."
raise RuntimeError(err)
# check if any geometry are overlapped
if isinstance(self.geometry, sp_geom.CompoundGeometry):
polygons = [sec_geom.geom for sec_geom in self.geometry.geoms]
overlapped_regions = sp_geom.check_geometry_overlaps(lop=polygons)
if overlapped_regions:
msg = "\nThe section geometry contains overlapping regions and the area"
msg += " of these overlapped regions will be 'double counted' in the"
msg += " plastic analysis which may result in incorrect values.\n"
msg += "If you do not intend for this double counting to occur, use"
msg += " a subtractive modelling approach to remove the overlapping"
msg += " region.\nPlease see "
msg += "https://sectionproperties.rtfd.io/en/stable/examples/geometry/"
msg += "advanced_geometry.html for more information."
warnings.warn(msg)
def calc_plastic(progress: Progress | None = None) -> None:
plastic_section = sp_plastic.PlasticSection(geometry=self.geometry)
plastic_section.calculate_plastic_properties(
section=self, verbose=verbose, progress=progress
)
if self.time_info:
progress = solver.create_progress()
progress_table = Table.grid()
panel = Panel.fit(
renderable=progress,
title="[bold]Plastic Analysis",
border_style="red",
padding=(1, 1),
)
progress_table.add_row(panel)
with Live(progress_table, refresh_per_second=10) as live:
calc_plastic(progress=progress)
panel.border_style = "green"
live.refresh()
print()
else:
calc_plastic()
[docs]
def calculate_stress(
self,
n: float = 0.0,
vx: float = 0.0,
vy: float = 0.0,
mxx: float = 0.0,
myy: float = 0.0,
m11: float = 0.0,
m22: float = 0.0,
mzz: float = 0.0,
) -> sp_stress_post.StressPost:
"""Calculates cross-section stresses.
Calculates the cross-section stress resulting from design actions and returns
a :class:`~sectionproperties.post.stress_post.StressPost` object allowing the
post-processing of the stress results.
Args:
n: Axial force
vx: Shear force acting in the x-direction
vy: Shear force acting in the y-direction
mxx: Bending moment about the centroidal xx-axis
myy: Bending moment about the centroidal yy-axis
m11: Bending moment about the centroidal 11-axis
m22: Bending moment about the centroidal 22-axis
mzz: Torsion moment about the centroidal zz-axis
Returns:
Object for post-processing cross-section stresses
Raises:
RuntimeError: If a geometric and warping analysis (if required) have not
been performed prior to calling this method
Note:
A geometric analysis must be performed prior to performing a stress
analysis. Further, if the shear force or torsion is non-zero a warping
analysis must also be performed.
"""
# check that a geometric analysis has been performed
if None in [
self.section_props.area,
self.section_props.ixx_c,
self.section_props.cx,
]:
msg = "Perform a geometric analysis before carrying out a stress analysis."
raise RuntimeError(msg)
# check for warping analysis requirement
if vx != 0 or vy != 0 or mzz != 0:
warping = True
if self.section_props.omega is None:
msg = "Perform a warping analysis before carrying out a stress "
msg += "analysis with non-zero shear forces or torsion moment."
raise RuntimeError(msg)
else:
warping = False
def calc_stress(progress: Progress | None = None) -> sp_stress_post.StressPost:
if progress:
task = progress.add_task(
description="[red]Calculating cross-section stresses",
total=len(self.elements),
)
else:
task = None
# create stress post object
stress_post = sp_stress_post.StressPost(section=self)
# get relevant section properties
cx, cy = self.get_c()
phi = self.get_phi()
if self.is_composite():
ea = self.get_ea()
ixx, iyy, ixy = self.get_eic()
i11, i22 = self.get_eip()
else:
ea = self.get_area()
ixx, iyy, ixy = self.get_ic()
i11, i22 = self.get_ip()
if warping:
if self.is_composite():
j = self.get_ej()
nu = self.get_nu_eff()
else:
j = self.get_j()
nu = 0.0
delta_s = self.section_props.delta_s
omega = self.section_props.omega
psi_shear = self.section_props.psi_shear
phi_shear = self.section_props.phi_shear
assert delta_s is not None
assert omega is not None
assert psi_shear is not None
assert phi_shear is not None
else:
j = 0.0
nu = 0.0
delta_s = 0.0
omega = np.zeros(shape=self.num_nodes)
psi_shear = np.zeros(shape=self.num_nodes)
phi_shear = np.zeros(shape=self.num_nodes)
# loop through all material groups
for group in stress_post.material_groups:
# allocate nodal weights vector for nodal averaging
nodal_weights = np.zeros(self.num_nodes)
sr = group.stress_result
# loop through all elements in the material group
for el in group.elements:
# get element omega and psi
omega_el = omega[el.node_ids]
psi_shear_el = psi_shear[el.node_ids]
phi_shear_el = phi_shear[el.node_ids]
(
sig_zz_n_el,
sig_zz_mxx_el,
sig_zz_myy_el,
sig_zz_m11_el,
sig_zz_m22_el,
sig_zx_mzz_el,
sig_zy_mzz_el,
sig_zx_vx_el,
sig_zy_vx_el,
sig_zx_vy_el,
sig_zy_vy_el,
weights,
) = el.element_stress(
n=n,
mxx=mxx,
myy=myy,
m11=m11,
m22=m22,
mzz=mzz,
vx=vx,
vy=vy,
ea=ea,
cx=cx,
cy=cy,
ixx=ixx,
iyy=iyy,
ixy=ixy,
i11=i11,
i22=i22,
phi=phi,
j=j,
nu=nu,
omega=omega_el,
psi_shear=psi_shear_el,
phi_shear=phi_shear_el,
delta_s=delta_s,
)
# add stresses to global vectors
sr.sig_zz_n[el.node_ids] += sig_zz_n_el * weights
sr.sig_zz_mxx[el.node_ids] += sig_zz_mxx_el * weights
sr.sig_zz_myy[el.node_ids] += sig_zz_myy_el * weights
sr.sig_zz_m11[el.node_ids] += sig_zz_m11_el * weights
sr.sig_zz_m22[el.node_ids] += sig_zz_m22_el * weights
sr.sig_zx_mzz[el.node_ids] += sig_zx_mzz_el * weights
sr.sig_zy_mzz[el.node_ids] += sig_zy_mzz_el * weights
sr.sig_zx_vx[el.node_ids] += sig_zx_vx_el * weights
sr.sig_zy_vx[el.node_ids] += sig_zy_vx_el * weights
sr.sig_zx_vy[el.node_ids] += sig_zx_vy_el * weights
sr.sig_zy_vy[el.node_ids] += sig_zy_vy_el * weights
# add nodal weights
nodal_weights[el.node_ids] += weights
if progress and task is not None:
progress.update(task_id=task, advance=1)
# nodal averaging
for i, weight in enumerate(nodal_weights):
if weight != 0:
sr.sig_zz_n[i] *= 1 / weight
sr.sig_zz_mxx[i] *= 1 / weight
sr.sig_zz_myy[i] *= 1 / weight
sr.sig_zz_m11[i] *= 1 / weight
sr.sig_zz_m22[i] *= 1 / weight
sr.sig_zx_mzz[i] *= 1 / weight
sr.sig_zy_mzz[i] *= 1 / weight
sr.sig_zx_vx[i] *= 1 / weight
sr.sig_zy_vx[i] *= 1 / weight
sr.sig_zx_vy[i] *= 1 / weight
sr.sig_zy_vy[i] *= 1 / weight
# calculate combined stresses
sr.calculate_combined_stresses()
if progress and task is not None:
msg = "[bold green]:white_check_mark: Stress analysis complete"
progress.update(task_id=task, description=msg)
return stress_post
if self.time_info:
progress = solver.create_progress()
progress_table = Table.grid()
panel = Panel.fit(
renderable=progress,
title="[bold]Stress Analysis",
border_style="red",
padding=(1, 1),
)
progress_table.add_row(panel)
with Live(progress_table, refresh_per_second=10) as live:
stress_post = calc_stress(progress=progress)
panel.border_style = "green"
live.refresh()
print()
else:
stress_post = calc_stress()
# return the stress_post object
return stress_post
[docs]
def assemble_torsion(
self,
progress: Progress | None = None,
task: TaskID | None = None,
) -> tuple[csc_matrix, npt.NDArray[np.float64]]:
"""Assembles the warping stiffness matrix.
Assembles stiffness matrices to be used for the computation of warping
properties and the torsion load vector (``f_torsion``). A Lagrangian multiplier
(``k_lg``) stiffness matrix is returned. The stiffness matrix is assembled using
the sparse COO format and returned in the sparse CSC format.
Args:
progress: Rich progress object
task: Rich task object
Returns:
Lagrangian multiplier stiffness matrix and torsion load vector (``k_lg``,
``f_torsion``)
"""
# initialise variables
n_size = self.num_nodes # size of matrix
row: list[float] = [] # list holding row indices
col: list[float] = [] # list holding column indices
data: list[float] = [] # list holding stiffness matrix entries
f_torsion = np.zeros(n_size) # force vector array
c_torsion = np.zeros(n_size) # constraint vector array
# loop through all elements in the mesh
for el in self.elements:
# determine number of nodes in the current element
n = len(el.node_ids)
# calculate the element stiffness matrix and torsion load vector
k_el, f_el, c_el = el.torsion_properties()
# assemble the torsion load vector
f_torsion[el.node_ids] += f_el
c_torsion[el.node_ids] += c_el
for node_id in el.node_ids:
row.extend([node_id] * n)
col.extend(list(el.node_ids) * n)
data.extend(k_el.flatten())
if progress and task:
progress.update(task_id=task, advance=1)
# construct Lagrangian multiplier matrix:
# column vector
row.extend(range(n_size))
col.extend([n_size] * n_size)
data.extend(c_torsion)
# row vector
row.extend([n_size] * n_size)
col.extend(range(n_size))
data.extend(c_torsion)
k_lg = coo_matrix(
(data, (row, col)), shape=(n_size + 1, n_size + 1), dtype=float
)
return csc_matrix(k_lg, dtype=float), f_torsion
[docs]
def plot_mesh(
self,
alpha: float = 0.5,
materials: bool = True,
mask: list[bool] | None = None,
title: str = "Finite Element Mesh",
**kwargs: Any,
) -> matplotlib.axes.Axes:
r"""Plots the finite element mesh.
Args:
alpha: Transparency of the mesh outlines: :math:`0 \leq \alpha \leq 1`
materials: If set to True shades the elements with the specified material
colors
mask: Mask array, of length ``num_nodes``, to mask out triangles
title: Plot title
kwargs: Passed to :func:`~sectionproperties.post.post.plotting_context`
Returns:
Matplotlib axes object
Example:
The following example plots the mesh for a rectanglular steel-timber
composite section:
.. plot::
:include-source: True
:caption: Composite section mesh
from sectionproperties.pre import Material
from sectionproperties.pre.library import rectangular_section
from sectionproperties.analysis import Section
steel = Material(
name="Steel",
elastic_modulus=200e3,
poissons_ratio=0.3,
density=7.85e-6,
yield_strength=250,
color="grey",
)
timber = Material(
name="Timber",
elastic_modulus=8e3,
poissons_ratio=0.35,
density=6.5e-7,
yield_strength=20,
color="burlywood",
)
geom_steel = rectangular_section(d=50, b=50, material=steel)
geom_timber = rectangular_section(d=50, b=50, material=timber)
geom = geom_timber.align_to(geom_steel, on="right") + geom_steel
geom.create_mesh(mesh_sizes=[10, 5])
Section(geometry=geom).plot_mesh()
"""
# create plot and setup the plot
with post.plotting_context(title=title, **kwargs) as (fig, ax):
assert ax
# create mesh triangulation
triang = tri.Triangulation(
self._mesh_nodes[:, 0],
self._mesh_nodes[:, 1],
self._mesh_elements[:, 0:3],
mask=mask,
)
# if the material colors are to be displayed
if materials and self.materials:
color_array = []
legend_labels = []
c = [] # Indices of elements for mapping colors
# create an array of finite element colors
for idx, element in enumerate(self.elements):
color_array.append(element.material.color)
c.append(idx)
# create a list of unique material legend entries
for i, material in enumerate(self.materials):
# if the material has not be entered yet
if i == 0 or material not in self.materials[0:i]:
# add the material color and name to the legend list
patch = mpatches.Patch(
color=material.color, label=material.name
)
legend_labels.append(patch)
cmap = ListedColormap(colors=color_array) # custom colormap
# plot the mesh colors
ax.tripcolor(
triang,
c,
cmap=cmap,
)
# display the legend
ax.legend(
loc="center left",
bbox_to_anchor=(1, 0.5),
handles=legend_labels,
)
# plot the mesh
ax.triplot(
triang,
lw=0.5,
color="black",
alpha=alpha,
)
return ax
[docs]
def plot_centroids(
self,
alpha: float = 0.5,
title: str = "Centroids",
**kwargs: Any,
) -> matplotlib.axes.Axes:
r"""Plots the calculated centroids over the mesh.
Plots the elastic centroid, the shear centre, the plastic centroids and the
principal axis, if they have been calculated, on top of the finite element mesh.
Args:
alpha: Transparency of the mesh outlines: :math:`0 \leq \alpha \leq 1`
title: Plot title
kwargs: Passed to :func:`~sectionproperties.post.post.plotting_context`
Returns:
Matplotlib axes object
Example:
The following example analyses a 200 PFC section and displays a plot of the
centroids:
.. plot::
:include-source: True
:caption: 200PFC centroids
from sectionproperties.pre.library import channel_section
from sectionproperties.analysis import Section
geom = channel_section(d=200, b=75, t_f=12, t_w=6, r=12, n_r=8)
geom.create_mesh(mesh_sizes=[20])
section = Section(geometry=geom)
section.calculate_geometric_properties()
section.calculate_warping_properties()
section.calculate_plastic_properties()
section.plot_centroids()
"""
# create plot and setup the plot
with post.plotting_context(title=title, **kwargs) as (fig, ax):
assert ax
# plot the finite element mesh
self.plot_mesh(alpha=alpha, **dict(kwargs, ax=ax))
# if the elastic centroid has been calculated
try:
cx, cy = self.get_c()
ax.scatter(
x=cx,
y=cy,
edgecolors="r",
facecolors="none",
marker="o",
s=100,
label="Elastic centroid",
)
except AssertionError:
pass
# if the shear centre has been calculated
try:
x_s, y_s = self.get_sc()
ax.scatter(x=x_s, y=y_s, c="r", marker="+", s=100, label="Shear centre")
except AssertionError:
pass
# if the global plastic centroid has been calculated
try:
x_pc, y_pc = self.get_pc()
ax.scatter(
x=x_pc,
y=y_pc,
c="r",
marker="x",
s=100,
label="Global plastic centroid",
)
except AssertionError:
pass
# if the principal plastic centroid has been calculated
try:
x11_pc, y22_pc = self.get_pc_p()
ax.scatter(
x11_pc,
y22_pc,
edgecolors="r",
facecolors="none",
marker="s",
s=100,
label="Principal plastic centroid",
)
except AssertionError:
pass
# if the principal axis has been calculated
try:
phi = self.get_phi()
cx, cy = self.get_c()
post.draw_principal_axis(ax=ax, phi=phi * np.pi / 180, cx=cx, cy=cy)
except AssertionError:
pass
# display the legend
ax.legend(loc="center left", bbox_to_anchor=(1, 0.5))
return ax
[docs]
def plot_warping_function(
self,
title: str = "Warping Function",
level: int = 20,
cmap: str = "viridis",
alpha: float = 0.2,
with_lines: bool = True,
**kwargs: Any,
) -> matplotlib.axes.Axes:
r"""Plots the warping function over the mesh.
Args:
title: Plot title
level: Number of contour levels
cmap: Colormap
with_lines: If set to True, contour lines are displayed
alpha: Transparency of the mesh outlines: :math:`0 \leq \alpha \leq 1`
kwargs: Passed to :func:`~sectionproperties.post.post.plotting_context`
Raises:
RuntimeError: If a warping analysis has not been performed
Returns:
Matplotlib axes object
"""
if self.section_props.omega is None:
raise RuntimeError(
"Perform a warping analysis before plotting the warping function."
)
# create plot and setup the plot
with post.plotting_context(title=title, **kwargs) as (fig, ax):
assert ax
# create triangulation
triang = tri.Triangulation(
self._mesh_nodes[:, 0],
self._mesh_nodes[:, 1],
self._mesh_elements[:, 0:3],
)
if with_lines:
ax.tricontour(
triang,
self.section_props.omega,
colors="k",
levels=level,
)
ax.tricontourf(
triang,
self.section_props.omega,
cmap=cmap,
levels=level,
)
# plot the finite element mesh
self.plot_mesh(alpha=alpha, materials=False, **dict(kwargs, ax=ax))
return ax
[docs]
def display_mesh_info(self) -> None:
"""Prints mesh statistics to the command line."""
console = Console()
console.print("Mesh Statistics:", style="bold underline magenta")
console.print(f"- {self.num_nodes} nodes")
console.print(f"- {len(self.elements)} elements")
regions = max(self._mesh_attributes) + 1
text = f"- {regions} region"
if regions > 1:
text += "s"
console.print(text)
console.print()
[docs]
def display_results(
self,
fmt: str = "8.6e",
) -> None:
"""Prints the results that have been calculated to the terminal.
Args:
fmt: Number formatting string, see
https://docs.python.org/3/library/string.html
"""
post.print_results(section=self, fmt=fmt)
[docs]
def is_composite(self) -> bool:
"""Returns whether or not a composite section is being analysed.
If the only material is the default material, a regular analysis is being
conducted. Otherwise, a composite analysis is being conducted.
Returns:
Whether or not a composite section is being analysed.
"""
return list(set(self.materials)) != [pre.DEFAULT_MATERIAL]
[docs]
def get_e_ref(
self,
e_ref: float | pre.Material,
) -> float:
"""Extract transformed elastic modulus from e_ref (float or material).
Args:
e_ref: Reference elastic modulus or material property by which to
transform the results
Returns:
Transformed elastic modulus
"""
if isinstance(e_ref, pre.Material):
return e_ref.elastic_modulus
else:
return e_ref
[docs]
def get_area(self) -> float:
"""Returns the cross-section area.
Returns:
Cross-section area
Raises:
AssertionError: If a geometric analysis has not been performed
"""
try:
assert self.section_props.area is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.area
[docs]
def get_perimeter(self) -> float:
"""Returns the cross-section perimeter.
Returns:
Cross-section perimeter
Raises:
AssertionError: If a geometric analysis has not been performed
"""
try:
assert self.section_props.perimeter is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.perimeter
[docs]
def get_mass(self) -> float:
"""Returns the cross-section mass.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
Returns:
Cross-section mass
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_area()."
raise RuntimeError(msg)
try:
assert self.section_props.mass is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.mass
[docs]
def get_ea(
self,
e_ref: float | pre.Material = 1,
) -> float:
"""Returns the cross-section axial rigidity.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted area (axial rigidity)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_area()."
raise RuntimeError(msg)
# obtain e_ref (reference elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.ea is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.ea / e_ref
[docs]
def get_q(self) -> tuple[float, float]:
"""Returns the cross-section first moments of area.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
First moments of area about the global axis (``qx``, ``qy``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a geometric analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_eq()."
raise RuntimeError(msg)
try:
assert self.section_props.qx is not None
assert self.section_props.qy is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.qx, self.section_props.qy
[docs]
def get_eq(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float]:
"""Returns the modulus-weighted cross-section first moments of area.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus-weighted first moments of area about the global axis (``e.qx``,
``e.qy``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_q()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.qx is not None
assert self.section_props.qy is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.qx / e_ref, self.section_props.qy / e_ref
[docs]
def get_ig(self) -> tuple[float, float, float]:
"""Returns the cross-section global second moments of area.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Second moments of area about the global axis (``ixx_g``, ``iyy_g``,
``ixy_g``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a geometric analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_eig()."
raise RuntimeError(msg)
try:
assert self.section_props.ixx_g is not None
assert self.section_props.ixy_g is not None
assert self.section_props.iyy_g is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.ixx_g,
self.section_props.iyy_g,
self.section_props.ixy_g,
)
[docs]
def get_eig(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float, float]:
"""Returns the modulus-weighted cross-section global second moments of area.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted second moments of area about the global axis (``e.ixx_g``,
``e.iyy_g``, ``e.ixy_g``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_ig()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.ixx_g is not None
assert self.section_props.ixy_g is not None
assert self.section_props.iyy_g is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.ixx_g / e_ref,
self.section_props.iyy_g / e_ref,
self.section_props.ixy_g / e_ref,
)
[docs]
def get_c(self) -> tuple[float, float]:
"""Returns the cross-section elastic centroid.
Returns:
Elastic centroid (``cx``, ``cy``)
Raises:
AssertionError: If a geometric analysis has not been performed
"""
try:
assert self.section_props.cx is not None
assert self.section_props.cy is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.cx, self.section_props.cy
[docs]
def get_ic(self) -> tuple[float, float, float]:
"""Returns the cross-section centroidal second moments of area.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Second moments of area about the centroidal axis (``ixx_c``, ``iyy_c``,
``ixy_c``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a geometric analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_eic()."
raise RuntimeError(msg)
try:
assert self.section_props.ixx_c is not None
assert self.section_props.ixy_c is not None
assert self.section_props.iyy_c is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.ixx_c,
self.section_props.iyy_c,
self.section_props.ixy_c,
)
[docs]
def get_eic(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float, float]:
"""Returns the modulus-weighted cross-section centroidal second moments of area.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted second moments of area about the centroidal axis
(``e.ixx_c``, ``e.iyy_c``, ``e.ixy_c``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_ic()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.ixx_c is not None
assert self.section_props.ixy_c is not None
assert self.section_props.iyy_c is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.ixx_c / e_ref,
self.section_props.iyy_c / e_ref,
self.section_props.ixy_c / e_ref,
)
[docs]
def get_z(self) -> tuple[float, float, float, float]:
"""Returns the cross-section centroidal elastic section moduli.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Elastic section moduli about the centroidal axis with respect to the top and
bottom fibres (``zxx_plus``, ``zxx_minus``, ``zyy_plus``, ``zyy_minus``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a geometric analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_ez()."
raise RuntimeError(msg)
try:
assert self.section_props.zxx_plus is not None
assert self.section_props.zxx_minus is not None
assert self.section_props.zyy_plus is not None
assert self.section_props.zyy_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.zxx_plus,
self.section_props.zxx_minus,
self.section_props.zyy_plus,
self.section_props.zyy_minus,
)
[docs]
def get_ez(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float, float, float]:
"""Returns the modulus-weighted cross-section centroidal elastic section moduli.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted elastic section moduli about the centroidal axis with
respect to the top and bottom fibres (``e.zxx_plus``, ``e.zxx_minus``,
``e.zyy_plus``, ``e.zyy_minus``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_z()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.zxx_plus is not None
assert self.section_props.zxx_minus is not None
assert self.section_props.zyy_plus is not None
assert self.section_props.zyy_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.zxx_plus / e_ref,
self.section_props.zxx_minus / e_ref,
self.section_props.zyy_plus / e_ref,
self.section_props.zyy_minus / e_ref,
)
[docs]
def get_rc(self) -> tuple[float, float]:
"""Returns the cross-section centroidal radii of gyration.
Returns:
Radii of gyration about the centroidal axis (``rx``, ``ry``)
Raises:
AssertionError: If a geometric analysis has not been performed
"""
try:
assert self.section_props.rx_c is not None
assert self.section_props.ry_c is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.rx_c, self.section_props.ry_c
[docs]
def get_ip(self) -> tuple[float, float]:
"""Returns the cross-section principal second moments of area.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Second moments of area about the principal axis (``i11_c``, ``i22_c``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a geometric analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_eip()."
raise RuntimeError(msg)
try:
assert self.section_props.i11_c is not None
assert self.section_props.i22_c is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.i11_c, self.section_props.i22_c
[docs]
def get_eip(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float]:
"""Returns the modulus-weighted cross-section principal second moments of area.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted second moments of area about the principal axis
(``e.i11_c``, ``e.i22_c``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_ip()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.i11_c is not None
assert self.section_props.i22_c is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.i11_c / e_ref, self.section_props.i22_c / e_ref
[docs]
def get_phi(self) -> float:
"""Returns the cross-section principal bending angle.
Returns:
Principal bending axis angle
Raises:
AssertionError: If a geometric analysis has not been performed
"""
try:
assert self.section_props.phi is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.phi
[docs]
def get_zp(self) -> tuple[float, float, float, float]:
"""Returns the cross-section principal elastic section moduli.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Elastic section moduli about the principal axis with respect to the top and
bottom fibres (``z11_plus``, ``z11_minus``, ``z22_plus``, ``z22_minus``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a geometric analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_ezp()."
raise RuntimeError(msg)
try:
assert self.section_props.z11_plus is not None
assert self.section_props.z11_minus is not None
assert self.section_props.z22_plus is not None
assert self.section_props.z22_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.z11_plus,
self.section_props.z11_minus,
self.section_props.z22_plus,
self.section_props.z22_minus,
)
[docs]
def get_ezp(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float, float, float]:
"""Returns the modulus-weighted cross-section principal elastic section moduli.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted elastic section moduli about the principal axis with
respect to the top and bottom fibres (``e.z11_plus``, ``e.z11_minus``,
``e.z22_plus``, ``e.z22_minus``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_zp()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.z11_plus is not None
assert self.section_props.z11_minus is not None
assert self.section_props.z22_plus is not None
assert self.section_props.z22_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return (
self.section_props.z11_plus / e_ref,
self.section_props.z11_minus / e_ref,
self.section_props.z22_plus / e_ref,
self.section_props.z22_minus / e_ref,
)
[docs]
def get_rp(self) -> tuple[float, float]:
"""Returns the cross-section principal radii of gyration.
Returns:
Radii of gyration about the principal axis (``r11``, ``r22``)
Raises:
AssertionError: If a geometric analysis has not been performed
"""
try:
assert self.section_props.r11_c is not None
assert self.section_props.r22_c is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.r11_c, self.section_props.r22_c
[docs]
def get_nu_eff(self) -> float:
"""Returns the cross-section effective Poisson's ratio.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
Returns:
Effective Poisson's ratio
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied)."
raise RuntimeError(msg)
try:
assert self.section_props.nu_eff is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.nu_eff
[docs]
def get_e_eff(self) -> float:
"""Returns the cross-section effective elastic modulus.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
Returns:
Effective elastic modulus based on area
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied)."
raise RuntimeError(msg)
try:
assert self.section_props.e_eff is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.e_eff
[docs]
def get_g_eff(self) -> float:
"""Returns the cross-section effective shear modulus.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
Returns:
Effective shear modulus based on area
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a geometric analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied)."
raise RuntimeError(msg)
try:
assert self.section_props.g_eff is not None
except AssertionError as e:
raise AssertionError("Conduct a geometric analysis.") from e
return self.section_props.g_eff
[docs]
def get_j(self) -> float:
"""Returns the cross-section St Venant torsion constant.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
St. Venant torsion constant
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a warping analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_ej()."
raise RuntimeError(msg)
try:
assert self.section_props.j is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.j
[docs]
def get_ej(
self,
e_ref: float | pre.Material = 1,
) -> float:
"""Returns the modulus-weighted cross-section St Venant torsion constant.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted St. Venant torsion constant
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a warping analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_j()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.j is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.j / e_ref
[docs]
def get_sc(self) -> tuple[float, float]:
"""Returns the cross-section centroidal shear centre (elasticity).
Returns:
Centroidal axis shear centre (elasticity approach) (``x_se``, ``y_se``)
Raises:
AssertionError: If a warping analysis has not been performed
"""
try:
assert self.section_props.x_se is not None
assert self.section_props.y_se is not None
assert self.section_props.cx is not None
assert self.section_props.cy is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
# add centroid location to move section back to original location
x_se = self.section_props.x_se + self.section_props.cx
y_se = self.section_props.y_se + self.section_props.cy
return x_se, y_se
[docs]
def get_sc_p(self) -> tuple[float, float]:
"""Returns the cross-section principal shear centre (elasticity).
Returns:
Principal axis shear centre (elasticity approach) (``x11_se``, ``y22_se``)
Raises:
AssertionError: If a warping analysis has not been performed
"""
try:
assert self.section_props.x11_se is not None
assert self.section_props.y22_se is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.x11_se, self.section_props.y22_se
[docs]
def get_sc_t(self) -> tuple[float, float]:
"""Returns the cross-section centroidal shear centre (Trefftz's method).
Returns:
Centroidal axis shear centre (Trefftz's method) (``x_st``, ``y_st``)
Raises:
AssertionError: If a warping analysis has not been performed
"""
try:
assert self.section_props.x_st is not None
assert self.section_props.y_st is not None
assert self.section_props.cx is not None
assert self.section_props.cy is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
# add centroid location to move section back to original location
x_st = self.section_props.x_st + self.section_props.cx
y_st = self.section_props.y_st + self.section_props.cy
return x_st, y_st
[docs]
def get_gamma(self) -> float:
"""Returns the cross-section warping constant.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Warping constant
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a warping analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using"
msg += " get_egamma()."
raise RuntimeError(msg)
try:
assert self.section_props.gamma is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.gamma
[docs]
def get_egamma(
self,
e_ref: float | pre.Material = 1,
) -> float:
"""Returns the modulus-weighted cross-section warping constant.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted warping constant
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a warping analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_gamma()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.gamma is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.gamma / e_ref
[docs]
def get_as(self) -> tuple[float, float]:
"""Returns the cross-section centroidal axis shear area.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Shear area for loading about the centroidal axis (``a_sx``, ``a_sy``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a warping analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_eas()."
raise RuntimeError(msg)
try:
assert self.section_props.a_sx is not None
assert self.section_props.a_sy is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.a_sx, self.section_props.a_sy
[docs]
def get_eas(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float]:
"""Returns modulus-weighted the cross-section centroidal axis shear area.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted shear area for loading about the centroidal axis
(``e.a_sx``, ``e.a_sy``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a warping analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_as()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.a_sx is not None
assert self.section_props.a_sy is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.a_sx / e_ref, self.section_props.a_sy / e_ref
[docs]
def get_as_p(self) -> tuple[float, float]:
"""Returns the cross-section princicpal axis shear area.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Shear area for loading about the princicpal bending axis (``a_s11``,
``a_s22``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a warping analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using"
msg += " get_eas_p()."
raise RuntimeError(msg)
try:
assert self.section_props.a_s11 is not None
assert self.section_props.a_s22 is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.a_s11, self.section_props.a_s22
[docs]
def get_eas_p(
self,
e_ref: float | pre.Material = 1,
) -> tuple[float, float]:
"""Returns the modulus-weighted cross-section princicpal axis shear area.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
This value can be transformed by providing a reference elastic modulus or a
material property.
Args:
e_ref: Reference elastic modulus or material property by which to transform
the results
Returns:
Modulus weighted shear area for loading about the princicpal bending axis
(``e.a_s11``, ``e.a_s22``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a warping analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_ig()."
raise RuntimeError(msg)
# calculate e_ref (transformed elastic modulus)
e_ref = self.get_e_ref(e_ref=e_ref)
try:
assert self.section_props.a_s11 is not None
assert self.section_props.a_s22 is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return self.section_props.a_s11 / e_ref, self.section_props.a_s22 / e_ref
[docs]
def get_beta(self) -> tuple[float, float, float, float]:
"""Returns the cross-section global monosymmetry constants.
Returns:
Monosymmetry constants for bending about both global axes (``beta_x_plus``,
``beta_x_minus``, ``beta_y_plus``, ``beta_y_minus``).
Note:
The ``plus`` value relates to the top flange in compression and the
``minus`` value relates to the bottom flange in compression.
Raises:
AssertionError: If a warping analysis has not been performed
"""
try:
assert self.section_props.beta_x_plus is not None
assert self.section_props.beta_x_minus is not None
assert self.section_props.beta_y_plus is not None
assert self.section_props.beta_y_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return (
self.section_props.beta_x_plus,
self.section_props.beta_x_minus,
self.section_props.beta_y_plus,
self.section_props.beta_y_minus,
)
[docs]
def get_beta_p(self) -> tuple[float, float, float, float]:
"""Returns the cross-section principal monosymmetry constants.
Returns:
Monosymmetry constant for bending about both principal axes
(``beta_11_plus``, ``beta_11_minus``, ``beta_22_plus``, ``beta_22_minus``)
Note:
The ``plus`` value relates to the top flange in compression and the
``minus`` value relates to the bottom flange in compression.
Raises:
AssertionError: If a warping analysis has not been performed
"""
try:
assert self.section_props.beta_11_plus is not None
assert self.section_props.beta_11_minus is not None
assert self.section_props.beta_22_plus is not None
assert self.section_props.beta_22_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a warping analysis.") from e
return (
self.section_props.beta_11_plus,
self.section_props.beta_11_minus,
self.section_props.beta_22_plus,
self.section_props.beta_22_minus,
)
[docs]
def get_pc(self) -> tuple[float, float]:
"""Returns the cross-section centroidal axis plastic centroid.
Returns:
Centroidal axis plastic centroid (``x_pc``, ``y_pc``)
Raises:
AssertionError: If a plastic analysis has not been performed
"""
try:
assert self.section_props.x_pc is not None
assert self.section_props.y_pc is not None
assert self.section_props.cx is not None
assert self.section_props.cy is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
# add centroid location to move section back to original location
x_pc = self.section_props.x_pc + self.section_props.cx
y_pc = self.section_props.y_pc + self.section_props.cy
return x_pc, y_pc
[docs]
def get_pc_p(self) -> tuple[float, float]:
"""Returns the cross-section principal axis plastic centroid.
Returns:
Principal bending axis plastic centroid (``x11_pc``, ``y22_pc``)
Raises:
AssertionError: If a plastic analysis has not been performed
"""
try:
assert self.section_props.phi is not None
assert self.section_props.x11_pc is not None
assert self.section_props.y22_pc is not None
assert self.section_props.cx is not None
assert self.section_props.cy is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
# determine the position of the plastic centroid in the global axis
x_pc, y_pc = fea.global_coordinate(
phi=self.section_props.phi,
x11=self.section_props.x11_pc,
y22=self.section_props.y22_pc,
)
# add centroid location to move section back to original location
return x_pc + self.section_props.cx, y_pc + self.section_props.cy
[docs]
def get_s(self) -> tuple[float, float]:
"""Returns the cross-section centroidal plastic section moduli.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Plastic section moduli about the centroidal axis (``sxx``, ``syy``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a plastic analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using get_mp()."
raise RuntimeError(msg)
try:
assert self.section_props.sxx is not None
assert self.section_props.syy is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
return self.section_props.sxx, self.section_props.syy
[docs]
def get_mp(self) -> tuple[float, float]:
"""Returns the cross-section plastic moment about the centroidal axis.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
Returns:
Plastic moment (:math:`M_p = f_y S`) about the centroidal axis (``mp_xx``,
``mp_yy``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a plastic analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_s()."
raise RuntimeError(msg)
try:
assert self.section_props.sxx is not None
assert self.section_props.syy is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
return self.section_props.sxx, self.section_props.syy
[docs]
def get_sp(self) -> tuple[float, float]:
"""Returns the cross-section principal axis plastic section moduli.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Plastic section moduli about the principal bending axis (``s11``, ``s22``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a plastic analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied). Consider using"
msg += " get_mp_p()."
raise RuntimeError(msg)
try:
assert self.section_props.s11 is not None
assert self.section_props.s22 is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
return self.section_props.s11, self.section_props.s22
[docs]
def get_mp_p(self) -> tuple[float, float]:
"""Returns the cross-section plastic moment about the principal axis.
This is a composite only property, as such this can only be returned if material
properties have been applied to the cross-section.
Returns:
Plastic moment (:math:`M_p = f_y S`) about the principal axis (``mp_11``,
``mp_22``)
Raises:
RuntimeError: If material properties have *not* been applied
AssertionError: If a plastic analysis has not been performed
"""
if not self.is_composite():
msg = "Attempting to get a composite only property for a geometric analysis"
msg += " (material properties have not been applied). Consider using"
msg += " get_sp()."
raise RuntimeError(msg)
try:
assert self.section_props.s11 is not None
assert self.section_props.s22 is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
return self.section_props.s11, self.section_props.s22
[docs]
def get_sf(self) -> tuple[float, float, float, float]:
"""Returns the cross-section centroidal axis shape factors.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Centroidal axis shape factors with respect to the top and bottom fibres
(``sf_xx_plus``, ``sf_xx_minus``, ``sf_yy_plus``, ``sf_yy_minus``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a plastic analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied)."
raise RuntimeError(msg)
try:
assert self.section_props.sf_xx_plus is not None
assert self.section_props.sf_xx_minus is not None
assert self.section_props.sf_yy_plus is not None
assert self.section_props.sf_yy_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
return (
self.section_props.sf_xx_plus,
self.section_props.sf_xx_minus,
self.section_props.sf_yy_plus,
self.section_props.sf_yy_minus,
)
[docs]
def get_sf_p(self) -> tuple[float, float, float, float]:
"""Returns the cross-section principal axis shape factors.
This is a geometric only property, as such this can only be returned if material
properties have *not* been applied to the cross-section.
Returns:
Principal bending axis shape factors with respect to the top and bottom
fibres (``sf_11_plus``, ``sf_11_minus``, ``sf_22_plus``, ``sf_22_minus``)
Raises:
RuntimeError: If material properties have been applied
AssertionError: If a plastic analysis has not been performed
"""
if self.is_composite():
msg = "Attempting to get a geometric only property for a composite analysis"
msg += " (material properties have been applied)."
raise RuntimeError(msg)
try:
assert self.section_props.sf_11_plus is not None
assert self.section_props.sf_11_minus is not None
assert self.section_props.sf_22_plus is not None
assert self.section_props.sf_22_minus is not None
except AssertionError as e:
raise AssertionError("Conduct a plastic analysis.") from e
return (
self.section_props.sf_11_plus,
self.section_props.sf_11_minus,
self.section_props.sf_22_plus,
self.section_props.sf_22_minus,
)
[docs]
def get_stress_at_points(
self,
pts: list[tuple[float, float]],
n: float = 0.0,
mxx: float = 0.0,
myy: float = 0.0,
m11: float = 0.0,
m22: float = 0.0,
mzz: float = 0.0,
vx: float = 0.0,
vy: float = 0.0,
agg_func: Callable[[list[float]], float] = np.average,
) -> list[tuple[float, float, float] | None]:
"""Calculates the stress for a list of points.
Calculates the stress at a set of points within an element for given design
actions and returns the global stress components for each point.
Args:
pts: A list of points ``[(x, y), ..., ]``
n: Axial force
mxx: Bending moment about the centroidal xx-axis
myy: Bending moment about the centroidal yy-axis
m11: Bending moment about the centroidal 11-axis
m22: Bending moment about the centroidal 22-axis
mzz: Torsion moment about the centroidal zz-axis
vx: Shear force acting in the x-direction
vy: Shear force acting in the y-direction
agg_func: A function that aggregates the stresses if the point is shared by
several elements. If the point, ``pt``, is shared by several elements
(e.g. if it is a node or on an edge), the stresses are retrieved from
each element and combined according to this function.
Raises:
RuntimeError: If a warping analysis has not been carried out and a shear
force or a torsion moment is supplied
Returns:
Resultant normal and shear stresses (``sigma_zz``, ``tau_xz``, ``tau_yz``)
for each point ``pt``. If a point is not in the section then ``None`` is
returned for that element in the list.
"""
# ensure warping analysis completed for shear and torsion
if vx != 0 or vy != 0 or mzz != 0:
warping = True
if self.section_props.omega is None:
msg = "Perform a warping analysis before carrying out a stress "
msg += "analysis with non-zero shear forces or torsion moment."
raise RuntimeError(msg)
else:
warping = False
action = {
"n": n,
"mxx": mxx,
"myy": myy,
"m11": m11,
"m22": m22,
"mzz": mzz,
"vx": vx,
"vy": vy,
}
# get relevant section properties
cx, cy = self.get_c()
if self.is_composite():
ixx, iyy, ixy = self.get_eic()
i11, i22 = self.get_eip()
else:
ixx, iyy, ixy = self.get_ic()
i11, i22 = self.get_ip()
if warping:
if self.is_composite():
j = self.get_ej()
nu = self.get_nu_eff()
else:
j = self.get_j()
nu = 0.0
delta_s = self.section_props.delta_s
omega = self.section_props.omega
psi_shear = self.section_props.psi_shear
phi_shear = self.section_props.phi_shear
assert delta_s is not None
assert omega is not None
assert psi_shear is not None
assert phi_shear is not None
else:
j = 0.0
nu = 0.0
delta_s = 0.0
omega = np.zeros(shape=self.num_nodes)
psi_shear = np.zeros(shape=self.num_nodes)
phi_shear = np.zeros(shape=self.num_nodes)
sect_prop = {
"ea": self.get_ea() if self.is_composite() else self.get_area(),
"cx": cx,
"cy": cy,
"ixx": ixx,
"iyy": iyy,
"ixy": ixy,
"i11": i11,
"i22": i22,
"phi": self.get_phi(),
"j": j,
"delta_s": delta_s,
"nu": nu,
}
stress_pts = []
for pt in pts:
query_geom = Point(pt)
tri_ids: list[int] = self.mesh_search_tree.query(
query_geom, predicate="intersects"
)
if len(tri_ids) == 0:
sig = None
elif len(tri_ids) == 1:
tri = self.elements[tri_ids[0]]
omega_el = omega[tri.node_ids]
psi_shear_el = psi_shear[tri.node_ids]
phi_shear_el = phi_shear[tri.node_ids]
sig = tri.local_element_stress(
p=pt,
**action,
**sect_prop,
omega=omega_el,
psi_shear=psi_shear_el,
phi_shear=phi_shear_el,
)
else:
sigs = []
for idx in tri_ids:
tri = self.elements[idx]
omega_el = omega[tri.node_ids]
psi_shear_el = psi_shear[tri.node_ids]
phi_shear_el = phi_shear[tri.node_ids]
sigs.append(
tri.local_element_stress(
p=pt,
**action,
**sect_prop,
omega=omega_el,
psi_shear=psi_shear_el,
phi_shear=phi_shear_el,
)
)
sig = (
float(agg_func([sig[0] for sig in sigs])),
float(agg_func([sig[1] for sig in sigs])),
float(agg_func([sig[2] for sig in sigs])),
)
stress_pts.append(sig)
return stress_pts
[docs]
@dataclass
class MaterialGroup:
"""Class for storing elements of different materials.
A MaterialGroup object contains the finite element objects for a specified
``material``. The ``stress_result`` variable provides storage for stresses related
each material.
Args:
num_nodes (int): Number of nodes in the finite element mesh
material (Material): Material object for the current MaterialGroup
stress_result (StressResult): A
:class:`~sectionproperties.post.stress_post.StressResult` object for saving
the stresses of the current material
elements (list[Tri6]): A list of finite element objects that are of the current
material type
el_ids (list[int]): A list of the element IDs of the elements that are of the
current material type
"""
num_nodes: int
material: pre.Material
stress_result: sp_stress_post.StressResult = field(init=False)
elements: list[fea.Tri6] = field(init=False, default_factory=list)
el_ids: list[int] = field(init=False, default_factory=list)
def __post_init__(self) -> None:
"""Allocates the stress_result parameter."""
self.stress_result = sp_stress_post.StressResult(self.num_nodes)
[docs]
def add_element(
self,
element: fea.Tri6,
) -> None:
"""Adds an element and its element ID to the MaterialGroup.
Args:
element: Element to add to the MaterialGroup
"""
# add Tri6 element to the list of elements
self.elements.append(element)
self.el_ids.append(element.el_id)