diff --git a/docs/source/rst/analysis.rst b/docs/source/rst/analysis.rst index c9c1314a..53071376 100644 --- a/docs/source/rst/analysis.rst +++ b/docs/source/rst/analysis.rst @@ -71,8 +71,9 @@ A stress analysis calculates the section stresses arising from a set of forces and moments. Executing this method returns a :class:`~sectionproperties.analysis.section.StressResult` object which stores the section stresses and provides stress plotting functions. -.. note:: A geometric and warping analysis must be performed on the Section - object before a stress analysis is carried out. +.. note:: A geometric analysis must be performed on the Section object before a stress + analysis is carried out. Further, if the shear force or twisting moment is non-zero + a warping analysis must also be performed. .. warning:: The stress analysis in *sectionproperties* is linear-elastic and does not account for the yielding of materials or any non-linearities. diff --git a/sectionproperties/analysis/section.py b/sectionproperties/analysis/section.py index e8643cb7..5a19615e 100644 --- a/sectionproperties/analysis/section.py +++ b/sectionproperties/analysis/section.py @@ -1043,8 +1043,8 @@ def calc_plastic(progress=None): calc_plastic() def calculate_stress(self, N=0, Vx=0, Vy=0, Mxx=0, Myy=0, M11=0, M22=0, Mzz=0): - """Calculates the cross-section stress resulting from design actions and returns a - :class:`~sectionproperties.analysis.section.StressPost` object allowing + """Calculates the cross-section stress resulting from design actions and returns + a :class:`~sectionproperties.analysis.section.StressPost` object allowing post-processing of the stress results. :param float N: Axial force @@ -1055,33 +1055,35 @@ def calculate_stress(self, N=0, Vx=0, Vy=0, Mxx=0, Myy=0, M11=0, M22=0, Mzz=0): :param float M11: Bending moment about the centroidal 11-axis :param float M22: Bending moment about the centroidal 22-axis :param float Mzz: Torsion moment about the centroidal zz-axis + :return: Object for post-processing cross-section stresses :rtype: :class:`~sectionproperties.analysis.section.StressPost` - Note that a geometric and warping analysis must be performed before a stress analysis is - carried out:: + Note that 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:: section = Section(geometry) section.calculate_geometric_properties() section.calculate_warping_properties() stress_post = section.calculate_stress(N=1e3, Vy=3e3, Mxx=1e6) - :raises RuntimeError: If a geometric and warping analysis have not been performed prior to - calling this method + :raises RuntimeError: If a geometric and warping analysis (if required) have not + been performed prior to calling this method """ # check that a geometric and warping analysis has been performed - if ( - None - in [ - self.section_props.area, - self.section_props.ixx_c, - self.section_props.cx, - self.section_props.j, - ] - and self.section_props.omega is None - ): - err = "Perform a geometric and warping analysis before carrying out a stress analysis." + if None in [ + self.section_props.area, + self.section_props.ixx_c, + self.section_props.cx, + ]: + err = "Perform a geometric analysis before carrying out a stress analysis." + raise RuntimeError(err) + + if self.section_props.omega is None and (Vx != 0 or Vy != 0 or Mzz != 0): + err = "Perform a warping analysis before carrying out a stress analysis " + err += "with non-zero shear forces or torsion moment." raise RuntimeError(err) def calc_stress(progress=None): @@ -1115,6 +1117,16 @@ def calc_stress(progress=None): # loop through all elements in the material group for el in group.elements: + # get element omega and psi + if self.section_props.omega is None: + omega_el = None + psi_shear_el = None + phi_shear_el = None + else: + omega_el = self.section_props.omega[el.node_ids] + psi_shear_el = self.section_props.psi_shear[el.node_ids] + phi_shear_el = self.section_props.phi_shear[el.node_ids] + ( sig_zz_n_el, sig_zz_mxx_el, @@ -1148,9 +1160,9 @@ def calc_stress(progress=None): phi, j, nu, - self.section_props.omega[el.node_ids], - self.section_props.psi_shear[el.node_ids], - self.section_props.phi_shear[el.node_ids], + omega_el, + psi_shear_el, + phi_shear_el, Delta_s, ) @@ -2214,6 +2226,12 @@ def get_stress_at_points( :rtype: List[Union[Tuple[float, float, float], None]] """ + # ensure warping analysis completed for shear and torsion + if self.section_props.omega is None and (Vx != 0 or Vy != 0 or Mzz != 0): + err = "Perform a warping analysis before carrying out a stress analysis " + err += "with non-zero shear forces or torsion moment." + raise RuntimeError(err) + action = { "N": N, "Mxx": Mxx, @@ -2251,26 +2269,46 @@ def get_stress_at_points( sig = None elif len(tri_ids) == 1: tri = self.elements[tri_ids[0]] + + if self.section_props.omega is None: + omega_el = None + psi_shear_el = None + phi_shear_el = None + else: + omega_el = self.section_props.omega[tri.node_ids] + psi_shear_el = self.section_props.psi_shear[tri.node_ids] + phi_shear_el = self.section_props.phi_shear[tri.node_ids] + sig = tri.local_element_stress( p=pt, **action, **sect_prop, - omega=self.section_props.omega[tri.node_ids], - psi_shear=self.section_props.psi_shear[tri.node_ids], - phi_shear=self.section_props.phi_shear[tri.node_ids], + omega=omega_el, + psi_shear=psi_shear_el, + phi_shear=phi_shear_el, ) else: sigs = [] for idx in tri_ids: tri = self.elements[idx] + + if self.section_props.omega is None: + omega_el = None + psi_shear_el = None + phi_shear_el = None + else: + omega_el = self.section_props.omega[tri.node_ids] + psi_shear_el = self.section_props.psi_shear[tri.node_ids] + phi_shear_el = self.section_props.phi_shear[tri.node_ids] + sigs.append( tri.local_element_stress( p=pt, **action, **sect_prop, - omega=self.section_props.omega[tri.node_ids], - psi_shear=self.section_props.psi_shear[tri.node_ids], - phi_shear=self.section_props.phi_shear[tri.node_ids], + omega=omega_el, + psi_shear=psi_shear_el, + phi_shear=phi_shear_el, ) ) sig = ( diff --git a/sectionproperties/tests/test_stress.py b/sectionproperties/tests/test_stress.py new file mode 100644 index 00000000..9d593974 --- /dev/null +++ b/sectionproperties/tests/test_stress.py @@ -0,0 +1,50 @@ +import pytest +import sectionproperties.pre.library.primitive_sections as primitive_sections +from sectionproperties.analysis.section import Section + +# define geometry +rect = primitive_sections.rectangular_section(b=50, d=100) +rect.create_mesh(mesh_sizes=0) # coarse mesh +sec = Section(rect) + + +# test stress runtime errors +def test_stress_runtime_errors(): + # check runtime error with no analysis performed + with pytest.raises(RuntimeError): + sec.calculate_stress() + + sec.calculate_geometric_properties() + + # check runtime errors with shear/torsion applied, no warping analysis + with pytest.raises(RuntimeError): + sec.calculate_stress(Vx=1) + sec.get_stress_at_points(pts=[[10, 10]], Vx=1) + + with pytest.raises(RuntimeError): + sec.calculate_stress(Vy=1) + sec.get_stress_at_points(pts=[[10, 10]], Vy=1) + + with pytest.raises(RuntimeError): + sec.calculate_stress(Mzz=1) + sec.get_stress_at_points(pts=[[10, 10]], Mzz=1) + + # check no runtime errors with no shear/torsion applied + sec.calculate_stress(N=1) + sec.calculate_stress(Mxx=1) + sec.calculate_stress(Myy=1) + sec.calculate_stress(M11=1) + sec.calculate_stress(M22=1) + + sec.calculate_warping_properties() + + # check no runtime errors after warping analysis + sec.calculate_stress(N=1) + sec.calculate_stress(Mxx=1) + sec.calculate_stress(Myy=1) + sec.calculate_stress(M11=1) + sec.calculate_stress(M22=1) + sec.calculate_stress(Vx=1) + sec.calculate_stress(Vy=1) + sec.calculate_stress(Mzz=1) + sec.get_stress_at_points(pts=[[10, 10]], Mzz=1) diff --git a/setup.cfg b/setup.cfg index 202f108f..38440504 100644 --- a/setup.cfg +++ b/setup.cfg @@ -51,4 +51,5 @@ rhino = test = rhino-shapley-interop>=0.0.4 + pytest pytest_check