Coverage for src/beamme/mesh_creation_functions/beam

1# The MIT License (MIT)

5# Permission is hereby granted, free of charge, to any person obtaining a copy

6# of this software and associated documentation files (the "Software"), to deal

7# in the Software without restriction, including without limitation the rights

8# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell

9# copies of the Software, and to permit persons to whom the Software is

10# furnished to do so, subject to the following conditions:

11#

12# The above copyright notice and this permission notice shall be included in

13# all copies or substantial portions of the Software.

14#

15# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR

16# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,

17# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE

18# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER

19# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,

20# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN

21# THE SOFTWARE.

22"""Generic function for beam creation."""

24from typing import Callable as _Callable

25from typing import Dict as _Dict

26from typing import List as _List

27from typing import Optional as _Optional

28from typing import Tuple as _Tuple

29from typing import Type as _Type

30from typing import Union as _Union

32import numpy as _np

34from beamme.core.conf import bme as _bme

35from beamme.core.element_beam import Beam as _Beam

36from beamme.core.geometry_set import GeometryName as _GeometryName

37from beamme.core.geometry_set import GeometrySet as _GeometrySet

38from beamme.core.material import MaterialBeamBase as _MaterialBeamBase

39from beamme.core.mesh import Mesh as _Mesh

40from beamme.core.node import NodeCosserat as _NodeCosserat

41from beamme.utils.nodes import get_single_node as _get_single_node

44def create_beam_mesh_generic(

45 mesh: _Mesh,

46 *,

47 beam_class: _Type[_Beam],

48 material: _MaterialBeamBase,

49 function_generator: _Callable,

50 interval: _Tuple[float, float],

51 n_el: _Optional[int] = None,

52 l_el: _Optional[float] = None,

53 interval_length: _Optional[float] = None,

54 set_nodal_arc_length: bool = False,

55 nodal_arc_length_offset: _Optional[float] = None,

56 node_positions_of_elements: _Optional[_List[float]] = None,

57 start_node: _Optional[_Union[_NodeCosserat, _GeometrySet]] = None,

58 end_node: _Optional[_Union[_NodeCosserat, _GeometrySet]] = None,

59 close_beam: bool = False,

60 vtk_cell_data: _Optional[_Dict[str, _Tuple]] = None,

61) -> _GeometryName:

62 """Generic beam creation function.

64 Remark for given start and/or end nodes:

65 If the rotation does not match, but the tangent vector is the same,

66 the created beams triads are rotated so the physical problem stays

67 the same (for axi-symmetric beam cross-sections) but the nodes can

68 be reused.

70 Args:

71 mesh:

72 Mesh that the created beam(s) should be added to.

73 beam_class:

74 Class of beam that will be used for this line.

75 material:

76 Material for this line.

77 function_generator:

78 The function_generator has to take two variables, point_a and

79 point_b (both within the interval) and return a function(xi) that

80 calculates the position and rotation along the beam, where

81 point_a -> xi = -1 and point_b -> xi = 1.

83 Usually, the Jacobian of the returned position field should be a unit

84 vector. Otherwise, the nodes may be spaced in an undesired way. All

85 standard mesh creation functions fulfill this property.

86 interval:

87 Start and end values for interval that will be used to create the

88 beam.

89 n_el:

90 Number of equally spaced beam elements along the line. Defaults to 1.

91 Mutually exclusive with l_el

92 l_el:

93 Desired length of beam elements. This requires the option interval_length

94 to be set. Mutually exclusive with n_el. Be aware, that this length

95 might not be achieved, if the elements are warped after they are

96 created.

97 interval_length:

98 Total length of the interval. Is required when the option l_el is given.

99 set_nodal_arc_length:

100 Flag if the arc length along the beam filament is set in the created

101 nodes. It is ensured that the arc length is consistent with possible

102 given start/end nodes.

103 nodal_arc_length_offset:

104 Offset of the stored nodal arc length w.r.t. to the one generated by

105 the function. Defaults to 0, the arc length set in the start node, or

106 the arc length in the end node minus total length (such that the arc

107 length at the end node matches).

108 node_positions_of_elements:

109 A list of normalized positions (within [0,1] and in ascending order)

110 that define the boundaries of beam elements along the created curve.

111 The given values will be mapped to the actual `interval` given as an

112 argument to this function. These values specify where elements start

113 and end, additional internal nodes (such as midpoints in higher-order

114 elements) may be placed automatically.

115 start_node:

116 Node to use as the first node for this line. Use this if the line

117 is connected to other lines (angles have to be the same, otherwise

118 connections should be used). If a geometry set is given, it can

119 contain one, and one node only.

120 end_node:

121 If this is a Node or GeometrySet, the last node of the created beam

122 is set to that node.

123 If it is True the created beam is closed within itself.

124 close_beam:

125 If it is True the created beam is closed within itself (mutually

126 exclusive with end_node).

127 vtk_cell_data:

128 With this argument, a vtk cell data can be set for the elements

129 created within this function. This can be used to check which

130 elements are created by which function.

131

132 Returns:

133 Geometry sets with the 'start' and 'end' node of the curve. Also a 'line' set

134 with all nodes of the curve.

135 """

136

137 # Check for mutually exclusive parameters

138 n_given_arguments = sum(

139 1

140 for argument in [n_el, l_el, node_positions_of_elements]

141 if argument is not None

142 )

143 if n_given_arguments == 0:

144 # No arguments were given, use a single element per default

145 n_el = 1

146 elif n_given_arguments > 1:

147 raise ValueError(

148 'The arguments "n_el", "l_el" and "node_positions_of_elements" are mutually exclusive'

149 )

150

151 if close_beam and end_node is not None:

152 raise ValueError(

153 'The arguments "close_beam" and "end_node" are mutually exclusive'

154 )

155

156 if set_nodal_arc_length:

157 if close_beam:

158 raise ValueError(

159 "The flags 'set_nodal_arc_length' and 'close_beam' are mutually exclusive."

160 )

161 elif nodal_arc_length_offset is not None:

162 raise ValueError(

163 'Providing the argument "nodal_arc_length_offset" without setting '

164 '"set_nodal_arc_length" to True does not make sense.'

165 )

166

167 # Cases where we have equally spaced elements

168 if n_el is not None or l_el is not None:

169 if l_el is not None:

170 # Calculate the number of elements in case a desired element length is provided

171 if interval_length is None:

172 raise ValueError(

173 'The parameter "l_el" requires "interval_length" to be set.'

174 )

175 n_el = max([1, round(interval_length / l_el)])

176 elif n_el is None:

177 raise ValueError("n_el should not be None at this point")

178

179 node_positions_of_elements = [i_node / n_el for i_node in range(n_el + 1)]

180 # A list for the element node positions was provided

181 elif node_positions_of_elements is not None:

182 # Check that the given positions are in ascending order and start with 0 and end with 1

183 for index, value, name in zip([0, -1], [0, 1], ["First", "Last"]):

184 if not _np.isclose(

185 value,

186 node_positions_of_elements[index],

187 atol=1e-12,

188 rtol=0.0,

189 ):

190 raise ValueError(

191 f"{name} entry of node_positions_of_elements must be {value}, got {node_positions_of_elements[index]}"

192 )

193 if not all(

194 x < y

195 for x, y in zip(node_positions_of_elements, node_positions_of_elements[1:])

196 ):

197 raise ValueError(

198 f"The given node_positions_of_elements must be in ascending order. Got {node_positions_of_elements}"

199 )

200

201 # Get the scale the node positions to the interval coordinates

202 interval_node_positions_of_elements = interval[0] + (

203 interval[1] - interval[0]

204 ) * _np.asarray(node_positions_of_elements)

205

206 # We need to make sure we have the number of elements for the case a given end node

207 n_el = len(interval_node_positions_of_elements) - 1

208

209 # Make sure the material is in the mesh.

210 mesh.add_material(material)

211

212 # List with nodes and elements that will be added in the creation of

213 # this beam.

214 elements = []

215 nodes = []

216

217 def check_given_node(node):

218 """Check that the given node is already in the mesh."""

219 if node not in mesh.nodes:

220 raise ValueError("The given node is not in the current mesh")

221

222 def get_relative_twist(rotation_node, rotation_function, name):

223 """Check if the rotation at a node and the one returned by the function

224 match.

225

226 If not, check if the first basis vector of the triads is the

227 same. If that is the case, a simple relative twist can be

228 applied to ensure that the triad field is continuous. This

229 relative twist can lead to issues if the beam cross-section is

230 not double symmetric.

231 """

232

233 if rotation_node == rotation_function:

234 return None

235 elif not _bme.allow_beam_rotation:

236 # The settings do not allow for a rotation of the beam

237 raise ValueError(

238 f"Given rotation of the {name} node does not match with given rotation function!"

239 )

240 else:

241 # Evaluate the relative rotation

242 # First check if the first basis vector is the same

243 relative_basis_1 = rotation_node.inv() * rotation_function * [1, 0, 0]

244 if _np.linalg.norm(relative_basis_1 - [1, 0, 0]) < _bme.eps_quaternion:

245 # Calculate the relative rotation

246 return rotation_function.inv() * rotation_node

247 else:

248 raise ValueError(

249 f"The tangent of the {name} node does not match with the given function!"

250 )

251

252 # Position and rotation at the start and end of the interval

253 function_over_whole_interval = function_generator(*interval)

254 relative_twist_start = None

255 relative_twist_end = None

256

257 # If a start node is given, set this as the first node for this beam.

258 if start_node is not None:

259 start_node = _get_single_node(start_node)

260 nodes = [start_node]

261 check_given_node(start_node)

262 _, start_rotation, _ = function_over_whole_interval(-1.0)

263 relative_twist_start = get_relative_twist(

264 start_node.rotation, start_rotation, "start"

265 )

266

267 # If an end node is given, check what behavior is wanted.

268 if end_node is not None:

269 end_node = _get_single_node(end_node)

270 check_given_node(end_node)

271 _, end_rotation, _ = function_over_whole_interval(1.0)

272 relative_twist_end = get_relative_twist(end_node.rotation, end_rotation, "end")

273

274 # Get the start value for the arc length functionality

275 if set_nodal_arc_length:

276 if nodal_arc_length_offset is not None:

277 # Let's use the given value, if it does not match with possible given

278 # start or end nodes, the check in the create beam function will detect

279 # that.

280 pass

281 elif start_node is not None and start_node.arc_length is not None:

282 nodal_arc_length_offset = start_node.arc_length

283 elif end_node is not None and end_node.arc_length is not None:

284 if interval_length is None:

285 raise ValueError(

286 "The end node has an arc length but the interval_length is not "

287 "given. This is not supported."

288 )

289 nodal_arc_length_offset = end_node.arc_length - interval_length

290 else:

291 # Default value

292 nodal_arc_length_offset = 0.0

293

294 # Check if a relative twist has to be applied

295 if relative_twist_start is not None and relative_twist_end is not None:

296 if relative_twist_start == relative_twist_end:

297 relative_twist = relative_twist_start

298 else:

299 raise ValueError(

300 "The relative twist required for the start and end node do not match"

301 )

302 elif relative_twist_start is not None:

303 relative_twist = relative_twist_start

304 elif relative_twist_end is not None:

305 relative_twist = relative_twist_end

306 else:

307 relative_twist = None

308

309 # Create the beams.

310 for i_el in range(n_el):

311 # If the beam is closed with itself, set the end node to be the

312 # first node of the beam. This is done when the second element is

313 # created, as the first node already exists here.

314 if i_el == 1 and close_beam:

315 end_node = nodes[0]

316

317 # Get the function to create this beam element.

318 function = function_generator(

319 interval_node_positions_of_elements[i_el],

320 interval_node_positions_of_elements[i_el + 1],

321 )

322

323 # Set the start node for the created beam.

324 if start_node is not None or i_el > 0:

325 first_node = nodes[-1]

326 else:

327 first_node = None

328

329 # If an end node is given, set this one for the last element.

330 if end_node is not None and i_el == n_el - 1:

331 last_node = end_node

332 else:

333 last_node = None

334

335 element = beam_class(material=material)

336 elements.append(element)

337 nodes.extend(

338 element.create_beam(

339 function,

340 start_node=first_node,

341 end_node=last_node,

342 relative_twist=relative_twist,

343 set_nodal_arc_length=set_nodal_arc_length,

344 nodal_arc_length_offset=nodal_arc_length_offset,

345 )

346 )

347

348 # Set vtk cell data on created elements.

349 if vtk_cell_data is not None:

350 for data_name, data_value in vtk_cell_data.items():

351 for element in elements:

352 if data_name in element.vtk_cell_data.keys():

353 raise KeyError(

354 'The cell data "{}" already exists!'.format(data_name)

355 )

356 element.vtk_cell_data[data_name] = data_value

357

358 # Add items to the mesh

359 mesh.elements.extend(elements)

360 if start_node is None:

361 mesh.nodes.extend(nodes)

362 else:

363 mesh.nodes.extend(nodes[1:])

364

365 # Set the last node of the beam.

366 if end_node is None:

367 end_node = nodes[-1]

368

369 # Set the nodes that are at the beginning and end of line (for search

370 # of overlapping points)

371 nodes[0].is_end_node = True

372 end_node.is_end_node = True

373

374 # Create geometry sets that will be returned.

375 return_set = _GeometryName()

376 return_set["start"] = _GeometrySet(nodes[0])

377 return_set["end"] = _GeometrySet(end_node)

378 return_set["line"] = _GeometrySet(elements)

379

380 return return_set

Coverage for src/beamme/mesh_creation_functions/beam_generic.py: 94%

125 statements