Elastic 3D beam structures

Elastic 3D beam structures #

Objectives

This tour explores the formulation of 3D beam-like elastic structures. One particularity of this tour is that the mesh topology is 1D (beams) but is embedded in a 3D ambient space. We will also show how to define a local frame containing the beam axis direction and two perpendicular directions for the definition of the cross-section geometrical properties.

See also

This tour bears similarities with the Linear elastic 2D truss and Linear shell model tours.

Download sources

Beam kinematics#

The variational formulation for 3D beams requires distinguishing between motion along the beam direction and along both perpendicular directions. These two directions often correspond to principal directions of inertia and distinguish between strong and weak bending directions. The user must therefore specify this local orientation frame throughout the whole structure. In this example, we will first compute the vector \(\boldsymbol{t}\) tangent to the beam direction, and two perpendicular directions \(\boldsymbol{a}_1\) and \(\boldsymbol{a}_2\). \(\boldsymbol{a}_1\) will always be perpendicular to \(t\) and the vertical direction, while \(\boldsymbol{a}_2\) will be such that \((\boldsymbol{t}, \boldsymbol{a}_1, \boldsymbol{a}_2)\) is a direct orthonormal frame.

Classical Euler-Bernoulli beam theory would typically require \(C^1\)-continuous interpolation of the displacement field, using Hermite polynomials for instance. Such a discretization scheme is unfortunately not available natively in dolfinx. Instead, we will consider a Timoshenko beam model with Saint-Venant uniform torsion theory for the torsional part. The beam kinematics will then be described by a 3D displacement field \(\boldsymbol{u}\) and an independent 3D rotation field \(\boldsymbol{\theta}\), that is 6 degrees of freedom at each point. \((u_t,u_1,u_2)\) will respectively denote the displacement components in the section local frame \((\boldsymbol{t}, \boldsymbol{a}_1, \boldsymbol{a}_2)\). Similarly, \(\theta_1\) and \(\theta_2\) will denote the section rotations about the two section axes whereas \(\theta_t\) is the twist angle.

The beam strains are given by (straight beam with constant local frame):

  • the normal strain: \(\delta = \dfrac{d u_t}{ds}\)

  • the bending curvatures: \(\chi_1 = \dfrac{d\theta_1}{ds}\) and \(\chi_2 = \dfrac{d\theta_2}{ds}\)

  • the shear strains: \(\gamma_1 = \dfrac{du_1}{ds}-\theta_2\) and \(\gamma_2 = \dfrac{du_2}{ds}+\theta_1\).

  • the torsional strain: \(\omega = \dfrac{d\theta_t}{ds}\).

where \(\dfrac{dv}{ds} = \nabla v \cdot \boldsymbol{t}\) is the tangential gradient. Associated with these generalized strains are the corresponding generalized stresses:

  • the normal force \(N\)

  • the bending moments \(M_1\) and \(M_2\)

  • the shear forces \(Q_1\) and \(Q_2\)

  • the torsional moment \(M_T\)

Constitutive equations#

The beam constitutive equations (assuming no normal force/bending moment coupling and that \(a_1\) and \(a_2\) are the principal axes of inertia) read as:

\[\begin{split}\begin{Bmatrix} N \\ Q_1 \\ Q_2 \\ M_T \\ M_1 \\ M_2 \end{Bmatrix} = \begin{bmatrix} ES & 0 & 0 & 0 & 0 & 0\\ 0 & GS_1 & 0 & 0 & 0 & 0\\ 0 & 0 & GS_2 & 0 & 0 & 0\\ 0 & 0 & 0 & GJ & 0 & 0\\ 0 & 0 & 0 & 0 & EI_1 & 0\\ 0 & 0 & 0 & 0 & 0 & EI_2 \end{bmatrix}\begin{Bmatrix} \delta \\ \gamma_1 \\ \gamma_2 \\ \omega \\ \chi_1 \\ \chi_2 \end{Bmatrix}\end{split}\]

where \(S\) is the cross-section area, \(E\) the material Young modulus and \(G\) the shear modulus, \(I_1=\int_S x_2^2 dS\) (resp. \(I_2=\int_S x_1^2 dS\)) are the bending second moment of inertia about axis \(a_1\) (resp. \(a_2\)), \(J\) is the torsion inertia, \(S_1\) (resp. \(S_2\)) is the shear area in direction \(a_1\) (resp. \(a_2\)).

Variational formulation#

The 3D beam variational formulation finally reads as: Find \((\boldsymbol{u},\boldsymbol{\theta})\in V\) such that:

\[\int_L (N\widehat{\delta}+Q_1\widehat{\gamma}_1+Q_2\widehat{\gamma}_2+M_T\widehat{\omega}+M_1\widehat{\chi}_1+M_2\widehat{\chi}_2)ds = \int_L \boldsymbol{f}\cdot\widehat{\boldsymbol{u}}ds \quad \forall (\widehat{\boldsymbol{u}},\widehat{\boldsymbol{\theta}})\in V\]

where we considered only a distributed loading \(\boldsymbol{f}\) and where \(\widehat{\delta},\ldots,\widehat{\chi}_2\) are the generalized strains associated with test functions \((\widehat{\boldsymbol{u}},\widehat{\boldsymbol{\theta}})\).

Implementation#

We start by loading the relevant modules.

import numpy as np
import gmsh
import ufl
import basix
from mpi4py import MPI
from dolfinx import fem, io
import dolfinx.fem.petsc

We generate the mesh (domain) using the Gmsh Python API by turning a circular dome into an elliptic one.

Hide code cell source

 
# Initialize gmsh
gmsh.initialize()

# Parameters
N = 20
h = 10
a = 20
M = 20
d = 1

# Creating points and lines
p0 = gmsh.model.occ.addPoint(a, 0, 0, d)
lines = []

for i in range(1, N + 1):
    p1 = gmsh.model.occ.addPoint(
        a * np.cos(0.99 * i / N * np.pi / 2.0),
        0,
        h * np.sin(i / N * np.pi / 2.0),
        d,
    )
    ll = gmsh.model.occ.addLine(p0, p1)
    lines.append(ll)
    p0 = p1
in_lines = [(1, line) for line in lines]
# Extruding lines
for j in range(1, M + 1):
    out = gmsh.model.occ.revolve(
        in_lines,
        0,
        0,
        0,
        0,
        0,
        1,
        angle=2 * np.pi / M,
    )
    in_lines = out[::4]


gmsh.model.occ.remove_all_duplicates()
gmsh.model.occ.synchronize()
lines = gmsh.model.get_entities(1)

gmsh.model.add_physical_group(1, [l[1] for l in lines], 1)
gmsh.model.mesh.generate(dim=1)

mesh_data = io.gmsh.model_to_mesh(gmsh.model, MPI.COMM_WORLD, 0)
domain = mesh_data.mesh
domain.geometry.x[:, 1] *= 2 / 3.0

# Finalize gmsh
gmsh.finalize()

Hide code cell output

 
Info    : [  0%] Fragments                                                                                  
Info    : [ 10%] Fragments                                                                                  

Info    : [ 20%] Fragments                                                                                  
Info    : [ 30%] Fragments                                                                                  

Info    : [ 40%] Fragments                                                                                  
Info    : [ 50%] Fragments                                                                                  
Info    : [ 60%] Fragments                                                                                  
Info    : [ 70%] Fragments - Filling splits of vertices                                                                                
Info    : [ 80%] Fragments - Splitting faces                                                                                

Info    : Meshing 1D...
Info    : [  0%] Meshing curve 23 (Line)
Info    : [ 10%] Meshing curve 25 (Line)
Info    : [ 10%] Meshing curve 27 (Line)
Info    : [ 10%] Meshing curve 29 (Line)
Info    : [ 10%] Meshing curve 31 (Line)
Info    : [ 10%] Meshing curve 33 (Line)
Info    : [ 10%] Meshing curve 35 (Line)
Info    : [ 10%] Meshing curve 37 (Line)
Info    : [ 10%] Meshing curve 39 (Line)
Info    : [ 10%] Meshing curve 41 (Line)
Info    : [ 10%] Meshing curve 43 (Line)
Info    : [ 10%] Meshing curve 45 (Line)
Info    : [ 10%] Meshing curve 47 (Line)
Info    : [ 10%] Meshing curve 49 (Line)
Info    : [ 10%] Meshing curve 51 (Line)
Info    : [ 10%] Meshing curve 53 (Line)
Info    : [ 10%] Meshing curve 55 (Line)
Info    : [ 10%] Meshing curve 57 (Line)
Info    : [ 10%] Meshing curve 59 (Line)
Info    : [ 10%] Meshing curve 61 (Line)
Info    : [ 10%] Meshing curve 62 (Circle)
Info    : [ 10%] Meshing curve 63 (Circle)
Info    : [ 10%] Meshing curve 64 (Line)
Info    : [ 10%] Meshing curve 65 (Circle)
Info    : [ 10%] Meshing curve 66 (Line)
Info    : [ 10%] Meshing curve 67 (Circle)
Info    : [ 10%] Meshing curve 68 (Line)
Info    : [ 10%] Meshing curve 69 (Circle)
Info    : [ 10%] Meshing curve 70 (Line)
Info    : [ 10%] Meshing curve 71 (Circle)
Info    : [ 10%] Meshing curve 72 (Line)
Info    : [ 10%] Meshing curve 73 (Circle)
Info    : [ 10%] Meshing curve 74 (Line)
Info    : [ 10%] Meshing curve 75 (Circle)
Info    : [ 10%] Meshing curve 76 (Line)
Info    : [ 10%] Meshing curve 77 (Circle)
Info    : [ 10%] Meshing curve 78 (Line)
Info    : [ 10%] Meshing curve 79 (Circle)
Info    : [ 10%] Meshing curve 80 (Line)
Info    : [ 10%] Meshing curve 81 (Circle)
Info    : [ 10%] Meshing curve 82 (Line)
Info    : [ 10%] Meshing curve 83 (Circle)
Info    : [ 10%] Meshing curve 84 (Line)
Info    : [ 10%] Meshing curve 85 (Circle)
Info    : [ 10%] Meshing curve 86 (Line)
Info    : [ 10%] Meshing curve 87 (Circle)
Info    : [ 10%] Meshing curve 88 (Line)
Info    : [ 10%] Meshing curve 89 (Circle)
Info    : [ 10%] Meshing curve 90 (Line)
Info    : [ 10%] Meshing curve 91 (Circle)
Info    : [ 10%] Meshing curve 92 (Line)
Info    : [ 10%] Meshing curve 93 (Circle)
Info    : [ 10%] Meshing curve 94 (Line)
Info    : [ 10%] Meshing curve 95 (Circle)
Info    : [ 10%] Meshing curve 96 (Line)
Info    : [ 10%] Meshing curve 97 (Circle)
Info    : [ 10%] Meshing curve 98 (Line)
Info    : [ 10%] Meshing curve 99 (Circle)
Info    : [ 10%] Meshing curve 100 (Line)
Info    : [ 10%] Meshing curve 101 (Circle)
Info    : [ 10%] Meshing curve 102 (Line)
Info    : [ 10%] Meshing curve 103 (Circle)
Info    : [ 10%] Meshing curve 104 (Circle)
Info    : [ 10%] Meshing curve 105 (Line)
Info    : [ 10%] Meshing curve 106 (Circle)
Info    : [ 10%] Meshing curve 107 (Line)
Info    : [ 10%] Meshing curve 108 (Circle)
Info    : [ 10%] Meshing curve 109 (Line)
Info    : [ 10%] Meshing curve 110 (Circle)
Info    : [ 10%] Meshing curve 111 (Line)
Info    : [ 10%] Meshing curve 112 (Circle)
Info    : [ 10%] Meshing curve 113 (Line)
Info    : [ 10%] Meshing curve 114 (Circle)
Info    : [ 10%] Meshing curve 115 (Line)
Info    : [ 10%] Meshing curve 116 (Circle)
Info    : [ 10%] Meshing curve 117 (Line)
Info    : [ 10%] Meshing curve 118 (Circle)
Info    : [ 10%] Meshing curve 119 (Line)
Info    : [ 10%] Meshing curve 120 (Circle)
Info    : [ 10%] Meshing curve 121 (Line)
Info    : [ 10%] Meshing curve 122 (Circle)
Info    : [ 10%] Meshing curve 123 (Line)
Info    : [ 20%] Meshing curve 124 (Circle)
Info    : [ 20%] Meshing curve 125 (Line)
Info    : [ 20%] Meshing curve 126 (Circle)
Info    : [ 20%] Meshing curve 127 (Line)
Info    : [ 20%] Meshing curve 128 (Circle)
Info    : [ 20%] Meshing curve 129 (Line)
Info    : [ 20%] Meshing curve 130 (Circle)
Info    : [ 20%] Meshing curve 131 (Line)
Info    : [ 20%] Meshing curve 132 (Circle)
Info    : [ 20%] Meshing curve 133 (Line)
Info    : [ 20%] Meshing curve 134 (Circle)
Info    : [ 20%] Meshing curve 135 (Line)
Info    : [ 20%] Meshing curve 136 (Circle)
Info    : [ 20%] Meshing curve 137 (Line)
Info    : [ 20%] Meshing curve 138 (Circle)
Info    : [ 20%] Meshing curve 139 (Line)
Info    : [ 20%] Meshing curve 140 (Circle)
Info    : [ 20%] Meshing curve 141 (Line)
Info    : [ 20%] Meshing curve 142 (Circle)
Info    : [ 20%] Meshing curve 143 (Line)
Info    : [ 20%] Meshing curve 144 (Circle)
Info    : [ 20%] Meshing curve 145 (Circle)
Info    : [ 20%] Meshing curve 146 (Line)
Info    : [ 20%] Meshing curve 147 (Circle)
Info    : [ 20%] Meshing curve 148 (Line)
Info    : [ 20%] Meshing curve 149 (Circle)
Info    : [ 20%] Meshing curve 150 (Line)
Info    : [ 20%] Meshing curve 151 (Circle)
Info    : [ 20%] Meshing curve 152 (Line)
Info    : [ 20%] Meshing curve 153 (Circle)
Info    : [ 20%] Meshing curve 154 (Line)
Info    : [ 20%] Meshing curve 155 (Circle)
Info    : [ 20%] Meshing curve 156 (Line)
Info    : [ 20%] Meshing curve 157 (Circle)
Info    : [ 20%] Meshing curve 158 (Line)
Info    : [ 20%] Meshing curve 159 (Circle)
Info    : [ 20%] Meshing curve 160 (Line)
Info    : [ 20%] Meshing curve 161 (Circle)
Info    : [ 20%] Meshing curve 162 (Line)
Info    : [ 20%] Meshing curve 163 (Circle)
Info    : [ 20%] Meshing curve 164 (Line)
Info    : [ 20%] Meshing curve 165 (Circle)
Info    : [ 20%] Meshing curve 166 (Line)
Info    : [ 20%] Meshing curve 167 (Circle)
Info    : [ 20%] Meshing curve 168 (Line)
Info    : [ 20%] Meshing curve 169 (Circle)
Info    : [ 20%] Meshing curve 170 (Line)
Info    : [ 20%] Meshing curve 171 (Circle)
Info    : [ 20%] Meshing curve 172 (Line)
Info    : [ 20%] Meshing curve 173 (Circle)
Info    : [ 20%] Meshing curve 174 (Line)
Info    : [ 20%] Meshing curve 175 (Circle)
Info    : [ 20%] Meshing curve 176 (Line)
Info    : [ 20%] Meshing curve 177 (Circle)
Info    : [ 20%] Meshing curve 178 (Line)
Info    : [ 20%] Meshing curve 179 (Circle)
Info    : [ 20%] Meshing curve 180 (Line)
Info    : [ 20%] Meshing curve 181 (Circle)
Info    : [ 20%] Meshing curve 182 (Line)
Info    : [ 20%] Meshing curve 183 (Circle)
Info    : [ 20%] Meshing curve 184 (Line)
Info    : [ 20%] Meshing curve 185 (Circle)
Info    : [ 20%] Meshing curve 186 (Circle)
Info    : [ 20%] Meshing curve 187 (Line)
Info    : [ 20%] Meshing curve 188 (Circle)
Info    : [ 20%] Meshing curve 189 (Line)
Info    : [ 20%] Meshing curve 190 (Circle)
Info    : [ 20%] Meshing curve 191 (Line)
Info    : [ 20%] Meshing curve 192 (Circle)
Info    : [ 20%] Meshing curve 193 (Line)
Info    : [ 20%] Meshing curve 194 (Circle)
Info    : [ 20%] Meshing curve 195 (Line)
Info    : [ 20%] Meshing curve 196 (Circle)
Info    : [ 20%] Meshing curve 197 (Line)
Info    : [ 20%] Meshing curve 198 (Circle)
Info    : [ 20%] Meshing curve 199 (Line)
Info    : [ 20%] Meshing curve 200 (Circle)
Info    : [ 20%] Meshing curve 201 (Line)
Info    : [ 20%] Meshing curve 202 (Circle)
Info    : [ 20%] Meshing curve 203 (Line)
Info    : [ 20%] Meshing curve 204 (Circle)
Info    : [ 20%] Meshing curve 205 (Line)
Info    : [ 30%] Meshing curve 206 (Circle)
Info    : [ 30%] Meshing curve 207 (Line)
Info    : [ 30%] Meshing curve 208 (Circle)
Info    : [ 30%] Meshing curve 209 (Line)
Info    : [ 30%] Meshing curve 210 (Circle)
Info    : [ 30%] Meshing curve 211 (Line)
Info    : [ 30%] Meshing curve 212 (Circle)
Info    : [ 30%] Meshing curve 213 (Line)
Info    : [ 30%] Meshing curve 214 (Circle)
Info    : [ 30%] Meshing curve 215 (Line)
Info    : [ 30%] Meshing curve 216 (Circle)
Info    : [ 30%] Meshing curve 217 (Line)
Info    : [ 30%] Meshing curve 218 (Circle)
Info    : [ 30%] Meshing curve 219 (Line)
Info    : [ 30%] Meshing curve 220 (Circle)
Info    : [ 30%] Meshing curve 221 (Line)
Info    : [ 30%] Meshing curve 222 (Circle)
Info    : [ 30%] Meshing curve 223 (Line)
Info    : [ 30%] Meshing curve 224 (Circle)
Info    : [ 30%] Meshing curve 225 (Line)
Info    : [ 30%] Meshing curve 226 (Circle)
Info    : [ 30%] Meshing curve 227 (Circle)
Info    : [ 30%] Meshing curve 228 (Line)
Info    : [ 30%] Meshing curve 229 (Circle)
Info    : [ 30%] Meshing curve 230 (Line)
Info    : [ 30%] Meshing curve 231 (Circle)
Info    : [ 30%] Meshing curve 232 (Line)
Info    : [ 30%] Meshing curve 233 (Circle)
Info    : [ 30%] Meshing curve 234 (Line)
Info    : [ 30%] Meshing curve 235 (Circle)
Info    : [ 30%] Meshing curve 236 (Line)
Info    : [ 30%] Meshing curve 237 (Circle)
Info    : [ 30%] Meshing curve 238 (Line)
Info    : [ 30%] Meshing curve 239 (Circle)
Info    : [ 30%] Meshing curve 240 (Line)
Info    : [ 30%] Meshing curve 241 (Circle)
Info    : [ 30%] Meshing curve 242 (Line)
Info    : [ 30%] Meshing curve 243 (Circle)
Info    : [ 30%] Meshing curve 244 (Line)
Info    : [ 30%] Meshing curve 245 (Circle)
Info    : [ 30%] Meshing curve 246 (Line)
Info    : [ 30%] Meshing curve 247 (Circle)
Info    : [ 30%] Meshing curve 248 (Line)
Info    : [ 30%] Meshing curve 249 (Circle)
Info    : [ 30%] Meshing curve 250 (Line)
Info    : [ 30%] Meshing curve 251 (Circle)
Info    : [ 30%] Meshing curve 252 (Line)
Info    : [ 30%] Meshing curve 253 (Circle)
Info    : [ 30%] Meshing curve 254 (Line)
Info    : [ 30%] Meshing curve 255 (Circle)
Info    : [ 30%] Meshing curve 256 (Line)
Info    : [ 30%] Meshing curve 257 (Circle)
Info    : [ 30%] Meshing curve 258 (Line)
Info    : [ 30%] Meshing curve 259 (Circle)
Info    : [ 30%] Meshing curve 260 (Line)
Info    : [ 30%] Meshing curve 261 (Circle)
Info    : [ 30%] Meshing curve 262 (Line)
Info    : [ 30%] Meshing curve 263 (Circle)
Info    : [ 30%] Meshing curve 264 (Line)
Info    : [ 30%] Meshing curve 265 (Circle)
Info    : [ 30%] Meshing curve 266 (Line)
Info    : [ 30%] Meshing curve 267 (Circle)
Info    : [ 30%] Meshing curve 268 (Circle)
Info    : [ 30%] Meshing curve 269 (Line)
Info    : [ 30%] Meshing curve 270 (Circle)
Info    : [ 30%] Meshing curve 271 (Line)
Info    : [ 30%] Meshing curve 272 (Circle)
Info    : [ 30%] Meshing curve 273 (Line)
Info    : [ 30%] Meshing curve 274 (Circle)
Info    : [ 30%] Meshing curve 275 (Line)
Info    : [ 30%] Meshing curve 276 (Circle)
Info    : [ 30%] Meshing curve 277 (Line)
Info    : [ 30%] Meshing curve 278 (Circle)
Info    : [ 30%] Meshing curve 279 (Line)
Info    : [ 30%] Meshing curve 280 (Circle)
Info    : [ 30%] Meshing curve 281 (Line)
Info    : [ 30%] Meshing curve 282 (Circle)
Info    : [ 30%] Meshing curve 283 (Line)
Info    : [ 30%] Meshing curve 284 (Circle)
Info    : [ 30%] Meshing curve 285 (Line)
Info    : [ 30%] Meshing curve 286 (Circle)
Info    : [ 30%] Meshing curve 287 (Line)
Info    : [ 40%] Meshing curve 288 (Circle)
Info    : [ 40%] Meshing curve 289 (Line)
Info    : [ 40%] Meshing curve 290 (Circle)
Info    : [ 40%] Meshing curve 291 (Line)
Info    : [ 40%] Meshing curve 292 (Circle)
Info    : [ 40%] Meshing curve 293 (Line)
Info    : [ 40%] Meshing curve 294 (Circle)
Info    : [ 40%] Meshing curve 295 (Line)
Info    : [ 40%] Meshing curve 296 (Circle)
Info    : [ 40%] Meshing curve 297 (Line)
Info    : [ 40%] Meshing curve 298 (Circle)
Info    : [ 40%] Meshing curve 299 (Line)
Info    : [ 40%] Meshing curve 300 (Circle)
Info    : [ 40%] Meshing curve 301 (Line)
Info    : [ 40%] Meshing curve 302 (Circle)
Info    : [ 40%] Meshing curve 303 (Line)
Info    : [ 40%] Meshing curve 304 (Circle)
Info    : [ 40%] Meshing curve 305 (Line)
Info    : [ 40%] Meshing curve 306 (Circle)
Info    : [ 40%] Meshing curve 307 (Line)
Info    : [ 40%] Meshing curve 308 (Circle)
Info    : [ 40%] Meshing curve 309 (Circle)
Info    : [ 40%] Meshing curve 310 (Line)
Info    : [ 40%] Meshing curve 311 (Circle)
Info    : [ 40%] Meshing curve 312 (Line)
Info    : [ 40%] Meshing curve 313 (Circle)
Info    : [ 40%] Meshing curve 314 (Line)
Info    : [ 40%] Meshing curve 315 (Circle)
Info    : [ 40%] Meshing curve 316 (Line)
Info    : [ 40%] Meshing curve 317 (Circle)
Info    : [ 40%] Meshing curve 318 (Line)
Info    : [ 40%] Meshing curve 319 (Circle)
Info    : [ 40%] Meshing curve 320 (Line)
Info    : [ 40%] Meshing curve 321 (Circle)
Info    : [ 40%] Meshing curve 322 (Line)
Info    : [ 40%] Meshing curve 323 (Circle)
Info    : [ 40%] Meshing curve 324 (Line)
Info    : [ 40%] Meshing curve 325 (Circle)
Info    : [ 40%] Meshing curve 326 (Line)
Info    : [ 40%] Meshing curve 327 (Circle)
Info    : [ 40%] Meshing curve 328 (Line)
Info    : [ 40%] Meshing curve 329 (Circle)
Info    : [ 40%] Meshing curve 330 (Line)
Info    : [ 40%] Meshing curve 331 (Circle)
Info    : [ 40%] Meshing curve 332 (Line)
Info    : [ 40%] Meshing curve 333 (Circle)
Info    : [ 40%] Meshing curve 334 (Line)
Info    : [ 40%] Meshing curve 335 (Circle)
Info    : [ 40%] Meshing curve 336 (Line)
Info    : [ 40%] Meshing curve 337 (Circle)
Info    : [ 40%] Meshing curve 338 (Line)
Info    : [ 40%] Meshing curve 339 (Circle)
Info    : [ 40%] Meshing curve 340 (Line)
Info    : [ 40%] Meshing curve 341 (Circle)
Info    : [ 40%] Meshing curve 342 (Line)
Info    : [ 40%] Meshing curve 343 (Circle)
Info    : [ 40%] Meshing curve 344 (Line)
Info    : [ 40%] Meshing curve 345 (Circle)
Info    : [ 40%] Meshing curve 346 (Line)
Info    : [ 40%] Meshing curve 347 (Circle)
Info    : [ 40%] Meshing curve 348 (Line)
Info    : [ 40%] Meshing curve 349 (Circle)
Info    : [ 40%] Meshing curve 350 (Circle)
Info    : [ 40%] Meshing curve 351 (Line)
Info    : [ 40%] Meshing curve 352 (Circle)
Info    : [ 40%] Meshing curve 353 (Line)
Info    : [ 40%] Meshing curve 354 (Circle)
Info    : [ 40%] Meshing curve 355 (Line)
Info    : [ 40%] Meshing curve 356 (Circle)
Info    : [ 40%] Meshing curve 357 (Line)
Info    : [ 40%] Meshing curve 358 (Circle)
Info    : [ 40%] Meshing curve 359 (Line)
Info    : [ 40%] Meshing curve 360 (Circle)
Info    : [ 40%] Meshing curve 361 (Line)
Info    : [ 40%] Meshing curve 362 (Circle)
Info    : [ 40%] Meshing curve 363 (Line)
Info    : [ 40%] Meshing curve 364 (Circle)
Info    : [ 40%] Meshing curve 365 (Line)
Info    : [ 40%] Meshing curve 366 (Circle)
Info    : [ 40%] Meshing curve 367 (Line)
Info    : [ 40%] Meshing curve 368 (Circle)
Info    : [ 40%] Meshing curve 369 (Line)
Info    : [ 50%] Meshing curve 370 (Circle)
Info    : [ 50%] Meshing curve 371 (Line)
Info    : [ 50%] Meshing curve 372 (Circle)
Info    : [ 50%] Meshing curve 373 (Line)
Info    : [ 50%] Meshing curve 374 (Circle)
Info    : [ 50%] Meshing curve 375 (Line)
Info    : [ 50%] Meshing curve 376 (Circle)
Info    : [ 50%] Meshing curve 377 (Line)
Info    : [ 50%] Meshing curve 378 (Circle)
Info    : [ 50%] Meshing curve 379 (Line)
Info    : [ 50%] Meshing curve 380 (Circle)
Info    : [ 50%] Meshing curve 381 (Line)
Info    : [ 50%] Meshing curve 382 (Circle)
Info    : [ 50%] Meshing curve 383 (Line)
Info    : [ 50%] Meshing curve 384 (Circle)
Info    : [ 50%] Meshing curve 385 (Line)
Info    : [ 50%] Meshing curve 386 (Circle)
Info    : [ 50%] Meshing curve 387 (Line)
Info    : [ 50%] Meshing curve 388 (Circle)
Info    : [ 50%] Meshing curve 389 (Line)
Info    : [ 50%] Meshing curve 390 (Circle)
Info    : [ 50%] Meshing curve 391 (Circle)
Info    : [ 50%] Meshing curve 392 (Line)
Info    : [ 50%] Meshing curve 393 (Circle)
Info    : [ 50%] Meshing curve 394 (Line)
Info    : [ 50%] Meshing curve 395 (Circle)
Info    : [ 50%] Meshing curve 396 (Line)
Info    : [ 50%] Meshing curve 397 (Circle)
Info    : [ 50%] Meshing curve 398 (Line)
Info    : [ 50%] Meshing curve 399 (Circle)
Info    : [ 50%] Meshing curve 400 (Line)
Info    : [ 50%] Meshing curve 401 (Circle)
Info    : [ 50%] Meshing curve 402 (Line)
Info    : [ 50%] Meshing curve 403 (Circle)
Info    : [ 50%] Meshing curve 404 (Line)
Info    : [ 50%] Meshing curve 405 (Circle)
Info    : [ 50%] Meshing curve 406 (Line)
Info    : [ 50%] Meshing curve 407 (Circle)
Info    : [ 50%] Meshing curve 408 (Line)
Info    : [ 50%] Meshing curve 409 (Circle)
Info    : [ 50%] Meshing curve 410 (Line)
Info    : [ 50%] Meshing curve 411 (Circle)
Info    : [ 50%] Meshing curve 412 (Line)
Info    : [ 50%] Meshing curve 413 (Circle)
Info    : [ 50%] Meshing curve 414 (Line)
Info    : [ 50%] Meshing curve 415 (Circle)
Info    : [ 50%] Meshing curve 416 (Line)
Info    : [ 50%] Meshing curve 417 (Circle)
Info    : [ 50%] Meshing curve 418 (Line)
Info    : [ 50%] Meshing curve 419 (Circle)
Info    : [ 50%] Meshing curve 420 (Line)
Info    : [ 50%] Meshing curve 421 (Circle)
Info    : [ 50%] Meshing curve 422 (Line)
Info    : [ 50%] Meshing curve 423 (Circle)
Info    : [ 50%] Meshing curve 424 (Line)
Info    : [ 50%] Meshing curve 425 (Circle)
Info    : [ 50%] Meshing curve 426 (Line)
Info    : [ 50%] Meshing curve 427 (Circle)
Info    : [ 50%] Meshing curve 428 (Line)
Info    : [ 50%] Meshing curve 429 (Circle)
Info    : [ 50%] Meshing curve 430 (Line)
Info    : [ 50%] Meshing curve 431 (Circle)
Info    : [ 50%] Meshing curve 432 (Circle)
Info    : [ 50%] Meshing curve 433 (Line)
Info    : [ 50%] Meshing curve 434 (Circle)
Info    : [ 50%] Meshing curve 435 (Line)
Info    : [ 50%] Meshing curve 436 (Circle)
Info    : [ 50%] Meshing curve 437 (Line)
Info    : [ 50%] Meshing curve 438 (Circle)
Info    : [ 50%] Meshing curve 439 (Line)
Info    : [ 50%] Meshing curve 440 (Circle)
Info    : [ 50%] Meshing curve 441 (Line)
Info    : [ 50%] Meshing curve 442 (Circle)
Info    : [ 50%] Meshing curve 443 (Line)
Info    : [ 50%] Meshing curve 444 (Circle)
Info    : [ 50%] Meshing curve 445 (Line)
Info    : [ 50%] Meshing curve 446 (Circle)
Info    : [ 50%] Meshing curve 447 (Line)
Info    : [ 50%] Meshing curve 448 (Circle)
Info    : [ 50%] Meshing curve 449 (Line)
Info    : [ 50%] Meshing curve 450 (Circle)
Info    : [ 50%] Meshing curve 451 (Line)
Info    : [ 60%] Meshing curve 452 (Circle)
Info    : [ 60%] Meshing curve 453 (Line)
Info    : [ 60%] Meshing curve 454 (Circle)
Info    : [ 60%] Meshing curve 455 (Line)
Info    : [ 60%] Meshing curve 456 (Circle)
Info    : [ 60%] Meshing curve 457 (Line)
Info    : [ 60%] Meshing curve 458 (Circle)
Info    : [ 60%] Meshing curve 459 (Line)
Info    : [ 60%] Meshing curve 460 (Circle)
Info    : [ 60%] Meshing curve 461 (Line)
Info    : [ 60%] Meshing curve 462 (Circle)
Info    : [ 60%] Meshing curve 463 (Line)
Info    : [ 60%] Meshing curve 464 (Circle)
Info    : [ 60%] Meshing curve 465 (Line)
Info    : [ 60%] Meshing curve 466 (Circle)
Info    : [ 60%] Meshing curve 467 (Line)
Info    : [ 60%] Meshing curve 468 (Circle)
Info    : [ 60%] Meshing curve 469 (Line)
Info    : [ 60%] Meshing curve 470 (Circle)
Info    : [ 60%] Meshing curve 471 (Line)
Info    : [ 60%] Meshing curve 472 (Circle)
Info    : [ 60%] Meshing curve 473 (Circle)
Info    : [ 60%] Meshing curve 474 (Line)
Info    : [ 60%] Meshing curve 475 (Circle)
Info    : [ 60%] Meshing curve 476 (Line)
Info    : [ 60%] Meshing curve 477 (Circle)
Info    : [ 60%] Meshing curve 478 (Line)
Info    : [ 60%] Meshing curve 479 (Circle)
Info    : [ 60%] Meshing curve 480 (Line)
Info    : [ 60%] Meshing curve 481 (Circle)
Info    : [ 60%] Meshing curve 482 (Line)
Info    : [ 60%] Meshing curve 483 (Circle)
Info    : [ 60%] Meshing curve 484 (Line)
Info    : [ 60%] Meshing curve 485 (Circle)
Info    : [ 60%] Meshing curve 486 (Line)
Info    : [ 60%] Meshing curve 487 (Circle)
Info    : [ 60%] Meshing curve 488 (Line)
Info    : [ 60%] Meshing curve 489 (Circle)
Info    : [ 60%] Meshing curve 490 (Line)
Info    : [ 60%] Meshing curve 491 (Circle)
Info    : [ 60%] Meshing curve 492 (Line)
Info    : [ 60%] Meshing curve 493 (Circle)
Info    : [ 60%] Meshing curve 494 (Line)
Info    : [ 60%] Meshing curve 495 (Circle)
Info    : [ 60%] Meshing curve 496 (Line)
Info    : [ 60%] Meshing curve 497 (Circle)
Info    : [ 60%] Meshing curve 498 (Line)
Info    : [ 60%] Meshing curve 499 (Circle)
Info    : [ 60%] Meshing curve 500 (Line)
Info    : [ 60%] Meshing curve 501 (Circle)
Info    : [ 60%] Meshing curve 502 (Line)
Info    : [ 60%] Meshing curve 503 (Circle)
Info    : [ 60%] Meshing curve 504 (Line)
Info    : [ 60%] Meshing curve 505 (Circle)
Info    : [ 60%] Meshing curve 506 (Line)
Info    : [ 60%] Meshing curve 507 (Circle)
Info    : [ 60%] Meshing curve 508 (Line)
Info    : [ 60%] Meshing curve 509 (Circle)
Info    : [ 60%] Meshing curve 510 (Line)
Info    : [ 60%] Meshing curve 511 (Circle)
Info    : [ 60%] Meshing curve 512 (Line)
Info    : [ 60%] Meshing curve 513 (Circle)
Info    : [ 60%] Meshing curve 514 (Circle)
Info    : [ 60%] Meshing curve 515 (Line)
Info    : [ 60%] Meshing curve 516 (Circle)
Info    : [ 60%] Meshing curve 517 (Line)
Info    : [ 60%] Meshing curve 518 (Circle)
Info    : [ 60%] Meshing curve 519 (Line)
Info    : [ 60%] Meshing curve 520 (Circle)
Info    : [ 60%] Meshing curve 521 (Line)
Info    : [ 60%] Meshing curve 522 (Circle)
Info    : [ 60%] Meshing curve 523 (Line)
Info    : [ 60%] Meshing curve 524 (Circle)
Info    : [ 60%] Meshing curve 525 (Line)
Info    : [ 60%] Meshing curve 526 (Circle)
Info    : [ 60%] Meshing curve 527 (Line)
Info    : [ 60%] Meshing curve 528 (Circle)
Info    : [ 60%] Meshing curve 529 (Line)
Info    : [ 60%] Meshing curve 530 (Circle)
Info    : [ 60%] Meshing curve 531 (Line)
Info    : [ 60%] Meshing curve 532 (Circle)
Info    : [ 60%] Meshing curve 533 (Line)
Info    : [ 70%] Meshing curve 534 (Circle)
Info    : [ 70%] Meshing curve 535 (Line)
Info    : [ 70%] Meshing curve 536 (Circle)
Info    : [ 70%] Meshing curve 537 (Line)
Info    : [ 70%] Meshing curve 538 (Circle)
Info    : [ 70%] Meshing curve 539 (Line)
Info    : [ 70%] Meshing curve 540 (Circle)
Info    : [ 70%] Meshing curve 541 (Line)
Info    : [ 70%] Meshing curve 542 (Circle)
Info    : [ 70%] Meshing curve 543 (Line)
Info    : [ 70%] Meshing curve 544 (Circle)
Info    : [ 70%] Meshing curve 545 (Line)
Info    : [ 70%] Meshing curve 546 (Circle)
Info    : [ 70%] Meshing curve 547 (Line)
Info    : [ 70%] Meshing curve 548 (Circle)
Info    : [ 70%] Meshing curve 549 (Line)
Info    : [ 70%] Meshing curve 550 (Circle)
Info    : [ 70%] Meshing curve 551 (Line)
Info    : [ 70%] Meshing curve 552 (Circle)
Info    : [ 70%] Meshing curve 553 (Line)
Info    : [ 70%] Meshing curve 554 (Circle)
Info    : [ 70%] Meshing curve 555 (Circle)
Info    : [ 70%] Meshing curve 556 (Line)
Info    : [ 70%] Meshing curve 557 (Circle)
Info    : [ 70%] Meshing curve 558 (Line)
Info    : [ 70%] Meshing curve 559 (Circle)
Info    : [ 70%] Meshing curve 560 (Line)
Info    : [ 70%] Meshing curve 561 (Circle)
Info    : [ 70%] Meshing curve 562 (Line)
Info    : [ 70%] Meshing curve 563 (Circle)
Info    : [ 70%] Meshing curve 564 (Line)
Info    : [ 70%] Meshing curve 565 (Circle)
Info    : [ 70%] Meshing curve 566 (Line)
Info    : [ 70%] Meshing curve 567 (Circle)
Info    : [ 70%] Meshing curve 568 (Line)
Info    : [ 70%] Meshing curve 569 (Circle)
Info    : [ 70%] Meshing curve 570 (Line)
Info    : [ 70%] Meshing curve 571 (Circle)
Info    : [ 70%] Meshing curve 572 (Line)
Info    : [ 70%] Meshing curve 573 (Circle)
Info    : [ 70%] Meshing curve 574 (Line)
Info    : [ 70%] Meshing curve 575 (Circle)
Info    : [ 70%] Meshing curve 576 (Line)
Info    : [ 70%] Meshing curve 577 (Circle)
Info    : [ 70%] Meshing curve 578 (Line)
Info    : [ 70%] Meshing curve 579 (Circle)
Info    : [ 70%] Meshing curve 580 (Line)
Info    : [ 70%] Meshing curve 581 (Circle)
Info    : [ 70%] Meshing curve 582 (Line)
Info    : [ 70%] Meshing curve 583 (Circle)
Info    : [ 70%] Meshing curve 584 (Line)
Info    : [ 70%] Meshing curve 585 (Circle)
Info    : [ 70%] Meshing curve 586 (Line)
Info    : [ 70%] Meshing curve 587 (Circle)
Info    : [ 70%] Meshing curve 588 (Line)
Info    : [ 70%] Meshing curve 589 (Circle)
Info    : [ 70%] Meshing curve 590 (Line)
Info    : [ 70%] Meshing curve 591 (Circle)
Info    : [ 70%] Meshing curve 592 (Line)
Info    : [ 70%] Meshing curve 593 (Circle)
Info    : [ 70%] Meshing curve 594 (Line)
Info    : [ 70%] Meshing curve 595 (Circle)
Info    : [ 70%] Meshing curve 596 (Circle)
Info    : [ 70%] Meshing curve 597 (Line)
Info    : [ 70%] Meshing curve 598 (Circle)
Info    : [ 70%] Meshing curve 599 (Line)
Info    : [ 70%] Meshing curve 600 (Circle)
Info    : [ 70%] Meshing curve 601 (Line)
Info    : [ 70%] Meshing curve 602 (Circle)
Info    : [ 70%] Meshing curve 603 (Line)
Info    : [ 70%] Meshing curve 604 (Circle)
Info    : [ 70%] Meshing curve 605 (Line)
Info    : [ 70%] Meshing curve 606 (Circle)
Info    : [ 70%] Meshing curve 607 (Line)
Info    : [ 70%] Meshing curve 608 (Circle)
Info    : [ 70%] Meshing curve 609 (Line)
Info    : [ 70%] Meshing curve 610 (Circle)
Info    : [ 70%] Meshing curve 611 (Line)
Info    : [ 70%] Meshing curve 612 (Circle)
Info    : [ 70%] Meshing curve 613 (Line)
Info    : [ 70%] Meshing curve 614 (Circle)
Info    : [ 70%] Meshing curve 615 (Line)
Info    : [ 80%] Meshing curve 616 (Circle)
Info    : [ 80%] Meshing curve 617 (Line)
Info    : [ 80%] Meshing curve 618 (Circle)
Info    : [ 80%] Meshing curve 619 (Line)
Info    : [ 80%] Meshing curve 620 (Circle)
Info    : [ 80%] Meshing curve 621 (Line)
Info    : [ 80%] Meshing curve 622 (Circle)
Info    : [ 80%] Meshing curve 623 (Line)
Info    : [ 80%] Meshing curve 624 (Circle)
Info    : [ 80%] Meshing curve 625 (Line)
Info    : [ 80%] Meshing curve 626 (Circle)
Info    : [ 80%] Meshing curve 627 (Line)
Info    : [ 80%] Meshing curve 628 (Circle)
Info    : [ 80%] Meshing curve 629 (Line)
Info    : [ 80%] Meshing curve 630 (Circle)
Info    : [ 80%] Meshing curve 631 (Line)
Info    : [ 80%] Meshing curve 632 (Circle)
Info    : [ 80%] Meshing curve 633 (Line)
Info    : [ 80%] Meshing curve 634 (Circle)
Info    : [ 80%] Meshing curve 635 (Line)
Info    : [ 80%] Meshing curve 636 (Circle)
Info    : [ 80%] Meshing curve 637 (Circle)
Info    : [ 80%] Meshing curve 638 (Line)
Info    : [ 80%] Meshing curve 639 (Circle)
Info    : [ 80%] Meshing curve 640 (Line)
Info    : [ 80%] Meshing curve 641 (Circle)
Info    : [ 80%] Meshing curve 642 (Line)
Info    : [ 80%] Meshing curve 643 (Circle)
Info    : [ 80%] Meshing curve 644 (Line)
Info    : [ 80%] Meshing curve 645 (Circle)
Info    : [ 80%] Meshing curve 646 (Line)
Info    : [ 80%] Meshing curve 647 (Circle)
Info    : [ 80%] Meshing curve 648 (Line)
Info    : [ 80%] Meshing curve 649 (Circle)
Info    : [ 80%] Meshing curve 650 (Line)
Info    : [ 80%] Meshing curve 651 (Circle)
Info    : [ 80%] Meshing curve 652 (Line)
Info    : [ 80%] Meshing curve 653 (Circle)
Info    : [ 80%] Meshing curve 654 (Line)
Info    : [ 80%] Meshing curve 655 (Circle)
Info    : [ 80%] Meshing curve 656 (Line)
Info    : [ 80%] Meshing curve 657 (Circle)
Info    : [ 80%] Meshing curve 658 (Line)
Info    : [ 80%] Meshing curve 659 (Circle)
Info    : [ 80%] Meshing curve 660 (Line)
Info    : [ 80%] Meshing curve 661 (Circle)
Info    : [ 80%] Meshing curve 662 (Line)
Info    : [ 80%] Meshing curve 663 (Circle)
Info    : [ 80%] Meshing curve 664 (Line)
Info    : [ 80%] Meshing curve 665 (Circle)
Info    : [ 80%] Meshing curve 666 (Line)
Info    : [ 80%] Meshing curve 667 (Circle)
Info    : [ 80%] Meshing curve 668 (Line)
Info    : [ 80%] Meshing curve 669 (Circle)
Info    : [ 80%] Meshing curve 670 (Line)
Info    : [ 80%] Meshing curve 671 (Circle)
Info    : [ 80%] Meshing curve 672 (Line)
Info    : [ 80%] Meshing curve 673 (Circle)
Info    : [ 80%] Meshing curve 674 (Line)
Info    : [ 80%] Meshing curve 675 (Circle)
Info    : [ 80%] Meshing curve 676 (Line)
Info    : [ 80%] Meshing curve 677 (Circle)
Info    : [ 80%] Meshing curve 678 (Circle)
Info    : [ 80%] Meshing curve 679 (Line)
Info    : [ 80%] Meshing curve 680 (Circle)
Info    : [ 80%] Meshing curve 681 (Line)
Info    : [ 80%] Meshing curve 682 (Circle)
Info    : [ 80%] Meshing curve 683 (Line)
Info    : [ 80%] Meshing curve 684 (Circle)
Info    : [ 80%] Meshing curve 685 (Line)
Info    : [ 80%] Meshing curve 686 (Circle)
Info    : [ 80%] Meshing curve 687 (Line)
Info    : [ 80%] Meshing curve 688 (Circle)
Info    : [ 80%] Meshing curve 689 (Line)
Info    : [ 80%] Meshing curve 690 (Circle)
Info    : [ 80%] Meshing curve 691 (Line)
Info    : [ 80%] Meshing curve 692 (Circle)
Info    : [ 80%] Meshing curve 693 (Line)
Info    : [ 80%] Meshing curve 694 (Circle)
Info    : [ 80%] Meshing curve 695 (Line)
Info    : [ 80%] Meshing curve 696 (Circle)
Info    : [ 80%] Meshing curve 697 (Line)
Info    : [ 90%] Meshing curve 698 (Circle)
Info    : [ 90%] Meshing curve 699 (Line)
Info    : [ 90%] Meshing curve 700 (Circle)
Info    : [ 90%] Meshing curve 701 (Line)
Info    : [ 90%] Meshing curve 702 (Circle)
Info    : [ 90%] Meshing curve 703 (Line)
Info    : [ 90%] Meshing curve 704 (Circle)
Info    : [ 90%] Meshing curve 705 (Line)
Info    : [ 90%] Meshing curve 706 (Circle)
Info    : [ 90%] Meshing curve 707 (Line)
Info    : [ 90%] Meshing curve 708 (Circle)
Info    : [ 90%] Meshing curve 709 (Line)
Info    : [ 90%] Meshing curve 710 (Circle)
Info    : [ 90%] Meshing curve 711 (Line)
Info    : [ 90%] Meshing curve 712 (Circle)
Info    : [ 90%] Meshing curve 713 (Line)
Info    : [ 90%] Meshing curve 714 (Circle)
Info    : [ 90%] Meshing curve 715 (Line)
Info    : [ 90%] Meshing curve 716 (Circle)
Info    : [ 90%] Meshing curve 717 (Line)
Info    : [ 90%] Meshing curve 718 (Circle)
Info    : [ 90%] Meshing curve 719 (Circle)
Info    : [ 90%] Meshing curve 720 (Line)
Info    : [ 90%] Meshing curve 721 (Circle)
Info    : [ 90%] Meshing curve 722 (Line)
Info    : [ 90%] Meshing curve 723 (Circle)
Info    : [ 90%] Meshing curve 724 (Line)
Info    : [ 90%] Meshing curve 725 (Circle)
Info    : [ 90%] Meshing curve 726 (Line)
Info    : [ 90%] Meshing curve 727 (Circle)
Info    : [ 90%] Meshing curve 728 (Line)
Info    : [ 90%] Meshing curve 729 (Circle)
Info    : [ 90%] Meshing curve 730 (Line)
Info    : [ 90%] Meshing curve 731 (Circle)
Info    : [ 90%] Meshing curve 732 (Line)
Info    : [ 90%] Meshing curve 733 (Circle)
Info    : [ 90%] Meshing curve 734 (Line)
Info    : [ 90%] Meshing curve 735 (Circle)
Info    : [ 90%] Meshing curve 736 (Line)
Info    : [ 90%] Meshing curve 737 (Circle)
Info    : [ 90%] Meshing curve 738 (Line)
Info    : [ 90%] Meshing curve 739 (Circle)
Info    : [ 90%] Meshing curve 740 (Line)
Info    : [ 90%] Meshing curve 741 (Circle)
Info    : [ 90%] Meshing curve 742 (Line)
Info    : [ 90%] Meshing curve 743 (Circle)
Info    : [ 90%] Meshing curve 744 (Line)
Info    : [ 90%] Meshing curve 745 (Circle)
Info    : [ 90%] Meshing curve 746 (Line)
Info    : [ 90%] Meshing curve 747 (Circle)
Info    : [ 90%] Meshing curve 748 (Line)
Info    : [ 90%] Meshing curve 749 (Circle)
Info    : [ 90%] Meshing curve 750 (Line)
Info    : [ 90%] Meshing curve 751 (Circle)
Info    : [ 90%] Meshing curve 752 (Line)
Info    : [ 90%] Meshing curve 753 (Circle)
Info    : [ 90%] Meshing curve 754 (Line)
Info    : [ 90%] Meshing curve 755 (Circle)
Info    : [ 90%] Meshing curve 756 (Line)
Info    : [ 90%] Meshing curve 757 (Circle)
Info    : [ 90%] Meshing curve 758 (Line)
Info    : [ 90%] Meshing curve 759 (Circle)
Info    : [ 90%] Meshing curve 760 (Circle)
Info    : [ 90%] Meshing curve 761 (Line)
Info    : [ 90%] Meshing curve 762 (Circle)
Info    : [ 90%] Meshing curve 763 (Line)
Info    : [ 90%] Meshing curve 764 (Circle)
Info    : [ 90%] Meshing curve 765 (Line)
Info    : [ 90%] Meshing curve 766 (Circle)
Info    : [ 90%] Meshing curve 767 (Line)
Info    : [ 90%] Meshing curve 768 (Circle)
Info    : [ 90%] Meshing curve 769 (Line)
Info    : [ 90%] Meshing curve 770 (Circle)
Info    : [ 90%] Meshing curve 771 (Line)
Info    : [ 90%] Meshing curve 772 (Circle)
Info    : [ 90%] Meshing curve 773 (Line)
Info    : [ 90%] Meshing curve 774 (Circle)
Info    : [ 90%] Meshing curve 775 (Line)
Info    : [ 90%] Meshing curve 776 (Circle)
Info    : [ 90%] Meshing curve 777 (Line)
Info    : [ 90%] Meshing curve 778 (Circle)
Info    : [ 90%] Meshing curve 779 (Line)
Info    : [100%] Meshing curve 780 (Circle)
Info    : [100%] Meshing curve 781 (Line)
Info    : [100%] Meshing curve 782 (Circle)
Info    : [100%] Meshing curve 783 (Line)
Info    : [100%] Meshing curve 784 (Circle)
Info    : [100%] Meshing curve 785 (Line)
Info    : [100%] Meshing curve 786 (Circle)
Info    : [100%] Meshing curve 787 (Line)
Info    : [100%] Meshing curve 788 (Circle)
Info    : [100%] Meshing curve 789 (Line)
Info    : [100%] Meshing curve 790 (Circle)
Info    : [100%] Meshing curve 791 (Line)
Info    : [100%] Meshing curve 792 (Circle)
Info    : [100%] Meshing curve 793 (Line)
Info    : [100%] Meshing curve 794 (Circle)
Info    : [100%] Meshing curve 795 (Line)
Info    : [100%] Meshing curve 796 (Circle)
Info    : [100%] Meshing curve 797 (Line)
Info    : [100%] Meshing curve 798 (Circle)
Info    : [100%] Meshing curve 799 (Line)
Info    : [100%] Meshing curve 800 (Circle)
Info    : [100%] Meshing curve 801 (Circle)
Info    : [100%] Meshing curve 802 (Line)
Info    : [100%] Meshing curve 803 (Circle)
Info    : [100%] Meshing curve 804 (Line)
Info    : [100%] Meshing curve 805 (Circle)
Info    : [100%] Meshing curve 806 (Line)
Info    : [100%] Meshing curve 807 (Circle)
Info    : [100%] Meshing curve 808 (Line)
Info    : [100%] Meshing curve 809 (Circle)
Info    : [100%] Meshing curve 810 (Line)
Info    : [100%] Meshing curve 811 (Circle)
Info    : [100%] Meshing curve 812 (Line)
Info    : [100%] Meshing curve 813 (Circle)
Info    : [100%] Meshing curve 814 (Line)
Info    : [100%] Meshing curve 815 (Circle)
Info    : [100%] Meshing curve 816 (Line)
Info    : [100%] Meshing curve 817 (Circle)
Info    : [100%] Meshing curve 818 (Line)
Info    : [100%] Meshing curve 819 (Circle)
Info    : [100%] Meshing curve 820 (Line)
Info    : [100%] Meshing curve 821 (Circle)
Info    : [100%] Meshing curve 822 (Line)
Info    : [100%] Meshing curve 823 (Circle)
Info    : [100%] Meshing curve 824 (Line)
Info    : [100%] Meshing curve 825 (Circle)
Info    : [100%] Meshing curve 826 (Line)
Info    : [100%] Meshing curve 827 (Circle)
Info    : [100%] Meshing curve 828 (Line)
Info    : [100%] Meshing curve 829 (Circle)
Info    : [100%] Meshing curve 830 (Line)
Info    : [100%] Meshing curve 831 (Circle)
Info    : [100%] Meshing curve 832 (Line)
Info    : [100%] Meshing curve 833 (Circle)
Info    : [100%] Meshing curve 834 (Line)
Info    : [100%] Meshing curve 835 (Circle)
Info    : [100%] Meshing curve 836 (Line)
Info    : [100%] Meshing curve 837 (Circle)
Info    : [100%] Meshing curve 838 (Line)
Info    : [100%] Meshing curve 839 (Circle)
Info    : [100%] Meshing curve 840 (Line)
Info    : [100%] Meshing curve 841 (Circle)
Info    : [100%] Meshing curve 842 (Circle)
Info    : [100%] Meshing curve 843 (Circle)
Info    : [100%] Meshing curve 844 (Circle)
Info    : [100%] Meshing curve 845 (Circle)
Info    : [100%] Meshing curve 846 (Circle)
Info    : [100%] Meshing curve 847 (Circle)
Info    : [100%] Meshing curve 848 (Circle)
Info    : [100%] Meshing curve 849 (Circle)
Info    : [100%] Meshing curve 850 (Circle)
Info    : [100%] Meshing curve 851 (Circle)
Info    : [100%] Meshing curve 852 (Circle)
Info    : [100%] Meshing curve 853 (Circle)
Info    : [100%] Meshing curve 854 (Circle)
Info    : [100%] Meshing curve 855 (Circle)
Info    : [100%] Meshing curve 856 (Circle)
Info    : [100%] Meshing curve 857 (Circle)
Info    : [100%] Meshing curve 858 (Circle)
Info    : [100%] Meshing curve 859 (Circle)
Info    : [100%] Meshing curve 860 (Circle)
Info    : [100%] Meshing curve 861 (Circle)
Info    : Done meshing 1D (Wall 0.0412191s, CPU 0.042784s)
Info    : 2200 nodes 3020 elements

We can check that the mesh is embedded in a 3D geometrical space while being of 1D topology.

gdim = domain.geometry.dim
tdim = domain.topology.dim

print(f"Geometrical dimension = {gdim}")
print(f"Topological dimension = {tdim}")

Geometrical dimension = 3
Topological dimension = 1

We use the ufl.Jacobian function to compute the transformation Jacobian between a reference element (interval here) and the current element. In our case, the Jacobian is of shape (3,1). Transforming it into a vector of unit length will give us the local tangent vector \(\boldsymbol{t}\).

dx_dX = ufl.Jacobian(domain)[:, 0]
t = dx_dX / ufl.sqrt(ufl.inner(dx_dX, dx_dX))

We now compute the section local axis. As mentioned earlier, \(\boldsymbol{a}_1\) will be perpendicular to \(\boldsymbol{t}\) and the vertical direction \(\boldsymbol{e}_z=(0,0,1)\). After normalization, \(\boldsymbol{a}_2\) is built by taking the cross product between \(\boldsymbol{t}\) and \(\boldsymbol{a}_1\), \(\boldsymbol{a}_2\) will therefore belong to the plane made by \(\boldsymbol{t}\) and the vertical direction.

ez = ufl.as_vector([0, 0, 1])
a1 = ufl.cross(t, ez)
a1 /= ufl.sqrt(ufl.dot(a1, a1))
a2 = ufl.cross(t, a1)
a2 /= ufl.sqrt(ufl.dot(a2, a2))

We now define the material and geometrical constants which will be used in the constitutive relation. We consider the case of a rectangular cross-section of width \(b\) and height \(h\) in directions \(\boldsymbol{a}_1\) and \(\boldsymbol{a}_2\). The bending inertia will therefore be \(I_1 = bh^3/12\) and \(I_2=hb^3/12\). The torsional inertia is \(J=\beta hb^3\) with \(\beta\approx 0.26\) for \(h=3b\). Finally, the shear areas are approximated by \(S_1=S_2=\kappa S\) with \(\kappa=5/6\).

thick = fem.Constant(domain, 0.3)
width = thick / 3
E = fem.Constant(domain, 70e3)
nu = fem.Constant(domain, 0.3)
G = E / 2 / (1 + nu)
rho = fem.Constant(domain, 2.7e-3)
g = fem.Constant(domain, 9.81)

S = thick * width
ES = E * S
EI1 = E * width * thick**3 / 12
EI2 = E * width**3 * thick / 12
GJ = G * 0.26 * thick * width**3
kappa = fem.Constant(domain, 5.0 / 6.0)
GS1 = kappa * G * S
GS2 = kappa * G * S

We now consider a mixed \(\mathbb{P}_1/\mathbb{P}_1\)-Lagrange interpolation for the displacement and rotation fields. The variational form is built using a function generalized_strains giving the vector of six generalized strains as well as a function generalized_stresses which computes the dot product of the strains with the above-mentioned constitutive matrix (diagonal here). Note that since the 1D beams are embedded in an ambient 3D space, we define a tangential gradient operator tgrad by taking the dot product with the local tangent vector \(t\).

Finally, similarly to Reissner-Mindlin plates, shear-locking issues might arise in the thin beam limit. To avoid this, reduced integration is performed on the shear part \(Q_1\widehat{\gamma}_1+Q_2\widehat{\gamma}_2\) of the variational form using a one-point rule.

Ue = basix.ufl.element("P", domain.basix_cell(), 1, shape=(gdim,))
W = fem.functionspace(domain, basix.ufl.mixed_element([Ue, Ue]))


u_, theta_ = ufl.TestFunctions(W)
du, dtheta = ufl.TrialFunctions(W)


def tgrad(u):
    return ufl.dot(ufl.grad(u), t)


def generalized_strains(u, theta):
    return ufl.as_vector(
        [
            ufl.dot(tgrad(u), t),
            ufl.dot(tgrad(u), a1) - ufl.dot(theta, a2),
            ufl.dot(tgrad(u), a2) + ufl.dot(theta, a1),
            ufl.dot(tgrad(theta), t),
            ufl.dot(tgrad(theta), a1),
            ufl.dot(tgrad(theta), a2),
        ]
    )


def generalized_stresses(u, theta):
    return ufl.dot(
        ufl.diag(ufl.as_vector([ES, GS1, GS2, GJ, EI1, EI2])),
        generalized_strains(u, theta),
    )


Sig = generalized_stresses(du, dtheta)
Eps_ = generalized_strains(u_, theta_)

dx_shear = ufl.dx(metadata={"quadrature_degree": 0})
k_form = (
    sum([Sig[i] * Eps_[i] * ufl.dx for i in [0, 3, 4, 5]])
    + (Sig[1] * Eps_[1] + Sig[2] * Eps_[2]) * dx_shear
)
l_form = -rho * S * g * u_[2] * ufl.dx

Clamped boundary conditions are considered at the bottom \(z=0\) level and the linear problem is finally solved.

def bottom(x):
    return np.isclose(x[2], 0.0)


Vu, _ = W.sub(0).collapse()
Vt, _ = W.sub(1).collapse()
u_dofs = fem.locate_dofs_geometrical((W.sub(0), Vu), bottom)
theta_dofs = fem.locate_dofs_geometrical((W.sub(1), Vt), bottom)
u0 = fem.Function(Vu)
theta0 = fem.Function(Vt)
bcs = [
    fem.dirichletbc(u0, u_dofs, W.sub(0)),
    fem.dirichletbc(theta0, theta_dofs, W.sub(1)),
]

w = fem.Function(W, name="Generalized_displacement")

problem = fem.petsc.LinearProblem(
    k_form,
    l_form,
    u=w,
    bcs=bcs,
    petsc_options={"ksp_type": "preonly", "pc_type": "lu"},
    petsc_options_prefix="beam_",
)
problem.solve()

Coefficient(FunctionSpace(Mesh(blocked element (Basix element (P, interval, 1, equispaced, unset, False, float64, []), (3,)), 0), mixed element (blocked element (Basix element (P, interval, 1, gll_warped, unset, False, float64, []), (3,)), blocked element (Basix element (P, interval, 1, gll_warped, unset, False, float64, []), (3,)))), 2)

We finally plot the deflected shape and the bending moment distribution using pyvista.

Hide code cell source

 
import pyvista
from dolfinx import plot

u = w.sub(0).collapse()

plotter = pyvista.Plotter()
topology, cell_types, geometry = plot.vtk_mesh(Vu)
grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
grid.point_data["Deflection"] = u.x.array.reshape(-1, 3)
warped = grid.warp_by_vector("Deflection", factor=500.0)
plotter.add_mesh(grid, show_edges=True, color="k", line_width=1, opacity=0.5)
plotter.add_mesh(warped, show_edges=True, line_width=5)
plotter.show()


M = fem.Function(Vu, name="Bending moments (M1,M2)")
Sig = generalized_stresses(w.sub(0), w.sub(1))
M_exp = fem.Expression(
    ufl.as_vector([Sig[4], Sig[5], 0]), Vu.element.interpolation_points
)
M.interpolate(M_exp)

grid.point_data["Bending_moments"] = M.x.array.reshape(-1, 3)
plotter = pyvista.Plotter()
plotter.add_mesh(grid, show_edges=True, line_width=5, scalars="Bending_moments")
plotter.show()

Attention

  • The \(\mathbb{P}_1/\mathbb{P}_1\) discretization choice made here might be too poor for achieving a decent error with a moderate mesh size. Increasing the order of interpolation, e.g. for the displacement \(\boldsymbol{u}\) is obviously possible but the integration degree for the reduced integration measure dx_shear must be adapted to this case.

  • To add concentrated forces, one can use the dS measure with appropriate tags or follow the strategy hinted in Linear elastic 2D truss .