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()