From 19f1bb9a14d22377e6335f2d281032c95087bd4d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Maur=C3=ADcio=20Bonatte?= Date: Sat, 26 Aug 2023 20:41:12 +0100 Subject: [PATCH] making random improvements --- secan/geometry.py | 167 +++++++++++++++++++---- secan/material.py | 331 +++++++++++++++++++++++++++++++++++++-------- secan/section.py | 5 +- secan/test/test.py | 6 +- 4 files changed, 423 insertions(+), 86 deletions(-) diff --git a/secan/geometry.py b/secan/geometry.py index c22ec99..fec8076 100644 --- a/secan/geometry.py +++ b/secan/geometry.py @@ -6,37 +6,120 @@ from abc import ABC, abstractmethod +from abc import ABC, abstractmethod +from typing import List, Tuple + class Geometry(ABC): + """ + Abstract base class representing a geometry. + + Attributes: + None + + Methods: + area -> float: Abstract property to get the area of the geometry. + boundary -> List[Tuple[float, float]]: Abstract property to get the boundary coordinates of the geometry. + get_normal_resistance(e0: float, k: float, center: float) -> np.ndarray: Abstract method to calculate normal resistance. + get_moment_resistance(e0: float, k: float, center: float) -> np.ndarray: Abstract method to calculate moment resistance. + get_stiffness(e0: float, k: float, center: float) -> np.ndarray: Abstract method to calculate stiffness. + set_material(material): Set the material for the geometry. + plot_geometry(graph): Method to plot the geometry (default implementation prints a message). + set_x_plot(graph, value): Method to set x-axis limits of a plot (adjusting with value). + """ @property @abstractmethod def area(self) -> float: + """ + Get the area of the geometry. + + Returns: + float: Area of the geometry. + """ pass @property @abstractmethod - def boundary(self) -> List[tuple[float, float]]: + def boundary(self) -> List[Tuple[float, float]]: + """ + Get the boundary coordinates of the geometry. + + Returns: + List[Tuple[float, float]]: List of boundary coordinates as tuples. + """ pass @abstractmethod def get_normal_resistance(self, e0: float, k: float, center: float) -> np.ndarray: + """ + Calculate normal resistance based on the given parameters. + + Args: + e0 (float): Normal strain value. + k (float): Curvature value. + center (float): Center value. + + Returns: + np.ndarray: Array of normal resistance values. + """ pass @abstractmethod def get_moment_resistance(self, e0: float, k: float, center: float) -> np.ndarray: + """ + Calculate moment resistance based on the given parameters. + + Args: + e0 (float): Normal strain value. + k (float): Curvature value. + center (float): Center value. + + Returns: + np.ndarray: Array of moment resistance values. + """ pass @abstractmethod def get_stiffness(self, e0: float, k: float, center: float) -> np.ndarray: + """ + Calculate stiffness based on the given parameters. + + Args: + e0 (float): Normal strain value. + k (float): Curvature value. + center (float): Center value. + + Returns: + np.ndarray: Array of stiffness values. + """ pass - def set_material(self, material): + def set_material(self, material) -> None: + """ + Set the material for the geometry. + + Args: + material: Material instance to be set for the geometry. + """ self.material = material - def plot_geometry(self, graph): + def plot_geometry(self, graph: object) -> None: + """ + Plot the geometry (default implementation prints a message). + + Args: + graph (object): Plotting object (e.g., Matplotlib axis). + """ print("Plot is not implemented for this section") - def set_x_plot(self, graph, value): + def set_x_plot(self, graph: object, value: float) -> None: + """ + Set x-axis limits of a plot (adjusting with value). + + Args: + graph (object): Plotting object (e.g., Matplotlib axis). + value (float): Value used to adjust x-axis limits. + """ xabs_max = abs(max(graph.get_xlim(), key=abs)) value *= 1.05 if value > xabs_max: @@ -45,10 +128,26 @@ def set_x_plot(self, graph, value): graph.set_xlim(xmin=-xabs_max, xmax=xabs_max) + class RectSection(Geometry): + """ + Rectangle section geometry class. + """ + def __init__(self, width: float, height: float, material, - center: tuple[float, float]=(0, 0), + center: tuple[float, float] = (0, 0), rotation=0, n_discret=200): + """ + Initialize the RectSection geometry. + + Args: + width (float): Width of the rectangle. + height (float): Height of the rectangle. + material: Material instance for the geometry. + center (tuple[float, float]): Center coordinates of the rectangle. + rotation: Rotation of the rectangle (not used). + n_discret (int): Number of discretization points. + """ self._width = width self._height = height self.center = np.array(center) @@ -93,50 +192,56 @@ def boundary(self) -> List[tuple[float, float]]: def get_area(self): return self.width*self.height - def get_strain(self, e0=0, k=0, pos=0): - return e0+k*(self.self.height/2-pos) - def get_strains(self, e0: float, k: float) -> np.ndarray: + """ + Calculate strains based on normal strain and curvature. + + Args: + e0 (float): Normal strain value. + k (float): Curvature value. + + Returns: + np.ndarray: Array of calculated strains. + """ return e0 + k * (self.height/2 - self.h_discret) - - def get_stress(self, e0: float, k: float, center: float) -> np.ndarray: - strains = self.get_strain(e0, k) - return self.material.get_stress(strains) def get_e0_sec(self, e0, k, center): + """ + Calculate normal strain of the section based on the normal strain, curvature and center of the main section. + + Args: + e0: Normal strain value. + k: Curvature value. + center: Center of the main section. + + Returns: + float: Calculated e0_sec. + """ return e0 + k * (center - self.center[1]) def get_normal_resistance_discrete(self, e0, k, center): e0_sec = self.get_e0_sec(e0, k, center) strains = self.get_strains(e0_sec, k) - normal = map(self.material.get_stress,strains) - normal = np.fromiter(normal, dtype=float) + normal = self.material.get_stress(strains) return normal * self.area_discret def get_normal_resistance(self, e0, k, center): - return sum(self.get_normal_resistance_discrete(e0, k, center)) + return np.sum(self.get_normal_resistance_discrete(e0, k, center)) def get_moment_resistance(self, e0, k, center): normal = self.get_normal_resistance_discrete(e0, k, center) - #dist = (center-self.center[1])+(self.center[1]-self.h_discret) - dist = (center-self.center[1])+(self._height/2-self.h_discret) - #print(self._height/2, self.center[1]) - return sum(normal * dist) + dist = (center - self.center[1]) + (self._height / 2 - self.h_discret) + return np.sum(normal * dist) def get_normal_stiff_discrete(self, e0, k, center): e0_sec = self.get_e0_sec(e0, k, center) strains = self.get_strains(e0_sec, k) - - normal = np.array([ - self.material.get_stiff(strain) for strain in strains - ]) + normal = self.material.get_stiff(strains) return normal * self.area_discret def get_stiffness(self, e0: float, k: float, - center: float) -> np.ndarray: - + center: float) -> np.ndarray: normal = self.get_normal_stiff_discrete(e0, k, center) - #dist = (center-self.center[1]) + (self.center[1]-self.h_discret) dist = (center-self.center[1])+(self._height/2-self.h_discret) a00 = normal.sum() @@ -175,6 +280,12 @@ def plot_strain(self, graph, e0, k, center): self.set_x_plot(graph, abs(max(strain, key=abs))) def plot_geometry(self, graph=None): + """ + Plot the rectangle geometry. + + Args: + graph: Matplotlib axis to plot the geometry. + """ if graph is None: fig, graph = plt.subplots(1, figsize=(10, 10)) graph.add_patch(Rectangle(self.boundary[0], @@ -186,7 +297,7 @@ def plot_geometry(self, graph=None): class Rebar(Geometry): def __init__(self, diameter, material, center=(0, 0)): - self._diameter = diameter + self.diameter = diameter self._area = pi * diameter**2 / 4 self.material = material self.center = np.array(center) @@ -211,7 +322,7 @@ def area(self): def area(self, new_area): if (new_area >= 0): self._area = new_area - self._diameter = (new_area*4/3.141592)**0.5 + self._diameter = (new_area*4/pi)**0.5 else: raise Exception("Area must be higher than 0") diff --git a/secan/material.py b/secan/material.py index 2459bb8..9035ead 100644 --- a/secan/material.py +++ b/secan/material.py @@ -4,37 +4,121 @@ class Material(ABC): + """ + Abstract base class representing a material. + + Attributes: + None + + Methods: + get_stiff(strain: float) -> float: Abstract method to calculate material stiffness. + get_stress(strain: float) -> float: Abstract method to calculate material stress. + plot(plot: object) -> None: Method to plot the material (default implementation prints a message). + """ @abstractmethod def get_stiff(self, strain: float) -> float: + """ + Calculate material stiffness based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stiffness. + """ pass @abstractmethod def get_stress(self, strain: float) -> float: + """ + Calculate material stress based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stress. + """ pass def plot(self, plot=None): + """ + Plot the material (default implementation prints a message). + + Args: + plot (object): Plotting object (optional). + """ print("Plot is not implemented for this material") class Linear(Material): + """ + Linear material class representing a material with constant Young's modulus. + + Attributes: + young (float): Young's modulus of the material. + + Methods: + get_stiff(strain: float) -> float: Calculate material stiffness based on the given strain. + get_stress(strain: float) -> float: Calculate material stress based on the given strain. + """ def __init__(self, young=0): + """ + Initialize the Linear material with a given Young's modulus. + + Args: + young (float): Young's modulus of the material. + """ self.young = young def get_stiff(self, strain=0): + """ + Calculate material stiffness based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stiffness. + """ return self.young def get_stress(self, strain=0): + """ + Calculate material stress based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stress. + """ return self.young*strain class Concrete(Material): # According to NBR6118 and EN1992 def __init__(self, fc=0): + """ + Initialize the Concrete material with a given fc value. + + Args: + fc (float): Compressive strength of concrete. + """ self.fc = fc @staticmethod def _get_constants(fc: float) -> tuple[float, float, float]: + """ + Compute material constants based on the given fc value. + + Args: + fc (float): Compressive strength of concrete. + + Returns: + tuple[float, float, float]: Tuple containing ec2, ecu, and n constants. + """ if 0e6 <= fc < 55e6: return -2 / 1e3, -3.5 / 1e3, 2 elif 55e6 <= fc <= 90e6: @@ -54,28 +138,68 @@ def fc(self, new_fc): self.ec2, self.ecu, self.n = self._get_constants(new_fc) def get_stiff(self, strain: float=0) -> float: - if (self.ec2 < strain <= 0): - return -1*(self.fc*self.n * (1-strain/self.ec2)**(self.n-1) / self.ec2) - elif (strain <= self.ec2): - return 0 - else: - return 0 + """ + Calculate material stiffness based on the given strain. + + Args: + strain (float): Strain value. + Returns: + float: Material stiffness. + """ + condition1 = (self.ec2 <= strain) & (strain <= 0) + condition2 = (strain <= self.ec2) + + result1 = -1*(self.fc*self.n * (1-strain/self.ec2)**(self.n-1) / self.ec2) + result2 = 0 + + stiff = np.select( + [condition1, condition2], + [result1, result2], + default=0 + ) + + return stiff + def get_stress(self, strain: float=0) -> float: - if (self.ec2 <= strain <= 0): - return -1*self.fc*(1-(1-strain/self.ec2)**self.n) - elif (self.ecu <= strain <= self.ec2): - return -self.fc - elif (strain <= self.ecu): - return 0 - else: - return 0 + """ + Calculate material stress based on the given strain. - def plot(self, graph=None): + Args: + strain (float): Strain value. + + Returns: + float: Material stress. + """ + condition1 = (self.ec2 <= strain) & (strain <= 0) + condition2 = (self.ecu <= strain) & (strain <= self.ec2) + condition3 = strain <= self.ecu + + result1 = -1 * self.fc * (1 - (1 - strain / self.ec2) ** self.n) + result2 = -self.fc + result3 = 0 + + stress = np.select( + [condition1, condition2, condition3], + [result1, result2, result3], + default=0 + ) + + return stress + + def plot(self, graph=None) -> None: + """ + Plot the stress-strain diagram for the Concrete material. + + Args: + graph: Matplotlib axis to plot the diagram. + """ if graph is None: fig, graph = plt.subplots(1, figsize=(10, 10)) + strain = np.arange(1.05*self.ecu, 1/5000, 1/50000) - stress = [self.get_stress(strain[i]) for i in range(len(strain))] + stress = self.get_stress(strain) + graph.set(xlabel='Strain') graph.set(ylabel='Stress') graph.set_title("Concrete Driagram") @@ -84,10 +208,31 @@ def plot(self, graph=None): class SteelIdeal(Material): + """ + SteelIdeal material class representing an idealized steel material. + + Attributes: + young (float): Young's modulus of the material. + fy (float): Yield strength of the material. + ultimate_strain (float): Ultimate strain of the material. + + Methods: + get_stiff(strain: float) -> float: Calculate material stiffness based on the given strain. + get_stress(strain: float) -> float: Calculate material stress based on the given strain. + plot(graph: object) -> None: Plot the stress-strain diagram for the SteelIdeal material. + """ def __init__(self, - young=0, - fy=0, - ultimate_strain=10e-3): + young: float = 0, + fy: float = 0, + ultimate_strain: float = 10e-3): + """ + Initialize the SteelIdeal material with given parameters. + + Args: + young (float): Young's modulus of the material. + fy (float): Yield strength of the material. + ultimate_strain (float): Ultimate strain of the material. + """ self.young = young self.fy = fy self.ultimate_strain = ultimate_strain @@ -99,82 +244,162 @@ def fy(self): @fy.setter def fy(self, new_fy): self._fy = new_fy - self.yeild_strain = new_fy / self.young + self.yield_strain = new_fy / self.young - def get_stiff(self, strain=0): - if (-self.yeild_strain <= strain <= self.yeild_strain): + def get_stiff(self, strain: float = 0) -> float: + """ + Calculate material stiffness based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stiffness. + """ + if (-self.yield_strain <= strain <= self.yield_strain): return self.young else: return 0 - def get_stress(self, strain=0): - if (-self.yeild_strain <= strain <= self.yeild_strain): + def get_stress(self, strain: float = 0) -> float: + """ + Calculate material stress based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stress. + """ + if (-self.yield_strain <= strain <= self.yield_strain): return self.young*strain elif (-self.ultimate_strain <= strain <= self.ultimate_strain): return self.fy * strain/abs(strain) else: return 0 - def plot(self, graph=None): + def plot(self, graph=None) -> None: + """ + Plot the stress-strain diagram for the SteelIdeal material. + + Args: + graph: Matplotlib axis to plot the diagram. + """ if graph is None: fig, graph = plt.subplots(1, figsize=(10, 10)) strain = np.arange(-self.ultimate_strain, self.ultimate_strain, self.ultimate_strain/100) stress = [self.get_stress(strain[i]) for i in range(len(strain))] - graph.set(xlabel='Strain') - graph.set(ylabel='Stress') - graph.set_title("Steel Driagram") + graph.set(xlabel='Strain', ylabel='Stress', title="Steel Diagram") graph.grid() graph.plot(strain, stress) class SteelHardening(SteelIdeal): + """ + SteelHardening material class representing an idealized steel material with hardening behavior. + + Attributes: + young (float): Young's modulus of the material. + fy (float): Yield strength of the material. + ft (float): Ultimate tensile strength of the material. + ultimate_strain (float): Ultimate strain of the material. + + Methods: + get_stiff(strain: float) -> float: Calculate material stiffness based on the given strain. + get_stress(strain: float) -> float: Calculate material stress based on the given strain. + plot(graph: object) -> None: Plot the stress-strain diagram for the SteelHardening material. + """ + def __init__(self, - young=0, - fy=0, - ft=0, - ultimate_strain=35e-3): + young: float = 0, + fy: float = 0, + ft: float = 0, + ultimate_strain: float = 35e-3): + """ + Initialize the SteelHardening material with given parameters. + + Args: + young (float): Young's modulus of the material. + fy (float): Yield strength of the material. + ft (float): Ultimate tensile strength of the material. + ultimate_strain (float): Ultimate strain of the material. + """ self.young = young self.fy = fy self.ultimate_strain = ultimate_strain self.ft = ft - + @property - def ft(self): - return self._ft - + def ft(self) -> float: + """ + Get the ultimate tensile strength property. + + Returns: + float: Ultimate tensile strength of the material. + """ + return self._ft + @ft.setter - def ft(self, new_ft): + def ft(self, new_ft: float) -> None: + """ + Set the ultimate tensile strength property and update related parameters. + + Args: + new_ft (float): New ultimate tensile strength of the material. + """ self._ft = new_ft - self.hardening_stiffness = (new_ft - self.fy) / (self.ultimate_strain - self.yeild_strain) + self.hardening_stiffness = (new_ft - self.fy) / (self.ultimate_strain - self.yield_strain) - def get_stiff(self, strain=0): - if (-self.yeild_strain <= strain <= self.yeild_strain): + def get_stiff(self, strain: float = 0) -> float: + """ + Calculate material stiffness based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stiffness. + """ + if (-self.yield_strain <= strain <= self.yield_strain): return self.young - elif self.yeild_strain < strain < self.ultimate_strain: + elif self.yield_strain < strain < self.ultimate_strain: return self.hardening_stiffness else: return 0 - def get_stress(self, strain=0): - if (-self.yeild_strain <= strain <= self.yeild_strain): - return self.young*strain - elif self.yeild_strain < strain < self.ultimate_strain: - return self.fy + (strain - self.yeild_strain) * self.hardening_stiffness - elif (-self.ultimate_strain <= strain <= -self.yeild_strain): - return self.fy * strain/abs(strain) + def get_stress(self, strain: float = 0) -> float: + """ + Calculate material stress based on the given strain. + + Args: + strain (float): Strain value. + + Returns: + float: Material stress. + """ + if (-self.yield_strain <= strain <= self.yield_strain): + return self.young * strain + elif self.yield_strain < strain < self.ultimate_strain: + return self.fy + (strain - self.yield_strain) * self.hardening_stiffness + elif (-self.ultimate_strain <= strain <= -self.yield_strain): + return self.fy * strain / abs(strain) else: return 0 - def plot(self, graph=None): + def plot(self, graph=None) -> None: + """ + Plot the stress-strain diagram for the SteelHardening material. + + Args: + graph: Matplotlib axis to plot the diagram. + """ if graph is None: fig, graph = plt.subplots(1, figsize=(10, 10)) strain = np.arange(-self.ultimate_strain, self.ultimate_strain, - self.ultimate_strain/100) + self.ultimate_strain / 100) stress = [self.get_stress(strain[i]) for i in range(len(strain))] - graph.set(xlabel='Strain') - graph.set(ylabel='Stress') - graph.set_title("Steel Driagram") + graph.set(xlabel='Strain', ylabel='Stress', title="Steel Diagram") graph.grid() graph.plot(strain, stress) diff --git a/secan/section.py b/secan/section.py index dfe232e..eff579b 100644 --- a/secan/section.py +++ b/secan/section.py @@ -58,7 +58,6 @@ def get_normal_res(self, e0: float, k: float) -> np.ndarray: def get_stiff(self, e0: float, k: float) -> np.ndarray: K = np.array([s.get_stiffness(e0, k, self.centroid[1]) for s in self.section]) return K.sum(axis=0) - def get_e0(self, k=0, e0=0, target_normal=0): for _ in range(self.n_ite_e0): @@ -88,17 +87,21 @@ def check_section(self, target_normal, target_moment, n_ite=10): residual = np.array([[target_normal-normal], [target_moment-moment]]) + if np.linalg.norm(residual) < self.tolerance_check_section: return e0, k stiff = self.get_stiff(e0, k) + if np.linalg.det(stiff) < self.tolerance_check_section: break # update section state inv_stiff = np.linalg.inv(stiff) + e0 += np.matmul(inv_stiff, residual )[0][0] k += np.matmul(inv_stiff, residual )[1][0] + return None, None def get_strain_base_top(self, f, inverted=False): diff --git a/secan/test/test.py b/secan/test/test.py index e4906aa..d4f2a1b 100644 --- a/secan/test/test.py +++ b/secan/test/test.py @@ -26,7 +26,6 @@ def test_steel_hardening(self): self.assertEqual(s.get_stress(25e-3), 957575757.5757575) self.assertEqual(s.get_stress(-10e-3), -400e6) - class TestGeometries(unittest.TestCase): def test_rectangle(self): @@ -46,7 +45,7 @@ def test_rectangle(self): self.assertEqual(rect.get_normal_resistance(e0=0.001, k=0.02, center=0), -1133313.0) self.assertEqual(rect.get_moment_resistance(e0=0.001, k=0.02, center=0.3), 172496.04929999996) - self.assertEqual(rect.get_moment_resistance(e0=0.001, k=0.02, center=0), 172496.04929999998) + self.assertEqual(rect.get_moment_resistance(e0=0.001, k=0.02, center=0), 172496.04929999996) def test_rebar(self): steel = SteelIdeal(200e9, 400e6) @@ -238,7 +237,7 @@ def test_example_3_4(self): self.assertAlmostEqual(self.conc.get_stress(-3.5e-3)/1e6, -18.21428, delta=self.delta) self.assertAlmostEqual(self.steel.get_stress(10e-3)/1e6, 434.7826, delta=self.delta) - self.assertAlmostEqual(self.steel.yeild_strain, 2.0704e-3) + self.assertAlmostEqual(self.steel.yield_strain, 2.0704e-3) self.assertAlmostEqual(self.s.area_rebar, 4.7123e-4, delta=self.delta) @@ -351,6 +350,5 @@ def test_example_4_5(self): self.assertAlmostEqual(rec_stiff[0][1]/1e9, -0.0102, delta=0.0004) ###### self.assertAlmostEqual(rec_stiff[1][1]/1e8, 0.0154, delta=0.0004) ###### - if __name__ == '__main__': unittest.main()