Skip to content

Documentation for mesh generation

ElasticTrussRebar

Bases: RebarInterface

This class can insert purely elastic rebar stiffnesses into the concrete matrix and the internal forces vector. Equations from http://what-when-how.com/the-finite-element-method/fem-for-trusses-finite-element-method-part-1/

Source code in reinforcement/rebar.py 
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
class ElasticTrussRebar(RebarInterface):
    """
    This class can insert purely elastic rebar stiffnesses into the concrete matrix and the internal forces vector.
    Equations from http://what-when-how.com/the-finite-element-method/fem-for-trusses-finite-element-method-part-1/

    """
    def __init__(self, concrete_mesh : dfx.mesh.Mesh, rebar_mesh : dfx.mesh.Mesh, function_space : dfx.fem.FunctionSpace, parameters: dict):
        """Initialize rebar class

        Args:
            concrete_mesh: The 3D mesh of the concrete structure.
            rebar_mesh: The line mesh of the steel rebar.
            function_space: The function space object of the concrete structure.
            parameters: A dictinary containing all needed parameters of the steel. Contains: 'A', 'E', 'rho'.

        """
        super().__init__(concrete_mesh, rebar_mesh, function_space, parameters)

    def apply_to_diagonal_mass(self, M : PETSc.Vec):
        """
        Adds nodal masses to a diagonal mass matrix.

        Args:
            M: The diagonal from of the mass matrix as a PETSc vector 

        """
        points = self.function_space.tabulate_dof_coordinates().flatten()
        diagonal_mass_1d = np.array([[1.0, 0.0], [0.0, 1.0]])
        T = np.zeros((2, 6))
        for dofs in self.dof_array.reshape(-1, 6):
            delta_x = points[dofs[3:]] - points[dofs[:3]]
            l_axial = np.linalg.norm(delta_x, 2)
            matrix_entries = delta_x / l_axial

            T[0, :3] = matrix_entries
            T[1, 3:] = matrix_entries

            mass_local = np.diag(
                T.T
                @ (
                    self.parameters["rho"]
                    * self.parameters["A"]
                    * l_axial
                    * diagonal_mass_1d
                )
                @ T
            )

            M.setValues(dofs, mass_local, addv=PETSc.InsertMode.ADD)
            M.assemble()

    def apply_to_stiffness(self, K : PETSc.Mat, u : PETSc.Vec):
        """
        Adds truss stiffness to global stiffness matrix,

        Args:
            K: The global stiffness matrix.
            u: The global displacements.
        """
        points = self.function_space.tabulate_dof_coordinates().flatten()
        K_1d = np.array([[1.0, -1.0], [-1.0, 1.0]])
        T = np.zeros((2, 6))
        for dofs in self.dof_array.reshape(-1, 6):
            delta_x = points[dofs[3:]] - points[dofs[:3]]
            l_axial = np.linalg.norm(delta_x, 2)

            matrix_entries = delta_x / l_axial
            T[0, :3] = matrix_entries
            T[1, 3:] = matrix_entries
            AEL = self.parameters["A"] * self.parameters["E"] / l_axial
            K_local = T.T @ (AEL * K_1d) @ T

            K.setOption(PETSc.Mat.Option.NEW_NONZERO_ALLOCATION_ERR, False)
            K.setValues(dofs, dofs, K_local.flat, addv=PETSc.InsertMode.ADD)
            K.assemble()

    def apply_to_forces(self, f_int : PETSc.Vec, u : PETSc.Vec, sign : float = 1.0):
        """
        Adds truss nodal internal forces to global forces vector.

        Args:
            f_int: The global forces vector.
            u: The global displacements.
            sign: Sign of the added values.
        """
        assert sign in [1.0, -1.0]
        points = self.function_space.tabulate_dof_coordinates().flatten()
        f_1d = np.array([-1.0, 1.0])
        T = np.zeros((2, 6))
        for dofs in self.dof_array.reshape(-1, 6):
            delta_x = points[dofs[3:]] - points[dofs[:3]]
            delta_u = u.array[dofs[3:]] - u.array[dofs[:3]]
            l_axial = np.linalg.norm(delta_x, 2)
            matrix_entries = delta_x / l_axial

            u_axial = np.inner(matrix_entries, delta_u)

            eps_axial = u_axial / l_axial

            T[0, :3] = matrix_entries
            T[1, 3:] = matrix_entries

            sigma_axial = self.parameters["E"] * eps_axial

            f_local = sign * T.T @ (self.parameters["A"] * sigma_axial * f_1d)

            f_int.setValues(dofs, f_local, addv=PETSc.InsertMode.ADD)
            f_int.assemble()

__init__(concrete_mesh, rebar_mesh, function_space, parameters)

Initialize rebar class

Parameters:

Name Type Description Default
concrete_mesh dfx.mesh.Mesh

The 3D mesh of the concrete structure.

 required
rebar_mesh dfx.mesh.Mesh

The line mesh of the steel rebar.

 required
function_space dfx.fem.FunctionSpace

The function space object of the concrete structure.

 required
parameters dict

A dictinary containing all needed parameters of the steel. Contains: 'A', 'E', 'rho'.

 required
Source code in reinforcement/rebar.py 
85
86
87
88
89
90
91
92
93
94
95
def __init__(self, concrete_mesh : dfx.mesh.Mesh, rebar_mesh : dfx.mesh.Mesh, function_space : dfx.fem.FunctionSpace, parameters: dict):
    """Initialize rebar class

    Args:
        concrete_mesh: The 3D mesh of the concrete structure.
        rebar_mesh: The line mesh of the steel rebar.
        function_space: The function space object of the concrete structure.
        parameters: A dictinary containing all needed parameters of the steel. Contains: 'A', 'E', 'rho'.

    """
    super().__init__(concrete_mesh, rebar_mesh, function_space, parameters)

apply_to_diagonal_mass(M)

Adds nodal masses to a diagonal mass matrix.

Parameters:

Name Type Description Default
M PETSc.Vec

The diagonal from of the mass matrix as a PETSc vector

 required
Source code in reinforcement/rebar.py 
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
def apply_to_diagonal_mass(self, M : PETSc.Vec):
    """
    Adds nodal masses to a diagonal mass matrix.

    Args:
        M: The diagonal from of the mass matrix as a PETSc vector 

    """
    points = self.function_space.tabulate_dof_coordinates().flatten()
    diagonal_mass_1d = np.array([[1.0, 0.0], [0.0, 1.0]])
    T = np.zeros((2, 6))
    for dofs in self.dof_array.reshape(-1, 6):
        delta_x = points[dofs[3:]] - points[dofs[:3]]
        l_axial = np.linalg.norm(delta_x, 2)
        matrix_entries = delta_x / l_axial

        T[0, :3] = matrix_entries
        T[1, 3:] = matrix_entries

        mass_local = np.diag(
            T.T
            @ (
                self.parameters["rho"]
                * self.parameters["A"]
                * l_axial
                * diagonal_mass_1d
            )
            @ T
        )

        M.setValues(dofs, mass_local, addv=PETSc.InsertMode.ADD)
        M.assemble()

apply_to_forces(f_int, u, sign=1.0)

Adds truss nodal internal forces to global forces vector.

Parameters:

Name Type Description Default
f_int PETSc.Vec

The global forces vector.

 required
u PETSc.Vec

The global displacements.

 required
sign float

Sign of the added values.

 1.0
Source code in reinforcement/rebar.py 
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
def apply_to_forces(self, f_int : PETSc.Vec, u : PETSc.Vec, sign : float = 1.0):
    """
    Adds truss nodal internal forces to global forces vector.

    Args:
        f_int: The global forces vector.
        u: The global displacements.
        sign: Sign of the added values.
    """
    assert sign in [1.0, -1.0]
    points = self.function_space.tabulate_dof_coordinates().flatten()
    f_1d = np.array([-1.0, 1.0])
    T = np.zeros((2, 6))
    for dofs in self.dof_array.reshape(-1, 6):
        delta_x = points[dofs[3:]] - points[dofs[:3]]
        delta_u = u.array[dofs[3:]] - u.array[dofs[:3]]
        l_axial = np.linalg.norm(delta_x, 2)
        matrix_entries = delta_x / l_axial

        u_axial = np.inner(matrix_entries, delta_u)

        eps_axial = u_axial / l_axial

        T[0, :3] = matrix_entries
        T[1, 3:] = matrix_entries

        sigma_axial = self.parameters["E"] * eps_axial

        f_local = sign * T.T @ (self.parameters["A"] * sigma_axial * f_1d)

        f_int.setValues(dofs, f_local, addv=PETSc.InsertMode.ADD)
        f_int.assemble()

apply_to_stiffness(K, u)

Adds truss stiffness to global stiffness matrix,

Parameters:

Name Type Description Default
K PETSc.Mat

The global stiffness matrix.

 required
u PETSc.Vec

The global displacements.

 required
Source code in reinforcement/rebar.py 
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
def apply_to_stiffness(self, K : PETSc.Mat, u : PETSc.Vec):
    """
    Adds truss stiffness to global stiffness matrix,

    Args:
        K: The global stiffness matrix.
        u: The global displacements.
    """
    points = self.function_space.tabulate_dof_coordinates().flatten()
    K_1d = np.array([[1.0, -1.0], [-1.0, 1.0]])
    T = np.zeros((2, 6))
    for dofs in self.dof_array.reshape(-1, 6):
        delta_x = points[dofs[3:]] - points[dofs[:3]]
        l_axial = np.linalg.norm(delta_x, 2)

        matrix_entries = delta_x / l_axial
        T[0, :3] = matrix_entries
        T[1, 3:] = matrix_entries
        AEL = self.parameters["A"] * self.parameters["E"] / l_axial
        K_local = T.T @ (AEL * K_1d) @ T

        K.setOption(PETSc.Mat.Option.NEW_NONZERO_ALLOCATION_ERR, False)
        K.setValues(dofs, dofs, K_local.flat, addv=PETSc.InsertMode.ADD)
        K.assemble()

RebarInterface

Bases: ABC

An interface for trusses. Contains methods to assign dofs.

Source code in reinforcement/rebar.py 
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
class RebarInterface(ABC):
    """
    An interface for trusses. Contains methods to assign dofs.
    """
    def __init__(self, concrete_mesh : dfx.mesh.Mesh, rebar_mesh : dfx.mesh.Mesh, function_space : dfx.fem.FunctionSpace, parameters: dict):
        """Initialize rebar class

        Args:
            concrete_mesh: The 3D mesh of the concrete structure.
            rebar_mesh: The line mesh of the steel rebar.
            function_space: The function space object of the concrete structure.
            parameters: A dictinary containing all needed parameters of the steel. Contains: 'A', 'E', 'rho'.

        """
        self.concrete_mesh = concrete_mesh
        self.rebar_mesh = rebar_mesh
        self.function_space = function_space
        self.parameters = parameters
        self.dof_array = np.array([], dtype=np.int32)
        self._assign_dofs()

    def _assign_dofs(self, tol=1e-6):
        # first all reinforcement dofs and their coordinates are found and saved in geometry_entities and points respectively
        fdim = self.rebar_mesh.topology.dim
        self.rebar_mesh.topology.create_connectivity(fdim, 0)
        num_lines_local = self.rebar_mesh.topology.index_map(fdim).size_local
        geometry_entities = dfx.cpp.mesh.entities_to_geometry(
            self.rebar_mesh, fdim, np.arange(num_lines_local, dtype=np.int32), False
        )
        dofs = []
        for line in geometry_entities:
            start = self.rebar_mesh.geometry.x[line][0]
            end = self.rebar_mesh.geometry.x[line][1]

            dofs_start = self._locate_concrete_dofs(start, tol)
            dofs_end = self._locate_concrete_dofs(end, tol)

            dofs.extend(dofs_start)
            dofs.extend(dofs_end)
        self.dof_array = np.array(dofs, dtype=np.int32).reshape(-1, 3)

    def _locate_concrete_dofs(self, point, tol):
        x, y, z = point

        def rebar_nodes(var):
            return np.logical_and(
                np.logical_and(np.abs(var[1] - y) < tol, np.abs(var[0] - x) < tol),
                np.abs(var[2] - z) < tol,
            )

        dofs_toappend = dfx.fem.locate_dofs_geometrical(
            self.function_space, rebar_nodes
        )
        try:
            assert len(dofs_toappend) == 1
        except AssertionError:
            raise Exception(
                f"{len(dofs_toappend)} dofs found at ({x},{y},{z}), expected 1. Try adjusting the tolerance"
            )

        return dofs_toappend * 3.0 + np.arange(3, dtype=np.float64)

    @abstractmethod
    def apply_to_forces(self, f_int, u):
        pass

    @abstractmethod
    def apply_to_stiffness(self, K, u):
        pass

__init__(concrete_mesh, rebar_mesh, function_space, parameters)

Initialize rebar class

Parameters:

Name Type Description Default
concrete_mesh dfx.mesh.Mesh

The 3D mesh of the concrete structure.

 required
rebar_mesh dfx.mesh.Mesh

The line mesh of the steel rebar.

 required
function_space dfx.fem.FunctionSpace

The function space object of the concrete structure.

 required
parameters dict

A dictinary containing all needed parameters of the steel. Contains: 'A', 'E', 'rho'.

 required
Source code in reinforcement/rebar.py 
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
def __init__(self, concrete_mesh : dfx.mesh.Mesh, rebar_mesh : dfx.mesh.Mesh, function_space : dfx.fem.FunctionSpace, parameters: dict):
    """Initialize rebar class

    Args:
        concrete_mesh: The 3D mesh of the concrete structure.
        rebar_mesh: The line mesh of the steel rebar.
        function_space: The function space object of the concrete structure.
        parameters: A dictinary containing all needed parameters of the steel. Contains: 'A', 'E', 'rho'.

    """
    self.concrete_mesh = concrete_mesh
    self.rebar_mesh = rebar_mesh
    self.function_space = function_space
    self.parameters = parameters
    self.dof_array = np.array([], dtype=np.int32)
    self._assign_dofs()