diff --git a/TEST.py b/TEST.py
deleted file mode 100644
index 060362c8658f6c88cc376a6b2be15cda1892eaf1..0000000000000000000000000000000000000000
--- a/TEST.py
+++ /dev/null
@@ -1,41 +0,0 @@
-
-import time
-import numpy as np
-import matplotlib.pyplot as plt
-
-
-def round_trip_connect(start, end):
- return [(i, i+1) for i in range(start, end)] + [(end, start)]
-
-def refine_grid(vertices, area):
- max_area = 1 ########################################################################
- return bool(area>max_area)
-
-def createMesh(vertices, edges) :
-
- from meshpy import triangle
- gridInfo = triangle.MeshInfo()
- gridInfo.set_points(vertices)
- gridInfo.set_facets(edges)
-
- mesh = triangle.build(gridInfo, refinement_func=refine_grid)
- mesh_points = np.array(mesh.points)
- mesh_elements = np.array(mesh.elements)
-
- return mesh_points, mesh_elements
-
-def plot(mesh_points, mesh_elements) :
- plt.triplot(mesh_points[:, 0], mesh_points[:, 1], mesh_elements)
- plt.gca().set_aspect('equal')
- for x,y in enumerate(mesh_points) :
- plt.text(y[0] - 0.01, y[1] + 0.01, x, ha='right')
- plt.show()
-
-
-vertices = [(-1,-1),(1,-1),(1,1),(-1,1)]
-
-edges = round_trip_connect(0, len(vertices)-1)
-
-mesh_points, mesh_elements = createMesh(vertices,edges)
-plot(mesh_points, mesh_elements)
-print(mesh_elements)
diff --git a/TEST_2.py b/TEST_2.py
deleted file mode 100644
index 01cbd36e9bc3394959b43e18c1478e222086edf0..0000000000000000000000000000000000000000
--- a/TEST_2.py
+++ /dev/null
@@ -1,32 +0,0 @@
-import numpy as np
-from fem_2d import Structure2DPlaneStress, Force2D, Constraint2D
-
-mesh_elements = np.array([[1, 2, 4],
- [3, 0, 4],
- [4, 2, 3],
- [0, 1, 4]], dtype='int32')
-
-mesh_points = np.array([[-1., -1.],
- [ 1., -1.],
- [ 1., 1.],
- [-1., 1.],
- [ 0., 0.]], dtype='float64')
-
-STR = Structure2DPlaneStress(mesh_elements, mesh_points,
- elastic_modulus = 10e6,
- poisson_ratio = 0.3,
- thickness = 0.1)
-
-
-STRUCT = list()
-STRUCT.append(STR)
-
-STR.node_list[0].constraint_2d = Constraint2D([1,1])
-STR.node_list[1].constraint_2d = Constraint2D([1,1])
-STR.node_list[2].force_2d = Force2D([1,1])
-STR.node_list[3].constraint_2d = Constraint2D([1,1])
-STR.node_list[4].constraint_2d = Constraint2D([0,1])
-
-STR.solve()
-STR.set_displacements()
-STR.print_details()
\ No newline at end of file
diff --git a/fem_2d.py b/fem_2d.py
index 220250c2b29422889833934175f55b9296ee1be7..2c64e8e27139fa6b80e3d27ea73c74afaf5e6e49 100644
--- a/fem_2d.py
+++ b/fem_2d.py
@@ -91,6 +91,7 @@ class Element2D:
self.node_3 = node_3
self.node_ids = node_ids
self.k_e = None
+ self.b_e = None
const_a = np.linalg.norm(self.node_1.position - self.node_2.position)
const_b = np.linalg.norm(self.node_2.position - self.node_3.position)
@@ -116,6 +117,7 @@ class Element2D:
[b_21, b_11, b_22, b_21, b_32, b_31]])
self.k_e = self.area * thickness * np.dot(np.dot(b_mat.T, c_mat), b_mat)
+ self.b_e = b_mat
return self.k_e
@@ -298,18 +300,11 @@ class Structure2DPlaneStress:
"""
print('\nStructure calculation results:')
print('\nDisplacements:')
- print('{:<10}{:>15}{:>15}{:>15}'.format('idx', 'u1', 'u2', 'u3'))
- for node_id in range(self.get_number_nodes()):
- displacement = self.get_node(node_id).get_displacement()
- print('{:<10}{:>15.6E}{:>15.6E}{:>15.6E}'.format(node_id,
- displacement[0],
- displacement[1],
- displacement[2]))
- print('\nElement forces:')
- print('{:<10}{:>15}'.format('idx', 'force'))
- for element in range(self.get_number_elements()):
- force = self.get_element(element).compute_force()
- print('{:<10}{:>15.6E}'.format(element, force))
+ print('{:<10}{:>15}{:>15}'.format('idx', 'u1', 'u2'))
+ for node_id in range(len(self.node_list)):
+ displacement = self.node_list[node_id].displacement
+ print('{:<10}{:>15.6E}{:>15.6E}'.format(node_id, displacement[0],
+ displacement[1]))
'''
class Visualizer:
@@ -463,4 +458,31 @@ class Visualizer:
else:
poly_set.color = 'red'
self.viewer.add_object(poly_set)
-'''
\ No newline at end of file
+'''
+def round_trip_connect(start, end):
+ return [(i, i+1) for i in range(start, end)] + [(end, start)]
+
+def refine_grid(vertices, area):
+ max_area = 1 ########################################################################
+ return bool(area>max_area)
+
+def createMesh(vertices, edges) :
+
+ from meshpy import triangle
+ gridInfo = triangle.MeshInfo()
+ gridInfo.set_points(vertices)
+ gridInfo.set_facets(edges)
+
+ mesh = triangle.build(gridInfo, refinement_func=refine_grid)
+ mesh_points = np.array(mesh.points)
+ mesh_elements = np.array(mesh.elements)
+
+ return mesh_points, mesh_elements
+
+def plot(mesh_points, mesh_elements):
+ import matplotlib.pyplot as plt
+ plt.triplot(mesh_points[:, 0], mesh_points[:, 1], mesh_elements)
+ plt.gca().set_aspect('equal')
+ for x,y in enumerate(mesh_points) :
+ plt.text(y[0] - 0.01, y[1] + 0.01, x, ha='right')
+ plt.show()
diff --git a/models_2d.py b/models_2d.py
new file mode 100644
index 0000000000000000000000000000000000000000..e4c1cbe1785fcd07077cf73159841a284185e317
--- /dev/null
+++ b/models_2d.py
@@ -0,0 +1,11 @@
+"""
+This is a missing docstring...
+"""
+from fem_2d import createMesh, round_trip_connect
+
+def simple_square():
+ vertices = [(-1,-1),(1,-1),(1,1),(-1,1)]
+ edges = round_trip_connect(0, len(vertices)-1)
+ return createMesh(vertices,edges)
+
+
diff --git a/tests_2d.py b/tests_2d.py
new file mode 100644
index 0000000000000000000000000000000000000000..3bf2e2725d3797e6d75775b25935ae7334b95a15
--- /dev/null
+++ b/tests_2d.py
@@ -0,0 +1,27 @@
+
+import numpy as np
+import matplotlib.pyplot as plt
+from fem_2d import Structure2DPlaneStress, Constraint2D, Force2D, plot
+from models_2d import simple_square
+
+# This test the simple square model
+SQ_mesh_points, SQ_mesh_elements = simple_square()
+
+STR = Structure2DPlaneStress(SQ_mesh_elements, SQ_mesh_points,
+ elastic_modulus = 10e6,
+ poisson_ratio = 0.3,
+ thickness = 0.1)
+
+
+STRUCT = list()
+STRUCT.append(STR)
+
+STR.node_list[0].constraint_2d = Constraint2D([1,1])
+STR.node_list[1].constraint_2d = Constraint2D([1,1])
+STR.node_list[2].force_2d = Force2D([1,1])
+STR.node_list[3].constraint_2d = Constraint2D([1,1])
+STR.node_list[4].constraint_2d = Constraint2D([1,1])
+
+STR.solve()
+STR.set_displacements()
+STR.print_results()
\ No newline at end of file