Coverage for src/sdynpy/core/sdynpy_geometry.py: 15%
2880 statements
« prev ^ index » next coverage.py v7.16.1, created at 2026-09-18 21:32 +0000
1"""
2Objects and procedures to handle operations on test or model geometries
4This module defines a Geometry object as well as all of the subcomponents of
5a geometry object: nodes, elements, tracelines and coordinate system. Geometry
6plotting is also handled in this module.
7"""
9"""
10Copyright 2022 National Technology & Engineering Solutions of Sandia,
11LLC (NTESS). Under the terms of Contract DE-NA0003525 with NTESS, the U.S.
12Government retains certain rights in this software.
14This program is free software: you can redistribute it and/or modify
15it under the terms of the GNU General Public License as published by
16the Free Software Foundation, either version 3 of the License, or
17(at your option) any later version.
19This program is distributed in the hope that it will be useful,
20but WITHOUT ANY WARRANTY; without even the implied warranty of
21MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22GNU General Public License for more details.
24You should have received a copy of the GNU General Public License
25along with this program. If not, see <https://www.gnu.org/licenses/>.
26"""
28import copy
29import numpy as np
30import pyvista as pv
31import pyvistaqt as pvqt
32import pyqtgraph as pg
33import vtk
34from qtpy import uic, QtCore
35from qtpy.QtGui import QKeySequence
36try:
37 from qtpy.QtGui import QAction
38except ImportError:
39 from qtpy.QtWidgets import QAction
40from qtpy.QtWidgets import QApplication, QMainWindow, QSizePolicy, QFileDialog, QMenu
41from .sdynpy_array import SdynpyArray
42from .sdynpy_colors import colormap, coord_colormap
43from .sdynpy_coordinate import CoordinateArray, coordinate_array, from_nodelist
44from ..signal_processing.sdynpy_rotation import R, lstsq_rigid_transform
45from ..signal_processing.sdynpy_camera import point_on_pixel
46from ..signal_processing.sdynpy_massprops import rigid_body_shape_scaling
47from .sdynpy_matrix import matrix
48import time
49import os
50from PIL import Image
51from scipy.spatial import Delaunay
52from scipy.linalg import block_diag
53import pandas as pd
54import openpyxl
55import warnings
56from scipy.io import savemat as scipy_savemat
58try:
59 repr(pv.GPUInfo())
60 IGNORE_PLOTS = False
61 # Set Qt Preferences
63 # Adding this to help out non-interactive consoles
64 os.environ["QT_AUTO_SCREEN_SCALE_FACTOR"] = "1"
65 if hasattr(QtCore.Qt, "AA_EnableHighDpiScaling"): # PyQt5 only
66 QApplication.setAttribute(QtCore.Qt.AA_EnableHighDpiScaling)
68 if hasattr(QtCore.Qt, "AA_UseHighDpiPixmaps"): # PyQt5 only
69 QApplication.setAttribute(QtCore.Qt.AA_UseHighDpiPixmaps)
71 QApplication.setHighDpiScaleFactorRoundingPolicy(
72 QtCore.Qt.HighDpiScaleFactorRoundingPolicy.PassThrough
73 )
74 if QApplication.instance() is None:
75 app = QApplication([''])
77except RuntimeError:
78 IGNORE_PLOTS = True
79 print('No GPU Found, Geometry plotting will not work!')
81# Enumerations
83_cs_type_strings = {0: 'Cartesian',
84 1: 'Polar',
85 2: 'Spherical'}
87_element_types = {
88 11: 'Rod',
89 21: 'Linear beam',
90 22: 'Tapered beam',
91 23: 'Curved beam',
92 24: 'Parabolic beam',
93 31: 'Straight pipe',
94 32: 'Curved pipe',
95 41: 'Plane Stress Linear Triangle',
96 42: 'Plane Stress Parabolic Triangle',
97 43: 'Plane Stress Cubic Triangle',
98 44: 'Plane Stress Linear Quadrilateral',
99 45: 'Plane Stress Parabolic Quadrilateral',
100 46: 'Plane Strain Cubic Quadrilateral',
101 51: 'Plane Strain Linear Triangle',
102 52: 'Plane Strain Parabolic Triangle',
103 53: 'Plane Strain Cubic Triangle',
104 54: 'Plane Strain Linear Quadrilateral',
105 55: 'Plane Strain Parabolic Quadrilateral',
106 56: 'Plane Strain Cubic Quadrilateral',
107 61: 'Plate Linear Triangle',
108 62: 'Plate Parabolic Triangle',
109 63: 'Plate Cubic Triangle',
110 64: 'Plate Linear Quadrilateral',
111 65: 'Plate Parabolic Quadrilateral',
112 66: 'Plate Cubic Quadrilateral',
113 71: 'Membrane Linear Quadrilateral',
114 72: 'Membrane Parabolic Triangle',
115 73: 'Membrane Cubic Triangle',
116 74: 'Membrane Linear Triangle',
117 75: 'Membrane Parabolic Quadrilateral',
118 76: 'Membrane Cubic Quadrilateral',
119 81: 'Axisymetric Solid Linear Triangle',
120 82: 'Axisymetric Solid Parabolic Triangle',
121 84: 'Axisymetric Solid Linear Quadrilateral',
122 85: 'Axisymetric Solid Parabolic Quadrilateral',
123 91: 'Thin Shell Linear Triangle',
124 92: 'Thin Shell Parabolic Triangle',
125 93: 'Thin Shell Cubic Triangle',
126 94: 'Thin Shell Linear Quadrilateral',
127 95: 'Thin Shell Parabolic Quadrilateral',
128 96: 'Thin Shell Cubic Quadrilateral',
129 101: 'Thick Shell Linear Wedge',
130 102: 'Thick Shell Parabolic Wedge',
131 103: 'Thick Shell Cubic Wedge',
132 104: 'Thick Shell Linear Brick',
133 105: 'Thick Shell Parabolic Brick',
134 106: 'Thick Shell Cubic Brick',
135 111: 'Solid Linear Tetrahedron',
136 112: 'Solid Linear Wedge',
137 113: 'Solid Parabolic Wedge',
138 114: 'Solid Cubic Wedge',
139 115: 'Solid Linear Brick',
140 116: 'Solid Parabolic Brick',
141 117: 'Solid Cubic Brick',
142 118: 'Solid Parabolic Tetrahedron',
143 121: 'Rigid Bar',
144 122: 'Rigid Element',
145 136: 'Node To Node Translational Spring',
146 137: 'Node To Node Rotational Spring',
147 138: 'Node To Ground Translational Spring',
148 139: 'Node To Ground Rotational Spring',
149 141: 'Node To Node Damper',
150 142: 'Node To Gound Damper',
151 151: 'Node To Node Gap',
152 152: 'Node To Ground Gap',
153 161: 'Lumped Mass',
154 171: 'Axisymetric Linear Shell',
155 172: 'Axisymetric Parabolic Shell',
156 181: 'Constraint',
157 191: 'Plastic Cold Runner',
158 192: 'Plastic Hot Runner',
159 193: 'Plastic Water Line',
160 194: 'Plastic Fountain',
161 195: 'Plastic Baffle',
162 196: 'Plastic Rod Heater',
163 201: 'Linear node-to-node interface',
164 202: 'Linear edge-to-edge interface',
165 203: 'Parabolic edge-to-edge interface',
166 204: 'Linear face-to-face interface',
167 208: 'Parabolic face-to-face interface',
168 212: 'Linear axisymmetric interface',
169 213: 'Parabolic axisymmetric interface',
170 221: 'Linear rigid surface',
171 222: 'Parabolic rigid surface',
172 231: 'Axisymetric linear rigid surface',
173 232: 'Axisymentric parabolic rigid surface'}
175_exodus_elem_type_map = {'triangle':41,
176 'quadrilateral':44,
177 'hexahedral':115,
178 'tet':111,
179 'tetrahedral':111,
180 'sphere':161,
181 'spring':136,
182 'bar':21,
183 'bar2':21,
184 'bar3':23,
185 'beam':21,
186 'beam2':21,
187 'beam3':23,
188 'truss':21,
189 'truss2':21,
190 'truss3':23,
191 'quad':94,
192 'quad4':94,
193 'quad5':94,
194 'quad8':95,
195 'quad9':95,
196 'shell':94,
197 'shell4':94,
198 'shell8':95,
199 'shell9':95,
200 'tri':91,
201 'tri3':91,
202 'tri6':92,
203 'tri7':92,
204 'trishell':91,
205 'trishell3':91,
206 'trishell6':92,
207 'trishell7':92,
208 'hex':115,
209 'hex8':115,
210 'hex9':116,
211 'hex20':116,
212 'hex27':117,
213 'tetra':111,
214 'tetra4':111,
215 'tetra8':118,
216 'tetra10':118,
217 'tetra14':118,
218 'wedge':112,
219 'hexshell':116,
220 # 'pyramid'
221 # 'pyramid5'
222 # 'pyramid8'
223 # 'pyramid13'
224 # 'pyramid18'
225 # Add new elements here
226 }
228_beam_elem_types = [11, 21, 22, 23, 24, 31, 32, 121]
229_face_element_types = [v for v in range(41, 107)]
230_solid_element_types = [v for v in range(111, 119)]
231_vtk_element_map = {
232 111: vtk.VTK_TETRA,
233 112: vtk.VTK_WEDGE,
234 113: vtk.VTK_QUADRATIC_WEDGE,
235 115: vtk.VTK_HEXAHEDRON,
236 116: vtk.VTK_QUADRATIC_HEXAHEDRON,
237 118: vtk.VTK_QUADRATIC_TETRA,
238}
240_vtk_connectivity_reorder = {
241 vtk.VTK_QUADRATIC_HEXAHEDRON: [0, 1, 2, 3, 4, 5, 6,
242 7, 8, 9, 10, 11, 16, 17, 18, 19, 12, 13, 14, 15]
243}
245MAX_NUMBER_REPR = 100
248class GeometryPlotter(pvqt.BackgroundPlotter):
249 """Class used to plot geometry
251 This class is essentially identical to PyVista's BackgroundPlotter;
252 however a small amount of additional functionality is added."""
254 def __init__(self, *args, **kwargs):
255 super().__init__(*args, **kwargs)
256 for action in self.main_menu.actions():
257 if action.menu() and action.text().lower() == 'view':
258 view_menu = action.menu()
259 break
260 for action in view_menu.actions():
261 if action.menu() and action.text().lower() == 'orientation marker':
262 orien_menu = action.menu()
263 break
264 orien_menu.addSeparator()
265 orien_menu.addAction('Show Origin', self.show_origin_marker)
266 orien_menu.addAction('Hide Origin', self.hide_origin_marker)
267 self._origin_marker = None
269 def show_origin_marker(self, arrow_scale=0.05):
270 if self._origin_marker is None:
271 bbox_diagonal = self.length
272 arrow_factor = bbox_diagonal * arrow_scale
273 self._origin_marker = pv.tools.create_axes_marker(label_size=np.array([0.25, 0.1])*arrow_scale)
274 self._origin_marker.SetTotalLength(arrow_factor, arrow_factor, arrow_factor)
275 self.renderer.AddActor(self._origin_marker)
276 self.renderer.Modified()
278 def hide_origin_marker(self):
279 if self._origin_marker is not None:
280 self.renderer.remove_actor(self._origin_marker)
281 self.renderer.Modified()
283 def save_rotation_animation(self, filename=None, frames=20, frame_rate=20,
284 individual_images=False):
285 """
286 Saves an animation of the rotating geometry to a file
288 Parameters
289 ----------
290 filename : str, optional
291 Path to the file being saved. The default is None, which returns
292 the images as an array rather than saves them to a file
293 frames : int, optional
294 Number of frames in the animation. The default is 20.
295 frame_rate : float, optional
296 Number of frames per second if the animation is saved to a GIF.
297 The default is 20.
298 individual_images : bool, optional
299 Boolean to specify whether the images are saved as individual PNG
300 files or a single GIF. The default is False.
302 Returns
303 -------
304 imgs : np.ndarray or None
305 Returns array of images if filename is None, otherwise returns None
307 """
308 imgs = []
309 phases = np.linspace(0, 2 * np.pi, frames, endpoint=False)
310 focus = np.array(self.camera.focal_point)
311 distance = np.array(self.camera.distance)
312 initial_view_from = np.array(self.camera.position) - focus
313 view_up = self.camera.up
314 for phase in phases:
315 # Compute a rotation matrix
316 rotmat = R(view_up, phase)
317 view_from = (rotmat @ initial_view_from[:, np.newaxis]).squeeze()
318 self.camera.position = np.array(focus) + distance * \
319 np.array(view_from) / np.linalg.norm(view_from)
320 self.camera.focal_point = focus
321 self.camera.up = view_up
322 self.render()
323 # Force GUI Update to ensure the OS doesn't delay events
324 QApplication.processEvents()
325 imgs.append(self.screenshot())
326 if filename is None:
327 return imgs
328 else:
329 if individual_images:
330 for i, img in enumerate(imgs):
331 fname = os.path.splitext(filename)
332 fname = fname[0] + '-{:d}'.format(i) + fname[1]
333 Image.fromarray(img).save(fname)
334 else:
335 imgs = [Image.fromarray(img).convert(
336 'P', palette=Image.ADAPTIVE, colors=256) for img in imgs]
337 imgs[0].save(fp=filename, format='GIF', append_images=imgs[1:],
338 save_all=True, duration=int(1 / frame_rate * 1000), loop=0)
341class TransientPlotter(GeometryPlotter):
342 """Class used to plot transient deformations"""
344 def __init__(self, geometry, displacement_data, displacement_scale=1.0,
345 frames_per_second=20,
346 undeformed_opacity=0.0, deformed_opacity=1.0, plot_kwargs={},
347 transformation_shapes=None, num_curves=50,
348 show: bool = True,
349 app=None,
350 window_size=None,
351 off_screen=None,
352 allow_quit_keypress=True,
353 toolbar=True,
354 menu_bar=True,
355 editor=False,
356 update_app_icon=None,
357 **kwargs):
358 """
359 Create a TransientPlotter object to plot displacements over time
361 Parameters
362 ----------
363 geometry : Geometry
364 Geometry on which the displacements will be plotted
365 displacement_data : TimeHistoryArray
366 Transient displacement data that will be applied
367 displacement_scale : float, optional
368 Scale factor applied to displacements. The default is 1.0.
369 frames_per_second : float, optional
370 Number of time steps to plot per second while the displacement is
371 animating. Default is 20.
372 undeformed_opacity : float, optional
373 Opacity of the undeformed geometry. The default is 0.0, or
374 completely transparent.
375 deformed_opacity : float, optional
376 Opacity of the deformed geometry. The default is 1.0, or completely
377 opaque.
378 plot_kwargs : dict, optional
379 Keyword arguments passed to the Geometry.plot function
380 transformation_shapes : ShapeArray
381 Shape matrix that will be used to expand the data. Must be the
382 same size as the `displacement_data`
383 num_curves : int, optional
384 Maximum number of curves to plot on the time selector. Default is
385 50
386 show : bool, optional
387 Show the plotting window. If ``False``, show this window by
388 running ``show()``. The default is True.
389 app : QApplication, optional
390 Creates a `QApplication` if left as `None`. The default is None.
391 window_size : list of int, optional
392 Window size in pixels. Defaults to ``[1024, 768]``
393 off_screen : TYPE, optional
394 Renders off screen when True. Useful for automated
395 screenshots or debug testing. The default is None.
396 allow_quit_keypress : bool, optional
397 Allow user to exit by pressing ``"q"``. The default is True.
398 toolbar : bool, optional
399 If True, display the default camera toolbar. Defaults to True.
400 menu_bar : bool, optional
401 If True, display the default main menu. Defaults to True.
402 editor : TYPE, optional
403 If True, display the VTK object editor. Defaults to False.
404 update_app_icon : bool, optional
405 If True, update_app_icon will be called automatically to update the
406 Qt app icon based on the current rendering output. If None, the
407 logo of PyVista will be used. If False, no icon will be set.
408 Defaults to None. The default is None.
409 title : str, optional
410 Title of plotting window.
411 multi_samples : int, optional
412 The number of multi-samples used to mitigate aliasing. 4 is a
413 good default but 8 will have better results with a potential
414 impact on performance.
415 line_smoothing : bool, optional
416 If True, enable line smothing
417 point_smoothing : bool, optional
418 If True, enable point smothing
419 polygon_smoothing : bool, optional
420 If True, enable polygon smothing
421 auto_update : float, bool, optional
422 Automatic update rate in seconds. Useful for automatically
423 updating the render window when actors are change without
424 being automatically ``Modified``. If set to ``True``, update
425 rate will be 1 second.
427 Returns
428 -------
429 None.
431 """
432 self.shape_select_toolbar: pvqt.plotting.QToolBar = None
433 self.animation_control_toolbar: pvqt.plotting.QToolBar = None
434 self.loop_action = None
435 self.plot_abscissa_action = None
436 self.timer: pvqt.plotting.QTimer = None
437 self.frames_per_second = frames_per_second
438 self.animation_speed = 1.0
439 self.last_time = time.time()
441 super(TransientPlotter, self).__init__(show, app, None if window_size is None else tuple(window_size),
442 off_screen,
443 allow_quit_keypress, toolbar,
444 menu_bar, editor, update_app_icon,
445 **kwargs)
447 self.text = self.add_text('', font_size=10)
448 self.loop_action.setChecked(False)
449 self.plot_abscissa_action.setChecked(True)
451 # Get the layout so we can add a new widget to it
452 layout = self.frame.layout()
453 self.time_selector_plot = pg.PlotWidget()
454 layout.addWidget(self.time_selector_plot)
455 sp = self.time_selector_plot.sizePolicy()
456 sp.setVerticalPolicy(QSizePolicy.Minimum)
457 self.time_selector_plot.setSizePolicy(sp)
458 self.time_selector_plot.sizeHint = lambda: QtCore.QSize(1006, 80)
459 # self.time_selector_plot.getPlotItem().hideAxis('left')
460 # self.time_selector_plot.getPlotItem().axes['left']['item'].setStyle(showValues=False)
462 # Now we need to plot the data in the plot
463 if num_curves < displacement_data.size:
464 indices = np.linspace(0, displacement_data.size - 1, num_curves).astype(int)
465 else:
466 indices = np.arange(displacement_data.size)
467 for data in displacement_data.flatten()[indices]:
468 self.time_selector_plot.plot(data.abscissa, data.ordinate)
469 # Plot the envelope
470 max_envelope = displacement_data.flatten().ordinate.max(axis=0)
471 min_envelope = displacement_data.flatten().ordinate.min(axis=0)
472 self.time_selector_plot.plot(displacement_data.flatten(
473 )[0].abscissa, max_envelope, pen=pg.mkPen(color='k', width=2))
474 self.time_selector_plot.plot(displacement_data.flatten(
475 )[0].abscissa, min_envelope, pen=pg.mkPen(color='k', width=2))
477 # Add an infinite line for a cursor
478 self.time_selector = pg.InfiniteLine(
479 0, movable=True, bounds=[displacement_data.abscissa.min(),
480 displacement_data.abscissa.max()])
481 self.time_selector_plot.addItem(self.time_selector)
482 self.time_selector.sigPositionChanged.connect(self.modify_abscissa)
484 # Set initial values
485 self.displacement_scale = displacement_scale
486 self.abscissa_index = 0
487 self.displacement_data = displacement_data.flatten()
488 self.transformation_shapes = None if transformation_shapes is None else transformation_shapes.flatten()
489 self.plot_geometry = geometry
490 self.num_abscissa = self.displacement_data.num_elements
491 self.abscissa = self.displacement_data[0].abscissa
492 self.abscissa_sig_figs = int(np.ceil(-np.log10(np.min(np.diff(self.abscissa))))) + 1
493 self.abscissa_format_spec = 'Time: {{:0.{:}f}}'.format(
494 self.abscissa_sig_figs if self.abscissa_sig_figs > 0 else 0)
496 # Plot the geometry
497 plotter, self.face_mesh_undeformed, self.point_mesh_undeformed, self.solid_mesh_undeformed = geometry.plot(plotter=self,
498 opacity=undeformed_opacity,
499 **plot_kwargs)
501 plotter, self.face_mesh_deformed, self.point_mesh_deformed, self.solid_mesh_deformed = geometry.plot(plotter=self,
502 opacity=deformed_opacity,
503 **plot_kwargs)
505 # Compute the node indices up front to make things faster when switching steps
506 if self.transformation_shapes is None:
507 coordinates = np.unique(self.displacement_data.coordinate)
508 else:
509 coordinates = np.unique(self.transformation_shapes.coordinate)
510 coordinates = coordinates[
511 np.isin(coordinates.node, self.plot_geometry.node.id)
512 & (np.abs(coordinates.direction) > 0)
513 & (np.abs(coordinates.direction) < 4)
514 ]
515 global_deflections = self.plot_geometry.global_deflection(coordinates).T
516 if self.transformation_shapes is None:
517 shape_array = self.displacement_data[coordinates[:, np.newaxis]].ordinate
518 else:
519 shape_array = self.transformation_shapes[coordinates].T
520 # Compute global displacements for each node
521 self.node_deflections = np.zeros((3 * self.plot_geometry.node.size, shape_array.shape[-1]))
522 for node_index, node in self.plot_geometry.node.ndenumerate():
523 coord_indices = coordinates.node == node.id
524 self.node_deflections[node_index[0] * 3:(node_index[0] + 1) *
525 3] = global_deflections[:, coord_indices] @ shape_array[coord_indices]
526 if self.transformation_shapes is None:
527 self.weighting = None
528 else:
529 self.weighting = self.displacement_data.ordinate
530 self.update_displacement()
532 def modify_abscissa(self):
533 abscissa = self.time_selector.value()
534 self.abscissa_index = np.argmin(np.abs(self.abscissa - abscissa))
535 self.update_displacement(update_selector=False)
537 def update_displacement(self, update_selector=True):
538 index = int(round(self.abscissa_index))
539 abscissa = self.abscissa[index]
540 if self.plot_abscissa_action.isChecked():
541 self.text.SetText(2, self.abscissa_format_spec.format(abscissa))
542 else:
543 self.text.SetText(2, '')
544 if self.weighting is None:
545 displacement = self.node_deflections[:, index] * self.displacement_scale
546 else:
547 displacement = self.node_deflections @ self.weighting[:,
548 np.newaxis, index] * self.displacement_scale
549 for mesh_deformed, mesh_undeformed in (
550 (self.face_mesh_deformed, self.face_mesh_undeformed),
551 (self.point_mesh_deformed, self.point_mesh_undeformed),
552 (self.solid_mesh_deformed, self.solid_mesh_undeformed)
553 ):
554 if mesh_deformed is not None:
555 mesh_deformed.points = (mesh_undeformed.points +
556 displacement.reshape(-1, 3) * self.displacement_scale
557 )
558 # Update the time selector
559 if update_selector:
560 self.time_selector.blockSignals(True)
561 self.time_selector.setValue(abscissa)
562 self.time_selector.blockSignals(False)
563 self.render()
565 def add_toolbars(self) -> None:
566 """
567 Adds toolbars to the BackgroundPlotter
569 """
570 super().add_toolbars()
571 self.shape_select_toolbar = self.app_window.addToolBar(
572 "Shape Selector"
573 )
575 self._add_action(
576 self.shape_select_toolbar, '|<', self.go_to_first_step
577 )
578 self._add_action(
579 self.shape_select_toolbar, '<', self.go_to_previous_step
580 )
581 self._add_action(
582 self.shape_select_toolbar, '>', self.go_to_next_step
583 )
584 self._add_action(
585 self.shape_select_toolbar, '>|', self.go_to_last_step
586 )
588 self.animation_control_toolbar = self.app_window.addToolBar(
589 "Animation Control"
590 )
592 self._add_action(
593 self.animation_control_toolbar, '< Play', self.play_animation_reverse
594 )
595 self._add_action(
596 self.animation_control_toolbar, 'Stop', self.stop_animation
597 )
598 self._add_action(
599 self.animation_control_toolbar, 'Play >', self.play_animation
600 )
602 def add_menu_bar(self) -> None:
603 """
604 Adds the menu bar to the BackgroundPlotter
606 """
607 super().add_menu_bar()
609 shape_menu = self.main_menu.addMenu("Shape")
610 scaling_menu = shape_menu.addMenu('Scaling')
611 scaling_menu.addAction('1/4x', self.select_scaling_0p25)
612 scaling_menu.addAction('1/2x', self.select_scaling_0p5)
613 self.scale_0p8_action = QAction('0.8x')
614 self.scale_0p8_action.triggered.connect(self.select_scaling_0p8)
615 self.scale_0p8_action.setShortcut(QKeySequence(','))
616 scaling_menu.addAction(self.scale_0p8_action)
617 self.scale_reset_action = QAction('Reset')
618 self.scale_reset_action.triggered.connect(self.select_scaling_1)
619 self.scale_reset_action.setShortcut(QKeySequence('/'))
620 scaling_menu.addAction(self.scale_reset_action)
621 self.scale_1p25_action = QAction('1.25x')
622 self.scale_1p25_action.triggered.connect(self.select_scaling_1p25)
623 self.scale_1p25_action.setShortcut(QKeySequence('.'))
624 scaling_menu.addAction(self.scale_1p25_action)
625 scaling_menu.addAction('2x', self.select_scaling_2p0)
626 scaling_menu.addAction('4x', self.select_scaling_4p0)
627 speed_menu = shape_menu.addMenu('Animation Speed')
628 self.speed_0p8_action = QAction('0.8x')
629 self.speed_0p8_action.triggered.connect(self.select_speed_0p8)
630 self.speed_0p8_action.setShortcut(QKeySequence(';'))
631 speed_menu.addAction(self.speed_0p8_action)
632 speed_menu.addAction('Reset', self.select_speed_1)
633 self.speed_1p25_action = QAction('1.25x')
634 self.speed_1p25_action.triggered.connect(self.select_speed_1p25)
635 self.speed_1p25_action.setShortcut(QKeySequence("'"))
636 speed_menu.addAction(self.speed_1p25_action)
637 shape_menu.addSeparator()
638 self.plot_abscissa_action = pvqt.plotting.QAction('Show Abscissa', checkable=True)
639 shape_menu.addAction(self.plot_abscissa_action)
640 self.loop_action = pvqt.plotting.QAction('Loop Animation', checkable=True)
641 shape_menu.addAction(self.loop_action)
643 def go_to_first_step(self):
644 self.stop_animation()
645 self.abscissa_index = 0
646 self.update_displacement()
648 def go_to_previous_step(self):
649 self.stop_animation()
650 self.abscissa_index -= 1
651 if self.abscissa_index < 0:
652 if self.loop_action.isChecked():
653 self.abscissa_index = self.num_abscissa - 1
654 else:
655 self.abscissa_index = 0
656 self.update_displacement()
658 def go_to_next_step(self):
659 self.stop_animation()
660 self.abscissa_index += 1
661 if self.abscissa_index > self.num_abscissa - 1:
662 if self.loop_action.isChecked():
663 self.abscissa_index = 0
664 else:
665 self.abscissa_index = self.num_abscissa - 1
666 self.update_displacement()
668 def go_to_last_step(self):
669 self.stop_animation()
670 self.abscissa_index = self.num_abscissa - 1
671 self.update_displacement()
673 def select_scaling_0p25(self):
674 """
675 Adjusts the current shape scaling by 0.25x
676 """
677 self.displacement_scale *= 0.25
678 self.last_time = time.time()
679 self.update_displacement()
681 def select_scaling_0p5(self):
682 """
683 Adjusts the current shape scaling by 0.5x
684 """
685 self.displacement_scale *= 0.5
686 self.last_time = time.time()
687 self.update_displacement()
689 def select_scaling_0p8(self):
690 """
691 Adjusts the current shape scaling by 0.8x
692 """
693 self.displacement_scale *= 0.8
694 self.last_time = time.time()
695 self.update_displacement()
697 def select_scaling_1(self):
698 """
699 Resets shape scaling to 1
700 """
701 self.displacement_scale = 1.0
702 self.last_time = time.time()
703 self.update_displacement()
705 def select_scaling_1p25(self):
706 """
707 Adjusts the current shape scaling by 1.25x
708 """
709 self.displacement_scale *= 1.25
710 self.last_time = time.time()
711 self.update_displacement()
713 def select_scaling_2p0(self):
714 """
715 Adjusts the current shape scaling by 2x
716 """
717 self.displacement_scale *= 2.0
718 self.last_time = time.time()
719 self.update_displacement()
721 def select_scaling_4p0(self):
722 """
723 Adjusts the current shape scaling by 4x
724 """
725 self.displacement_scale *= 4.0
726 self.last_time = time.time()
727 self.update_displacement()
729 def select_speed_1(self):
730 """
731 Resets animation speed to default quantities
732 """
733 self.animation_speed = 1.0
735 def select_speed_0p8(self):
736 """
737 Adjusts the current animation speed by 0.8
738 """
739 self.animation_speed *= 0.8
741 def select_speed_1p25(self):
742 """
743 Adjusts the current animation speed by 1.25
744 """
745 self.animation_speed *= 1.25
747 def play_animation(self):
748 self.timer = pvqt.plotting.QTimer()
749 self.timer.timeout.connect(self.update_abscissa)
750 self.timer.start(10)
751 self.last_time = time.time()
753 def play_animation_reverse(self):
754 self.timer = pvqt.plotting.QTimer()
755 self.timer.timeout.connect(self.update_abscissa_reverse)
756 self.timer.start(10)
757 self.last_time = time.time()
759 def stop_animation(self):
760 if self.timer is not None:
761 self.timer.stop()
762 self.timer = None
764 def update_abscissa(self):
765 now = time.time()
766 time_elapsed = now - self.last_time
767 index_change = self.frames_per_second * self.animation_speed * time_elapsed
768 self.last_time = now
769 self.abscissa_index += index_change
770 if self.abscissa_index > self.num_abscissa - 1:
771 if self.loop_action.isChecked():
772 self.abscissa_index = 0
773 else:
774 self.stop_animation()
775 self.abscissa_index = self.num_abscissa - 1
776 self.update_displacement()
778 def update_abscissa_reverse(self):
779 now = time.time()
780 time_elapsed = now - self.last_time
781 index_change = self.frames_per_second * self.animation_speed * time_elapsed
782 self.last_time = now
783 self.abscissa_index -= index_change
784 if self.abscissa_index < 0:
785 if self.loop_action.isChecked():
786 self.abscissa_index = self.num_abscissa - 1
787 else:
788 self.stop_animation()
789 self.abscissa_index = 0
790 self.update_displacement()
792 def save_animation(self, filename=None, frame_rate=20,
793 individual_images=False, start_time=None, stop_time=None,
794 step=None):
795 """
796 Saves the current shape animation to a file
798 Parameters
799 ----------
800 filename : str, optional
801 Path to the file being saved. The default is None, which returns
802 the images as an array rather than saves them to a file
803 frame_rate : float, optional
804 Number of frames per second if the animation is saved to a GIF.
805 The default is 20.
806 individual_images : bool, optional
807 Boolean to specify whether the images are saved as individual PNG
808 files or a single GIF. The default is False.
809 start_time : float, optional
810 Time value at which the animation will start. Default is first
811 time step
812 stop_time : float, optional
813 Time value at which the animation will end. Default is last time
814 step
815 step : int, optional
816 Only render every ``step`` time steps. Default is to render all
817 timesteps
818 Returns
819 -------
820 imgs : np.ndarray or None
821 Returns array of images if filename is None, otherwise returns None
823 """
824 self.stop_animation()
825 imgs = []
826 indices = np.arange(self.abscissa.size)
827 if start_time is not None:
828 indices = indices[self.abscissa[indices] >= start_time]
829 if stop_time is not None:
830 indices = indices[self.abscissa[indices] <= stop_time]
831 if step is not None:
832 indices = indices[::step]
833 for index in indices:
834 self.abscissa_index = index
835 self.update_displacement()
836 # Force GUI Update to ensure the OS doesn't delay events
837 QApplication.processEvents()
838 imgs.append(self.screenshot())
839 if filename is None:
840 return imgs
841 else:
842 if individual_images:
843 for i, img in enumerate(imgs):
844 fname = os.path.splitext(filename)
845 fname = fname[0] + '-{:d}'.format(i) + fname[1]
846 Image.fromarray(img).save(fname)
847 else:
848 imgs = [Image.fromarray(img).convert(
849 'P', palette=Image.ADAPTIVE, colors=256) for img in imgs]
850 imgs[0].save(fp=filename, format='GIF', append_images=imgs[1:],
851 save_all=True, duration=int(1 / frame_rate * 1000), loop=0)
853 def _close(self):
854 self.stop_animation()
855 super(TransientPlotter, self)._close()
858class ShapePlotter(GeometryPlotter):
859 """Class used to plot animated shapes"""
861 def __init__(self, geometry, shapes, plot_kwargs,
862 background_plotter_kwargs={'auto_update': 0.01},
863 undeformed_opacity=0.25, starting_scale=1.0,
864 deformed_opacity=1.0, shape_name='Mode',
865 show_damping=True):
866 """
867 Create a Shape Plotter object to plot shapes
869 Parameters
870 ----------
871 geometry : Geometry
872 Geometry on which the shapes will be plotted
873 shapes : ShapeArray
874 Shapes to plot on the geometry
875 plot_kwargs : dict
876 Keyword arguments passed to the Geometry.plot function
877 background_plotter_kwargs : dict, optional
878 Keyword arguments passed to the BackgroundPlotter constructor.
879 The default is {'auto_update':0.01}.
880 undeformed_opacity : float, optional
881 Opacity of the undeformed geometry. Set to zero for no undeformed
882 geometry. The default is 0.25.
883 starting_scale : float, optional
884 Starting scale of the shapes on the plot. The default is 1.0.
885 shape_name : str, optional
886 Name that will be used in the plotter to describe the shape. The
887 default is "Mode"
888 show_damping : bool, optional
889 Boolean that specifies whether the damping should be displayed in
890 the comment.
891 """
892 self.shape_select_toolbar: pvqt.plotting.QToolBar = None
893 self.animation_control_toolbar: pvqt.plotting.QToolBar = None
894 self.signed_amplitude_action = None
895 self.complex_action = None
896 self.real_action = None
897 self.imag_action = None
898 self.loop_action = None
899 self.plot_comments_action = None
900 self.save_animation_action = None
901 self.display_undeformed = None
902 self.geometry = geometry
903 self.shapes = shapes.ravel()
904 self.node_displacements = None
905 self.current_shape = 0
906 self.displacement_scale = 1.0
907 self.displacement_scale_modification = starting_scale
908 self.animation_speed = 1.0
909 self.last_time = time.time()
910 self.phase = 0.0
911 self.timer: pvqt.plotting.QTimer = None
912 self.shape_name = shape_name
913 self.show_damping = show_damping
914 self.shape_comment_table = None
916 super().__init__(**background_plotter_kwargs)
918 self.text = self.add_text('', font_size=10)
920 self.complex_action.setChecked(True)
921 self.loop_action.setChecked(True)
922 self.plot_comments_action.setChecked(True)
924 plotter, self.face_mesh_undeformed, self.point_mesh_undeformed, self.solid_mesh_undeformed = geometry.plot(plotter=self,
925 opacity=undeformed_opacity,
926 **plot_kwargs)
928 plotter, self.face_mesh_deformed, self.point_mesh_deformed, self.solid_mesh_deformed = geometry.plot(plotter=self,
929 opacity=deformed_opacity,
930 **plot_kwargs)
932 # Compute the node indices up front to make things faster when switching shapes
933 coordinates = self.shapes[self.current_shape].coordinate
934 nodes = coordinates.node
935 # Only do coordinates that are in the geometry node array
936 self.indices = np.isin(coordinates.node, self.geometry.node.id)
937 coordinates = coordinates[self.indices]
938 nodes = nodes[self.indices]
939 direction = coordinates.direction
940 self.coordinate_node_index_map = {}
941 for coordinate_index, node_id in np.ndenumerate(nodes):
942 if direction[coordinate_index] == 0 or abs(direction[coordinate_index]) > 3:
943 continue
944 if node_id not in self.coordinate_node_index_map:
945 self.coordinate_node_index_map[node_id] = []
946 self.coordinate_node_index_map[node_id].append(coordinate_index[0])
947 local_deformations = coordinates.local_direction()
948 ordered_nodes = self.geometry.node(nodes)
949 coordinate_systems = self.geometry.coordinate_system(ordered_nodes.disp_cs)
950 points = global_coord(self.geometry.coordinate_system(
951 ordered_nodes.def_cs), ordered_nodes.coordinate)
952 self.global_deflections = global_deflection(coordinate_systems, local_deformations, points)
954 self.compute_displacements()
955 self.show_comment()
956 self.render()
957 self.update_shape()
958 self.update_shape_mode(0)
959 QApplication.processEvents()
961 def compute_displacements(self, compute_scale=True) -> np.ndarray:
962 """
963 Computes displacements to apply to the geometry
965 Parameters
966 ----------
967 compute_scale : bool, optional
968 Renormalize the displacement scaling. The default is True.
970 """
971 shape_displacements = self.shapes[self.current_shape].shape_matrix[..., self.indices]
972 node_displacements = np.zeros(self.geometry.node.shape + (3,),
973 dtype=shape_displacements.dtype)
974 for node_index, node in self.geometry.node.ndenumerate():
975 try:
976 node_indices = self.coordinate_node_index_map[node.id]
977 except KeyError:
978 print('Node {:} not found in shape array'.format(node.id))
979 node_displacements[node_index] = 0
980 continue
981 node_deflection_scales = shape_displacements[node_indices]
982 node_deflection_directions = self.global_deflections[node_indices]
983 node_displacements[node_index] += np.sum(node_deflection_scales[:, np.newaxis] *
984 node_deflection_directions, axis=0)
985 self.node_displacements = node_displacements.reshape(-1, 3)
986 # Compute maximum displacement
987 if compute_scale:
988 max_displacement = np.max(np.linalg.norm(self.node_displacements, axis=-1))
989 global_coords = self.geometry.global_node_coordinate()
990 bbox_diagonal = np.linalg.norm(
991 np.max(global_coords, axis=0) - np.min(global_coords, axis=0))
992 self.displacement_scale = 0.05 * bbox_diagonal / \
993 max_displacement # Set to 5% of bounding box dimension
995 def update_shape_mode(self, phase=None):
996 """
997 Updates the mode that is being plotted
999 Parameters
1000 ----------
1001 phase : float, optional
1002 Sets the current phase of the shape. The default is None, which
1003 computes the new phase based on the time elapsed.
1005 """
1006 if phase is None:
1007 now = time.time()
1008 time_elapsed = now - self.last_time
1009 phase_change = 2 * np.pi / self.animation_speed * time_elapsed
1010 self.last_time = now
1011 self.phase += phase_change
1012 else:
1013 self.phase = phase
1014 if self.phase > 2 * np.pi:
1015 self.phase -= 2 * np.pi
1016 if not self.loop_action.isChecked():
1017 self.stop_animation()
1018 # print('Updating Phase: {:}'.format(self.phase))
1019 if self.complex_action.isChecked():
1020 # print('Complex')
1021 deformation = np.real(np.exp(self.phase * 1j) * self.node_displacements)
1022 elif self.real_action.isChecked():
1023 # print('Real')
1024 deformation = np.cos(self.phase) * np.real(self.node_displacements)
1025 elif self.imag_action.isChecked():
1026 # print('Imag')
1027 deformation = np.cos(self.phase) * np.imag(self.node_displacements)
1028 elif self.signed_amplitude_action.isChecked():
1029 # print('Signed Amplitude')
1030 imag = np.imag(self.node_displacements).flatten()[:, np.newaxis]
1031 real = np.real(self.node_displacements).flatten()[:, np.newaxis]
1032 best_fit = np.linalg.lstsq(real, imag)[0][0, 0]
1033 best_fit_phase = np.arctan(best_fit)
1034 rectified_displacements = self.node_displacements / np.exp(1j * best_fit_phase)
1035 deformation = np.cos(self.phase) * np.real(rectified_displacements)
1036 for mesh_deformed, mesh_undeformed in (
1037 (self.face_mesh_deformed, self.face_mesh_undeformed),
1038 (self.point_mesh_deformed, self.point_mesh_undeformed),
1039 (self.solid_mesh_deformed, self.solid_mesh_undeformed)
1040 ):
1041 if mesh_deformed is not None:
1042 mesh_deformed.points = (mesh_undeformed.points +
1043 deformation * self.displacement_scale * self.displacement_scale_modification
1044 )
1045 self.render()
1047 def update_shape(self):
1048 """
1049 Updates the shape that is being plotted
1050 """
1051 self.update_shape_mode()
1053 def save_animation_from_action(self):
1054 filename, file_filter = QFileDialog.getSaveFileName(
1055 self, 'Select File to Save Animation', filter='GIF (*.gif)')
1056 if filename == '':
1057 return
1058 self.save_animation(filename)
1060 def save_animation(self, filename=None, frames=20, frame_rate=20,
1061 individual_images=False):
1062 """
1063 Saves the current shape animation to a file
1065 Parameters
1066 ----------
1067 filename : str, optional
1068 Path to the file being saved. The default is None, which returns
1069 the images as an array rather than saves them to a file
1070 frames : int, optional
1071 Number of frames in the animation. The default is 20.
1072 frame_rate : float, optional
1073 Number of frames per second if the animation is saved to a GIF.
1074 The default is 20.
1075 individual_images : bool, optional
1076 Boolean to specify whether the images are saved as individual PNG
1077 files or a single GIF. The default is False.
1079 Returns
1080 -------
1081 imgs : np.ndarray or None
1082 Returns array of images if filename is None, otherwise returns None
1084 """
1085 self.stop_animation()
1086 imgs = []
1087 phases = np.linspace(0, 2 * np.pi, frames, endpoint=False)
1088 for phase in phases:
1089 self.update_shape_mode(phase)
1090 # Force GUI Update to ensure the OS doesn't delay events
1091 QApplication.processEvents()
1092 imgs.append(self.screenshot())
1093 if filename is None:
1094 return imgs
1095 else:
1096 if individual_images:
1097 for i, img in enumerate(imgs):
1098 fname = os.path.splitext(filename)
1099 fname = fname[0] + '-{:d}'.format(i) + fname[1]
1100 Image.fromarray(img).save(fname)
1101 else:
1102 imgs = [Image.fromarray(img).convert(
1103 'P', palette=Image.ADAPTIVE, colors=256) for img in imgs]
1104 imgs[0].save(fp=filename, format='GIF', append_images=imgs[1:],
1105 save_all=True, duration=int(1 / frame_rate * 1000), loop=0)
1107 def save_animation_all_shapes(self, filename_base='Shape_{:}.gif',
1108 frames=20, frame_rate=20,
1109 individual_images=False):
1110 """
1111 Helper function to easily save animations for all shapes
1113 Parameters
1114 ----------
1115 filename_base : str, optional
1116 Filename that will be passed to the format function to produce the
1117 actual file name for each shape. The default is 'Shape_{:}.gif'.
1118 frames : int, optional
1119 Number of frames in the animation. The default is 20.
1120 frame_rate : float, optional
1121 Number of frames per second if the animation is saved to a GIF.
1122 The default is 20.
1123 individual_images : bool, optional
1124 Boolean to specify whether the images are saved as individual PNG
1125 files or a single GIF. The default is False.
1127 """
1128 for shape in range(self.shapes.size):
1129 self.current_shape = shape
1130 self.compute_displacements()
1131 self.show_comment()
1132 self.save_animation(filename_base.format(shape + 1),
1133 frames, frame_rate, individual_images)
1135 def add_menu_bar(self) -> None:
1136 """
1137 Adds the menu bar to the BackgroundPlotter
1139 """
1140 super().add_menu_bar()
1142 # Get the file menu
1143 file_menu = [child for child in self.main_menu.children()
1144 if child.__class__.__name__ == 'QMenu'
1145 and child.title() == 'File'][0]
1146 self.save_animation_action = pvqt.plotting.QAction('Save Animation')
1147 self.save_animation_action.triggered.connect(self.save_animation_from_action)
1148 file_menu.insertAction(file_menu.actions()[1], self.save_animation_action)
1150 shape_menu = self.main_menu.addMenu("Shape")
1151 complex_display_menu = shape_menu.addMenu("Complex Display")
1152 self.signed_amplitude_action = pvqt.plotting.QAction('Signed Amplitude', checkable=True)
1153 self.signed_amplitude_action.triggered.connect(self.select_complex)
1154 complex_display_menu.addAction(self.signed_amplitude_action)
1155 self.complex_action = pvqt.plotting.QAction('Complex', checkable=True)
1156 self.complex_action.triggered.connect(self.select_complex)
1157 complex_display_menu.addAction(self.complex_action)
1158 self.real_action = pvqt.plotting.QAction('Real', checkable=True)
1159 self.real_action.triggered.connect(self.select_complex)
1160 complex_display_menu.addAction(self.real_action)
1161 self.imag_action = pvqt.plotting.QAction('Imaginary', checkable=True)
1162 self.imag_action.triggered.connect(self.select_complex)
1163 complex_display_menu.addAction(self.imag_action)
1164 scaling_menu = shape_menu.addMenu('Scaling')
1165 scaling_menu.addAction('1/4x', self.select_scaling_0p25)
1166 scaling_menu.addAction('1/2x', self.select_scaling_0p5)
1167 self.scale_0p8_action = QAction('0.8x')
1168 self.scale_0p8_action.triggered.connect(self.select_scaling_0p8)
1169 self.scale_0p8_action.setShortcut(QKeySequence(','))
1170 scaling_menu.addAction(self.scale_0p8_action)
1171 self.scale_reset_action = QAction('Reset')
1172 self.scale_reset_action.triggered.connect(self.select_scaling_1)
1173 self.scale_reset_action.setShortcut(QKeySequence('/'))
1174 scaling_menu.addAction(self.scale_reset_action)
1175 self.scale_1p25_action = QAction('1.25x')
1176 self.scale_1p25_action.triggered.connect(self.select_scaling_1p25)
1177 self.scale_1p25_action.setShortcut(QKeySequence('.'))
1178 scaling_menu.addAction(self.scale_1p25_action)
1179 scaling_menu.addAction('2x', self.select_scaling_2p0)
1180 scaling_menu.addAction('4x', self.select_scaling_4p0)
1181 speed_menu = shape_menu.addMenu('Animation Speed')
1182 self.speed_0p8_action = QAction('0.8x')
1183 self.speed_0p8_action.triggered.connect(self.select_speed_0p8)
1184 self.speed_0p8_action.setShortcut(QKeySequence(';'))
1185 speed_menu.addAction(self.speed_0p8_action)
1186 speed_menu.addAction('Reset', self.select_speed_1)
1187 self.speed_1p25_action = QAction('1.25x')
1188 self.speed_1p25_action.triggered.connect(self.select_speed_1p25)
1189 self.speed_1p25_action.setShortcut(QKeySequence("'"))
1190 speed_menu.addAction(self.speed_1p25_action)
1191 shape_menu.addSeparator()
1192 self.plot_comments_action = pvqt.plotting.QAction('Show Comments', checkable=True)
1193 self.plot_comments_action.triggered.connect(self.show_comment)
1194 shape_menu.addAction(self.plot_comments_action)
1195 self.loop_action = pvqt.plotting.QAction('Loop Animation', checkable=True)
1196 self.loop_action.triggered.connect(self.select_loop)
1197 shape_menu.addAction(self.loop_action)
1199 display_menu = self.main_menu.addMenu('Display')
1200 self.display_undeformed = pvqt.plotting.QAction('Show Undeformed', checkable=True)
1201 self.display_undeformed.triggered.connect(self.toggle_undeformed)
1202 display_menu.addAction(self.display_undeformed)
1204 def add_toolbars(self) -> None:
1205 """
1206 Adds toolbars to the BackgroundPlotter
1208 """
1209 super().add_toolbars()
1210 self.shape_select_toolbar = self.app_window.addToolBar(
1211 "Shape Selector"
1212 )
1214 self._add_action(
1215 self.shape_select_toolbar, '<<', self.prev_shape
1216 )
1217 self._add_action(
1218 self.shape_select_toolbar, '?', self.select_shape
1219 )
1220 self._add_action(
1221 self.shape_select_toolbar, '>>', self.next_shape
1222 )
1224 self.animation_control_toolbar = self.app_window.addToolBar(
1225 "Animation Control"
1226 )
1228 self._add_action(
1229 self.animation_control_toolbar, 'Play', self.play_animation
1230 )
1231 self._add_action(
1232 self.animation_control_toolbar, 'Stop', self.stop_animation
1233 )
1235 def play_animation(self):
1236 """Starts the animation playing"""
1237 self.timer = pvqt.plotting.QTimer()
1238 self.timer.timeout.connect(self.update_shape)
1239 self.timer.start(10)
1240 self.phase = 0.0
1241 self.last_time = time.time()
1243 def stop_animation(self):
1244 """Stops the animation from playing"""
1245 if self.timer is not None:
1246 self.timer.stop()
1247 self.timer = None
1249 def toggle_undeformed(self) -> None:
1250 pass
1252 def select_scaling_0p25(self):
1253 """
1254 Adjusts the current shape scaling by 0.25x
1255 """
1256 self.displacement_scale_modification *= 0.25
1257 self.phase = 0.0
1258 self.last_time = time.time()
1259 self.update_shape()
1261 def select_scaling_0p5(self):
1262 """
1263 Adjusts the current shape scaling by 0.5x
1264 """
1265 self.displacement_scale_modification *= 0.5
1266 self.phase = 0.0
1267 self.last_time = time.time()
1268 self.update_shape()
1270 def select_scaling_0p8(self):
1271 """
1272 Adjusts the current shape scaling by 0.8x
1273 """
1274 self.displacement_scale_modification *= 0.8
1275 self.phase = 0.0
1276 self.last_time = time.time()
1277 self.update_shape()
1279 def select_scaling_1(self):
1280 """
1281 Resets shape scaling to 1
1282 """
1283 self.displacement_scale_modification = 1.0
1284 self.phase = 0.0
1285 self.last_time = time.time()
1286 self.update_shape()
1288 def select_scaling_1p25(self):
1289 """
1290 Adjusts the current shape scaling by 1.25x
1291 """
1292 self.displacement_scale_modification *= 1.25
1293 self.phase = 0.0
1294 self.last_time = time.time()
1295 self.update_shape()
1297 def select_scaling_2p0(self):
1298 """
1299 Adjusts the current shape scaling by 2x
1300 """
1301 self.displacement_scale_modification *= 2.0
1302 self.phase = 0.0
1303 self.last_time = time.time()
1304 self.update_shape()
1306 def select_scaling_4p0(self):
1307 """
1308 Adjusts the current shape scaling by 4x
1309 """
1310 self.displacement_scale_modification *= 4.0
1311 self.phase = 0.0
1312 self.last_time = time.time()
1313 self.update_shape()
1315 def select_speed_1(self):
1316 """
1317 Resets animation speed to default quantities
1318 """
1319 self.animation_speed = 1.0
1321 def select_speed_0p8(self):
1322 """
1323 Adjusts the current animation speed by 0.8
1324 """
1325 self.animation_speed /= 0.8
1327 def select_speed_1p25(self):
1328 """
1329 Adjusts the current animation speed by 1.25
1330 """
1331 self.animation_speed /= 1.25
1333 def select_complex(self) -> None:
1334 """
1335 Adjusts how complex shapes are plotted
1336 """
1337 complex_actions = [
1338 self.signed_amplitude_action,
1339 self.complex_action,
1340 self.real_action,
1341 self.imag_action
1342 ]
1343 sender = self.sender()
1344 sender_index = complex_actions.index(sender)
1345 complex_actions.pop(sender_index)
1346 if not sender.isChecked():
1347 sender.setChecked(True)
1348 else:
1349 for action in complex_actions:
1350 action.setChecked(False)
1351 self.phase = 0.0
1352 self.last_time = time.time()
1353 self.update_shape()
1355 def select_loop(self) -> None:
1356 pass
1358 def next_shape(self) -> None:
1359 """
1360 Increases the index of the shape being plotted
1361 """
1362 self.current_shape += 1
1363 if self.current_shape >= self.shapes.size:
1364 self.current_shape = 0
1365 self.compute_displacements()
1366 self.phase = 0.0
1367 self.last_time = time.time()
1368 self.update_shape()
1369 self.show_comment()
1371 def prev_shape(self) -> None:
1372 """
1373 Decreases the index of the shape being plotted
1374 """
1375 self.current_shape -= 1
1376 if self.current_shape < 0:
1377 self.current_shape = self.shapes.size - 1
1378 self.compute_displacements()
1379 self.phase = 0.0
1380 self.last_time = time.time()
1381 self.update_shape()
1382 self.show_comment()
1384 def reset_shape(self) -> None:
1385 """
1386 Resets a shape if the value of the current shape is changed
1387 """
1388 self.compute_displacements()
1389 self.phase = 0.0
1390 self.last_time = time.time()
1391 self.update_shape()
1392 self.show_comment()
1394 def select_shape(self) -> None:
1395 from .sdynpy_shape import ShapeCommentTable
1396 self.shape_comment_table = ShapeCommentTable(self.shapes, self)
1398 def show_comment(self):
1399 """
1400 Shows or hides the shape comment
1401 """
1402 if self.plot_comments_action.isChecked():
1403 self.text.SetText(2, '\n'.join([
1404 '{:} {:}'.format(self.shape_name, self.current_shape + 1),
1405 ' Frequency: {:0.2f}'.format(self.shapes[self.current_shape].frequency)] +
1406 ([' Damping: {:.2f}%'.format(self.shapes[self.current_shape].damping * 100)] if self.show_damping else []) +
1407 [' ' + self.shapes[self.current_shape].comment1,
1408 ' ' + self.shapes[self.current_shape].comment2,
1409 ' ' + self.shapes[self.current_shape].comment3,
1410 ' ' + self.shapes[self.current_shape].comment4,
1411 ' ' + self.shapes[self.current_shape].comment5]))
1412 else:
1413 self.text.SetText(2, '')
1415 def _close(self):
1416 self.stop_animation()
1417 super(ShapePlotter, self)._close()
1420class DeflectionShapePlotter(ShapePlotter):
1421 """Class used to plot animated deflection shapes from spectra"""
1423 def __init__(self, geometry, deflection_shape_data, plot_kwargs,
1424 background_plotter_kwargs={'auto_update': 0.01},
1425 undeformed_opacity=0.25, starting_scale=1.0,
1426 deformed_opacity=1.0, num_curves=50):
1427 """
1428 Create a DeflectionShapePlotter object to plot shapes
1430 Parameters
1431 ----------
1432 geometry : Geometry
1433 Geometry on which the shapes will be plotted
1434 deflection_shape_data : NDDataArray
1435 Data array containing the deflection shapes to plot
1436 plot_kwargs : dict
1437 Keyword arguments passed to the Geometry.plot function
1438 background_plotter_kwargs : dict, optional
1439 Keyword arguments passed to the BackgroundPlotter constructor.
1440 The default is {'auto_update':0.01}.
1441 undeformed_opacity : float, optional
1442 Opacity of the undeformed geometry. Set to zero for no undeformed
1443 geometry. The default is 0.25.
1444 starting_scale : float, optional
1445 Starting scale of the shapes on the plot. The default is 1.0.
1446 num_curves : int, optional
1447 Maximum number of curves to plot on the frequency selector.
1448 Default is 50
1449 """
1450 # Convert function data into shapes
1451 from .sdynpy_shape import shape_array
1452 flattened_data = deflection_shape_data.flatten()
1453 response_coordinates = flattened_data.response_coordinate
1454 if not np.unique(response_coordinates).size == response_coordinates.size:
1455 raise ValueError('all response coordinates for operating data must be unique')
1456 shapes = shape_array(response_coordinates, flattened_data.ordinate.T,
1457 flattened_data[0].abscissa)
1458 super().__init__(geometry, shapes, plot_kwargs, background_plotter_kwargs,
1459 undeformed_opacity, starting_scale, deformed_opacity, shape_name='Shape',
1460 show_damping=False)
1462 # Get the layout so we can add a new widget to it
1463 layout = self.frame.layout()
1464 self.frequency_selector_plot = pg.PlotWidget()
1465 layout.addWidget(self.frequency_selector_plot)
1466 sp = self.frequency_selector_plot.sizePolicy()
1467 sp.setVerticalPolicy(QSizePolicy.Minimum)
1468 self.frequency_selector_plot.setSizePolicy(sp)
1469 self.frequency_selector_plot.sizeHint = lambda: QtCore.QSize(1006, 120)
1471 # Now we need to plot the data in the plot
1472 if num_curves < flattened_data.size:
1473 indices = np.linspace(0, flattened_data.size - 1, num_curves).astype(int)
1474 else:
1475 indices = np.arange(flattened_data.size)
1476 for data in flattened_data[indices]:
1477 self.frequency_selector_plot.plot(data.abscissa, np.abs(data.ordinate))
1478 # Plot the envelope
1479 max_envelope = np.linalg.svd(
1480 flattened_data.ordinate.T[..., np.newaxis], False, False, False).squeeze()
1481 self.frequency_selector_plot.plot(
1482 flattened_data[0].abscissa, max_envelope, pen=pg.mkPen(color='k', width=2))
1484 # Add an infinite line for a cursor
1485 self.frequency_selector = pg.InfiniteLine(
1486 0, movable=True, bounds=[flattened_data.abscissa.min(),
1487 flattened_data.abscissa.max()])
1488 self.frequency_selector_plot.addItem(self.frequency_selector)
1489 self.frequency_selector.sigPositionChanged.connect(self.modify_abscissa)
1491 self.frequency_selector_plot.setLogMode(False, True)
1492 self.frequency_selector_plot.setRange(
1493 yRange=(np.log10(max_envelope.min()), np.log10(max_envelope.max())))
1495 def modify_abscissa(self):
1496 frequency = self.frequency_selector.value()
1497 self.current_shape = np.argmin(np.abs(self.shapes.frequency - frequency))
1498 self.compute_displacements()
1499 self.update_shape()
1500 self.show_comment()
1502 def save_multiline_animation(self, filename=None, frame_rate=20,
1503 phase_change_per_frame=0.087,
1504 frequency_change_per_frame=0.5,
1505 start_phase=0.0,
1506 start_frequency=None, end_frequency=None,
1507 individual_images=False):
1508 self.stop_animation()
1509 imgs = []
1510 if start_frequency is None:
1511 start_frequency = self.shapes.frequency[0]
1512 if end_frequency is None:
1513 end_frequency = self.shapes.frequency[-1]
1514 frequency_frames = np.arange(start_frequency, end_frequency, frequency_change_per_frame)
1515 phase_frames = start_phase + (np.arange(frequency_frames.size) * phase_change_per_frame)
1516 for frequency, phase in zip(frequency_frames, phase_frames):
1517 self.current_shape = np.argmin(np.abs(self.shapes.frequency - frequency))
1518 self.compute_displacements()
1519 self.update_shape_mode(phase)
1520 self.show_comment()
1521 # Force GUI Update to ensure the OS doesn't delay events
1522 QApplication.processEvents()
1523 imgs.append(self.screenshot())
1524 if filename is None:
1525 return imgs
1526 else:
1527 if individual_images:
1528 for i, img in enumerate(imgs):
1529 fname = os.path.splitext(filename)
1530 fname = fname[0] + '-{:d}'.format(i) + fname[1]
1531 Image.fromarray(img).save(fname)
1532 else:
1533 imgs = [Image.fromarray(img).convert(
1534 'P', palette=Image.ADAPTIVE, colors=256) for img in imgs]
1535 imgs[0].save(fp=filename, format='GIF', append_images=imgs[1:],
1536 save_all=True, duration=int(1 / frame_rate * 1000), loop=0)
1539def split_list(seq, value):
1540 """
1541 Splits a list by a value that exists in the list
1543 Parameters
1544 ----------
1545 seq : iterable
1546 The sequence to split up
1547 value :
1548 Value within the list to use to split up the list
1550 Yields
1551 ------
1552 group : list
1553 Sublist of seq separated into groups.
1555 """
1556 group = []
1557 for val in seq:
1558 if val != value:
1559 group.append(val)
1560 elif group:
1561 yield group
1562 group = []
1563 yield group
1566def global_coord(coordinate_system, local_coord):
1567 '''Compute global coordinates from local coordinates
1569 Compute the global coordinates of a node from its local coordinates in the
1570 corresponding coordinate system.
1572 Parameters
1573 ----------
1575 coordinate_system : CoordinateSystemArray
1576 An array of coordinate systems of shape [...]
1577 local_coord : ndarray
1578 An array of coordinates of shape [...,3]
1580 Returns
1581 ------
1582 global_coord : ndarray
1583 An array of coordinates of shape [...,3]
1584 '''
1585 # Make sure we don't change local_coord
1586 local_coord = np.array(copy.deepcopy(local_coord))
1588 # First convert to cartesian if necessary
1589 cylindrical_css = coordinate_system.cs_type == 1
1590 cyl_coord = local_coord[cylindrical_css]
1591 cyl_coord[..., 0], cyl_coord[..., 1] = (cyl_coord[..., 0] * np.cos(cyl_coord[..., 1] * np.pi / 180),
1592 cyl_coord[..., 0] * np.sin(cyl_coord[..., 1] * np.pi / 180))
1593 local_coord[cylindrical_css] = cyl_coord
1595 spherical_css = coordinate_system.cs_type == 2
1596 sph_coord = local_coord[spherical_css]
1597 sph_coord[..., 0], sph_coord[..., 1], sph_coord[..., 2] = (
1598 sph_coord[..., 0] * np.sin(sph_coord[..., 1] * np.pi / 180) *
1599 np.cos(sph_coord[..., 2] * np.pi / 180),
1600 sph_coord[..., 0] * np.sin(sph_coord[..., 1] * np.pi / 180) *
1601 np.sin(sph_coord[..., 2] * np.pi / 180),
1602 sph_coord[..., 0] * np.cos(sph_coord[..., 1] * np.pi / 180))
1603 local_coord[spherical_css] = sph_coord
1604 # Convert to global cartesian from local cartesian
1605 return np.einsum('...ij,...i->...j', coordinate_system.matrix[..., :3, :], local_coord) + coordinate_system.matrix[..., 3, :]
1608def local_coord(coordinate_system, global_coord):
1609 """Compute local coordinates from global coordinates
1611 Compute the local coordinates of a node in the
1612 corresponding coordinate system from its global coordinates.
1614 Parameters
1615 ----------
1617 coordinate_system : CoordinateSystemArray
1618 An array of coordinate systems of shape [...]
1619 global_coord : ndarray
1620 An array of coordinates of shape [...,3]
1622 Returns
1623 ------
1624 local_coord : ndarray
1625 An array of coordinates of shape [...,3]
1626 """
1627 # Make sure we don't change global_coord
1628 global_coord = np.array(copy.deepcopy(global_coord))
1629 # Make sure we're working with a numpy array
1630 # Convert to the local cartesian
1631 local_coordinate = np.einsum(
1632 '...ji,...i->...j', coordinate_system.matrix[..., :3, :], global_coord - coordinate_system.matrix[..., 3, :])
1634 # Convert to cylindrical or spherical if necessary.
1635 cylindrical_css = coordinate_system.cs_type == 1
1636 cyl_coord = local_coordinate[cylindrical_css]
1637 cyl_coord[..., 0], cyl_coord[..., 1] = (np.sqrt(cyl_coord[..., 0]**2 + cyl_coord[..., 1]**2),
1638 np.arctan2(cyl_coord[..., 1], cyl_coord[..., 0]) * 180 / np.pi)
1639 local_coordinate[cylindrical_css] = cyl_coord
1641 spherical_css = coordinate_system.cs_type == 2
1642 sph_coord = local_coordinate[spherical_css]
1643 sph_coord[..., 0], sph_coord[..., 1] = (
1644 np.sqrt(sph_coord[..., 0]**2 + sph_coord[..., 1]**2 + sph_coord[..., 2]**2),
1645 np.arctan2(sph_coord[..., 1], sph_coord[..., 0]) * 180 / np.pi)
1646 # At this point local_coordinate[2,:] is still z and local_coordinate[0,:] is
1647 # now r instead of x
1648 sph_coord[..., 2] = np.arccos(sph_coord[..., 2] / sph_coord[..., 0]) * 180 / np.pi
1649 # Flip the 1 and 2 rows to be consistent with IMAT
1650 sph_coord = sph_coord[..., [0, 2, 1]]
1651 local_coordinate[spherical_css] = sph_coord
1652 return local_coordinate
1655def global_deflection(coordinate_system, local_deflection, global_point=None):
1656 """
1657 Compute deflections in the global coordinate system
1659 Parameters
1660 ----------
1661 coordinate_system : CoordinateSystemArray
1662 An array of coordinate systems of shape [...]
1663 local_deflection : ndarray
1664 An array of deflections of shape [...,3]
1665 global_point : np.ndarray, optional
1666 An array of coordinates in the global coordinate system of shape
1667 [...,3]. Must be specified if there are spherical or cylindrical
1668 coordinate systems. The default is None.
1670 Raises
1671 ------
1672 ValueError
1673 If cylindrical or spherical coordinate systems exist and global_point
1674 is not specified
1676 Returns
1677 -------
1678 global_deflection : ndarray
1679 An array of coordinates of shape [...,3] specifying the deflection
1680 directions
1682 """
1683 # Make sure we don't change local_deflection or point
1684 local_deflection = np.array(copy.deepcopy(local_deflection))
1685 if global_point is not None:
1686 global_point = np.array(copy.deepcopy(global_point))
1687 # Convert cylindrical and spherical to a local cartesian AT THE POINT
1688 # OF INTEREST (i.e. rotated in the coordinate system)
1689 local_rotation = np.zeros(coordinate_system.shape + (3, 3))
1690 for key, cs in coordinate_system.ndenumerate():
1691 if cs.cs_type == 0:
1692 this_local_rotation = np.eye(3)
1693 elif cs.cs_type == 1:
1694 if global_point is None:
1695 raise ValueError('point must be specified for Cylindrical coordinate systems')
1696 local_point = local_coord(cs, global_point[key])
1697 this_local_rotation = R(2, local_point[1], degrees=True)
1698 elif cs.cs_type == 2:
1699 if global_point is None:
1700 raise ValueError('point must be specified for Spherical coordinate systems')
1701 local_point = local_coord(cs, global_point[key])
1702 this_local_rotation = R(
1703 2, local_point[2], degrees=True) @ R(1, local_point[1], degrees=True)
1704 else:
1705 raise ValueError('invalid coordinate system type {:}'.format(cs.cs_type))
1706 local_rotation[key] = this_local_rotation
1707 return np.einsum('...ij,...j->...i', np.einsum('...ij,...ik->...jk', coordinate_system.matrix[..., :3, :], local_rotation), local_deflection)
1710class CoordinateSystemArray(SdynpyArray):
1711 """Coordinate System Information specifying position and directions.
1713 Use the coordinate_system_array helper function to create the array.
1714 """
1715 data_dtype = [('id', 'uint64'), ('name', '<U40'), ('color', 'uint16'),
1716 ('cs_type', 'uint16'), ('matrix', 'float64', (4, 3))]
1718 def __new__(subtype, shape, buffer=None, offset=0,
1719 strides=None, order=None):
1720 # Create the ndarray instance of our type, given the usual
1721 # ndarray input arguments. This will call the standard
1722 # ndarray constructor, but return an object of our type.
1723 # It also triggers a call to InfoArray.__array_finalize__
1724 obj = super(CoordinateSystemArray, subtype).__new__(subtype, shape,
1725 CoordinateSystemArray.data_dtype, buffer, offset, strides, order)
1726 # Finally, we must return the newly created object:
1727 return obj
1729 def __repr__(self):
1730 string_out = '{:>8s}, {:>6s}, {:>20s}, {:>5s}, {:>10s}\n'.format(
1731 'Index', 'ID', 'Name', 'Color', 'Type')
1732 if self.size == 0:
1733 string_out += '----------- Empty -------------\n'
1734 for i, (key, val) in enumerate(self.ndenumerate()):
1735 if i >= MAX_NUMBER_REPR:
1736 string_out += ' .\n .\n .\n'
1737 break
1738 string_out += '{:>8}, {:>6}, {:>20}, {:>5}, {:>10}\n'.format(
1739 str(key), val.id, val.name, val.color, _cs_type_strings[val.cs_type[()]])
1740 return string_out
1742 def __call__(self, index_by_id):
1743 """
1744 Select coordinate system by id rather than by index
1746 Parameters
1747 ----------
1748 index_by_id : int or np.ndarray
1749 ID number(s) to use to select coordinate systems.
1751 Raises
1752 ------
1753 ValueError
1754 If specified ID(s) not found in array.
1756 Returns
1757 -------
1758 output : CoordinateSystemArray
1759 Subset of CoordinateSystemArray with the specified IDs.
1761 """
1762 index_dict = {node.id[()]: index for index, node in self.ndenumerate()}
1763 ids = np.array(index_by_id)
1764 output = np.empty(ids.shape, dtype=self.dtype).view(self.__class__)
1765 for key, val in np.ndenumerate(ids):
1766 try:
1767 index = index_dict[val]
1768 except KeyError:
1769 raise ValueError('ID {:} not found in array'.format(val))
1770 output[key] = self[index]
1771 return output
1773 @staticmethod
1774 def from_unv(unv_data_dict, combine=True):
1775 """
1776 Load CoordinateSystemArrays from universal file data from read_unv
1778 Parameters
1779 ----------
1780 unv_data_dict : dict
1781 Dictionary containing data from read_unv
1782 combine : bool, optional
1783 If True, return as a single CoordinateSytemArray
1785 Returns
1786 -------
1787 output_arrays : CoordinateSystemArray
1788 Coordinate Systems read from unv
1790 """
1791 try:
1792 datasets = unv_data_dict[2420]
1793 except KeyError:
1794 if combine:
1795 return CoordinateSystemArray(0)
1796 else:
1797 return []
1798 output_arrays = []
1799 for dataset in datasets:
1800 output_arrays.append(coordinate_system_array(dataset.cs_labels,
1801 dataset.cs_names,
1802 dataset.cs_colors,
1803 dataset.cs_types,
1804 dataset.cs_matrices))
1805 if combine:
1806 output_arrays = np.concatenate(output_arrays)
1807 return output_arrays
1810def coordinate_system_array(id=1, name='', color=1, cs_type=0, matrix=np.concatenate((np.eye(3), np.zeros((1, 3))), axis=0),
1811 structured_array=None):
1812 """
1813 Creates an array that specifies coordinate systems in the geometry
1815 Creates an array of coordinate systems that specify directions of sensors
1816 in a test or analysis. Coordinate system arrays can be created using a
1817 numpy structured array or individual arrays for each attribute.
1818 Multidimensional arrays can be used.
1820 Parameters
1821 ----------
1822 id : ndarray
1823 Integer array corresponding to the id of the coordinate systems. Input
1824 will be cast to an integer (i.e. 2.0 -> 2, 1.9 -> 1)
1825 name : ndarray
1826 name of the coordinate systems
1827 color : ndarray
1828 color of the coordinate systems as integers
1829 cs_type : ndarray
1830 Coordinate system types (0 = Cartesian, 1 = Polar, 2 = Spherical)
1831 matrix : ndarray
1832 Coordinate system transformation matrix with shape [... x 4 x 3]
1833 structured_array : ndarray (structured)
1834 Alternatively to the individual attributes, a single numpy structured
1835 array can be passed, which should have the same name as the inputs to
1836 the function listed above.
1838 Returns
1839 -------
1840 coordinate_system_array : CoordinateSystemArray
1841 """
1842 if structured_array is not None:
1843 try:
1844 id = structured_array['id']
1845 name = structured_array['name']
1846 color = structured_array['color']
1847 cs_type = structured_array['cs_type']
1848 matrix = structured_array['matrix']
1849 except (ValueError, TypeError):
1850 raise ValueError(
1851 'structured_array must be numpy.ndarray with dtype names "id", "name", "color", "cs_type" and "matrix"')
1852 else:
1853 id = np.array(id)
1854 name = np.array(name)
1855 color = np.array(color)
1856 cs_type = np.array(cs_type)
1857 matrix = np.array(matrix)
1859 # Don't check shapes because we want them to be broadcastable, but id have to
1860 # be unique so we will use that for the shape
1861 coord_sys_array = CoordinateSystemArray(id.shape)
1862 coord_sys_array.id = id
1863 coord_sys_array.name = name
1864 coord_sys_array.color = color
1865 coord_sys_array.cs_type = cs_type
1866 coord_sys_array.matrix = matrix
1868 return coord_sys_array
1871class ElementArray(SdynpyArray):
1872 """Element information array
1874 Use the element_array helper function to create the array.
1875 """
1877 data_dtype = [('id', 'uint64'), ('type', 'uint8'), ('color', 'uint16'),
1878 ('connectivity', 'object')]
1880 def __new__(subtype, shape, buffer=None, offset=0,
1881 strides=None, order=None):
1882 # Create the ndarray instance of our type, given the usual
1883 # ndarray input arguments. This will call the standard
1884 # ndarray constructor, but return an object of our type.
1885 # It also triggers a call to InfoArray.__array_finalize__
1887 obj = super(ElementArray, subtype).__new__(subtype, shape,
1888 ElementArray.data_dtype, buffer, offset, strides, order)
1889 # Finally, we must return the newly created object:
1890 return obj
1892 def __repr__(self):
1893 string_out = '{:>8s}, {:>6s}, {:>4s}, {:>5s}, {:>7s}\n'.format(
1894 'Index', 'ID', 'Type', 'Color', '# Nodes')
1895 if self.size == 0:
1896 string_out += '----------- Empty -------------\n'
1897 for i, (key, val) in enumerate(self.ndenumerate()):
1898 if i >= MAX_NUMBER_REPR:
1899 string_out += ' .\n .\n .\n'
1900 break
1901 string_out += '{:>8}, {:>6}, {:>4}, {:>5}, {:>7}\n'.format(
1902 str(key), val.id, val.type, val.color, len(val.connectivity))
1903 return string_out
1905 def __call__(self, index_by_id):
1906 """
1907 Select elements by id rather than by index
1909 Parameters
1910 ----------
1911 index_by_id : int or np.ndarray
1912 ID number(s) to use to select elements.
1914 Raises
1915 ------
1916 ValueError
1917 If specified ID(s) not found in array.
1919 Returns
1920 -------
1921 output : ElementArray
1922 Subset of ElementArray with the specified IDs.
1924 """
1925 index_dict = {node.id[()]: index for index, node in self.ndenumerate()}
1926 ids = np.array(index_by_id)
1927 output = np.empty(ids.shape, dtype=self.dtype).view(self.__class__)
1928 for key, val in np.ndenumerate(ids):
1929 try:
1930 index = index_dict[val]
1931 except KeyError:
1932 raise ValueError('ID {:} not found in array'.format(val))
1933 output[key] = self[index]
1934 return output
1936 def reduce(self, node_list):
1937 """
1938 Keep only elements that have all nodes contained in node_list
1940 Parameters
1941 ----------
1942 node_list : iterable
1943 Iterable containing nodes to keep.
1945 Returns
1946 -------
1947 ElementArray
1948 ElementArray containing only elements with all nodes in node_list.
1950 """
1951 output_list = []
1952 for key, element in self.ndenumerate():
1953 if np.all(np.isin(element.connectivity, node_list)):
1954 output_list.append(element)
1955 return np.array(output_list, self.dtype).view(ElementArray)
1957 def remove(self, element_list):
1958 """
1959 Removes elements with id numbers in traceline_list
1961 Parameters
1962 ----------
1963 traceline_list : interable
1964 Iterable containing tracelines to discard.
1966 Returns
1967 -------
1968 TracelineArray
1969 TracelineArray containing tracelines with id numbers not it
1970 traceline_list.
1972 """
1973 output_list = []
1974 for key,element in self.ndenumerate():
1975 if element.id not in element_list:
1976 output_list.append(element)
1977 return np.array(output_list, self.dtype).view(ElementArray)
1979 @staticmethod
1980 def from_unv(unv_data_dict, combine=True):
1981 """
1982 Load ElementArrays from universal file data from read_unv
1984 Parameters
1985 ----------
1986 unv_data_dict : dict
1987 Dictionary containing data from read_unv
1988 combine : bool, optional
1989 If True, return as a single ElementArray
1991 Returns
1992 -------
1993 output_arrays : ElementArray
1994 Elements read from unv
1996 """
1997 try:
1998 datasets = unv_data_dict[2412]
1999 except KeyError:
2000 if combine:
2001 return ElementArray(0)
2002 else:
2003 return []
2004 output_arrays = []
2005 for dataset in datasets:
2006 output_arrays.append(element_array(dataset.element_labels,
2007 dataset.fe_descriptor_ids,
2008 dataset.colors,
2009 dataset.connectivities))
2010 if combine:
2011 output_arrays = np.concatenate(output_arrays)
2012 return output_arrays
2015def element_array(id=1, type=None, color=1,
2016 connectivity=None, structured_array=None):
2017 """
2018 Creates an array that specifies elements in the geometry
2020 Creates an array of elements that specify connectivity of sensors
2021 in a test or analysis. Element arrays can be created using a
2022 numpy structured array or individual arrays for each attribute.
2023 Multidimensional arrays can be used.
2025 Parameters
2026 ----------
2027 id : ndarray
2028 Integer array corresponding to the id of the node. Input
2029 will be cast to an integer (i.e. 2.0 -> 2, 1.9 -> 1)
2030 type : ndarray
2031 Element types. See notes for type mapping
2032 color : ndarray
2033 color of the elements as integers
2034 connectivity : ndarray
2035 An object array of iterables defining the connectivity for each element.
2036 structured_array : ndarray (structured)
2037 Alternatively to the individual attributes, a single numpy structured
2038 array can be passed, which should have the same name as the inputs to
2039 the function listed above.
2042 Returns
2043 -------
2044 element_array : ElementArray
2047 Notes
2048 -----
2049 Here is a list of element types
2051 11 : 'Rod',
2052 21 : 'Linear beam',
2053 22 : 'Tapered beam',
2054 23 : 'Curved beam',
2055 24 : 'Parabolic beam',
2056 31 : 'Straight pipe',
2057 32 : 'Curved pipe',
2058 41 : 'Plane Stress Linear Triangle',
2059 42 : 'Plane Stress Parabolic Triangle',
2060 43 : 'Plane Stress Cubic Triangle',
2061 44 : 'Plane Stress Linear Quadrilateral',
2062 45 : 'Plane Stress Parabolic Quadrilateral',
2063 46 : 'Plane Strain Cubic Quadrilateral',
2064 51 : 'Plane Strain Linear Triangle',
2065 52 : 'Plane Strain Parabolic Triangle',
2066 53 : 'Plane Strain Cubic Triangle',
2067 54 : 'Plane Strain Linear Quadrilateral',
2068 55 : 'Plane Strain Parabolic Quadrilateral',
2069 56 : 'Plane Strain Cubic Quadrilateral',
2070 61 : 'Plate Linear Triangle',
2071 62 : 'Plate Parabolic Triangle',
2072 63 : 'Plate Cubic Triangle',
2073 64 : 'Plate Linear Quadrilateral',
2074 65 : 'Plate Parabolic Quadrilateral',
2075 66 : 'Plate Cubic Quadrilateral',
2076 71 : 'Membrane Linear Quadrilateral',
2077 72 : 'Membrane Parabolic Triangle',
2078 73 : 'Membrane Cubic Triangle',
2079 74 : 'Membrane Linear Triangle',
2080 75 : 'Membrane Parabolic Quadrilateral',
2081 76 : 'Membrane Cubic Quadrilateral',
2082 81 : 'Axisymetric Solid Linear Triangle',
2083 82 : 'Axisymetric Solid Parabolic Triangle',
2084 84 : 'Axisymetric Solid Linear Quadrilateral',
2085 85 : 'Axisymetric Solid Parabolic Quadrilateral',
2086 91 : 'Thin Shell Linear Triangle',
2087 92 : 'Thin Shell Parabolic Triangle',
2088 93 : 'Thin Shell Cubic Triangle',
2089 94 : 'Thin Shell Linear Quadrilateral',
2090 95 : 'Thin Shell Parabolic Quadrilateral',
2091 96 : 'Thin Shell Cubic Quadrilateral',
2092 101: 'Thick Shell Linear Wedge',
2093 102: 'Thick Shell Parabolic Wedge',
2094 103: 'Thick Shell Cubic Wedge',
2095 104: 'Thick Shell Linear Brick',
2096 105: 'Thick Shell Parabolic Brick',
2097 106: 'Thick Shell Cubic Brick',
2098 111: 'Solid Linear Tetrahedron',
2099 112: 'Solid Linear Wedge',
2100 113: 'Solid Parabolic Wedge',
2101 114: 'Solid Cubic Wedge',
2102 115: 'Solid Linear Brick',
2103 116: 'Solid Parabolic Brick',
2104 117: 'Solid Cubic Brick',
2105 118: 'Solid Parabolic Tetrahedron',
2106 121: 'Rigid Bar',
2107 122: 'Rigid Element',
2108 136: 'Node To Node Translational Spring',
2109 137: 'Node To Node Rotational Spring',
2110 138: 'Node To Ground Translational Spring',
2111 139: 'Node To Ground Rotational Spring',
2112 141: 'Node To Node Damper',
2113 142: 'Node To Gound Damper',
2114 151: 'Node To Node Gap',
2115 152: 'Node To Ground Gap',
2116 161: 'Lumped Mass',
2117 171: 'Axisymetric Linear Shell',
2118 172: 'Axisymetric Parabolic Shell',
2119 181: 'Constraint',
2120 191: 'Plastic Cold Runner',
2121 192: 'Plastic Hot Runner',
2122 193: 'Plastic Water Line',
2123 194: 'Plastic Fountain',
2124 195: 'Plastic Baffle',
2125 196: 'Plastic Rod Heater',
2126 201: 'Linear node-to-node interface',
2127 202: 'Linear edge-to-edge interface',
2128 203: 'Parabolic edge-to-edge interface',
2129 204: 'Linear face-to-face interface',
2130 208: 'Parabolic face-to-face interface',
2131 212: 'Linear axisymmetric interface',
2132 213: 'Parabolic axisymmetric interface',
2133 221: 'Linear rigid surface',
2134 222: 'Parabolic rigid surface',
2135 231: 'Axisymetric linear rigid surface',
2136 232: 'Axisymentric parabolic rigid surface'
2137 """
2138 if structured_array is not None:
2139 try:
2140 id = structured_array['id']
2141 type = structured_array['type']
2142 color = structured_array['color']
2143 connectivity = structured_array['connectivity']
2144 except (ValueError, TypeError):
2145 raise ValueError(
2146 'structured_array must be numpy.ndarray with dtype names "id", "type", "color", and "connectivity"')
2147 else:
2148 id = np.array(id)
2149 type = np.array(type)
2150 color = np.array(color)
2151 # Need to make connectivity an object array so don't do this!
2152 # connectivity = np.array(connectivity,dtype=object)
2154 # Don't check shapes because we want them to be broadcastable, but id have to
2155 # be unique so we will use that for the shape
2156 earray = ElementArray(id.shape)
2157 earray.id = id
2158 earray.type = type
2159 earray.color = color
2160 if earray.ndim == 0:
2161 # This is gross, but it works!
2162 np.ndarray.__getitem__(earray, 'connectivity')[()] = np.array(connectivity, dtype='uint32')
2163 else:
2164 connectivity = np.array(connectivity, dtype=object)
2165 for key, val in np.ndenumerate(id):
2166 earray.connectivity[key] = np.array(connectivity[key], dtype='uint32')
2167 return earray
2170class NodeArray(SdynpyArray):
2171 """Node information array
2173 Use the node_array helper function to create the array.
2174 """
2176 data_dtype = [('id', 'uint64'), ('coordinate', 'float64', (3,)), ('color', 'uint16'),
2177 ('def_cs', 'uint64'), ('disp_cs', 'uint64')]
2179 def __new__(subtype, shape, buffer=None, offset=0,
2180 strides=None, order=None):
2181 # Create the ndarray instance of our type, given the usual
2182 # ndarray input arguments. This will call the standard
2183 # ndarray constructor, but return an object of our type.
2184 # It also triggers a call to InfoArray.__array_finalize__
2186 obj = super(NodeArray, subtype).__new__(subtype, shape,
2187 NodeArray.data_dtype, buffer, offset, strides, order)
2188 # Finally, we must return the newly created object:
2189 return obj
2191 def __repr__(self):
2192 string_out = '{:>8s}, {:>6s}, {:>8s}, {:>8s}, {:>8s}, {:>5s}, {:>5s}\n'.format(
2193 'Index', 'ID', 'X', 'Y', 'Z', 'DefCS', 'DisCS')
2194 if self.size == 0:
2195 string_out += '----------- Empty -------------\n'
2196 for i, (key, val) in enumerate(self.ndenumerate()):
2197 if i >= MAX_NUMBER_REPR:
2198 string_out += ' .\n .\n .\n'
2199 break
2200 string_out += '{:>8}, {:>6}, {:>8.3f}, {:>8.3f}, {:>8.3f}, {:>5}, {:>5}\n'.format(
2201 str(key), val.id, *val.coordinate, val.def_cs, val.disp_cs)
2202 return string_out
2204 def __call__(self, index_by_id):
2205 """
2206 Select nodes by id rather than by index
2208 Parameters
2209 ----------
2210 index_by_id : int or np.ndarray
2211 ID number(s) to use to select nodes.
2213 Raises
2214 ------
2215 ValueError
2216 If specified ID(s) not found in array.
2218 Returns
2219 -------
2220 output : NodeArray
2221 Subset of NodeArray with the specified IDs.
2223 """
2224 index_dict = {node.id[()]: index for index, node in self.ndenumerate()}
2225 ids = np.array(index_by_id)
2226 output = np.empty(ids.shape, dtype=self.dtype).view(self.__class__)
2227 for key, val in np.ndenumerate(ids):
2228 try:
2229 index = index_dict[val]
2230 except KeyError:
2231 raise ValueError('ID {:} not found in array'.format(val))
2232 output[key] = self[index]
2233 return output
2235 def reduce(self, node_list):
2236 """
2237 Keep only nodes that are contained in node_list
2239 Parameters
2240 ----------
2241 node_list : iterable
2242 Iterable containing nodes to keep.
2244 Returns
2245 -------
2246 NodeArray
2247 NodeArray containing only nodes in node_list.
2249 """
2250 return self(node_list)
2252 @staticmethod
2253 def from_unv(unv_data_dict, combine=True):
2254 """
2255 Load NodeArrays from universal file data from read_unv
2257 Parameters
2258 ----------
2259 unv_data_dict : dict
2260 Dictionary containing data from read_unv
2261 combine : bool, optional
2262 If True, return as a single NodeArray
2264 Returns
2265 -------
2266 output_arrays : NodeArray
2267 Nodes read from unv
2269 """
2270 try:
2271 datasets = unv_data_dict[2411]
2272 except KeyError:
2273 if combine:
2274 return NodeArray(0)
2275 else:
2276 return []
2277 output_arrays = []
2278 for dataset in datasets:
2279 output_arrays.append(node_array(dataset.node_labels,
2280 dataset.coordinates,
2281 dataset.colors,
2282 dataset.export_coordinate_systems,
2283 dataset.displacement_coordinate_systems))
2284 if combine:
2285 output_arrays = np.concatenate(output_arrays)
2286 return output_arrays
2288 def by_position(self, position_array):
2289 """
2290 Select node by closest position
2292 Parameters
2293 ----------
2294 position_array : np.ndarray
2295 A (...,3) shape array containing positions of nodes to keep
2297 Returns
2298 -------
2299 NodeArray
2300 NodArray containing nodes that were closest to the positions in
2301 position_array.
2303 """
2304 node_differences = np.linalg.norm(
2305 self.coordinate - position_array[(Ellipsis,) + (np.newaxis,) * self.ndim + (slice(None),)], axis=-1)
2306 node_differences = node_differences.reshape(*position_array.shape[:-1], -1)
2307 node_indices = np.argmin(
2308 node_differences,
2309 axis=-1)
2310 return self.flatten()[node_indices]
2312 @staticmethod
2313 def project_to_minimum_plane(coordinates, return_3D=True):
2314 """
2315 Projects coordinates to a single best-fit plane
2317 Parameters
2318 ----------
2319 coordinates : np.ndarray
2320 A (...,3) coordinate array.
2321 return_3D : bool, optional
2322 If True, return the 3D coordinates of the projected points.
2323 Otherwise return the projected 2D coordinate. The default is True.
2325 Returns
2326 -------
2327 np.ndarray
2328 Points projected to a best-fit plane.
2330 """
2331 coords = coordinates.T.copy()
2332 mean = np.mean(coords, axis=-1, keepdims=True)
2333 coords -= mean
2334 U, S, V = np.linalg.svd(coords)
2335 coords = U.T @ coords
2336 coords[2, :] = 0
2337 if return_3D:
2338 return (U @ coords + mean).T
2339 else:
2340 return coords[:2].T
2342 def global_coordinate(self, cs_array):
2343 """
2344 Get the global coordinate of each node
2346 Parameters
2347 ----------
2348 cs_array : CoordinateSystemArray
2349 CoordinateSystemArray consisting of the local coordinate systems for
2350 each node
2352 Returns
2353 -------
2354 points : np.ndarray
2355 Coordinates of the nodes in the global coordinate system.
2357 """
2358 points = global_coord(cs_array(self.def_cs), self.coordinate)
2359 return points
2361 def triangulate(self, cs_array, projection_function=None,
2362 return_element_array=True, element_color=1,
2363 condition_threshold=None):
2364 """
2365 Creates elements for a node set using Delaunay Triangulation
2367 Parameters
2368 ----------
2369 cs_array : CoordinateSystemArray
2370 CoordinateSystemArray containing coordinate systems for each node.
2371 projection_function : function, optional
2372 Function to use to project 3D coordinates to 2D coordinates for
2373 triangulation. The default is None.
2374 return_element_array : bool, optional
2375 Returns an ElementArray if True, otherwise it simply returns the
2376 triangle simplices. The default is True.
2377 element_color : np.ndarray or int, optional
2378 Integers representing colors applied to the elements.
2379 The default is 1.
2380 condition_threshold : float, optional
2381 Condition number threshold used to remove triangles that are poorly
2382 shaped. The default is None.
2384 Returns
2385 -------
2386 ElementArray or np.ndarray
2387 ElementArray containing elements or np.ndarray containing triangle
2388 simplices.
2390 """
2391 if projection_function is None:
2392 def projection_function(coords): return self.project_to_minimum_plane(coords, False)
2393 # Compute global positions of the nodes
2394 global_positions = self.global_coordinate(cs_array)
2395 # Now compute the projection function
2396 projected_positions = projection_function(global_positions)
2397 # Now triangulate
2398 tri = Delaunay(projected_positions)
2399 simplices = tri.simplices
2400 if condition_threshold is not None:
2401 # Now remove small triangles via condition number
2402 tri_points = projected_positions[tri.simplices, :]
2403 tri_points -= np.mean(tri_points, axis=1, keepdims=True)
2404 simplices = simplices[np.linalg.cond(tri_points) < condition_threshold]
2405 if not return_element_array:
2406 return simplices
2407 else:
2408 # Create an element array
2409 connectivity = self.id[simplices]
2410 return element_array(id=np.arange(connectivity.shape[0]) + 1,
2411 type=61, color=element_color,
2412 connectivity=connectivity)
2414 def by_grid(self, grid_spacing, x_min=None, y_min=None, z_min=None,
2415 x_max=None, y_max=None, z_max=None):
2416 """
2417 Selects nodes in a grid
2419 Parameters
2420 ----------
2421 grid_spacing : float
2422 Approximate grid spacing between selected nodes
2423 x_min : float, optional
2424 Minimum X-value of the grid.
2425 The default is the minimum node x-coordinate.
2426 y_min : float, optional
2427 Minimun Y-value of the grid.
2428 The default is the minimum node y-coordinate.
2429 z_min : float, optional
2430 Minimum Z-value of the grid.
2431 The default is the minimum node z-coordinate.
2432 x_max : float, optional
2433 Maximum X-value of the grid.
2434 The default is the maximum node x-coordinate.
2435 y_max : float, optional
2436 Maximum Y-value of the grid.
2437 The default is the maximum node y-coordinate.
2438 z_max : float, optional
2439 Maximum Z-value of the grid.
2440 The default is the maximum node z-coordinate.
2442 Returns
2443 -------
2444 NodeArray
2445 NodeArray containing the closest nodes to the grid points
2446 """
2447 if x_min is None:
2448 x_min = self.coordinate[:, 0].min()
2449 if y_min is None:
2450 y_min = self.coordinate[:, 1].min()
2451 if z_min is None:
2452 z_min = self.coordinate[:, 2].min()
2453 if x_max is None:
2454 x_max = self.coordinate[:, 0].max()
2455 if y_max is None:
2456 y_max = self.coordinate[:, 1].max()
2457 if z_max is None:
2458 z_max = self.coordinate[:, 2].max()
2459 min_coords = np.array((x_min, y_min, z_min))
2460 max_coords = np.array((x_max, y_max, z_max))
2461 range_coords = max_coords - min_coords
2462 num_grids = np.ceil(range_coords / grid_spacing).astype(int)
2463 # Now we will create the grid in each dimension using linspace
2464 grid_arrays = [np.linspace(min_coord, max_coord, num_coord)
2465 for min_coord, max_coord, num_coord
2466 in zip(min_coords, max_coords, num_grids)]
2467 # We will then use meshgrid to assemble the array grid. We will flatten the
2468 # point dimension and transpose so the array has shape (n_points x 3)
2469 grid_points = np.array(np.meshgrid(*grid_arrays, indexing='ij')).reshape(3, -1).T
2470 # Now we will select nodes that are closest to the points in the grid
2471 candidate_nodes = self.by_position(grid_points)
2472 # Reduce to only nodes that are within one grid spacing of their requested point
2473 # this reduces the nodes that are selected by points far from the model
2474 candidate_nodes = candidate_nodes[
2475 np.linalg.norm(candidate_nodes.coordinate - grid_points, axis=-1) < grid_spacing]
2476 # Remove any duplicates
2477 candidate_node_ids = np.unique(candidate_nodes.id)
2478 return candidate_nodes(candidate_node_ids)
2481def node_array(id=1, coordinate=np.array((0, 0, 0)), color=1, def_cs=1, disp_cs=1,
2482 structured_array=None):
2483 """
2484 Creates an array that specifies nodes in the geometry
2486 Creates an array of nodes that specify positions of sensors
2487 in a test or analysis. Node arrays can be created using a
2488 numpy structured array or individual arrays for each attribute.
2489 Multidimensional arrays can be used.
2491 Parameters
2492 ----------
2493 id : ndarray
2494 Integer array corresponding to the id of the node. Input
2495 will be cast to an integer (i.e. 2.0 -> 2, 1.9 -> 1)
2496 coordinate : ndarray
2497 Positions of the nodes in space
2498 color : ndarray
2499 color of the nodes as integers
2500 def_cs : ndarray
2501 Coordinate system where the node's position is defined
2502 disp_cs : ndarray
2503 Coordinate system where the node's displacements are defined
2504 structured_array : ndarray (structured)
2505 Alternatively to the individual attributes, a single numpy structured
2506 array can be passed, which should have the same name as the inputs to
2507 the function listed above.
2510 Returns
2511 -------
2512 node_array : NodeArray
2513 """
2514 if structured_array is not None:
2515 try:
2516 id = structured_array['id']
2517 coordinate = structured_array['coordinate']
2518 color = structured_array['color']
2519 def_cs = structured_array['def_cs']
2520 disp_cs = structured_array['disp_cs']
2521 except (ValueError, TypeError):
2522 raise ValueError(
2523 'structured_array must be numpy.ndarray with dtype names "id", "coordinate", "color", "def_cs" and "disp_cs"')
2524 else:
2525 id = np.array(id)
2526 coordinate = np.array(coordinate)
2527 color = np.array(color)
2528 def_cs = np.array(def_cs)
2529 disp_cs = np.array(disp_cs)
2531 # Don't check shapes because we want them to be broadcastable, but id have to
2532 # be unique so we will use that for the shape
2533 narray = NodeArray(id.shape)
2534 narray.id = id
2535 narray.coordinate = coordinate
2536 narray.color = color
2537 narray.def_cs = def_cs
2538 narray.disp_cs = disp_cs
2540 return narray
2543class TracelineArray(SdynpyArray):
2544 """Traceline information array
2546 Use the traceline_array helper function to create the array.
2547 """
2548 data_dtype = [('id', 'uint64'), ('color', 'uint16'), ('description', '<U40'),
2549 ('connectivity', 'object')]
2551 def __new__(subtype, shape, buffer=None, offset=0,
2552 strides=None, order=None):
2553 # Create the ndarray instance of our type, given the usual
2554 # ndarray input arguments. This will call the standard
2555 # ndarray constructor, but return an object of our type.
2556 # It also triggers a call to InfoArray.__array_finalize__
2558 obj = super(TracelineArray, subtype).__new__(subtype, shape,
2559 TracelineArray.data_dtype, buffer, offset, strides, order)
2560 # Finally, we must return the newly created object:
2561 return obj
2563 def __repr__(self):
2564 string_out = '{:>8s}, {:>6s}, {:>20s}, {:>5s}, {:>7s}\n'.format(
2565 'Index', 'ID', 'Description', 'Color', '# Nodes')
2566 if self.size == 0:
2567 string_out += '----------- Empty -------------\n'
2568 for i, (key, val) in enumerate(self.ndenumerate()):
2569 if i >= MAX_NUMBER_REPR:
2570 string_out += ' .\n .\n .\n'
2571 break
2572 string_out += '{:>8s}, {:>6d}, {:>20s}, {:>5d}, {:>7d}\n'.format(
2573 str(key), val.id, val.description, val.color, len(val.connectivity))
2574 return string_out
2576 def __call__(self, index_by_id):
2577 """
2578 Select nodes by id rather than by index
2580 Parameters
2581 ----------
2582 index_by_id : int or np.ndarray
2583 ID number(s) to use to select nodes.
2585 Raises
2586 ------
2587 ValueError
2588 If specified ID(s) not found in array.
2590 Returns
2591 -------
2592 output : TracelineArray
2593 Subset of TracelineArray with the specified IDs.
2595 """
2596 index_dict = {node.id[()]: index for index, node in self.ndenumerate()}
2597 ids = np.array(index_by_id)
2598 output = np.empty(ids.shape, dtype=self.dtype).view(self.__class__)
2599 for key, val in np.ndenumerate(ids):
2600 try:
2601 index = index_dict[val]
2602 except KeyError:
2603 raise ValueError('ID {:} not found in array'.format(val))
2604 output[key] = self[index]
2605 return output
2607 def reduce(self, node_list):
2608 """
2609 Keep only tracelines fully contain nodes in node_list
2611 Parameters
2612 ----------
2613 node_list : iterable
2614 Iterable containing nodes to keep.
2616 Returns
2617 -------
2618 TracelineArray
2619 TracelineArray containing lines that contain only nodes in node_list.
2621 """
2622 output_list = []
2623 for key, traceline in self.ndenumerate():
2624 if np.all(np.isin(traceline.connectivity, node_list) | (traceline.connectivity == 0)):
2625 output_list.append(traceline)
2626 return np.array(output_list, self.dtype).view(TracelineArray)
2628 def remove(self, traceline_list):
2629 """
2630 Removes tracelines with id numbers in traceline_list
2632 Parameters
2633 ----------
2634 traceline_list : interable
2635 Iterable containing tracelines to discard.
2637 Returns
2638 -------
2639 TracelineArray
2640 TracelineArray containing tracelines with id numbers not it
2641 traceline_list.
2643 """
2644 output_list = []
2645 for key,traceline in self.ndenumerate():
2646 if traceline.id not in traceline_list:
2647 output_list.append(traceline)
2648 return np.array(output_list, self.dtype).view(TracelineArray)
2650 @staticmethod
2651 def from_unv(unv_data_dict, combine=True):
2652 """
2653 Load TracelineArrays from universal file data from read_unv
2655 Parameters
2656 ----------
2657 unv_data_dict : dict
2658 Dictionary containing data from read_unv
2659 combine : bool, optional
2660 If True, return as a single TracelineArray object
2662 Returns
2663 -------
2664 output_arrays : TracelineArray
2665 Tracelines read from unv
2667 """
2668 try:
2669 datasets = unv_data_dict[82]
2670 except KeyError:
2671 if combine:
2672 return TracelineArray(0)
2673 else:
2674 return []
2675 output_arrays = []
2676 for dataset in datasets:
2677 output_arrays.append(traceline_array(dataset.traceline_number,
2678 dataset.identification,
2679 dataset.color,
2680 dataset.nodes)[np.newaxis])
2681 if combine:
2682 output_arrays = np.concatenate(output_arrays)
2683 return output_arrays
2686def traceline_array(id=1, description='', color=1,
2687 connectivity=None, structured_array=None):
2688 """
2689 Creates an array that specifies tracelines in the geometry
2691 Creates an array of tracelines that specify connectivity of sensors
2692 in a test or analysis. Traceline arrays can be created using a
2693 numpy structured array or individual arrays for each attribute.
2694 Multidimensional arrays can be used.
2696 Parameters
2697 ----------
2698 id : ndarray
2699 Integer array corresponding to the id of the node. Input
2700 will be cast to an integer (i.e. 2.0 -> 2, 1.9 -> 1)
2701 Description : ndarray
2702 Traceline name
2703 color : ndarray
2704 color of the elements as integers
2705 connectivity : ndarray
2706 An object array of iterables defining the connectivity for each traceline.
2707 structured_array : ndarray (structured)
2708 Alternatively to the individual attributes, a single numpy structured
2709 array can be passed, which should have the same name as the inputs to
2710 the function listed above.
2713 Returns
2714 -------
2715 traceline_array : TracelineArray
2717 """
2718 if structured_array is not None:
2719 try:
2720 id = structured_array['id']
2721 description = structured_array['description']
2722 color = structured_array['color']
2723 connectivity = structured_array['connectivity']
2724 except (ValueError, TypeError):
2725 raise ValueError(
2726 'structured_array must be numpy.ndarray with dtype names "id", "color", "description", and "connectivity"')
2727 else:
2728 id = np.array(id)
2729 description = np.array(description)
2730 color = np.array(color)
2731 # Don't do this, need to keep as an object dtype
2732 # connectivity = np.array(connectivity)
2734 # Don't check shapes because we want them to be broadcastable, but id have to
2735 # be unique so we will use that for the shape
2736 tlarray = TracelineArray(id.shape)
2737 tlarray.id = id
2738 tlarray.description = description
2739 tlarray.color = color
2740 if tlarray.ndim == 0:
2741 np.ndarray.__getitem__(tlarray, 'connectivity')[()] = np.array(connectivity, dtype='uint32')
2742 else:
2743 connectivity = np.array(connectivity, dtype=object)
2744 for key, val in np.ndenumerate(id):
2745 tlarray.connectivity[key] = np.array(connectivity[key], dtype='uint32')
2747 return tlarray
2750from ..fem.sdynpy_exodus import Exodus, ExodusInMemory # noqa: E402
2751from .sdynpy_shape import shape_array # noqa: E402
2754class Geometry:
2755 """Container for nodes, coordinate systems, tracelines, and elements
2757 Geometry is the class that is most useful for working with the positioning
2758 and spatial visualization of test or analysis data. It contains functions
2759 to plot and animate 3D geometry"""
2761 def __init__(self, node=NodeArray((0,)), coordinate_system=CoordinateSystemArray((0,)),
2762 traceline=TracelineArray((0,)), element=ElementArray((0,))):
2763 """
2764 Initialize a geometry object with nodes, coordinate systems, tracelines,
2765 and elements.
2767 All input arguments will be flattened when passed to the Geometry, as
2768 the geometry does not support multi-dimensional object arrays.
2770 Parameters
2771 ----------
2772 node : NodeArray, optional
2773 The set of nodes in the geometry. The default is NodeArray((0,)),
2774 an empty NodeArray.
2775 coordinate_system : CoordinateSystemArray, optional
2776 The set of coordinate systems in the geometry. The default is
2777 CoordinateSystemArray((0,)), an empty CoordinateSystemArray
2778 traceline : TracelineArray, optional
2779 The set of tracelines defined in the geometry. The default is
2780 TracelineArray((0,)), an empty TracelineArray
2781 element : ElementArray, optional
2782 The set of elements defined in the geometry. The default is
2783 ElementArray((0,)), an empty element array.
2785 Returns
2786 -------
2787 None.
2789 """
2790 self.node = node.flatten()
2791 self.coordinate_system = coordinate_system.flatten()
2792 self.traceline = traceline.flatten()
2793 self.element = element.flatten()
2795 def __repr__(self):
2796 string_out = '\n'.join(val.capitalize(
2797 ) + '\n' + repr(self.__dict__[val]) for val in ['node', 'coordinate_system', 'traceline', 'element'])
2798 return string_out
2800 def reduce(self, node_list):
2801 """
2802 Reduce the geometry to only contain nodes in `node_list`
2804 Elements and tracelines will only be kept if all nodes in each element
2805 or traceline are in `node_list`. Coordinate systems will only be kept
2806 if they are required by a node in `node_list`
2808 Parameters
2809 ----------
2810 node_list : iterable
2811 An iterable of integer node id numbers.
2813 Returns
2814 -------
2815 Geometry
2816 A geometry only containing the nodes in `node_list`
2818 """
2819 node = self.node.reduce(node_list).copy()
2820 traceline = self.traceline.reduce(node_list).copy()
2821 element = self.element.reduce(node_list).copy()
2822 css = np.unique(np.concatenate((np.unique(node.def_cs), node.disp_cs)))
2823 coordinate_system = self.coordinate_system(css).copy()
2824 return Geometry(node, coordinate_system, traceline, element)
2826 def add_traceline(self, connectivity, id=None, description='', color=1):
2827 """
2828 Adds a traceline to the geometry
2830 Parameters
2831 ----------
2832 connectivity : iterable
2833 Iterable containing the node_ids to connect with the traceline
2834 id : TYPE, optional
2835 The id number of the new traceline. The default is None, which
2836 results in the new traceline having an id of one greater than the
2837 current maximum traceline it.
2838 description : str, optional
2839 A string description of the new traceline. The default is ''.
2840 color : int, optional
2841 An integer corresponding to the color of the new traceline.
2842 The default is 1.
2844 Returns
2845 -------
2846 None. Modifications are made in-place to the current geometry.
2848 """
2849 if id is None:
2850 if self.traceline.size > 0:
2851 id = np.max(self.traceline.id) + 1
2852 else:
2853 id = 1
2854 new_traceline = traceline_array((id,), description, color,
2855 [np.array(connectivity)])
2856 self.traceline = np.concatenate((self.traceline, new_traceline), axis=0)
2858 def add_element(self, elem_type, connectivity, id=None, color=1):
2859 if id is None:
2860 if self.element.size > 0:
2861 id = np.max(self.element.id) + 1
2862 else:
2863 id = 1
2864 new_element = element_array((id,), elem_type, color, [np.array(connectivity)])
2865 self.element = np.concatenate((self.element, new_element))
2867 def add_node(self, coordinate, id=None, color=1, def_cs=1, disp_cs=1):
2868 """
2869 Adds a node to the geometry
2871 Parameters
2872 ----------
2873 coordinate : iterable
2874 A length-3 iterable containing the node's position in its
2875 definition coordinate system.
2876 id : int, optional
2877 The id number of the new node. The default is None, which results
2878 in the new node having an id of one greater than the current
2879 maximum node id.
2880 color : int, optional
2881 An integer corresponding to the color of the new node.
2882 The default is 1.
2883 def_cs : int, optional
2884 The id of the coordinate system in which the node's position is
2885 defined. The default is 1.
2886 disp_cs : int, optional
2887 The id of the coordinate system in which the node's displacements
2888 are defined. The default is 1.
2890 Returns
2891 -------
2892 None. Modifications are made in-place to the current geometry.
2894 """
2895 if id is None:
2896 if self.node.size > 0:
2897 id = np.max(self.node.id) + 1
2898 else:
2899 id = 1
2900 new_node = node_array((id,), np.array(coordinate, dtype='float64'),
2901 color, def_cs, disp_cs)
2902 self.node = np.concatenate((self.node, new_node), axis=0)
2904 def add_coordinate_system(self, matrix=None, id=None, name='', color=1, cs_type=0):
2905 """
2906 Adds a coordinate system to the geometry
2908 Parameters
2909 ----------
2910 matrix : ndarray, optional
2911 A 4x3 transformation matrix (rows 0-2 are the rotation matrix, row
2912 3 is the origin). The default is None, which results in a global
2913 Cartesian coordinate system (identity rotation, zero origin).
2914 id : int, optional
2915 The id number of the new coordinate system. The default is None,
2916 which results in the new coordinate system having an id of one
2917 greater than the current maximum coordinate system id.
2918 name : str, optional
2919 A string name for the new coordinate system. The default is ''.
2920 color : int, optional
2921 An integer corresponding to the color of the new coordinate system.
2922 The default is 1.
2923 cs_type : int, optional
2924 The coordinate system type (0 = Cartesian, 1 = Polar,
2925 2 = Spherical). The default is 0.
2927 Returns
2928 -------
2929 None. Modifications are made in-place to the current geometry.
2931 """
2932 if id is None:
2933 if self.coordinate_system.size > 0:
2934 id = np.max(self.coordinate_system.id) + 1
2935 else:
2936 id = 1
2937 if matrix is None:
2938 matrix = np.concatenate((np.eye(3), np.zeros((1, 3))), axis=0)
2939 new_cs = coordinate_system_array((id,), name, color, cs_type,
2940 np.array(matrix, dtype='float64'))
2941 self.coordinate_system = np.concatenate((self.coordinate_system, new_cs), axis=0)
2943 def edit(self=None, node_size=10, line_width=2):
2944 """
2945 Open an interactive Qt GUI for editing this geometry in place.
2947 Can also be called on the class itself — ``new_geom =
2948 Geometry.edit()`` — to create a new empty geometry and build it from
2949 scratch in the editor. In that case the new Geometry object is
2950 returned (populated in place as you edit; if the call blocks in the Qt
2951 event loop, it returns once the window is closed).
2953 A window opens with a 3-D view of the geometry (left, resizable via the
2954 splitter) beside editable tables for the nodes, coordinate systems,
2955 tracelines, and elements (right). Each coordinate system is drawn as a
2956 labeled X/Y/Z arrow triad at its origin (labels read e.g. '1X+' for
2957 coordinate system 1's X+ axis). The active tab determines what
2958 clicking in the
2959 3-D view selects — nodes, coordinate systems, tracelines, or elements —
2960 and what the buttons operate on. Hold Shift to select multiple items;
2961 press Escape to deselect everything; changing tabs deselects everything.
2963 On the Nodes and Coord Systems tabs, "Add" immediately appends a new
2964 item. On the Tracelines tab, "Add" is a toggle: while active, click the
2965 first node and Shift-click the rest in order, then "Create" commits
2966 them as a new traceline. On the Elements tab, "Add" is replaced by
2967 three mutually exclusive toggles — "Add Beams" (type 21), "Add
2968 Triangles" (type 41), and "Add Quadrilaterals" (type 44). While one is
2969 active, clicking in the 3-D view picks nodes instead of elements
2970 (plain click restarts the pick list, Shift-click extends it), and each
2971 element is committed automatically as soon as enough nodes have been
2972 picked (2, 3, or 4), ready for the next element of the same type.
2974 All changes are applied to this Geometry object in place, so the object
2975 reflects the edits once the window is closed. This method depends only
2976 on ``qtpy``/``pyvistaqt`` (already required by :meth:`plot`) and on the
2977 ``add_node``, ``add_coordinate_system``, ``add_traceline``, and
2978 ``add_element`` methods, so it is self-contained and can be shared as a
2979 single method.
2981 Parameters
2982 ----------
2983 node_size : int, optional
2984 Pixel size used to display nodes in the 3-D view. The default is 10.
2985 line_width : int, optional
2986 Pixel width used to display tracelines and element edges. The
2987 default is 2.
2989 Returns
2990 -------
2991 window or geometry
2992 When called on an instance (``g.edit()``), the editor QMainWindow
2993 is returned (also stored on ``self._editor_window``). When called
2994 on the class (``Geometry.edit()``), the newly created Geometry is
2995 returned instead.
2996 """
2997 import sys
2998 import pyvista as pv
2999 from qtpy import QtCore
3000 from qtpy.QtCore import QItemSelectionModel
3001 from qtpy.QtWidgets import (
3002 QApplication, QMainWindow, QWidget, QSplitter, QVBoxLayout,
3003 QHBoxLayout, QTabWidget, QTableWidget, QTableWidgetItem, QPushButton,
3004 QLabel, QMessageBox, QAbstractItemView, QComboBox, QStyledItemDelegate,
3005 QInputDialog, QDialog, QFormLayout, QDoubleSpinBox, QDialogButtonBox)
3006 try:
3007 from qtpy.QtWidgets import QShortcut
3008 except ImportError:
3009 from qtpy.QtGui import QShortcut
3010 from pyvistaqt import QtInteractor
3012 # Called on the class (Geometry.edit()) rather than an instance:
3013 # create a new empty geometry to build in the editor and return it
3014 created_new = self is None
3015 if created_new:
3016 self = Geometry()
3017 geometry = self
3018 # A node's global position requires its definition coordinate system to
3019 # exist, so guarantee at least a global Cartesian system is present.
3020 if geometry.coordinate_system.size == 0:
3021 geometry.add_coordinate_system()
3023 app = QApplication.instance()
3024 created_app = app is None
3025 if created_app:
3026 app = QApplication(sys.argv)
3028 # ------------------------------------------------------------------ #
3029 # Editor state and small helpers
3030 # ------------------------------------------------------------------ #
3031 # add_mode: None (selecting), 'traceline', or an element type (21/41/44)
3032 # picked: ordered node ids accumulated while an add mode is active
3033 # selected: ids of selected items of the active tab's kind
3034 state = {'add_mode': None, 'picked': [], 'selected': [],
3035 'show_dofs': False}
3036 flags = {'populating': False}
3037 elem_node_counts = {21: 2, 41: 3, 44: 4}
3039 def _fmt(v):
3040 if isinstance(v, (float, np.floating)):
3041 return repr(float(v))
3042 return str(int(v))
3044 def _conn_to_str(connectivity):
3045 return ', '.join(str(int(c)) for c in np.atleast_1d(connectivity))
3047 def _str_to_conn(text):
3048 return np.array([int(t) for t in text.replace(',', ' ').split()],
3049 dtype='uint32')
3051 # Color names in the order of the sdynpy (I-deas) rendering colormap:
3052 # index 0 = Black ... index 15 = White
3053 color_names = ['Black', 'Blue', 'Gray Blue', 'Light Blue', 'Cyan',
3054 'Dark Olive', 'Dark Green', 'Green', 'Yellow',
3055 'Golden Orange', 'Orange', 'Red', 'Magenta',
3056 'Light Magenta', 'Pink', 'White']
3058 def _color_name(value):
3059 value = int(value)
3060 return color_names[value] if 0 <= value < 16 else str(value)
3062 def _color_value(text):
3063 text = text.strip()
3064 try:
3065 return int(text)
3066 except ValueError:
3067 pass
3068 for value, name in enumerate(color_names):
3069 if name.lower() == text.lower():
3070 return value
3071 raise ValueError('Unknown color: {}'.format(text))
3073 cs_type_names = ['Cartesian', 'Cylindrical', 'Spherical']
3075 def _cs_type_name(value):
3076 value = int(value)
3077 return cs_type_names[value] if 0 <= value < 3 else str(value)
3079 def _cs_type_value(text):
3080 text = text.strip()
3081 try:
3082 return int(text)
3083 except ValueError:
3084 pass
3085 for value, name in enumerate(cs_type_names):
3086 if name.lower() == text.lower():
3087 return value
3088 if text.lower() == 'sypherical': # sdynpy's Excel template spelling
3089 return 2
3090 raise ValueError('Unknown coordinate system type: {}'.format(text))
3092 # ------------------------------------------------------------------ #
3093 # Window scaffolding
3094 # ------------------------------------------------------------------ #
3095 window = QMainWindow()
3096 window.setWindowTitle('SDynPy Geometry Editor')
3097 window.resize(1250, 780)
3098 central = QWidget()
3099 window.setCentralWidget(central)
3100 splitter = QSplitter(QtCore.Qt.Horizontal)
3101 splitter.setHandleWidth(8)
3102 splitter.setChildrenCollapsible(False)
3103 QHBoxLayout(central).addWidget(splitter)
3105 # multi_samples=0: multisampling causes partially-rendered/faded lines
3106 # on some macOS/VTK combinations, and it must be disabled at
3107 # construction (before the GL context is realized) to take effect
3108 interactor = QtInteractor(window, multi_samples=0)
3109 # plot() calls plotter.show(); the interactor is already visible via its
3110 # parent window, so neutralize its show() (headless idiom).
3111 interactor.show = lambda *a, **k: None
3112 interactor.add_axes()
3113 splitter.addWidget(interactor)
3115 panel = QWidget()
3116 panel_layout = QVBoxLayout(panel)
3117 splitter.addWidget(panel)
3118 splitter.setStretchFactor(0, 2)
3119 splitter.setStretchFactor(1, 1)
3120 splitter.setSizes([2 * window.width() // 3, window.width() // 3])
3122 tabs = QTabWidget()
3123 panel_layout.addWidget(tabs)
3124 # Transient feedback (pick progress, errors) lives below the table
3125 status = QLabel('')
3126 panel_layout.addWidget(status)
3128 node_cols = ['ID', 'X', 'Y', 'Z', 'Color', 'Def CS', 'Disp CS']
3129 cs_cols = ['ID', 'Name', 'Color', 'Type', 'Origin X', 'Origin Y', 'Origin Z']
3130 tl_cols = ['ID', 'Color', 'Description', 'Connectivity']
3131 elem_cols = ['ID', 'Type', 'Color', 'Connectivity']
3133 def _make_table(headers):
3134 t = QTableWidget(0, len(headers))
3135 t.setHorizontalHeaderLabels(headers)
3136 t.setSelectionBehavior(QAbstractItemView.SelectRows)
3137 # Excel-style: plain click selects one row (deselecting the rest),
3138 # Shift-click selects the range in between, Cmd-click (Ctrl on
3139 # Windows) adds/removes individual rows
3140 t.setSelectionMode(QAbstractItemView.ExtendedSelection)
3141 t.horizontalHeader().setStretchLastSection(True)
3142 return t
3144 node_table = _make_table(node_cols)
3145 cs_table = _make_table(cs_cols)
3146 tl_table = _make_table(tl_cols)
3147 elem_table = _make_table(elem_cols)
3148 tabs.addTab(node_table, 'Nodes')
3149 tabs.addTab(cs_table, 'Coord Systems')
3150 tabs.addTab(tl_table, 'Tracelines')
3151 tabs.addTab(elem_table, 'Elements')
3152 all_tables = [node_table, cs_table, tl_table, elem_table]
3154 def _ids_for_tab(tab):
3155 return (geometry.node.id, geometry.coordinate_system.id,
3156 geometry.traceline.id, geometry.element.id)[tab]
3158 class _ChoiceDelegate(QStyledItemDelegate):
3159 # Dropdown of named choices for columns stored as integers
3160 def __init__(self, parent, items):
3161 super().__init__(parent)
3162 self._items = items
3164 def createEditor(self, parent, option, index):
3165 combo = QComboBox(parent)
3166 combo.addItems(self._items)
3167 return combo
3169 def setEditorData(self, editor, index):
3170 position = editor.findText(str(index.data() or ''))
3171 if position >= 0:
3172 editor.setCurrentIndex(position)
3174 def setModelData(self, editor, model, index):
3175 model.setData(index, editor.currentText())
3177 # Named-choice columns per table: {column: list of names}.
3178 # Nodes/Tracelines/Elements have a Color column; Coord Systems also
3179 # has the Type column.
3180 choice_columns_by_table = ({4: color_names},
3181 {2: color_names, 3: cs_type_names},
3182 {1: color_names},
3183 {2: color_names})
3184 for table, choice_columns in zip(all_tables, choice_columns_by_table):
3185 for column, items in choice_columns.items():
3186 table.setItemDelegateForColumn(
3187 column, _ChoiceDelegate(table, items))
3189 # ------------------------------------------------------------------ #
3190 # Buttons
3191 # ------------------------------------------------------------------ #
3192 buttons = QHBoxLayout()
3193 panel_layout.addLayout(buttons)
3194 btn_add = QPushButton('Add')
3195 btn_beam = QPushButton('Add Beams')
3196 btn_tri = QPushButton('Add Triangles')
3197 btn_quad = QPushButton('Add Quadrilaterals')
3198 for b in (btn_beam, btn_tri, btn_quad):
3199 b.setCheckable(True)
3200 b.setVisible(False)
3201 btn_create = QPushButton('Create')
3202 btn_create.setVisible(False)
3203 btn_delete = QPushButton('Delete')
3204 btn_dofs = QPushButton('Show DOF Arrows')
3205 btn_dofs.setCheckable(True)
3206 btn_rotate = QPushButton('Rotate Coordinate Systems')
3207 btn_rotate.setVisible(False)
3208 for w in (btn_add, btn_beam, btn_tri, btn_quad, btn_create, btn_delete,
3209 btn_dofs, btn_rotate):
3210 buttons.addWidget(w)
3212 # ------------------------------------------------------------------ #
3213 # Geometry queries for picking
3214 # ------------------------------------------------------------------ #
3215 def _scene_scale():
3216 if geometry.node.size > 1:
3217 extent = np.max(np.ptp(geometry.global_node_coordinate(), axis=0))
3218 if extent > 0:
3219 return extent
3220 return 1.0
3222 def _traceline_segments(connectivity):
3223 # Split a traceline's connectivity on 0 (pen-up) entries
3224 segments, current = [], []
3225 for c in np.atleast_1d(connectivity):
3226 if int(c) == 0:
3227 if len(current) >= 2:
3228 segments.append(np.array(current))
3229 current = []
3230 else:
3231 current.append(int(c))
3232 if len(current) >= 2:
3233 segments.append(np.array(current))
3234 return segments
3236 def _point_segment_distance(p, a, b):
3237 ab = b - a
3238 denom = np.dot(ab, ab)
3239 t = 0.0 if denom == 0 else np.clip(np.dot(p - a, ab) / denom, 0, 1)
3240 return np.linalg.norm(p - (a + t * ab))
3242 def _polyline_distance(p, node_ids, close=False):
3243 ids = [int(i) for i in np.atleast_1d(node_ids)
3244 if np.any(geometry.node.id == i)]
3245 if len(ids) == 0:
3246 return np.inf
3247 points = geometry.global_node_coordinate(np.array(ids))
3248 if len(ids) == 1:
3249 return np.linalg.norm(points[0] - p)
3250 if close and len(points) > 2:
3251 points = np.concatenate((points, points[:1]))
3252 return min(_point_segment_distance(p, points[k], points[k + 1])
3253 for k in range(len(points) - 1))
3255 def _nearest_node(p):
3256 if geometry.node.size == 0:
3257 return None
3258 coords = geometry.global_node_coordinate()
3259 return int(geometry.node.id[
3260 int(np.argmin(np.linalg.norm(coords - p[None, :], axis=1)))])
3262 def _nearest_cs(p):
3263 if geometry.coordinate_system.size == 0:
3264 return None
3265 origins = geometry.coordinate_system.matrix[:, 3, :]
3266 return int(geometry.coordinate_system.id[
3267 int(np.argmin(np.linalg.norm(origins - p[None, :], axis=1)))])
3269 def _nearest_traceline(p):
3270 best, best_distance = None, np.inf
3271 for tl in geometry.traceline:
3272 segments = _traceline_segments(tl.connectivity)
3273 distance = min((_polyline_distance(p, s) for s in segments),
3274 default=np.inf)
3275 if distance < best_distance:
3276 best, best_distance = int(tl.id), distance
3277 return best
3279 def _nearest_element(p):
3280 best, best_distance = None, np.inf
3281 for elem in geometry.element:
3282 connectivity = np.atleast_1d(elem.connectivity)
3283 distance = _polyline_distance(p, connectivity,
3284 close=len(connectivity) > 2)
3285 if distance < best_distance:
3286 best, best_distance = int(elem.id), distance
3287 return best
3289 # ------------------------------------------------------------------ #
3290 # Table <-> geometry synchronization
3291 # ------------------------------------------------------------------ #
3292 def refresh_tables():
3293 flags['populating'] = True
3294 node_table.setRowCount(geometry.node.size)
3295 for r in range(geometry.node.size):
3296 n = geometry.node[r]
3297 for c, v in enumerate((_fmt(n.id), _fmt(n.coordinate[0]),
3298 _fmt(n.coordinate[1]), _fmt(n.coordinate[2]),
3299 _color_name(n.color), _fmt(n.def_cs),
3300 _fmt(n.disp_cs))):
3301 node_table.setItem(r, c, QTableWidgetItem(v))
3303 cs_table.setRowCount(geometry.coordinate_system.size)
3304 for r in range(geometry.coordinate_system.size):
3305 cs = geometry.coordinate_system[r]
3306 cs_table.setItem(r, 0, QTableWidgetItem(_fmt(cs.id)))
3307 cs_table.setItem(r, 1, QTableWidgetItem(str(cs.name)))
3308 cs_table.setItem(r, 2, QTableWidgetItem(_color_name(cs.color)))
3309 cs_table.setItem(r, 3, QTableWidgetItem(_cs_type_name(cs.cs_type)))
3310 for c in range(3):
3311 cs_table.setItem(r, 4 + c, QTableWidgetItem(_fmt(cs.matrix[3, c])))
3313 tl_table.setRowCount(geometry.traceline.size)
3314 for r in range(geometry.traceline.size):
3315 tl = geometry.traceline[r]
3316 tl_table.setItem(r, 0, QTableWidgetItem(_fmt(tl.id)))
3317 tl_table.setItem(r, 1, QTableWidgetItem(_color_name(tl.color)))
3318 tl_table.setItem(r, 2, QTableWidgetItem(str(tl.description)))
3319 tl_table.setItem(r, 3, QTableWidgetItem(_conn_to_str(tl.connectivity)))
3321 elem_table.setRowCount(geometry.element.size)
3322 for r in range(geometry.element.size):
3323 el = geometry.element[r]
3324 elem_table.setItem(r, 0, QTableWidgetItem(_fmt(el.id)))
3325 elem_table.setItem(r, 1, QTableWidgetItem(_fmt(el.type)))
3326 elem_table.setItem(r, 2, QTableWidgetItem(_color_name(el.color)))
3327 elem_table.setItem(r, 3, QTableWidgetItem(_conn_to_str(el.connectivity)))
3328 flags['populating'] = False
3330 def _apply_node(row, col, txt):
3331 if col == 0:
3332 geometry.node.id[row] = int(txt)
3333 elif col in (1, 2, 3):
3334 geometry.node.coordinate[row, col - 1] = float(txt)
3335 elif col == 4:
3336 geometry.node.color[row] = _color_value(txt)
3337 elif col == 5:
3338 geometry.node.def_cs[row] = int(txt)
3339 elif col == 6:
3340 geometry.node.disp_cs[row] = int(txt)
3342 def _apply_cs(row, col, txt):
3343 if col == 0:
3344 geometry.coordinate_system.id[row] = int(txt)
3345 elif col == 1:
3346 geometry.coordinate_system.name[row] = txt
3347 elif col == 2:
3348 geometry.coordinate_system.color[row] = _color_value(txt)
3349 elif col == 3:
3350 geometry.coordinate_system.cs_type[row] = _cs_type_value(txt)
3351 else:
3352 geometry.coordinate_system.matrix[row, 3, col - 4] = float(txt)
3354 def _apply_tl(row, col, txt):
3355 if col == 0:
3356 geometry.traceline.id[row] = int(txt)
3357 elif col == 1:
3358 geometry.traceline.color[row] = _color_value(txt)
3359 elif col == 2:
3360 geometry.traceline.description[row] = txt
3361 elif col == 3:
3362 geometry.traceline.connectivity[row] = _str_to_conn(txt)
3364 def _apply_elem(row, col, txt):
3365 if col == 0:
3366 geometry.element.id[row] = int(txt)
3367 elif col == 1:
3368 geometry.element.type[row] = int(txt)
3369 elif col == 2:
3370 geometry.element.color[row] = _color_value(txt)
3371 elif col == 3:
3372 geometry.element.connectivity[row] = _str_to_conn(txt)
3374 def _remap_node_ids(mapping):
3375 # Renumber node references in traceline/element connectivity after
3376 # node ids change. Applied as one simultaneous mapping so
3377 # overlapping renames (e.g. offsetting 1,2,3 by +1) don't chain.
3378 if not mapping:
3379 return
3380 for idx in range(geometry.traceline.size):
3381 conn = np.atleast_1d(geometry.traceline.connectivity[idx])
3382 geometry.traceline.connectivity[idx] = np.array(
3383 [mapping.get(int(c), int(c)) if int(c) != 0 else 0
3384 for c in conn], dtype='uint64')
3385 for idx in range(geometry.element.size):
3386 conn = np.atleast_1d(geometry.element.connectivity[idx])
3387 geometry.element.connectivity[idx] = np.array(
3388 [mapping.get(int(c), int(c)) for c in conn], dtype='uint64')
3390 def _remap_cs_ids(mapping):
3391 # Renumber coordinate-system references on nodes after cs ids change
3392 if not mapping or geometry.node.size == 0:
3393 return
3394 for field in (geometry.node.def_cs, geometry.node.disp_cs):
3395 for i in range(len(field)):
3396 field[i] = mapping.get(int(field[i]), int(field[i]))
3398 def _change_ids(tab_index, new_id_by_row):
3399 # Assign new ids for the given rows and remap all references.
3400 # Validated up front so a bad value changes nothing.
3401 id_array = _ids_for_tab(tab_index)
3402 mapping = {}
3403 for r, new_id in new_id_by_row.items():
3404 new_id = int(new_id)
3405 if new_id < 0:
3406 raise ValueError('ids must be non-negative')
3407 old_id = int(id_array[r])
3408 if old_id != new_id:
3409 mapping[old_id] = new_id
3410 for r, new_id in new_id_by_row.items():
3411 id_array[r] = int(new_id)
3412 if tab_index == 0:
3413 _remap_node_ids(mapping)
3414 elif tab_index == 1:
3415 _remap_cs_ids(mapping)
3417 def _finish_table_edit(table, tab_index, selected):
3418 # Repopulate (normalizes display, reverts bad input) and restore
3419 # the row selection it clears
3420 refresh_tables()
3421 flags['populating'] = True
3422 selection_model = table.selectionModel()
3423 for r in selected:
3424 if r < table.rowCount():
3425 selection_model.select(
3426 table.model().index(r, 0),
3427 QItemSelectionModel.Select | QItemSelectionModel.Rows)
3428 flags['populating'] = False
3429 ids_arr = _ids_for_tab(tab_index)
3430 state['selected'] = [int(ids_arr[r]) for r in selected
3431 if r < len(ids_arr)]
3432 render()
3434 def _apply_edit(table, tab_index, apply_func, row, col, txt, selected):
3435 # Editing a cell in a multi-row selection applies the new value to
3436 # that column for every selected row
3437 if len(selected) > 1 and row in selected:
3438 apply_rows = selected
3439 else:
3440 apply_rows = [row]
3441 try:
3442 if col == 0:
3443 # id edits also renumber everything referencing those ids
3444 _change_ids(tab_index,
3445 {r: int(txt) for r in apply_rows})
3446 else:
3447 for r in apply_rows:
3448 apply_func(r, col, txt)
3449 except (ValueError, TypeError, OverflowError):
3450 pass
3451 _finish_table_edit(table, tab_index, selected)
3453 def _make_cell_changed(table, tab_index, apply_func):
3454 def handler(row, col):
3455 if flags['populating']:
3456 return
3457 txt = table.item(row, col).text()
3458 selected = sorted({i.row() for i in table.selectedIndexes()})
3459 _apply_edit(table, tab_index, apply_func, row, col, txt, selected)
3460 return handler
3462 def _make_context_edit(table, tab_index, apply_func, choice_columns):
3463 # Right-clicking a cell edits it WITHOUT disturbing the current
3464 # selection (double-click editing collapses the selection first),
3465 # so this is the way to edit a column for many rows at once
3466 def show_dialog(row, col, current, selected):
3467 n_rows = len(selected) if (len(selected) > 1 and row in selected) else 1
3468 header = table.horizontalHeaderItem(col).text()
3469 label = '{} for {} selected row{}:'.format(
3470 header, n_rows, 's' if n_rows > 1 else '')
3471 if col == 0 and n_rows > 1:
3472 # IDs must stay distinct, so a multi-row ID edit adds an
3473 # offset to each selected ID instead of setting them equal
3474 offset, ok = QInputDialog.getInt(
3475 window, 'Add ID Offset',
3476 'Offset to add to {} selected IDs:'.format(n_rows), 0)
3477 if not ok or offset == 0:
3478 return
3479 ids_before = _ids_for_tab(tab_index).copy()
3480 try:
3481 _change_ids(tab_index,
3482 {r: int(ids_before[r]) + offset
3483 for r in selected})
3484 except (ValueError, TypeError, OverflowError):
3485 pass
3486 _finish_table_edit(table, tab_index, selected)
3487 return
3488 if col in choice_columns:
3489 items = choice_columns[col]
3490 start = items.index(current) if current in items else 0
3491 txt, ok = QInputDialog.getItem(
3492 window, 'Edit ' + header, label, items, start, False)
3493 else:
3494 txt, ok = QInputDialog.getText(
3495 window, 'Edit ' + header, label, text=current)
3496 if not ok:
3497 return
3498 _apply_edit(table, tab_index, apply_func, row, col, txt, selected)
3500 def handler(pos):
3501 item = table.itemAt(pos)
3502 if item is None:
3503 return
3504 row, col = item.row(), item.column()
3505 current = item.text()
3506 selected = sorted({i.row() for i in table.selectedIndexes()})
3508 # Defer the modal dialog until the right-click's mouse release
3509 # has actually been processed -- opening it while the button is
3510 # still down leaves the table holding a mouse grab that eats
3511 # every subsequent click in the window (tabs stop responding).
3512 # A 0 ms defer can still fire before the release arrives, so
3513 # poll until no mouse button is pressed.
3514 def open_when_released():
3515 if QApplication.mouseButtons() != QtCore.Qt.NoButton:
3516 QtCore.QTimer.singleShot(50, open_when_released)
3517 return
3518 show_dialog(row, col, current, selected)
3519 QtCore.QTimer.singleShot(0, open_when_released)
3520 return handler
3522 table_wiring = ((node_table, 0, _apply_node),
3523 (cs_table, 1, _apply_cs),
3524 (tl_table, 2, _apply_tl),
3525 (elem_table, 3, _apply_elem))
3526 for (table, tab_index, apply_func), choice_columns in zip(
3527 table_wiring, choice_columns_by_table):
3528 table.cellChanged.connect(
3529 _make_cell_changed(table, tab_index, apply_func))
3530 table.setContextMenuPolicy(QtCore.Qt.CustomContextMenu)
3531 table.customContextMenuRequested.connect(
3532 _make_context_edit(table, tab_index, apply_func, choice_columns))
3534 def _clear_table_selections():
3535 flags['populating'] = True
3536 for t in all_tables:
3537 t.clearSelection()
3538 flags['populating'] = False
3540 def _sync_current_table_selection():
3541 tab = tabs.currentIndex()
3542 table = all_tables[tab]
3543 ids_arr = _ids_for_tab(tab)
3544 flags['populating'] = True
3545 table.clearSelection()
3546 selection_model = table.selectionModel()
3547 for row in range(len(ids_arr)):
3548 if int(ids_arr[row]) in state['selected']:
3549 selection_model.select(
3550 table.model().index(row, 0),
3551 QItemSelectionModel.Select | QItemSelectionModel.Rows)
3552 flags['populating'] = False
3554 def table_selection_changed():
3555 if flags['populating']:
3556 return
3557 tab = tabs.currentIndex()
3558 table = all_tables[tab]
3559 ids_arr = _ids_for_tab(tab)
3560 rows = sorted({i.row() for i in table.selectedIndexes()})
3561 state['selected'] = [int(ids_arr[r]) for r in rows if r < len(ids_arr)]
3562 render()
3563 for t in all_tables:
3564 t.itemSelectionChanged.connect(table_selection_changed)
3566 # ------------------------------------------------------------------ #
3567 # 3-D rendering
3568 # ------------------------------------------------------------------ #
3569 def _draw_red_polyline(connectivity):
3570 for segment in _traceline_segments(connectivity):
3571 valid = [i for i in segment if np.any(geometry.node.id == i)]
3572 if len(valid) < 2:
3573 continue
3574 points = geometry.global_node_coordinate(np.array(valid))
3575 line = pv.PolyData(points,
3576 lines=[len(points)] + list(range(len(points))))
3577 interactor.add_mesh(line, color='red', line_width=line_width + 3,
3578 reset_camera=False)
3580 def _draw_red_points(points):
3581 interactor.add_mesh(pv.PolyData(np.atleast_2d(points)), color='red',
3582 point_size=node_size + 8, reset_camera=False,
3583 render_points_as_spheres=True)
3585 def _draw_all_lines():
3586 # Workaround for a VTK bug on some platforms where scalar-mapped
3587 # (cmap-colored) line cells are only partially rendered: re-draw
3588 # every traceline and beam-type element as direct-colored meshes
3589 # (no scalar mapping), grouped by color. Lines that plot() DID
3590 # draw are simply overdrawn in the identical color.
3591 groups = {}
3592 for tl in geometry.traceline:
3593 for segment in _traceline_segments(tl.connectivity):
3594 groups.setdefault(int(tl.color), []).append(segment)
3595 for elem in geometry.element:
3596 if elem.type in _beam_elem_types:
3597 groups.setdefault(int(elem.color), []).append(
3598 np.atleast_1d(elem.connectivity))
3599 for color_value, segments in groups.items():
3600 all_points, lines, offset = [], [], 0
3601 for segment in segments:
3602 valid = [int(i) for i in segment
3603 if np.any(geometry.node.id == i)]
3604 if len(valid) < 2:
3605 continue
3606 points = geometry.global_node_coordinate(np.array(valid))
3607 all_points.append(points)
3608 lines.extend([len(valid)] + list(range(offset, offset + len(valid))))
3609 offset += len(valid)
3610 if not all_points:
3611 continue
3612 mesh = pv.PolyData(np.concatenate(all_points), lines=lines)
3613 rgb = colormap(np.clip(color_value, 0, 15) / 15)[:3]
3614 interactor.add_mesh(mesh, color=rgb, line_width=line_width,
3615 reset_camera=False)
3617 def _draw_node_dofs():
3618 # Draw labeled DOF arrows (translation directions in each node's
3619 # displacement coordinate system) at every node. A temporary
3620 # geometry with only the nodes and coordinate systems is used so
3621 # plot_coordinate adds just the arrows and labels -- the geometry
3622 # itself is already drawn.
3623 if geometry.node.size == 0:
3624 return
3625 arrow_geometry = Geometry(geometry.node.copy(),
3626 geometry.coordinate_system.copy())
3627 arrow_geometry.plot_coordinate(
3628 label_dofs=True,
3629 arrow_scale=0.1 * _scene_scale(),
3630 arrow_scale_type='absolute',
3631 plot_kwargs={'plotter': interactor, 'node_size': 0,
3632 'line_width': 0})
3634 def _draw_coordinate_systems():
3635 # Draw one labeled X/Y/Z arrow triad per coordinate system at its
3636 # origin by calling plot_coordinate on a temporary geometry that
3637 # has one node per coordinate system (node id = cs id, placed at
3638 # the cs origin, displaced in that cs). Labels read e.g. '1X+'
3639 # for coordinate system 1's X+ axis.
3640 if geometry.coordinate_system.size == 0:
3641 return
3642 css = geometry.coordinate_system
3643 identity_cs_id = int(np.max(css.id)) + 1
3644 triad_nodes = node_array(
3645 css.id.copy(), np.atleast_2d(css.matrix[:, 3, :]).copy(),
3646 1, identity_cs_id, css.id.copy())
3647 triad_geometry = Geometry(
3648 triad_nodes,
3649 np.concatenate((css.copy(),
3650 coordinate_system_array((identity_cs_id,)))))
3651 triad_geometry.plot_coordinate(
3652 label_dofs=True,
3653 arrow_scale=0.1 * _scene_scale(),
3654 arrow_scale_type='absolute',
3655 plot_kwargs={'plotter': interactor, 'node_size': 0,
3656 'line_width': 0})
3658 def _draw_highlights():
3659 tab = tabs.currentIndex()
3660 highlight_nodes = list(state['picked'])
3661 if state['add_mode'] is None:
3662 if tab == 0:
3663 highlight_nodes.extend(state['selected'])
3664 elif tab == 1:
3665 rows = [int(np.argmax(geometry.coordinate_system.id == i))
3666 for i in state['selected']
3667 if np.any(geometry.coordinate_system.id == i)]
3668 if rows:
3669 _draw_red_points(
3670 geometry.coordinate_system.matrix[rows, 3, :])
3671 elif tab == 2:
3672 for i in state['selected']:
3673 match = geometry.traceline.id == i
3674 if np.any(match):
3675 _draw_red_polyline(
3676 geometry.traceline.connectivity[np.argmax(match)])
3677 elif tab == 3:
3678 for i in state['selected']:
3679 match = geometry.element.id == i
3680 if not np.any(match):
3681 continue
3682 row = int(np.argmax(match))
3683 connectivity = np.atleast_1d(
3684 geometry.element.connectivity[row])
3685 valid = [c for c in connectivity
3686 if np.any(geometry.node.id == c)]
3687 if len(valid) < 3 or geometry.element.type[row] in _beam_elem_types:
3688 _draw_red_polyline(connectivity)
3689 else:
3690 points = geometry.global_node_coordinate(np.array(valid))
3691 face = pv.PolyData(
3692 points, faces=[len(points)] + list(range(len(points))))
3693 # Translucent so it draws after (on top of) the
3694 # coincident opaque element face
3695 interactor.add_mesh(face, color='red', opacity=0.7,
3696 reset_camera=False)
3697 valid = [n for n in dict.fromkeys(highlight_nodes)
3698 if np.any(geometry.node.id == n)]
3699 if valid:
3700 _draw_red_points(geometry.global_node_coordinate(np.array(valid)))
3702 def render():
3703 try:
3704 cam = interactor.camera_position
3705 except Exception:
3706 cam = None
3707 interactor.clear()
3708 try:
3709 if geometry.node.size > 0:
3710 geometry.plot(plotter=interactor, node_size=node_size,
3711 line_width=line_width, label_nodes=True)
3712 _draw_all_lines()
3713 # Coordinate-system triads only show on the Coord Systems tab
3714 if tabs.currentIndex() == 1:
3715 _draw_coordinate_systems()
3716 # Per-node DOF arrows only show on the Nodes tab when toggled
3717 if tabs.currentIndex() == 0 and state['show_dofs']:
3718 _draw_node_dofs()
3719 _draw_highlights()
3720 except Exception as e:
3721 status.setText('Render error: {}'.format(e))
3722 try:
3723 if cam is not None:
3724 interactor.camera_position = cam
3725 interactor.render()
3726 except Exception as e:
3727 status.setText('Render error: {}'.format(e))
3729 # ------------------------------------------------------------------ #
3730 # 3-D picking (tab-dependent selection / add-mode node picking)
3731 # ------------------------------------------------------------------ #
3732 def on_pick(point, *args):
3733 if point is None:
3734 return
3735 p = np.asarray(point, dtype=float).ravel()[:3]
3736 # Shift or Cmd (Ctrl on Windows; Qt maps Command to ControlModifier
3737 # on macOS) extends the selection/pick instead of restarting it
3738 shift = bool(QApplication.keyboardModifiers()
3739 & (QtCore.Qt.ShiftModifier | QtCore.Qt.ControlModifier))
3740 if state['add_mode'] is not None:
3741 nid = _nearest_node(p)
3742 if nid is None:
3743 return
3744 # Same selection rule as everywhere else: a plain click starts
3745 # over with just this node; Shift-click extends the pick list
3746 if shift:
3747 if not state['picked'] or state['picked'][-1] != nid:
3748 state['picked'].append(nid)
3749 else:
3750 state['picked'] = [nid]
3751 mode = state['add_mode']
3752 if mode in elem_node_counts and \
3753 len(state['picked']) == elem_node_counts[mode]:
3754 geometry.add_element(mode, state['picked'])
3755 state['picked'] = []
3756 status.setText('Added element {} (type {})'.format(
3757 int(geometry.element.id[-1]), mode))
3758 refresh_tables()
3759 else:
3760 status.setText('Picked nodes: {}'.format(state['picked']))
3761 render()
3762 return
3763 tab = tabs.currentIndex()
3764 item_id = (_nearest_node, _nearest_cs,
3765 _nearest_traceline, _nearest_element)[tab](p)
3766 if item_id is None:
3767 return
3768 if shift:
3769 if item_id in state['selected']:
3770 state['selected'].remove(item_id)
3771 else:
3772 state['selected'].append(item_id)
3773 else:
3774 state['selected'] = [item_id]
3775 _sync_current_table_selection()
3776 render()
3778 interactor.enable_point_picking(callback=on_pick, show_message=False,
3779 left_clicking=True, picker='cell',
3780 show_point=False)
3782 # ------------------------------------------------------------------ #
3783 # Button and mode handlers
3784 # ------------------------------------------------------------------ #
3785 def deselect_all():
3786 state['selected'] = []
3787 state['picked'] = []
3788 _clear_table_selections()
3789 status.setText('')
3790 render()
3792 def add_clicked(checked=False):
3793 tab = tabs.currentIndex()
3794 if tab == 0:
3795 if geometry.node.size > 0:
3796 coord = geometry.global_node_coordinate().mean(axis=0)
3797 else:
3798 coord = np.zeros(3)
3799 geometry.add_node(coord)
3800 refresh_tables()
3801 render()
3802 elif tab == 1:
3803 geometry.add_coordinate_system()
3804 refresh_tables()
3805 render()
3806 elif tab == 2:
3807 if checked:
3808 state['add_mode'] = 'traceline'
3809 state['picked'] = []
3810 state['selected'] = []
3811 _clear_table_selections()
3812 status.setText('Click the first node, Shift-click the rest '
3813 'in order, then press Create.')
3814 else:
3815 state['add_mode'] = None
3816 state['picked'] = []
3817 status.setText('')
3818 render()
3820 def elem_mode_handler(button, elem_type):
3821 def handler(checked):
3822 if checked:
3823 for b in (btn_beam, btn_tri, btn_quad):
3824 if b is not button:
3825 b.blockSignals(True)
3826 b.setChecked(False)
3827 b.blockSignals(False)
3828 state['add_mode'] = elem_type
3829 state['picked'] = []
3830 state['selected'] = []
3831 _clear_table_selections()
3832 status.setText(
3833 'Click the first node, then Shift-click {} more per '
3834 'element (type {}).'.format(
3835 elem_node_counts[elem_type] - 1, elem_type))
3836 else:
3837 state['add_mode'] = None
3838 state['picked'] = []
3839 status.setText('')
3840 render()
3841 return handler
3843 def create_clicked():
3844 if state['add_mode'] == 'traceline' and len(state['picked']) >= 2:
3845 geometry.add_traceline(state['picked'])
3846 state['picked'] = []
3847 status.setText('Traceline added; keep clicking to build another.')
3848 refresh_tables()
3849 render()
3851 def tab_changed(index):
3852 state['add_mode'] = None
3853 state['picked'] = []
3854 state['selected'] = []
3855 for b in (btn_beam, btn_tri, btn_quad):
3856 b.blockSignals(True)
3857 b.setChecked(False)
3858 b.blockSignals(False)
3859 b.setVisible(index == 3)
3860 btn_add.blockSignals(True)
3861 btn_add.setChecked(False)
3862 btn_add.blockSignals(False)
3863 btn_add.setCheckable(index == 2)
3864 btn_add.setVisible(index != 3)
3865 btn_create.setVisible(index == 2)
3866 btn_dofs.setVisible(index == 0)
3867 btn_rotate.setVisible(index == 1)
3868 _clear_table_selections()
3869 status.setText('')
3870 render()
3872 def _prune_connectivity(remove_ids):
3873 remove = set(int(x) for x in np.atleast_1d(remove_ids))
3874 if geometry.traceline.size > 0:
3875 keep = [not (set(int(c) for c in tl.connectivity) & remove)
3876 for tl in geometry.traceline]
3877 geometry.traceline = geometry.traceline[np.array(keep, dtype=bool)]
3878 if geometry.element.size > 0:
3879 keep = [not (set(int(c) for c in el.connectivity) & remove)
3880 for el in geometry.element]
3881 geometry.element = geometry.element[np.array(keep, dtype=bool)]
3883 def _selected_rows(table):
3884 return sorted({i.row() for i in table.selectedIndexes()})
3886 def delete_selected():
3887 tab = tabs.currentIndex()
3888 if tab == 0:
3889 rows = _selected_rows(node_table)
3890 if not rows:
3891 return
3892 remove_ids = geometry.node.id[rows]
3893 keep = np.ones(geometry.node.size, dtype=bool)
3894 keep[rows] = False
3895 geometry.node = geometry.node[keep]
3896 _prune_connectivity(remove_ids)
3897 elif tab == 1:
3898 rows = _selected_rows(cs_table)
3899 if not rows:
3900 return
3901 cs_ids = geometry.coordinate_system.id[rows]
3902 used = np.concatenate((geometry.node.def_cs, geometry.node.disp_cs)) \
3903 if geometry.node.size > 0 else np.array([], dtype='uint64')
3904 in_use = np.isin(cs_ids, used)
3905 if in_use.any():
3906 QMessageBox.warning(
3907 window, 'Coordinate system in use',
3908 'Cannot delete coordinate system(s) still referenced by '
3909 'nodes: {}'.format(cs_ids[in_use].tolist()))
3910 return
3911 keep = np.ones(geometry.coordinate_system.size, dtype=bool)
3912 keep[rows] = False
3913 geometry.coordinate_system = geometry.coordinate_system[keep]
3914 elif tab == 2:
3915 rows = _selected_rows(tl_table)
3916 if not rows:
3917 return
3918 keep = np.ones(geometry.traceline.size, dtype=bool)
3919 keep[rows] = False
3920 geometry.traceline = geometry.traceline[keep]
3921 elif tab == 3:
3922 rows = _selected_rows(elem_table)
3923 if not rows:
3924 return
3925 keep = np.ones(geometry.element.size, dtype=bool)
3926 keep[rows] = False
3927 geometry.element = geometry.element[keep]
3928 state['selected'] = []
3929 state['picked'] = []
3930 refresh_tables()
3931 render()
3933 btn_add.clicked.connect(add_clicked)
3934 btn_beam.toggled.connect(elem_mode_handler(btn_beam, 21))
3935 btn_tri.toggled.connect(elem_mode_handler(btn_tri, 41))
3936 btn_quad.toggled.connect(elem_mode_handler(btn_quad, 44))
3937 def toggle_dofs(checked):
3938 state['show_dofs'] = checked
3939 render()
3941 def rotate_cs():
3942 rows = _selected_rows(cs_table)
3943 if not rows:
3944 status.setText('Select coordinate system(s) to rotate first.')
3945 return
3946 dialog = QDialog(window)
3947 dialog.setWindowTitle('Rotate Coordinate Systems')
3948 form = QFormLayout(dialog)
3949 angle_box = QDoubleSpinBox()
3950 angle_box.setRange(-360.0, 360.0)
3951 angle_box.setDecimals(3)
3952 angle_box.setSuffix(' deg')
3953 axis_combo = QComboBox()
3954 axis_combo.addItems(['X', 'Y', 'Z'])
3955 form.addRow('Rotation angle:', angle_box)
3956 form.addRow('Rotation axis:', axis_combo)
3957 buttons_box = QDialogButtonBox(
3958 QDialogButtonBox.Ok | QDialogButtonBox.Cancel)
3959 buttons_box.accepted.connect(dialog.accept)
3960 buttons_box.rejected.connect(dialog.reject)
3961 form.addRow(buttons_box)
3962 if not dialog.exec_():
3963 return
3964 angle = angle_box.value()
3965 axis_index = axis_combo.currentIndex()
3966 # Same convention as the Excel-template reader:
3967 # new rotation = R(axis, angle, degrees=True).T @ existing
3968 rotation = R(axis_index, angle, degrees=True).T
3969 for r in rows:
3970 geometry.coordinate_system.matrix[r, :3, :] = \
3971 rotation @ geometry.coordinate_system.matrix[r, :3, :]
3972 status.setText('Rotated {} coordinate system(s) {} deg about {}.'.format(
3973 len(rows), angle, axis_combo.currentText()))
3974 _finish_table_edit(cs_table, 1, rows)
3976 btn_create.clicked.connect(create_clicked)
3977 btn_delete.clicked.connect(delete_selected)
3978 btn_dofs.toggled.connect(toggle_dofs)
3979 btn_rotate.clicked.connect(rotate_cs)
3980 tabs.currentChanged.connect(tab_changed)
3981 escape_shortcut = QShortcut(QKeySequence(QtCore.Qt.Key_Escape), window)
3982 escape_shortcut.activated.connect(deselect_all)
3983 try:
3984 interactor.add_key_event('Escape', deselect_all)
3985 except Exception:
3986 pass
3988 # ------------------------------------------------------------------ #
3989 # Launch
3990 # ------------------------------------------------------------------ #
3991 # Shut the VTK render window down cleanly when the editor closes --
3992 # otherwise its render timer keeps firing into a destroyed native
3993 # window, which spams wglMakeCurrent/'handle is invalid' errors on
3994 # Windows the next time the editor is opened
3995 original_close_event = window.closeEvent
3997 def close_event(event):
3998 try:
3999 interactor.close()
4000 except Exception:
4001 pass
4002 original_close_event(event)
4003 window.closeEvent = close_event
4005 refresh_tables()
4006 window.show()
4007 window.raise_()
4008 render()
4009 interactor.reset_camera()
4010 self._editor_window = window
4011 if created_app:
4012 app.exec_()
4013 return geometry if created_new else window
4015 @classmethod
4016 def from_unv(cls, unv_dict):
4017 """
4018 Create a geometry from universal file format data from readunv
4020 Parameters
4021 ----------
4022 unv_dict : dict
4023 Dictionary containing data from read_unv
4025 Returns
4026 -------
4027 Geometry
4028 Geometry object read from the unv data
4030 """
4031 try:
4032 coordinate_systems = CoordinateSystemArray.from_unv(unv_dict, combine=True)
4033 except KeyError:
4034 coordinate_systems = CoordinateSystemArray((0,))
4035 # Nodes
4036 try:
4037 nodes = NodeArray.from_unv(unv_dict, combine=True)
4038 cs_indices = coordinate_systems.id
4039 for node_index, node in nodes.ndenumerate():
4040 cs_index = np.where(cs_indices == node.def_cs)
4041 matching_cs = coordinate_systems[cs_index]
4042 if matching_cs.size == 0:
4043 print('Warning: Coordinate System {:} not found in universal file structure. Adding Global Cartesian with ID {:}'.format(
4044 node.def_cs, node.def_cs))
4045 coordinate_systems = np.concatenate(
4046 (coordinate_systems, coordinate_system_array((node.def_cs,))))
4047 cs_indices = coordinate_systems.id
4048 # cs_index = np.where(cs_indices == node.def_cs)
4049 # node.coordinate = local_coord(coordinate_systems[cs_index][0],node.coordinate)
4050 except KeyError:
4051 nodes = NodeArray((0,))
4052 try:
4053 elements = ElementArray.from_unv(unv_dict, combine=True)
4054 except KeyError:
4055 elements = ElementArray((0,))
4056 try:
4057 tracelines = TracelineArray.from_unv(unv_dict, combine=True)
4058 except KeyError:
4059 tracelines = TracelineArray((0,))
4061 return cls(nodes, coordinate_systems, tracelines, elements)
4063 from_uff = from_unv
4065 @staticmethod
4066 def write_excel_template(path_to_xlsx):
4067 """
4068 Writes an Excel File Template for Creating Geometry
4070 Parameters
4071 ----------
4072 path_to_xlsx : string
4073 Path to write xlsx Excel file
4075 Returns
4076 -------
4077 Nothing
4079 Notes
4080 -----
4081 See documentation for from_excel_template for instructions on filling
4082 out the template to create a geometry.
4083 """
4085 # Create a new workbook
4086 workbook = openpyxl.Workbook()
4088 # Add multiple sheets to the workbook
4089 sheet_names = ['Coordinate Systems', 'Nodes', 'Elements', 'Trace Lines']
4090 for sheet_name in sheet_names:
4091 workbook.create_sheet(sheet_name)
4093 # Remove the default sheet created at the start
4094 del workbook['Sheet']
4096 # Define column names for each sheet
4097 columns = {
4098 'Coordinate Systems': ['ID', 'Name', 'Color', 'Type', 'X Location', 'Y Location', 'Z Location', 'Rotation 1 Angle', 'Rotation 1 Axis', 'Rotation 2 Angle', 'Rotation 2 Axis', 'Rotation 3 Angle', 'Rotation 3 Axis'],
4099 'Nodes': ['ID', 'Color', 'X Location', 'Y Location', 'Z Location', 'Displacement CS', 'Definition CS'],
4100 'Elements': ['ID', 'Color', 'Type', 'Node 1', 'Node 2', 'Node 3', 'Node 4', 'Node 5'],
4101 'Trace Lines': ['ID', 'Description', 'Color', 'Node 1', 'Node 2', 'Node 3', 'Node 4', 'Node 5']
4102 }
4104 # Create data validations
4105 data_validation_ID = openpyxl.worksheet.datavalidation.DataValidation(type="whole", operator='greaterThan', formula1='-1',allow_blank=True,error='Must be interger greater than -1',errorStyle='stop',showErrorMessage=True)
4106 data_validation_Color = openpyxl.worksheet.datavalidation.DataValidation(type="list", formula1='"Black,Blue,Gray Blue,Light Blue,Cyan,Dark Olive,Dark Green,Green,Yellow,Golden Orange,Orange,Red,Magenta,Light Magenta,Pink,White"',allow_blank=True,error='Must be value from list',errorStyle='stop',showErrorMessage=True)
4107 data_validation_CSType = openpyxl.worksheet.datavalidation.DataValidation(type="list", formula1='"Cartesian,Cylindrical,Sypherical"',allow_blank=True,error='Must be value from list',errorStyle='stop',showErrorMessage=True)
4108 data_validation_Axis = openpyxl.worksheet.datavalidation.DataValidation(type="list", formula1='"X,Y,Z"',allow_blank=True,error='Must be value from list',errorStyle='stop',showErrorMessage=True)
4109 data_validation_ElemType = openpyxl.worksheet.datavalidation.DataValidation(type="list", formula1='"beam,triangle,quadrilateral,tetrahedral"',allow_blank=True,error='Must be value from list',errorStyle='stop',showErrorMessage=True)
4110 data_validation_Node = openpyxl.worksheet.datavalidation.DataValidation(type="whole", operator='greaterThan', formula1='-1',allow_blank=True,error='Must be value from list',errorStyle='stop',showErrorMessage=True)
4111 data_validation_Decimal = openpyxl.worksheet.datavalidation.DataValidation(type="decimal",allow_blank=True,error='Must be decimal',errorStyle='stop',showErrorMessage=True)
4113 # Specify Columns for data validation
4114 table_height = 10
4115 data_validation_ID_A = copy.deepcopy(data_validation_ID)
4116 data_validation_ID_A.add('A2:A'+str(table_height+1))
4117 data_validation_Color_B = copy.deepcopy(data_validation_Color)
4118 data_validation_Color_B.add('B2:B'+str(table_height+1))
4119 data_validation_Color_C = copy.deepcopy(data_validation_Color)
4120 data_validation_Color_C.add('C2:C'+str(table_height+1))
4121 data_validation_CSType.add('D2:D'+str(table_height+1))
4122 data_validation_Axis.add('I2:I'+str(table_height+1))
4123 data_validation_Axis.add('K2:K'+str(table_height+1))
4124 data_validation_Axis.add('M2:M'+str(table_height+1))
4125 data_validation_Decimal_Coordinate_System = data_validation_Decimal
4126 data_validation_Decimal_Coordinate_System.add('E2:H'+str(table_height+1))
4127 data_validation_Decimal_Coordinate_System.add('J2:J'+str(table_height+1))
4128 data_validation_Decimal_Coordinate_System.add('L2:L'+str(table_height+1))
4129 data_validation_Decimal_Nodes = copy.deepcopy(data_validation_Decimal)
4130 data_validation_Decimal_Nodes.add('C2:E'+str(table_height+1))
4131 data_validation_ID_Nodes = copy.deepcopy(data_validation_ID)
4132 data_validation_ID_Nodes.add('F2:G'+str(table_height+1))
4133 data_validation_ElemType.add('C2:C'+str(table_height+1))
4134 data_validation_Node.add('D2:H'+str(table_height+1))
4136 # Add data and tables to each sheet
4137 for sheet_name in sheet_names:
4138 sheet = workbook[sheet_name]
4140 # Define the column headers
4141 headers = columns[sheet_name]
4143 # Write the headers into the first row
4144 for col_num, header in enumerate(headers, start=1):
4145 sheet.cell(row=1, column=col_num, value=header)
4147 table_range = f'A1:{chr(65 + len(headers) - 1)}{table_height + 1}'
4149 table = openpyxl.worksheet.table.Table(displayName=sheet_name.replace(' ','_'), ref=table_range)
4151 # Apply table style
4152 style = openpyxl.worksheet.table.TableStyleInfo(name="TableStyleMedium9", showFirstColumn=False,showLastColumn=False, showRowStripes=True, showColumnStripes=False)
4153 table.tableStyleInfo = style
4155 # Add table to the worksheet
4156 sheet.add_table(table)
4158 # Apply Data Validations for "Coordinate Systems" Sheet
4159 sheet = workbook['Coordinate Systems']
4160 sheet.add_data_validation(data_validation_ID_A)
4161 sheet.add_data_validation(data_validation_Color_C)
4162 sheet.add_data_validation(data_validation_CSType)
4163 sheet.add_data_validation(data_validation_Decimal_Coordinate_System)
4164 sheet.add_data_validation(data_validation_Axis)
4166 # Apply Data Validations for "Nodes" Sheet
4167 sheet = workbook['Nodes']
4168 sheet.add_data_validation(data_validation_ID_A)
4169 sheet.add_data_validation(data_validation_ID_Nodes)
4170 sheet.add_data_validation(data_validation_Color_B)
4171 sheet.add_data_validation(data_validation_Decimal_Nodes)
4173 # Apply Data Validations for "Elements" Sheet
4174 sheet = workbook['Elements']
4175 sheet.add_data_validation(data_validation_ID_A)
4176 sheet.add_data_validation(data_validation_Color_B)
4177 sheet.add_data_validation(data_validation_ElemType)
4178 sheet.add_data_validation(data_validation_Node)
4180 # Apply Data Validations for "Trace Lines" Sheet
4181 sheet = workbook['Trace Lines']
4182 sheet.add_data_validation(data_validation_ID_A)
4183 sheet.add_data_validation(data_validation_Color_C)
4184 sheet.add_data_validation(data_validation_Node)
4186 # Save the workbook
4187 workbook.save(path_to_xlsx)
4189 @classmethod
4190 def from_excel_template(cls, path_to_xlsx):
4191 """
4192 Create a geometry from Excel file template
4194 Parameters
4195 ----------
4196 path_to_xlsx : string
4197 Path to xlsx Excel file containing geometry information
4199 Returns
4200 -------
4201 Geometry
4202 Geometry object created from the Excel file
4204 Notes
4205 -----
4206 To use this function, first save out an excel template file using the
4207 `write_excel_template` function. This will construct an excel
4208 workbook with four worksheets on which the different portions of the
4209 Geometry are defined.
4211 On the Coordinate Systems tab, users will define the various global and
4212 local coordinate systems in their geometry. Each geometry requires
4213 an ID number. A Name can optionally be given. The Color should be
4214 selectred from the dropdown list of colors in the Excel template. The
4215 type of the coordinate system should be selected from the dropdown list of
4216 coordinate system types.
4218 The origin of the coordinate system can be specified with the X Location,
4219 Y Location, Z Location columns. Then rotations of the coordinate system
4220 can be specified using rotations about axes. Up to three axes and angles
4221 can be specified to create arbitrary compound rotations. The rotation
4222 axes should be X, Y, or Z, and the rotation angles are in degrees.
4224 On the Nodes tab, all tabs must be filled out. ID numbers must be unique.
4225 Colors should be selected from the dropdown list in the Excel template. The
4226 position of the node is specified using the X Location, Y Location, Z
4227 Location columns. Each node has a displacement coordinate system in which
4228 its position is defined, and a a definition coordinate system in which
4229 its displacements are defined. These columns should consist of integers
4230 corresponding to ID numbers from the Coordinate Systems tab.
4232 On the Elements tab, elements connecting nodes are defined. The ID
4233 column must consist of unique integer identifiers. The Color tab should
4234 be selected from the dropdown list. The type can be selected from the dropdown list.
4235 If the type column
4236 is empty, an element type based on the number of connections given will
4237 be used. Defult element types for connection length of 2 is "Type 21 - Linear Beam",
4238 for a connection length of 3 is "Type 41 - Plane Stress Linear Triangle",
4239 and for a connection length of 4 is "Type 44 - Plane Stress Linear Quadrilateral".
4240 The columns Node 1 through Node 4 contain the nodes in each element. Only the required
4241 number of nodes must be filled out (e.g. a tri element would only contain
4242 3 nodes, so only columns Node 1, Node 2, and Node 3 would be filled).
4244 On the Trace Lines tab, lines connecting the nodes are defined. The
4245 ID column must consist of unique integer identifiers. The Description
4246 column contains a string description of the line. The Color should be
4247 selected from the dropdown list. The
4248 Node 1 through Node 5 columns should contain the nodes for each line.
4249 Only the number of nodes in the line must be filled out, so if a line
4250 only connects 5 nodes, only Node 1 though Node 5 must be filled. More columns can
4251 can be added for additional Nodes if longer tracelines are needed.
4254 """
4256 # % Read Geometry Information from Excel File
4257 df_coordinate_systems = pd.read_excel(io = path_to_xlsx,sheet_name='Coordinate Systems')
4258 df_nodes = pd.read_excel(io = path_to_xlsx,sheet_name='Nodes')
4259 df_elements = pd.read_excel(io = path_to_xlsx,sheet_name='Elements')
4260 df_trace_lines = pd.read_excel(io = path_to_xlsx,sheet_name='Trace Lines')
4262 color_number_dict = {
4263 'Black':1,
4264 'Blue':2,
4265 'Gray Blue':3,
4266 'Light Blue':4,
4267 'Cyan':5,
4268 'Dark Olive':6,
4269 'Dark Green':7,
4270 'Green':8,
4271 'Yellow':9,
4272 'Golden Orange':10,
4273 'Orange':11,
4274 'Red':12,
4275 'Magenta':13,
4276 'Light Magenta':14,
4277 'Pink':15,
4278 'White':16,
4279 }
4281 cs_type_dict = {
4282 'Cartesian':0,
4283 'Cylindrical':1,
4284 'Sypherical':2,
4285 }
4287 df_coordinate_systems = df_coordinate_systems.replace({'Color':color_number_dict})
4288 df_coordinate_systems = df_coordinate_systems.replace({'Type':cs_type_dict})
4289 df_nodes = df_nodes.replace({'Color':color_number_dict})
4290 df_elements = df_elements.replace({'Color':color_number_dict})
4291 df_trace_lines = df_trace_lines.replace({'Color':color_number_dict})
4293 # % Add Coordinate Systems
4294 # Start with an empty coordinate system
4295 coordinate_systems = coordinate_system_array([])
4296 for row_index, df_coordinate_system in df_coordinate_systems.iterrows():
4297 # Create Transformation Matrix
4298 T = np.eye(3)
4299 for rotation in range(1,(len(df_coordinate_system)-7)//2+1):
4300 rotation_angle = df_coordinate_system['Rotation ' + str(rotation) + ' Angle']
4301 rotation_axis = df_coordinate_system['Rotation ' + str(rotation) + ' Axis']
4302 rotation_axis_dict = {'X':0,'Y':1,'Z':2}
4303 if pd.notnull(rotation_angle) and pd.notnull(rotation_axis):
4304 T = R(rotation_axis_dict[rotation_axis],rotation_angle,degrees=True).T@T
4306 # Add Coordinate System to Geometry
4307 origin = np.array((df_coordinate_system['X Location'],df_coordinate_system['Y Location'],df_coordinate_system['Z Location']))[np.newaxis,:]
4308 coordinate_systems = np.concatenate((coordinate_systems,coordinate_system_array(id=df_coordinate_system['ID'], name=df_coordinate_system['Name'], color=df_coordinate_system['Color']-1, cs_type=df_coordinate_system['Type'], matrix = np.concatenate((T,origin),axis=0))),axis=None)
4310 # % Add Nodes
4311 # print(df_nodes['Color'])
4312 nodes = node_array(
4313 id=np.array(df_nodes['ID']),
4314 coordinate=np.concatenate(
4315 [np.array(df_nodes['X Location'])[:,np.newaxis],np.array(df_nodes['Y Location'])[:,np.newaxis],np.array(df_nodes['Z Location'])[:,np.newaxis]],axis=-1),
4316 color=np.array(df_nodes['Color'])-1,
4317 def_cs=np.array(df_nodes['Definition CS']),
4318 disp_cs=np.array(df_nodes['Displacement CS']))
4320 # % Add Elements
4321 elements = element_array([])
4322 element_connection_length_dict = {
4323 2:21, # 21 : 'Linear Beam'
4324 3:41, # 41 : 'Plane Stress Linear Triangle'
4325 4:44, # 44 : 'Plane Stress Linear Quadrilateral'
4326 }
4327 for row_index, df_element in df_elements.iterrows():
4328 connectivity = df_element[3:].values
4329 connectivity = connectivity[pd.notnull(connectivity)]
4330 if isinstance(df_element['Type'],str):
4331 element_type = _exodus_elem_type_map[df_element['Type'].lower()]
4332 elif isinstance(df_element['Type'],int):
4333 element_type = df_element['Type']
4334 elif pd.isnull(df_element['Type']):
4335 element_type = element_connection_length_dict[len(connectivity)]
4336 new_element = element_array(id=df_element['ID'], type=element_type, color=df_element['Color']-1,connectivity=connectivity)
4337 elements = np.concatenate((elements,new_element),axis=None)
4339 # % Add Tracelines
4340 tracelines = traceline_array([])
4341 for row_index, df_trace_line in df_trace_lines.iterrows():
4342 connectivity=df_trace_line[3:].values
4343 connectivity = connectivity[pd.notnull(connectivity)]
4344 new_traceline = traceline_array(connectivity=connectivity,id=df_trace_line['ID'],description=df_trace_line['Description'],color=df_trace_line['Color']-1)
4345 tracelines = np.concatenate((tracelines,new_traceline),axis=None)
4347 return cls(nodes, coordinate_systems, tracelines, elements)
4349 @classmethod
4350 def from_exodus(cls, exo: Exodus, blocks=None, local=False,
4351 preferred_local_orientation=np.array((0.0, 0.0, 1.0)),
4352 secondary_preferred_local_orientation=np.array((1.0, 0.0, 0.0)),
4353 local_nodes=None):
4354 """
4355 Generate a geometry from exodus file data
4357 Parameters
4358 ----------
4359 exo : Exodus or ExodusInMemory
4360 The exodus data from which geometry will be created.
4361 blocks : iterable, optional
4362 Iterable containing a set of block ids to import when creating the
4363 geometry. The default is None, which uses all blocks in the model.
4364 local : bool, optional
4365 Flag to specify whether or not to create local coordinate systems
4366 for each node in the model. This can be useful when creating
4367 instrumnetation positions from a finite element model, where the
4368 sensor will be oriented perpendicular to the surface it is mounted
4369 on. The default is False, which returns all data in a single
4370 global coordinate system.
4371 preferred_local_orientation : np.ndarray, optional
4372 A preferred direction for the local coordinate system. The first
4373 constraint is that the local Z+ axis is perpendicular to the
4374 surface. The coordinate system will then try to align itself as
4375 much as it can with this direction. The default is
4376 np.array((0.0,0.0,1.0)), which points the local coordinate systems
4377 along the global Z+ axis.
4378 secondary_preferred_local_orientation : np.ndarray, optional
4379 A secondary preferred direction is only used if the surface normal
4380 direction is parallel to the primary preferred_local_orientation.
4381 The default is np.array((1.0,0.0,0.0)), which points the local
4382 coordinate system along the local Z+ axis.
4383 local_nodes : np.ndarray, optional
4384 If specified, only create local coordinate systems at the specified
4385 nodes.
4387 Returns
4388 -------
4389 Geometry
4390 A geometry consisting of the finite element nodes and element
4391 connectivity.
4393 """
4394 if isinstance(exo, Exodus):
4395 exo = exo.load_into_memory(close=False, variables=[], timesteps=[], blocks=blocks)
4396 # Get nodes
4397 coordinates = exo.nodes.coordinates
4398 node_map = np.arange(
4399 coordinates.shape[0]) + 1 if exo.nodes.node_num_map is None else exo.nodes.node_num_map
4400 nodes = node_array(node_map, coordinates)
4401 tl_connectivity = []
4402 elem_connectivity = []
4403 tl_color = []
4404 elem_color = []
4405 elem_types = []
4406 if blocks is None:
4407 blocks = exo.blocks
4408 else:
4409 blocks = [block for block in exo.blocks if block.id in blocks]
4410 for i, block in enumerate(blocks):
4411 color_index = (i % 14) + 1
4412 try:
4413 connectivity = node_map[block.connectivity]
4414 except AttributeError:
4415 continue
4416 elem_type = block.elem_type
4417 if connectivity.shape[-1] <= 1:
4418 # Skip conmasses
4419 continue
4420 for element in connectivity:
4421 elem_connectivity.append(element)
4422 elem_color.append(color_index)
4423 elem_types.append(_exodus_elem_type_map[elem_type.lower()])
4424 tls = traceline_array(np.arange(len(tl_connectivity)) + 1,
4425 color=tl_color,
4426 connectivity=tl_connectivity)
4427 elems = element_array(np.arange(len(elem_connectivity)) + 1,
4428 elem_types, elem_color, elem_connectivity)
4429 if local:
4430 if local_nodes is not None:
4431 node_map_dict = {val: index for index, val in enumerate(node_map)}
4432 local_node_indices = np.array([node_map_dict[node] for node in local_nodes])
4433 local_node_index_map = {node_index: index for index, node_index in enumerate(local_node_indices)}
4434 preferred_local_orientation /= np.linalg.norm(preferred_local_orientation)
4435 secondary_preferred_local_orientation /= np.linalg.norm(
4436 secondary_preferred_local_orientation)
4437 # Go through each node, find the block that it's in, get its elements
4438 normal_sum = np.zeros(
4439 (coordinates.shape[0] if local_nodes is None else len(local_nodes),
4440 3))
4441 for block in blocks:
4442 if block.connectivity.shape[-1] <= 1:
4443 # Skip conmasses
4444 continue
4445 full_block_node_positions = coordinates[block.connectivity, :]
4446 normal_vectors = np.cross(full_block_node_positions[:, 2, :] - full_block_node_positions[:, 0, :],
4447 full_block_node_positions[:, 1, :] - full_block_node_positions[:, 0, :])
4448 normal_vectors /= np.linalg.norm(normal_vectors, axis=-1, keepdims=True)
4449 block_nodes = np.unique(block.connectivity)
4450 if local_nodes is not None:
4451 block_nodes = local_node_indices[np.isin(local_node_indices, block_nodes)]
4452 for index in block_nodes:
4453 node_elements = np.any(block.connectivity == index, axis=-1)
4454 node_normals = normal_vectors[node_elements]
4455 if local_nodes is None:
4456 normal_sum[index] += np.sum(node_normals, axis=0)
4457 else:
4458 normal_sum[local_node_index_map[index]] += np.sum(node_normals, axis=0)
4459 # Put small nodes to reasonable values to avoid divide by zero
4460 normal_sum[np.linalg.norm(normal_sum, axis=-1) < 1e-14] = preferred_local_orientation
4461 node_rotations = np.zeros((3,) + normal_sum.shape)
4462 node_rotations[2] = -normal_sum / np.linalg.norm(normal_sum, axis=-1, keepdims=True)
4463 # Find the nodes that are basically in the direction of the preferred direction
4464 nodes_in_preferred_direction = np.abs(
4465 np.einsum('ij,j->i', node_rotations[2], preferred_local_orientation)) > 0.99
4466 # Compute cross products if we are not parallel with the preferred orientation
4467 node_rotations[1, ~nodes_in_preferred_direction] = np.cross(
4468 node_rotations[2, ~nodes_in_preferred_direction], preferred_local_orientation)
4469 node_rotations[1, nodes_in_preferred_direction] = np.cross(
4470 node_rotations[2, nodes_in_preferred_direction], secondary_preferred_local_orientation)
4471 node_rotations[1] /= np.linalg.norm(node_rotations[1], axis=-1, keepdims=True)
4472 node_rotations[0] = np.cross(node_rotations[1], node_rotations[2])
4473 node_rotations[0] /= np.linalg.norm(node_rotations[0], axis=-1, keepdims=True)
4474 node_rotations = np.moveaxis(node_rotations, 0, 1)
4475 # Create the coordinate system matrix
4476 cs_matrix = np.concatenate((node_rotations, np.zeros(
4477 (node_rotations.shape[0], 1, node_rotations.shape[-1]))), axis=1)
4478 # Put a global coordinate system at the end
4479 cs_matrix = np.concatenate((cs_matrix, np.eye(4, 3)[np.newaxis]), axis=0)
4480 if local_nodes is None:
4481 cs_ids = np.concatenate((node_map, (node_map.max() + 1,)), axis=0)
4482 cs_descriptions = ['Node {:} Disp'.format(id) for id in node_map] + ['Global']
4483 css = coordinate_system_array(cs_ids, cs_descriptions, 1, 0, cs_matrix)
4484 nodes.def_cs = cs_ids[-1]
4485 nodes.disp_cs = cs_ids[:-1]
4486 else:
4487 cs_ids = np.concatenate((local_nodes, (node_map.max() + 1,)), axis=0)
4488 cs_descriptions = ['Node {:} Disp'.format(id) for id in local_nodes] + ['Global']
4489 css = coordinate_system_array(cs_ids, cs_descriptions, 1, 0, cs_matrix)
4490 nodes.def_cs = cs_ids[-1]
4491 nodes.disp_cs = cs_ids[-1]
4492 nodes.disp_cs[local_node_indices] = cs_ids[:-1]
4493 else:
4494 css = coordinate_system_array((1,))
4495 return Geometry(nodes, css, tls, elems)
4497 @classmethod
4498 def camera_visualization(cls, K, RT, image_size, size=1, colors=1):
4499 """
4500 Create a geometry used to visualize a camera specified with K and RT
4502 Parameters
4503 ----------
4504 K : ndarray
4505 A (...,3,3) matrix containing the intrinsic parameters of the cameras
4506 RT : ndarray
4507 A (...,3,4) matrix containing the extrinsic parameters of the cameras
4508 image_size : ndarray
4509 A (...,2) matrix containing the width,height of each camera
4510 size : float
4511 The distance that the rays will project from the pinhole of the camera
4512 colors : int or ndarray
4513 The colors assigned to each camera
4515 Returns
4516 -------
4517 geometry : cls
4518 A geoemtry containing the cameras that can be used for visualization
4520 """
4521 if K.ndim == 2:
4522 K = K[np.newaxis, ...]
4523 if K.ndim > 3:
4524 K = K.reshape(-1, 3, 3)
4525 if RT.ndim == 2:
4526 RT = RT[np.newaxis, ...]
4527 if RT.ndim > 3:
4528 RT = RT.reshape(-1, 3, 4)
4529 if image_size.ndim == 1:
4530 image_size = image_size[np.newaxis, ...]
4531 if image_size.ndim > 2:
4532 image_size = image_size.reshape(-1, 2)
4533 if np.array(colors).ndim == 0:
4534 colors = colors * np.ones(K.shape[0])
4535 if np.array(colors).ndim > 1:
4536 colors = np.array(colors).flatten()
4538 corners = np.array([(0, 0),
4539 (1, 0),
4540 (1, 1),
4541 (0, 1)])
4543 geometry = cls(coordinate_system=coordinate_system_array((1,)))
4545 for i, (thisK, thisRT, this_size, this_color) in enumerate(zip(K, RT, image_size, colors)):
4546 pixel_corners = corners * this_size
4548 # Compute the positions
4549 image_corners = point_on_pixel(thisK, thisRT, pixel_corners, size)
4550 pinhole = -thisRT[:3, :3].T @ thisRT[..., :3, 3:]
4551 # Create a geometry
4552 pinhole_node = node_array(((i + 1) * 10,), pinhole.T, disp_cs=i + 2, color=this_color)
4553 pixel_nodes = node_array((i + 1) * 10 + np.arange(4) + 1,
4554 image_corners, disp_cs=i + 2, color=this_color)
4555 cs = coordinate_system_array((i + 2,), color=this_color)
4556 cs.matrix[..., :3, :3] = thisRT[:3, :3]
4557 geometry.node = np.concatenate((geometry.node, pinhole_node, pixel_nodes))
4558 geometry.coordinate_system = np.concatenate((geometry.coordinate_system, cs))
4559 geometry.add_traceline(np.array([0, 1, 2, 3, 4, 1]) + 10 * (i + 1), color=this_color)
4560 geometry.add_traceline(np.array([0, 2]) + 10 * (i + 1), color=this_color)
4561 geometry.add_traceline(np.array([0, 3]) + 10 * (i + 1), color=this_color)
4562 geometry.add_traceline(np.array([0, 4]) + 10 * (i + 1), color=this_color)
4563 return geometry
4565 def compress_ids(self, compress_nodes=True, compress_elements=True,
4566 compress_tracelines=True, compress_cs=True,
4567 return_maps=False):
4568 """
4569 Compresses ID numbers to make node, element, etc. contiguous
4571 Parameters
4572 ----------
4573 compress_nodes : bool, optional
4574 If True, compress the node ids. The default is True.
4575 compress_elements : bool, optional
4576 If True, compress the element ids. The default is True.
4577 compress_tracelines : bool, optional
4578 If True, compress the traceline ids. The default is True.
4579 compress_cs : bool, optional
4580 If True, compress the coordinate system ids. The default is True.
4581 return_maps : bool, optional
4582 If True, return id maps for nodes, elements, tracelines, and
4583 coordinate systems. Maps will be equal to `None` if no compression
4584 was done on that field. The default is False.
4586 Returns
4587 -------
4588 mapped_geometry : Geometry
4589 Geometry with contiguous id numbers for the selected fields.
4590 node_map : id_map or `None`
4591 Mapping from the old set of nodes to the new set of nodes, only
4592 returned if `return_maps` is True. If `compress_nodes` is False,
4593 this will be `None`.
4594 traceline_map : id_map or `None`
4595 Mapping from the old set of tracelines to the new set of tracelines, only
4596 returned if `return_maps` is True. If `compress_tracelines` is False,
4597 this will be `None`.
4598 element_map : id_map or `None`
4599 Mapping from the old set of elements to the new set of elements, only
4600 returned if `return_maps` is True. If `compress_elements` is False,
4601 this will be `None`.
4602 cs_map : id_map or `None`
4603 Mapping from the old set of coordinate systems to the new set of
4604 coordinate systems, only
4605 returned if `return_maps` is True. If `compress_cs` is False,
4606 this will be `None`.
4608 """
4609 if compress_nodes:
4610 to_nodes = np.arange(self.node.size)+1
4611 from_nodes = self.node.id.flatten()
4612 node_map = id_map(from_nodes, to_nodes)
4613 else:
4614 node_map = None
4615 if compress_cs:
4616 to_cs = np.arange(self.coordinate_system.size)+1
4617 from_cs = self.coordinate_system.id.flatten()
4618 cs_map = id_map(from_cs, to_cs)
4619 else:
4620 cs_map = None
4621 if compress_elements:
4622 to_elements = np.arange(self.element.size)+1
4623 from_elements = self.coordinate_system.id.flatten()
4624 element_map = id_map(from_elements, to_elements)
4625 else:
4626 element_map = None
4627 if compress_tracelines:
4628 to_tls = np.arange(self.traceline.size)+1
4629 from_tls = self.coordinate_system.id.flatten()
4630 traceline_map = id_map(from_tls, to_tls)
4631 else:
4632 traceline_map = None
4633 mapped_geometry = self.map_ids(node_map, traceline_map, element_map, cs_map)
4634 if return_maps:
4635 return mapped_geometry, node_map, traceline_map, element_map, cs_map
4636 else:
4637 return mapped_geometry
4639 @classmethod
4640 def from_imat_struct(cls, imat_fem_struct):
4641 """
4642 Constructs a Geometry from an imat_fem class saved to a Matlab structure
4644 In IMAT, a structure can be created from an `imat_fem` by using the get()
4645 function. This can then be saved to a .mat file and loaded using
4646 `scipy.io.loadmat`. The output from loadmat can be passed into this function
4648 Parameters
4649 ----------
4650 imat_fem_struct : np.ndarray
4651 structure from loadmat containing data from an imat_fem
4653 Returns
4654 -------
4655 Geometry
4656 Geometry constructed from the data in the imat structure
4658 """
4659 node_struct = imat_fem_struct['node'][0, 0]
4660 cs_array = node_struct['cs'][0, 0].reshape(-1, 3)
4661 nodes = node_array(id=node_struct['id'][0, 0].flatten(),
4662 coordinate=node_struct['coord'][0, 0].reshape(-1, 3),
4663 color=node_struct['color'][0, 0].flatten(),
4664 def_cs=cs_array[:, 0], disp_cs=cs_array[:, 1])
4665 cs_struct = imat_fem_struct['cs'][0, 0]
4666 cs_ids = cs_struct['id'][0, 0].flatten()
4667 if cs_ids.size > 0:
4668 cs_names = np.concatenate(cs_struct['name'][0, 0].flatten()).flatten()
4669 else:
4670 cs_names = []
4671 css = coordinate_system_array(id=cs_ids,
4672 name=cs_names,
4673 color=cs_struct['color'][0, 0].flatten(),
4674 cs_type=cs_struct['type'][0, 0].flatten(),
4675 matrix=cs_struct['matrix'][0, 0].reshape(4, 3, -1).transpose(2, 0, 1))
4676 elem_struct = imat_fem_struct['elem'][0, 0]
4677 elems = element_array(id=elem_struct['id'][0, 0].flatten(),
4678 type=elem_struct['type'][0, 0].flatten(),
4679 color=elem_struct['color'][0, 0].flatten(),
4680 connectivity=[val.flatten() for val in elem_struct['conn'][0, 0].flatten()])
4681 tl_struct = imat_fem_struct['tl'][0, 0]
4682 tl_ids = tl_struct['id'][0, 0].flatten()
4683 if tl_ids.size > 0:
4684 tl_names = np.concatenate(tl_struct['desc'][0, 0].flatten()).flatten()
4685 else:
4686 tl_names = []
4687 tls = traceline_array(id=tl_ids,
4688 description=tl_names,
4689 color=tl_struct['color'][0, 0].flatten(),
4690 connectivity=[val.flatten() for val in tl_struct['conn'][0, 0].flatten()])
4691 return cls(nodes, css, tls, elems)
4693 def global_node_coordinate(self, node_ids=None):
4694 """
4695 Position of the Geometry's nodes in the global coordinate system
4697 Parameters
4698 ----------
4699 node_ids : np.ndarray
4700 An array of node id numbers to keep in the output coordinate. If
4701 not specified, all nodes will be returned in the order they are in
4702 the model.
4705 Returns
4706 -------
4707 np.ndarray
4708 numpy array with shape (n,3) where n is the number of nodes in the
4709 geometry. This array contains the 3d position of each node in the
4710 global coordinate system.
4712 """
4713 if node_ids is None:
4714 node = self.node
4715 else:
4716 node = self.node(node_ids)
4717 cs_array = self.coordinate_system(node.def_cs)
4718 return node.global_coordinate(cs_array)
4720 def node_by_global_position(self, global_position_array):
4721 """
4722 Select node by closest position
4724 Parameters
4725 ----------
4726 global_position_array : np.ndarray
4727 A (...,3) shape array containing positions of nodes to keep
4729 Returns
4730 -------
4731 NodeArray
4732 NodArray containing nodes that were closest to the positions in
4733 position_array.
4735 """
4736 node_differences = np.linalg.norm(
4737 self.global_node_coordinate() - global_position_array[(Ellipsis,) + (np.newaxis,) * self.node.ndim + (slice(None),)], axis=-1)
4738 node_differences = node_differences.reshape(*global_position_array.shape[:-1], -1)
4739 node_indices = np.argmin(
4740 node_differences,
4741 axis=-1)
4742 return self.node.flatten()[node_indices]
4744 def node_by_global_grid(self,grid_spacing, x_min = None, y_min = None,
4745 z_min = None, x_max = None, y_max = None,
4746 z_max = None):
4747 """
4748 Selects nodes in a grid
4750 Parameters
4751 ----------
4752 grid_spacing : float
4753 Approximate grid spacing between selected nodes
4754 x_min : float, optional
4755 Minimum X-value of the grid.
4756 The default is the minimum node x-coordinate.
4757 y_min : float, optional
4758 Minimun Y-value of the grid.
4759 The default is the minimum node y-coordinate.
4760 z_min : float, optional
4761 Minimum Z-value of the grid.
4762 The default is the minimum node z-coordinate.
4763 x_max : float, optional
4764 Maximum X-value of the grid.
4765 The default is the maximum node x-coordinate.
4766 y_max : float, optional
4767 Maximum Y-value of the grid.
4768 The default is the maximum node y-coordinate.
4769 z_max : float, optional
4770 Maximum Z-value of the grid.
4771 The default is the maximum node z-coordinate.
4773 Returns
4774 -------
4775 NodeArray
4776 NodeArray containing the closest nodes to the grid points
4777 """
4778 coordinate = self.global_node_coordinate()
4779 if x_min is None:
4780 x_min = coordinate[:, 0].min()
4781 if y_min is None:
4782 y_min = coordinate[:, 1].min()
4783 if z_min is None:
4784 z_min = coordinate[:, 2].min()
4785 if x_max is None:
4786 x_max = coordinate[:, 0].max()
4787 if y_max is None:
4788 y_max = coordinate[:, 1].max()
4789 if z_max is None:
4790 z_max = coordinate[:, 2].max()
4791 min_coords = np.array((x_min, y_min, z_min))
4792 max_coords = np.array((x_max, y_max, z_max))
4793 range_coords = max_coords - min_coords
4794 num_grids = np.ceil(range_coords / grid_spacing).astype(int)
4795 # Now we will create the grid in each dimension using linspace
4796 grid_arrays = [np.linspace(min_coord, max_coord, num_coord)
4797 for min_coord, max_coord, num_coord
4798 in zip(min_coords, max_coords, num_grids)]
4799 # We will then use meshgrid to assemble the array grid. We will flatten the
4800 # point dimension and transpose so the array has shape (n_points x 3)
4801 grid_points = np.array(np.meshgrid(*grid_arrays, indexing='ij')).reshape(3, -1).T
4802 # Now we will select nodes that are closest to the points in the grid
4803 candidate_nodes = self.node_by_global_position(grid_points)
4804 # Reduce to only nodes that are within one grid spacing of their requested point
4805 # this reduces the nodes that are selected by points far from the model
4806 candidate_nodes = candidate_nodes[
4807 np.linalg.norm(self.global_node_coordinate(candidate_nodes.id) - grid_points, axis=-1) < grid_spacing]
4808 # Remove any duplicates
4809 candidate_node_ids = np.unique(candidate_nodes.id)
4810 return candidate_nodes(candidate_node_ids)
4812 def global_deflection(self, coordinate_array):
4813 """
4814 Direction of local deflection in the global coordinate system
4816 Parameters
4817 ----------
4818 coordinate_array : CoordinateArray
4819 A list of coordinates for which the global deformation direction
4820 will be computed
4822 Returns
4823 -------
4824 global_deflections : np.ndarray
4825 numpy array with shape (n,3) where n is the number of coordinates
4826 in the specified `coordinate_array`. This array contains the 3d
4827 direction of motion for each coordinate in the global coordinate
4828 system.
4830 """
4831 local_deformations = coordinate_array.local_direction()
4832 ordered_nodes = self.node(coordinate_array.node)
4833 coordinate_systems = self.coordinate_system(ordered_nodes.disp_cs)
4834 points = global_coord(self.coordinate_system(
4835 ordered_nodes.def_cs), ordered_nodes.coordinate)
4836 global_deflections = global_deflection(coordinate_systems, local_deformations, points)
4837 return global_deflections
4839 @staticmethod
4840 def overlay_geometries(geometries, color_override=None, return_node_id_offset=False):
4841 """
4842 Combines several geometries, offsetting the id numbers to avoid conflicts
4844 Parameters
4845 ----------
4846 geometries : iterable
4847 An iterable of geometry objects that will be combined into a single
4848 geometry
4849 color_override : iterable, optional
4850 An iterble of integers specifying colors, which will override the
4851 existing geometry colors. This should have the same length as the
4852 `geometries` input. The default is None, which keeps the original
4853 geometry colors.
4854 return_node_id_offset : bool, optional
4855 Specifies whether or not the applied node offset should be returned,
4856 which is useful if the new node numbers are to be mapped to the old
4857 node numbers. The default is False, which only returns the combined
4858 geometry.
4860 Returns
4861 -------
4862 Geometry
4863 A geometry consisting of a combination of the specified geometries
4864 node_offset
4865 An integer specifying the node id offset applied to avoid conflicts
4867 """
4868 node_ids = np.concatenate([geometry.node.id for geometry in geometries])
4869 cs_ids = np.concatenate([geometry.coordinate_system.id for geometry in geometries])
4870 tl_ids = np.concatenate([geometry.traceline.id for geometry in geometries])
4871 elem_ids = np.concatenate([geometry.element.id for geometry in geometries])
4872 try:
4873 max_node_id = np.max(node_ids)
4874 except ValueError:
4875 max_node_id = 1
4876 try:
4877 max_cs_id = np.max(cs_ids)
4878 except ValueError:
4879 max_cs_id = 1
4880 try:
4881 max_tl_id = np.max(tl_ids)
4882 except ValueError:
4883 max_tl_id = 1
4884 try:
4885 max_elem_id = np.max(elem_ids)
4886 except ValueError:
4887 max_elem_id = 1
4888 node_length = int(np.floor(np.log10(max_node_id)) + 1)
4889 tl_length = int(np.floor(np.log10(max_tl_id)) + 1)
4890 elem_length = int(np.floor(np.log10(max_elem_id)) + 1)
4891 cs_length = int(np.floor(np.log10(max_cs_id)) + 1)
4892 new_geometries = [geometry.modify_ids(
4893 node_change=10**(node_length) * (i + 1),
4894 traceline_change=10**(tl_length) * (i + 1),
4895 element_change=10**(elem_length) * (i + 1),
4896 coordinate_system_change=10**(cs_length) * (i + 1))
4897 for i, geometry in enumerate(geometries)]
4898 if color_override is not None:
4899 for i, color in enumerate(color_override):
4900 new_geometries[i].node.color = color
4901 new_geometries[i].traceline.color = color
4902 new_geometries[i].element.color = color
4903 new_geometries[i].coordinate_system.color = color
4905 new_geometry = new_geometries[0]
4906 for i in range(1, len(new_geometries)):
4907 new_geometry = new_geometry + new_geometries[i]
4909 if return_node_id_offset:
4910 return new_geometry, 10**(node_length)
4911 else:
4912 return new_geometry
4914 def map_coordinate(self,dof_list,to_geometry,node_map=None,
4915 plot_maps = False):
4916 """
4917 Map degrees of freedom from one geometry to the closest aligned on a
4918 second geometry.
4920 Parameters
4921 ----------
4922 dof_list : CoordinateArray
4923 An array of degrees of freedom on the geometry to map.
4924 to_geometry : Geometry
4925 The geometry to which the degrees of freedom will be mapped.
4926 node_map : None, str, or id_map, optional
4927 If specified, it will map the node ids from the present geometry to
4928 the mapping geometry. If the string 'position' is specified, then
4929 the closest node by position on the mapping geometry to the present
4930 geometry will be chosen. To explicitly map between nodes, pass a
4931 sdpy.id_map. If not specified, it is assumed that the node numbers
4932 in the present geometry map to the nodes in the mapping geometry.
4933 plot_maps : bool
4934 If True, it will plot the CoordinateArray on the geometries to
4935 visualize the two mapped degrees of freedom.
4937 Raises
4938 ------
4939 ValueError
4940 If an invalid mapping is specified.
4942 Returns
4943 -------
4944 output_dof_list : CoordinateArray
4945 A CoordinateArray of the same size and shape of dof_list with each
4946 entry in dof_list corresponding to the same position in output_dof_list.
4948 """
4949 dof_list = dof_list.copy()
4950 if node_map == 'position':
4951 from_nodes = np.unique(dof_list.node)
4952 to_nodes = to_geometry.node_by_global_position(self.global_node_coordinate(from_nodes)).id
4953 to_geometry = to_geometry.reduce(to_nodes)
4954 node_map = id_map(from_nodes,to_nodes)
4955 self = self.reduce(from_nodes).map_ids(node_id_map = node_map)
4956 dof_list.node = node_map(dof_list.node)
4957 elif isinstance(node_map,id_map):
4958 self = self.reduce(node_map.from_ids).map_ids(node_id_map = node_map)
4959 to_geometry = to_geometry.reduce(node_map.to_nodes)
4960 dof_list.node = node_map(dof_list.node)
4961 elif node_map is None:
4962 self = self.reduce(np.unique(dof_list.node))
4963 to_geometry = to_geometry.reduce(np.unique(dof_list.node))
4964 else:
4965 raise ValueError('Invalid node_map. Must be an id_map, "position", or None.')
4967 output_dof_list = dof_list.copy()
4968 global_directions = self.global_deflection(dof_list)
4969 to_directions = {id_num:to_geometry.global_deflection(coordinate_array(id_num,[1,2,3]))
4970 for id_num in to_geometry.node.id}
4971 for index,value in dof_list.ndenumerate():
4972 dirs = to_directions[value.node]
4973 dots = np.einsum('ij,j',dirs,global_directions[index])
4974 direction = np.argmax(abs(dots))
4975 sign = np.sign(dots[direction])
4976 output_dof_list[index].direction = sign*(direction+1)
4977 if plot_maps:
4978 plotter = self.plot_coordinate(dof_list[index])
4979 to_geometry.plot_coordinate(output_dof_list[index],plot_kwargs={'plotter':plotter})
4981 return output_dof_list
4983 def convert_elements_to_tracelines(self,keep_ids = False, start_id = None,
4984 in_place = False):
4985 """
4986 Transforms all elements in the geometry into tracelines
4988 Parameters
4989 ----------
4990 keep_ids : bool, optional
4991 If True, the element ids will be used as the new traceline ids.
4992 The default is False.
4993 start_id : int, optional
4994 The starting id number for the tracelines converted from elements.
4995 Only used if keep_ids is False. If not specified, each subsequent
4996 element will use the current highest traceline id number + 1.
4997 in_place : bool, optional
4998 If True, modify the geometry in place. If False, a copy will be
4999 created and returned. The default is False.
5001 Raises
5002 ------
5003 ValueError
5004 If conflicts are expected to occur between the newly created
5005 traceline IDs and the current existing traceline IDs.
5006 NotImplementedError
5007 If solid element types are converted.
5009 Returns
5010 -------
5011 geometry : Geometry
5012 Only returned if in_place is False.
5014 """
5015 if in_place:
5016 geometry = self
5017 else:
5018 geometry = self.copy()
5019 if keep_ids:
5020 if np.any(np.isin(geometry.element.id,geometry.traceline.id)):
5021 raise ValueError('Keeping element IDs would result in a conflict with existing traceline IDs')
5022 if start_id is not None:
5023 if np.any(np.isin(np.arange(geometry.element.size)+start_id,geometry.traceline.id)):
5024 raise ValueError('Setting new traceline ID starting at {:} would result in a conflict with existing traceline IDs'.format(start_id))
5025 if not keep_ids and start_id is None:
5026 if len(geometry.traceline) > 0:
5027 start_id = np.max(geometry.traceline.id) + 1
5028 else:
5029 start_id = 1
5030 for idx, element in enumerate(geometry.element.flatten()):
5031 if keep_ids:
5032 new_id = element.id
5033 else:
5034 new_id = start_id + idx
5035 if element.type <= 96:
5036 conn = element.connectivity
5037 conn = np.concatenate((conn,[conn[0]]))
5038 else:
5039 raise NotImplementedError('Solid Elements are not implemented yet!')
5040 geometry.add_traceline(conn,new_id,color=element.color)
5041 geometry.element = ElementArray((0,))
5042 if not in_place:
5043 return geometry
5045 def split_tracelines_into_segments(self, in_place = False):
5046 """
5047 Splits long tracelines into many length-2 tracelines
5049 Parameters
5050 ----------
5051 in_place : bool, optional
5052 If True, the current geometry will be modified in place. The
5053 default is False.
5055 Returns
5056 -------
5057 geometry : Geometry
5058 Only returned if in_place is False.
5060 """
5061 temp_geometry = Geometry()
5062 for traceline in self.traceline:
5063 connectivity = np.array((traceline.connectivity[:-1],traceline.connectivity[1:])).T
5064 connectivity = connectivity[np.all(connectivity != 0,axis=-1)]
5065 for conn in connectivity:
5066 temp_geometry.add_traceline(conn, traceline.id, traceline.description, traceline.color)
5067 if not in_place:
5068 geometry = self.copy()
5069 geometry.traceline = temp_geometry.traceline
5070 return geometry
5071 else:
5072 self.traceline = temp_geometry.traceline
5074 def remove_duplicate_tracelines(self, in_place = False):
5075 """
5076 Removes tracelines that are effectively equivalent to each other.
5078 If two tracelines have the same connectivity (or reversed connectivity),
5079 then the first will be kept and the second will be discarded. This will
5080 not take into account other fields such as description or color.
5082 Parameters
5083 ----------
5084 in_place : bool, optional
5085 If True, the current geometry will be modified in place. The
5086 default is False.
5088 Returns
5089 -------
5090 geometry : Geometry
5091 Only returned if in_place is False.
5093 """
5094 conn_set = set([])
5095 keep_array = []
5096 for traceline in self.traceline:
5097 conn = tuple(traceline.connectivity)
5098 if conn in conn_set or conn[::-1] in conn_set:
5099 keep_array.append(False)
5100 else:
5101 keep_array.append(True)
5102 conn_set.add(conn)
5103 if not in_place:
5104 geometry = self.copy()
5105 geometry.traceline = geometry.traceline[keep_array]
5106 return geometry
5107 else:
5108 self.traceline = self.traceline[keep_array]
5110 def plot(self, node_size: int = 5, line_width: int = 1, opacity=1.0, view_up=None, view_from=None, plotter=None,
5111 show_edges=False, label_nodes = False, label_tracelines = False,
5112 label_elements = False, label_font_size = 16,
5113 plot_individual_items = False):
5114 """
5115 Plots the geometry in an interactive 3D window.
5117 Parameters
5118 ----------
5119 node_size : int, optional
5120 Size to display the nodes in pixels. Set to 0 to not display nodes.
5121 The default is 5.
5122 line_width : int, optional
5123 Width to display tracelines and element edges in pixels. Set to 0
5124 to not show tracelines or edges. The default is 1.
5125 opacity : float, optional
5126 A float between 0 and 1 to specify the transparency of the geometry.
5127 Set to 1 for completely opaque, and 0 for completely transparent
5128 (invisible). The default is 1.0, no transparency.
5129 view_up : np.ndarray, optional
5130 Set the "up" direction in the plot by passing in a size-3 numpy
5131 array. The default is None.
5132 view_from : np.ndarray, optional
5133 Specify the direction from which the geometry is viewed. The
5134 default is None.
5135 plotter : BackgroundPlotter, optional
5136 A plotter can be specified to plot the geometry in an existing
5137 plot. The default is None, which creates a new window and plot.
5138 show_edges : bool, optional
5139 Specify whether or not to draw edges on elements. The default is False.
5140 label_nodes : bool or iteratble, optional
5141 Specify whether or not to label the nodes in the geometry. If True,
5142 all nodes will be plotted. If label_nodes is an array or list, it
5143 should contain the node numbers to label. The default is False,
5144 which does not label nodes.
5145 label_tracelines : bool or iteratble, optional
5146 Specify whether or not to label the tracelines in the geometry. If True,
5147 all tracelines will be plotted. If label_tracelines is an array or list, it
5148 should contain the traceline numbers to label. The default is False,
5149 which does not label tracelines.
5150 label_elements : bool or iteratble, optional
5151 Specify whether or not to label the elements in the geometry. If True,
5152 all elements will be plotted. If label_elements is an array or list, it
5153 should contain the element numbers to label. The default is False,
5154 which does not label elements.
5155 label_font_size : int, optional
5156 Specify the font size of the node labels. The default is 16.
5157 plot_individual_items : bool, optional
5158 If True, then each item will be plotted one at a time. This will
5159 make the plotting much slower, but it is required for certain
5160 features (like exporting to 3D PDF content) where having all data
5161 in one mesh is not feasible. Note that volume elements will NOT
5162 be plotted in this case.
5164 Raises
5165 ------
5166 KeyError
5167 If referenced id numbers are not found in the corresponding object,
5168 for example if a traceline references node 11 but there is no node
5169 11 in the NodeArray
5170 ValueError
5171 If an invalid or unknown element type is used
5173 Returns
5174 -------
5175 plotter : BackgroundPlotter
5176 A reference to the plotter object that the geometry was plotted in
5177 face_mesh : TYPE
5178 A reference to the mesh used to plot surface elements
5179 point_mesh : TYPE
5180 A reference to the mesh used to plot nodes and tracelines
5181 solid_mesh : TYPE
5182 A reference to the mesh used to plot volume elements
5184 """
5185 if IGNORE_PLOTS:
5186 return None
5187 # Get part information
5188 nodes = self.node.flatten()
5189 css = self.coordinate_system.flatten()
5190 elems = self.element.flatten()
5191 tls = self.traceline.flatten()
5193 num_labels = ((1 if (label_nodes is True or (label_nodes is not False and len(label_nodes) > 0)) else 0)
5194 + (1 if (label_tracelines is True or (label_tracelines is not False and len(label_tracelines) > 0)) else 0)
5195 + (1 if (label_elements is True or (label_elements is not False and len(label_elements) > 0)) else 0))
5197 coordinate_system = CoordinateSystemArray(nodes.shape)
5198 cs_map = {cs['id']:idx for idx,cs in np.ndenumerate(self.coordinate_system)}
5199 for key, node in nodes.ndenumerate():
5200 coordinate_system[key] = css[cs_map[node.def_cs]]
5201 global_node_positions = global_coord(coordinate_system, nodes.coordinate)
5202 node_index_dict = {node.id[()]: index[0] for index, node in nodes.ndenumerate()}
5203 node_index_map = np.vectorize(node_index_dict.__getitem__)
5204 tl_index_dict = {tl.id[()]: index[0] for index, tl in tls.ndenumerate()}
5205 tl_index_map = np.vectorize(tl_index_dict.__getitem__)
5206 elem_index_dict = {elem.id[()]: index[0] for index, elem in elems.ndenumerate()}
5207 elem_index_map = np.vectorize(elem_index_dict.__getitem__)
5209 # Now go through and get the element and line connectivity
5210 if not plot_individual_items:
5211 face_element_connectivity = []
5212 face_element_colors = []
5213 solid_element_connectivity = []
5214 solid_element_colors = []
5215 solid_element_types = []
5216 node_colors = []
5217 line_connectivity = []
5218 line_colors = []
5219 for index, node in nodes.ndenumerate():
5220 # element_connectivity.append(1)
5221 # element_connectivity.append(index[0])
5222 # element_colors.append(node.color)
5223 node_colors.append(node.color)
5224 for index, element in elems.ndenumerate():
5225 # Check which type of element it is
5226 if element.type in _beam_elem_types: # Beamlike element, use a line
5227 line_connectivity.append(len(element.connectivity))
5228 try:
5229 line_connectivity.extend(node_index_map(element.connectivity))
5230 except KeyError:
5231 raise KeyError(
5232 'Element {:} contains a node id not found in the node array'.format(element.id))
5233 line_colors.append(element.color)
5234 elif element.type in _face_element_types:
5235 face_element_connectivity.append(len(element.connectivity))
5236 try:
5237 face_element_connectivity.extend(node_index_map(element.connectivity))
5238 except KeyError:
5239 raise KeyError(
5240 'Element {:} contains a node id not found in the node array'.format(element.id))
5241 face_element_colors.append(element.color)
5242 elif element.type in _solid_element_types:
5243 try:
5244 solid_element_types.append(_vtk_element_map[element.type])
5245 except KeyError:
5246 raise ValueError('Do not know equivalent VTK element type for {:}: {:}'.format(
5247 element.type, _element_types[element.type]))
5248 solid_element_connectivity.append(len(element.connectivity))
5249 try:
5250 if solid_element_types[-1] in _vtk_connectivity_reorder:
5251 solid_element_connectivity.extend(node_index_map(
5252 element.connectivity[..., _vtk_connectivity_reorder[solid_element_types[-1]]]))
5253 else:
5254 solid_element_connectivity.extend(node_index_map(element.connectivity))
5255 except KeyError:
5256 raise KeyError(
5257 'Element {:} contains a node id not found in the node array'.format(element.id))
5258 solid_element_colors.append(element.color)
5259 else:
5260 raise ValueError('Unknown element type {:}'.format(element.type))
5262 for index, tl in tls.ndenumerate():
5263 for conn_group in split_list(tl.connectivity, 0):
5264 if len(conn_group) == 0:
5265 continue
5266 line_connectivity.append(len(conn_group))
5267 try:
5268 line_connectivity.extend(node_index_map(conn_group))
5269 except KeyError:
5270 raise KeyError(
5271 'Traceline {:} contains a node id not found in the node array'.format(tl.id))
5272 line_colors.append(tl.color)
5273 else:
5274 face_element_connectivity = []
5275 face_element_colors = []
5276 node_colors = []
5277 line_connectivity = []
5278 line_colors = []
5279 for index, node in nodes.ndenumerate():
5280 # element_connectivity.append(1)
5281 # element_connectivity.append(index[0])
5282 # element_colors.append(node.color)
5283 node_colors.append(node.color)
5284 for index, element in elems.ndenumerate():
5285 # Check which type of element it is
5286 if element.type in _beam_elem_types: # Beamlike element, use a line
5287 try:
5288 line_connectivity.append(node_index_map(element.connectivity))
5289 except KeyError:
5290 raise KeyError(
5291 'Element {:} contains a node id not found in the node array'.format(element.id))
5292 line_colors.append(element.color)
5293 elif element.type in _face_element_types:
5294 try:
5295 if len(element.connectivity) == 3:
5296 face_element_connectivity.append(node_index_map(element.connectivity))
5297 face_element_colors.append(element.color)
5298 else:
5299 face_element_connectivity.append(node_index_map(element.connectivity[[0,1,3]]))
5300 face_element_colors.append(element.color)
5301 face_element_connectivity.append(node_index_map(element.connectivity[[1,2,3]]))
5302 face_element_colors.append(element.color)
5303 except KeyError:
5304 raise KeyError(
5305 'Element {:} contains a node id not found in the node array'.format(element.id))
5306 elif element.type in _solid_element_types:
5307 warnings.warn('Solid Elements are currently not supported and will be skipped!')
5308 else:
5309 raise ValueError('Unknown element type {:}'.format(element.type))
5311 for index, tl in tls.ndenumerate():
5312 for conn_group in split_list(tl.connectivity, 0):
5313 if len(conn_group) == 0:
5314 continue
5315 try:
5316 mapped_conn_group = node_index_map(conn_group)
5317 except KeyError:
5318 raise KeyError(
5319 'Traceline {:} contains a node id not found in the node array'.format(tl.id))
5320 for indices in zip(mapped_conn_group[:-1],mapped_conn_group[1:]):
5321 line_connectivity.append(indices)
5322 line_colors.append(tl.color)
5324 # Now we start to plot
5325 if plotter is None:
5326 plotter = GeometryPlotter(editor=False)
5328 if not plot_individual_items:
5329 # Need to split up between point mesh, face/line mesh, and solid mesh
5330 if len(face_element_connectivity) == 0 and len(line_connectivity) == 0:
5331 face_mesh = None
5332 else:
5333 face_mesh = pv.PolyData(global_node_positions,
5334 faces=None if len(
5335 face_element_connectivity) == 0 else face_element_connectivity,
5336 lines=None if len(line_connectivity) == 0 else line_connectivity)
5337 face_mesh.cell_data['color'] = line_colors + face_element_colors
5339 plotter.add_mesh(face_mesh, scalars='color', cmap=colormap, clim=[0, 15],
5340 show_edges=show_edges, # True if line_width > 0 else False,
5341 show_scalar_bar=False, line_width=line_width,
5342 opacity=opacity)
5343 if len(solid_element_connectivity) == 0:
5344 solid_mesh = None
5345 else:
5346 solid_mesh = pv.UnstructuredGrid(np.array(solid_element_connectivity),
5347 np.array(solid_element_types, dtype='uint8'),
5348 np.array(global_node_positions))
5349 solid_mesh.cell_data['color'] = solid_element_colors
5351 plotter.add_mesh(solid_mesh, scalars='color', cmap=colormap, clim=[0, 15],
5352 show_edges=show_edges, # True if line_width > 0 else False,
5353 show_scalar_bar=False, line_width=line_width,
5354 opacity=opacity)
5356 if node_size > 0:
5357 point_mesh = pv.PolyData(global_node_positions)
5358 point_mesh.cell_data['color'] = node_colors
5359 plotter.add_mesh(point_mesh, scalars=node_colors, cmap=colormap, clim=[0, 15],
5360 show_edges=show_edges, # True if line_width > 0 else False,
5361 show_scalar_bar=False, point_size=node_size,
5362 opacity=opacity)
5363 else:
5364 point_mesh = None
5366 if (label_nodes is True or (label_nodes is not False and len(label_nodes) > 0)):
5367 try:
5368 label_node_ids = [id_num for id_num in label_nodes]
5369 label_node_indices = node_index_map(label_node_ids)
5370 except TypeError:
5371 label_node_ids = [id_num for id_num in self.node.id]
5372 label_node_indices = node_index_map(label_node_ids)
5374 node_label_mesh = pv.PolyData(global_node_positions[label_node_indices])
5375 node_label_mesh['NODE_ID'] = [
5376 ('N' if num_labels > 1 else '') + str(val) for val in label_node_ids]
5377 plotter.add_point_labels(node_label_mesh, 'NODE_ID', tolerance=0.0, shape=None,
5378 show_points=False, always_visible=True, font_size=label_font_size)
5380 if (label_tracelines is True or (label_tracelines is not False and len(label_tracelines) > 0)):
5381 try:
5382 label_traceline_ids = [id_num for id_num in label_tracelines]
5383 label_traceline_indices = tl_index_map(label_traceline_ids)
5384 except TypeError:
5385 label_traceline_ids = [id_num for id_num in self.traceline.id]
5386 label_traceline_indices = tl_index_map(label_traceline_ids)
5388 # Need to get the point to plot for each traceline
5389 traceline_locations = []
5390 for idx in label_traceline_indices:
5391 tl = tls[idx]
5392 connectivity = tl.connectivity[tl.connectivity != 0]
5393 # Find the middle of the traceline
5394 center_nodes = [connectivity[len(connectivity)//2-1],connectivity[len(connectivity)//2]]
5395 position = np.mean(global_node_positions[node_index_map(center_nodes)],axis=0)
5396 traceline_locations.append(position)
5398 tl_label_mesh = pv.PolyData(np.array(traceline_locations))
5399 tl_label_mesh['TL_ID'] = [
5400 ('L' if num_labels > 1 else '') + str(val) for val in label_traceline_ids]
5401 plotter.add_point_labels(tl_label_mesh, 'TL_ID', tolerance=0.0, shape=None,
5402 show_points=False, always_visible=True, font_size=label_font_size)
5404 if (label_elements is True or (label_elements is not False and len(label_elements) > 0)):
5405 try:
5406 label_element_ids = [id_num for id_num in label_elements]
5407 label_element_indices = elem_index_map(label_element_ids)
5408 except TypeError:
5409 label_element_ids = [id_num for id_num in self.element.id]
5410 label_element_indices = elem_index_map(label_element_ids)
5412 # Need to get the point to plot for each traceline
5413 element_locations = []
5414 for idx in label_element_indices:
5415 elem = elems[idx]
5416 connectivity = elem.connectivity
5417 # Find the middle of the traceline
5418 position = np.mean(global_node_positions[node_index_map(connectivity)],axis=0)
5419 element_locations.append(position)
5421 elem_label_mesh = pv.PolyData(np.array(element_locations))
5422 elem_label_mesh['ELEM_ID'] = [
5423 ('E' if num_labels > 1 else '') + str(val) for val in label_element_ids]
5424 plotter.add_point_labels(elem_label_mesh, 'ELEM_ID', tolerance=0.0, shape=None,
5425 show_points=False, always_visible=True, font_size=label_font_size)
5427 else:
5428 face_map = []
5429 face_mesh = []
5430 for conn,color in zip(face_element_connectivity,face_element_colors):
5431 node_indices,inverse = np.unique(conn,return_inverse=True)
5432 node_positions = global_node_positions[node_indices][inverse]
5433 node_connectivity = np.arange(node_positions.shape[0])
5434 mesh = pv.PolyData(node_positions,faces=[len(node_connectivity)]+[c for c in node_connectivity])
5435 mesh.cell_data['color'] = color
5436 mesh.SetObjectName('Elem: {:} {:} {:}'.format(*self.node.id[node_indices[inverse]]))
5437 plotter.add_mesh(mesh,scalars='color',cmap=colormap, clim = [0,15],
5438 show_edges=show_edges, show_scalar_bar = False,
5439 line_width=line_width,opacity=opacity,label='Elem: {:} {:} {:}'.format(*self.node.id[node_indices[inverse]]))
5440 face_map.append(self.node.id[node_indices[inverse]])
5441 face_mesh.append(mesh)
5443 line_map = []
5444 line_mesh = []
5445 for conn,color in zip(line_connectivity,line_colors):
5446 node_indices,inverse = np.unique(conn,return_inverse=True)
5447 node_positions = global_node_positions[node_indices][inverse]
5448 node_connectivity = np.arange(node_positions.shape[0])
5449 mesh = pv.PolyData(node_positions,lines=[len(node_connectivity)]+[c for c in node_connectivity])
5450 mesh.cell_data['color'] = color
5451 mesh.SetObjectName('Line: {:} {:}'.format(*self.node.id[node_indices[inverse]]))
5452 plotter.add_mesh(mesh,scalars='color',cmap=colormap, clim = [0,15],
5453 show_edges=show_edges, show_scalar_bar = False,
5454 line_width=line_width,opacity=opacity,label='Line: {:} {:}'.format(*self.node.id[node_indices[inverse]]))
5455 line_map.append(self.node.id[node_indices[inverse]])
5456 line_mesh.append(mesh)
5458 point_map = []
5459 point_mesh = []
5460 for node,position,color in zip(self.node.id,global_node_positions,self.node.color):
5461 if node_size > 0:
5462 mesh = pv.PolyData(position)
5463 mesh.cell_data['color'] = color
5464 plotter.add_mesh(mesh, scalars = color, cmap=colormap, clim=[0,15],
5465 show_edges=show_edges, show_scalar_bar=False, point_size=node_size,
5466 opacity=opacity)
5467 point_map.append(node)
5468 point_mesh.append(mesh)
5470 if label_nodes:
5471 try:
5472 label_node_ids = [id_num for id_num in label_nodes]
5473 label_node_indices = node_index_map(label_node_ids)
5474 except TypeError:
5475 label_node_ids = [id_num for id_num in self.node.id]
5476 label_node_indices = node_index_map(label_node_ids)
5478 for index,label in zip(label_node_indices,label_node_ids):
5479 node_label_mesh = pv.PolyData(global_node_positions[index])
5480 node_label_mesh['NODE_ID'] = [str(label)]
5481 plotter.add_point_labels(node_label_mesh, 'NODE_ID', tolerance=0.0, shape=None,
5482 show_points=False, always_visible=True, font_size=label_font_size)
5484 if view_from is not None:
5485 focus = plotter.camera.focal_point
5486 distance = plotter.camera.distance
5487 plotter.camera.position = np.array(
5488 focus) + distance * np.array(view_from) / np.linalg.norm(view_from)
5489 plotter.camera.focal_point = focus
5490 if view_up is not None:
5491 plotter.camera.up = view_up
5492 plotter.show()
5493 plotter.render()
5494 QApplication.processEvents()
5495 if plot_individual_items:
5496 return plotter, (face_mesh,face_map), (line_mesh,line_map), (point_mesh,point_map)
5497 else:
5498 return plotter, face_mesh, point_mesh, solid_mesh
5500 def global_geometry_table(self,coordinate_array_or_nodelist = None):
5501 if coordinate_array_or_nodelist is None:
5502 coord = coordinate_array(self.node.id,[1,2,3],force_broadcast=True)
5503 elif isinstance(coordinate_array_or_nodelist,CoordinateArray):
5504 coord = coordinate_array_or_nodelist.flatten()
5505 else:
5506 coord = coordinate_array(coordinate_array_or_nodelist,[1,2,3],force_broadcast=True)
5507 global_positions = self.global_node_coordinate(coord.node)
5508 global_directions = self.global_deflection(coord)
5509 coord_string = coord.string_array()
5510 df_dict = {}
5511 df_dict['Degree of Freedom'] = coord_string
5512 for label,column in zip('XYZ',global_positions.T):
5513 df_dict[f'Coordinate {label}'] = column
5514 for label,column in zip('XYZ',global_directions.T):
5515 df_dict[f'Direction {label}'] = column
5516 df = pd.DataFrame(df_dict)
5517 return df
5519 def plot_coordinate(self, coordinates: CoordinateArray = None,
5520 arrow_scale=0.1,
5521 arrow_scale_type='bbox', label_dofs=False,
5522 label_font_size=16, opacity=1.0,
5523 arrow_ends_on_node=False, plot_kwargs={},
5524 plot_individual_items = False):
5525 """
5526 Plots coordinate arrows on the geometry
5528 Parameters
5529 ----------
5530 coordinates : CoordinateArray, optional
5531 Coordinates to draw on the geometry. If no coordinates are specified,
5532 all translation degrees of freedom at each node will be plotted.
5533 arrow_scale : float | np.ndarray, optional
5534 Size of the arrows in proportion to the length of the diagonal of
5535 the bounding box of the Geometry if `arrow_scale_type` is 'bbox',
5536 otherwise the raw length of the arrow. If an ndarray is provided,
5537 it should match the shape of `coordinates`. The default is 0.1.
5538 arrow_scale_type : str, optional
5539 Specifies how to compute the size of the arrows. If 'bbox', then
5540 the arrow is scaled based on the size of the geometry. Otherwise,
5541 the arrow size is the specified length. The default is 'bbox'.
5542 label_dofs : bool, optional
5543 Specify whether or not to label the coordinates with strings.
5544 The default is False.
5545 label_font_size : int, optional
5546 Specifies the font size for the node labels.
5547 Default is 16.
5548 opacity : float, optional
5549 A float between 0 and 1 to specify the transparency of the geometry.
5550 Set to 1 for completely opaque, and 0 for completely transparent
5551 (invisible). The default is 1.0, no transparency.
5552 arrow_ends_on_node : bool, optional
5553 If True, arrow tip ends at the node, otherwise the arrow begins at
5554 node.
5555 Defualt is False
5556 plot_kwargs : dict, optional
5557 Any additional keywords that should be passed to the Geometry.plot
5558 function. The default is {}.
5559 plot_individual_items : bool, optional
5560 If True, then each item will be plotted one at a time. This will
5561 make the plotting much slower, but it is required for certain
5562 features (like exporting to 3D PDF content) where having all data
5563 in one mesh is not feasible. Note that volume elements will NOT
5564 be plotted in this case.
5566 Returns
5567 -------
5568 plotter : BackgroundPlotter
5569 A reference to the window the geometry was plotted in.
5571 """
5572 if IGNORE_PLOTS:
5573 return None
5574 plot_kwargs = plot_kwargs.copy()
5575 plot_kwargs['plot_individual_items'] = plot_individual_items
5576 plotter, face_mesh, point_mesh, solid_mesh = self.plot(opacity=opacity, **plot_kwargs)
5577 # Add deflection directions to each node
5578 # Get local deformation for each coordinate direction
5580 if coordinates is None:
5581 coordinates = from_nodelist(self.node.id)
5583 def build_coord_mesh(coordinates, geom, scale_factors=None):
5584 nodes = coordinates.node
5585 indices = np.isin(coordinates.node, geom.node.id)
5586 coordinates = coordinates[indices]
5587 nodes = nodes[indices]
5588 local_deformations = coordinates.local_direction()
5589 coordinate_systems = geom.coordinate_system(geom.node(nodes).disp_cs)
5590 points = global_coord(geom.coordinate_system(
5591 geom.node(nodes).def_cs), geom.node(nodes).coordinate)
5592 global_deflections = global_deflection(coordinate_systems, local_deformations, points)
5593 # Now add the point array to the mesh
5594 coord_mesh = pv.PolyData(points)
5595 coord_mesh.point_data['Coordinates'] = global_deflections
5596 coord_mesh.point_data['Direction'] = abs(coordinates.direction)
5597 coord_mesh.point_data['Node ID'] = nodes
5598 if scale_factors is not None:
5599 coord_mesh.point_data['Scale'] = scale_factors
5600 return coord_mesh, points, global_deflections
5602 coordinates = coordinates.flatten()
5604 if isinstance(arrow_scale, int | float):
5605 scale = False
5606 arrow_scale_array = np.ones(len(coordinates))
5607 global_arrow_scale = arrow_scale
5608 elif isinstance(arrow_scale, np.ndarray):
5609 scale = 'Scale'
5610 arrow_scale_array = arrow_scale
5611 global_arrow_scale = 1
5613 if arrow_scale_type == 'bbox':
5614 node_coords = self.global_node_coordinate()
5615 bbox_diagonal = np.linalg.norm(
5616 np.max(node_coords, axis=0) - np.min(node_coords, axis=0))
5617 arrow_factor = bbox_diagonal * global_arrow_scale
5618 else:
5619 arrow_factor = global_arrow_scale
5621 # Create Straight Arrows and Labels
5622 coordinates_straight = coordinates[np.abs(coordinates.direction) < 4]
5623 if arrow_scale_array is not None:
5624 straight_scaling = arrow_scale_array[np.abs(coordinates.direction) < 4]
5625 else:
5626 straight_scaling = None
5627 if coordinates_straight.size:
5628 if plot_individual_items:
5629 for coordinate_straight in coordinates_straight:
5630 [coord_mesh_straight, points_straight, global_deflections_straight] = build_coord_mesh(
5631 coordinate_straight[np.newaxis], self, straight_scaling)
5633 if arrow_ends_on_node:
5634 arrow_start = (-1.0, 0.0, 0.0)
5635 else:
5636 arrow_start = (0.0, 0.0, 0.0)
5637 coord_arrows_stright = coord_mesh_straight.glyph(
5638 orient='Coordinates', scale=scale, factor=arrow_factor, geom=pv.Arrow(start=arrow_start))
5639 color = ((1.0,0.0,0.0) if abs(coordinate_straight.direction) == 1 else
5640 (0.0,1.0,0.0) if abs(coordinate_straight.direction) == 2 else
5641 (0.0,0.0,1.0))
5642 plotter.add_mesh(coord_arrows_stright, color = color, show_scalar_bar=False)
5643 else:
5644 [coord_mesh_straight, points_straight, global_deflections_straight] = build_coord_mesh(
5645 coordinates_straight, self, straight_scaling)
5647 if arrow_ends_on_node:
5648 arrow_start = (-1.0, 0.0, 0.0)
5649 else:
5650 arrow_start = (0.0, 0.0, 0.0)
5651 coord_arrows_stright = coord_mesh_straight.glyph(
5652 orient='Coordinates', scale=scale, factor=arrow_factor, geom=pv.Arrow(start=arrow_start))
5653 plotter.add_mesh(coord_arrows_stright, scalars='Direction',
5654 cmap=coord_colormap, clim=[1, 6], show_scalar_bar=False)
5656 if label_dofs:
5657 if arrow_ends_on_node:
5658 dof_label_mesh = pv.PolyData(
5659 points_straight - global_deflections_straight * arrow_factor * straight_scaling[:,np.newaxis])
5660 else:
5661 dof_label_mesh = pv.PolyData(
5662 points_straight + global_deflections_straight * arrow_factor * straight_scaling[:,np.newaxis])
5663 dof_label_mesh['DOF'] = [val for val in coordinates_straight.string_array()]
5664 plotter.add_point_labels(dof_label_mesh, 'DOF', tolerance=0.0, shape=None,
5665 show_points=False, always_visible=True, font_size=label_font_size)
5667 # Create Rotational Arrows and Labels
5668 coordinates_rotation = coordinates[np.abs(coordinates.direction) > 3]
5669 if arrow_scale_array is not None:
5670 rotation_scaling = arrow_scale_array[np.abs(coordinates.direction) > 3]
5671 else:
5672 rotation_scaling = None
5673 if coordinates_rotation.size:
5674 if plot_individual_items:
5675 for coordinate_rotation in coordinates_rotation:
5676 [coord_mesh_rotation, points_rotation, global_deflections_rotation] = build_coord_mesh(
5677 coordinate_rotation[np.newaxis], self, rotation_scaling)
5679 arrow_index = np.arange(4)
5680 head_location_angles = 1 / 8 * np.pi + 1 / 2 * np.pi * arrow_index
5681 tail_location_angles = 3 / 8 * np.pi + 1 / 2 * np.pi * arrow_index
5683 r = 1
5684 arrow_head_locations = np.array([np.zeros(
5685 head_location_angles.size), r * np.sin(head_location_angles), r * np.cos(head_location_angles)]).T
5686 arrow_tail_locations = np.array([np.zeros(
5687 tail_location_angles.size), r * np.sin(tail_location_angles), r * np.cos(tail_location_angles)]).T
5688 cone_vectors = np.array([np.zeros(head_location_angles.size), np.sin(
5689 head_location_angles - np.pi / 2), np.cos(head_location_angles - np.pi / 2)]).T
5691 arc = pv.merge([pv.CircularArc(pointa=start, pointb=stop, center=(0, 0, 0)).tube(
5692 radius=0.05) for start, stop in zip(arrow_head_locations, arrow_tail_locations)])
5693 cone = pv.merge([pv.Cone(center=cone_center, direction=cone_vector, height=0.3, radius=0.1, resolution=20)
5694 for cone_center, cone_vector in zip(arrow_head_locations, cone_vectors)])
5696 curved_arrow = pv.merge([arc, cone])
5697 coord_arrows_rotation = coord_mesh_rotation.glyph(
5698 orient='Coordinates', scale=scale, factor=arrow_factor, geom=curved_arrow)
5699 color = ((1.0,0.0,0.0) if abs(coordinate_straight.direction) == 1 else
5700 (0.0,1.0,0.0) if abs(coordinate_straight.direction) == 2 else
5701 (0.0,0.0,1.0))
5702 plotter.add_mesh(coord_arrows_rotation, color=color, show_scalar_bar=False)
5703 else:
5704 [coord_mesh_rotation, points_rotation, global_deflections_rotation] = build_coord_mesh(
5705 coordinates_rotation, self, rotation_scaling)
5707 arrow_index = np.arange(4)
5708 head_location_angles = 1 / 8 * np.pi + 1 / 2 * np.pi * arrow_index
5709 tail_location_angles = 3 / 8 * np.pi + 1 / 2 * np.pi * arrow_index
5711 r = 1
5712 arrow_head_locations = np.array([np.zeros(
5713 head_location_angles.size), r * np.sin(head_location_angles), r * np.cos(head_location_angles)]).T
5714 arrow_tail_locations = np.array([np.zeros(
5715 tail_location_angles.size), r * np.sin(tail_location_angles), r * np.cos(tail_location_angles)]).T
5716 cone_vectors = np.array([np.zeros(head_location_angles.size), np.sin(
5717 head_location_angles - np.pi / 2), np.cos(head_location_angles - np.pi / 2)]).T
5719 arc = pv.merge([pv.CircularArc(pointa=start, pointb=stop, center=(0, 0, 0)).tube(
5720 radius=0.05) for start, stop in zip(arrow_head_locations, arrow_tail_locations)])
5721 cone = pv.merge([pv.Cone(center=cone_center, direction=cone_vector, height=0.3, radius=0.1, resolution=20)
5722 for cone_center, cone_vector in zip(arrow_head_locations, cone_vectors)])
5724 curved_arrow = pv.merge([arc, cone])
5725 coord_arrows_rotation = coord_mesh_rotation.glyph(
5726 orient='Coordinates', scale=scale, factor=arrow_factor, geom=curved_arrow)
5727 # Make Colors of Rotations same as Straight
5728 coord_arrows_rotation['Direction'] = coord_arrows_rotation['Direction'] - 3
5730 plotter.add_mesh(coord_arrows_rotation, scalars='Direction',
5731 cmap=coord_colormap, clim=[1, 6], show_scalar_bar=False)
5733 if label_dofs:
5734 nodes = coordinates_rotation.node
5735 indices = np.isin(coordinates_rotation.node, self.node.id)
5736 nodes = nodes[indices]
5737 local_deformations = coordinates_rotation.local_direction()
5738 local_deformation_transformation = np.array(([arrow_head_locations[0, 0], arrow_head_locations[0, 1], arrow_head_locations[0, 2]],
5739 [arrow_head_locations[0, 2], arrow_head_locations[0,
5740 0], arrow_head_locations[0, 1]],
5741 [arrow_head_locations[0, 1], arrow_head_locations[0, 2], arrow_head_locations[0, 0]]))
5742 local_deformations_new = np.transpose(
5743 local_deformation_transformation.T @ local_deformations.T)
5744 coordinate_systems = self.coordinate_system(self.node(nodes).disp_cs)
5745 points = global_coord(self.coordinate_system(
5746 self.node(nodes).def_cs), self.node(nodes).coordinate)
5747 global_deflections_rotation_new = global_deflection(
5748 coordinate_systems, local_deformations_new, points)
5750 dof_label_mesh = pv.PolyData(
5751 points_rotation + global_deflections_rotation_new * arrow_factor * rotation_scaling[:,np.newaxis])
5753 dof_label_mesh['DOF'] = [val for val in coordinates_rotation.string_array()]
5754 plotter.add_point_labels(dof_label_mesh, 'DOF', tolerance=0.0, shape=None,
5755 show_points=False, always_visible=True, font_size=label_font_size)
5756 return plotter
5758 def plot_shape(self, shape, plot_kwargs={}, background_plotter_kwargs={'editor': False}, undeformed_opacity=0.25, deformed_opacity=1.0, starting_scale=1.0):
5759 """
5760 Plot mode shapes on the geometry
5762 Parameters
5763 ----------
5764 shape : ShapeArray
5765 The set of shapes to plot
5766 plot_kwargs : dict, optional
5767 Any additional keywords that should be passed to the Geometry.plot
5768 function. The default is {}.
5769 background_plotter_kwargs : dict, optional
5770 Any additional arguments that should be passed to the
5771 BackgroundPlotter initializer. The default is {'editor':False}.
5772 undeformed_opacity : float, optional
5773 A float between 0 and 1 to specify the transparency of the undeformed
5774 geometry. Set to 1 for completely opaque, and 0 for completely
5775 transparent (invisible). The default is 0.25.
5776 deformed_opacity : float, optional
5777 A float between 0 and 1 to specify the transparency of the deformed
5778 geometry. Set to 1 for completely opaque, and 0 for completely
5779 transparent (invisible). The default is 1.0.
5780 starting_scale : float, optional
5781 The starting scale factor of the animation. The default is 1.0.
5783 Returns
5784 -------
5785 ShapePlotter
5786 A reference to the ShapePlotter class that is created to plot the
5787 animated shapes
5789 """
5790 if IGNORE_PLOTS:
5791 return None
5792 return ShapePlotter(self, shape, plot_kwargs, background_plotter_kwargs,
5793 undeformed_opacity, deformed_opacity=deformed_opacity, starting_scale=starting_scale)
5795 def plot_deflection_shape(self, deflection_shape_data, plot_kwargs={},
5796 background_plotter_kwargs={'editor': False},
5797 undeformed_opacity=0.25, deformed_opacity=1.0,
5798 starting_scale=1.0):
5799 """
5800 Plot deflection shapes shapes on the geometry
5802 Parameters
5803 ----------
5804 deflection_shape_data : NDDataArray
5805 Data array containing the deflection shapes to plot
5806 plot_kwargs : dict, optional
5807 Any additional keywords that should be passed to the Geometry.plot
5808 function. The default is {}.
5809 background_plotter_kwargs : dict, optional
5810 Any additional arguments that should be passed to the
5811 BackgroundPlotter initializer. The default is {'editor':False}.
5812 undeformed_opacity : float, optional
5813 A float between 0 and 1 to specify the transparency of the undeformed
5814 geometry. Set to 1 for completely opaque, and 0 for completely
5815 transparent (invisible). The default is 0.25.
5816 deformed_opacity : float, optional
5817 A float between 0 and 1 to specify the transparency of the deformed
5818 geometry. Set to 1 for completely opaque, and 0 for completely
5819 transparent (invisible). The default is 1.0.
5820 starting_scale : float, optional
5821 The starting scale factor of the animation. The default is 1.0.
5823 Returns
5824 -------
5825 deflection_shape_data
5826 A reference to the deflection_shape_data object that is created to
5827 plot the animated shapes
5829 """
5830 if IGNORE_PLOTS:
5831 return None
5832 return DeflectionShapePlotter(self, deflection_shape_data, plot_kwargs, background_plotter_kwargs,
5833 undeformed_opacity, deformed_opacity=deformed_opacity, starting_scale=starting_scale)
5835 def plot_transient(self, displacement_data, displacement_scale=1.0,
5836 frames_per_second=20,
5837 undeformed_opacity=0.0, deformed_opacity=1.0,
5838 plot_kwargs={},
5839 transformation_shapes=None,
5840 num_curves=50,
5841 show: bool = True,
5842 app=None,
5843 window_size=None,
5844 off_screen=None,
5845 allow_quit_keypress=True,
5846 toolbar=True,
5847 menu_bar=True,
5848 editor=False,
5849 update_app_icon=None,
5850 **kwargs):
5851 """
5852 Create a TransientPlotter object to plot displacements over time
5854 Parameters
5855 ----------
5856 displacement_data : TimeHistoryArray
5857 Transient displacement data that will be applied
5858 displacement_scale : float, optional
5859 Scale factor applied to displacements. The default is 1.0.
5860 frames_per_second : float, optional
5861 Number of time steps to plot per second while the displacement is
5862 animating. Default is 20.
5863 undeformed_opacity : float, optional
5864 Opacity of the undeformed geometry. The default is 0.0, or
5865 completely transparent.
5866 deformed_opacity : float, optional
5867 Opacity of the deformed geometry. The default is 1.0, or completely
5868 opaque.
5869 plot_kwargs : dict, optional
5870 Keyword arguments passed to the Geometry.plot function
5871 transformation_shapes : ShapeArray
5872 Shape matrix that will be used to expand the data. Must be the
5873 same size as the `displacement_data`
5874 num_curves : int, optional
5875 Maximum number of curves to plot on the time selector. Default is
5876 50.
5877 show : bool, optional
5878 Show the plotting window. If ``False``, show this window by
5879 running ``show()``. The default is True.
5880 app : QApplication, optional
5881 Creates a `QApplication` if left as `None`. The default is None.
5882 window_size : list of int, optional
5883 Window size in pixels. Defaults to ``[1024, 768]``
5884 off_screen : TYPE, optional
5885 Renders off screen when True. Useful for automated
5886 screenshots or debug testing. The default is None.
5887 allow_quit_keypress : bool, optional
5888 Allow user to exit by pressing ``"q"``. The default is True.
5889 toolbar : bool, optional
5890 If True, display the default camera toolbar. Defaults to True.
5891 menu_bar : bool, optional
5892 If True, display the default main menu. Defaults to True.
5893 editor : TYPE, optional
5894 If True, display the VTK object editor. Defaults to False.
5895 update_app_icon : bool, optional
5896 If True, update_app_icon will be called automatically to update the
5897 Qt app icon based on the current rendering output. If None, the
5898 logo of PyVista will be used. If False, no icon will be set.
5899 Defaults to None. The default is None.
5900 title : str, optional
5901 Title of plotting window.
5902 multi_samples : int, optional
5903 The number of multi-samples used to mitigate aliasing. 4 is a
5904 good default but 8 will have better results with a potential
5905 impact on performance.
5906 line_smoothing : bool, optional
5907 If True, enable line smothing
5908 point_smoothing : bool, optional
5909 If True, enable point smothing
5910 polygon_smoothing : bool, optional
5911 If True, enable polygon smothing
5912 auto_update : float, bool, optional
5913 Automatic update rate in seconds. Useful for automatically
5914 updating the render window when actors are change without
5915 being automatically ``Modified``. If set to ``True``, update
5916 rate will be 1 second.
5918 Returns
5919 -------
5920 TransientPlotter
5921 """
5922 if IGNORE_PLOTS:
5923 return None
5924 return TransientPlotter(self, displacement_data, displacement_scale,
5925 frames_per_second,
5926 undeformed_opacity, deformed_opacity,
5927 plot_kwargs,
5928 transformation_shapes,
5929 num_curves,
5930 show,
5931 app,
5932 window_size,
5933 off_screen,
5934 allow_quit_keypress,
5935 toolbar,
5936 menu_bar,
5937 editor,
5938 update_app_icon,
5939 **kwargs)
5941 def __add__(self, geometry):
5942 if not isinstance(self, geometry.__class__):
5943 return NotImplemented
5944 final_fields = {}
5945 for field in ['node', 'traceline', 'element', 'coordinate_system']:
5946 # Look through and see which are in common
5947 common_ids = np.intersect1d(getattr(self, field).id, getattr(geometry, field).id)
5948 if common_ids.size > 0:
5949 equal_ids = getattr(self, field)(common_ids) == getattr(geometry, field)(common_ids)
5950 if not all(equal_ids):
5951 raise ValueError('Both geometries contain {:} with ID {:} but they are not equivalent'.format(
5952 field, common_ids[~equal_ids]))
5953 self_ids = np.concatenate(
5954 (np.setdiff1d(getattr(self, field).id, getattr(geometry, field).id), common_ids))
5955 geometry_ids = np.setdiff1d(getattr(geometry, field).id, getattr(self, field).id)
5956 final_fields[field] = np.concatenate(
5957 (getattr(self, field)(self_ids), getattr(geometry, field)(geometry_ids)))
5958 return Geometry(**final_fields)
5960 def modify_ids(self, node_change=0, traceline_change=0, element_change=0, coordinate_system_change=0):
5961 """
5962 Shifts the id numbers in the geometry
5964 Parameters
5965 ----------
5966 node_change : int, optional
5967 The amount to shift the node ids. The default is 0.
5968 traceline_change : int, optional
5969 The amount to shift the traceline ids.. The default is 0.
5970 element_change : int, optional
5971 The amount to shift the element ids. The default is 0.
5972 coordinate_system_change : int, optional
5973 The amount to shift the coordinate system ids. The default is 0.
5975 Returns
5976 -------
5977 geom_out : Geometry
5978 A copy of the original geometry with id numbers modified
5980 """
5981 geom_out = self.copy()
5982 geom_out.node.id += node_change
5983 geom_out.traceline.id += traceline_change
5984 geom_out.element.id += element_change
5985 geom_out.coordinate_system.id += coordinate_system_change
5986 geom_out.node.def_cs += coordinate_system_change
5987 geom_out.node.disp_cs += coordinate_system_change
5988 geom_out.traceline.connectivity += node_change
5989 geom_out.element.connectivity += node_change
5990 # Set values of tracelines back to zero if they just equal node_change,
5991 # because that means they were zero and should represent picking up the
5992 # pen
5993 for connectivity_array in geom_out.traceline.connectivity:
5994 connectivity_array[connectivity_array == node_change] = 0
5995 return geom_out
5997 def map_ids(self, node_id_map=None, traceline_id_map=None, element_id_map=None, coordinate_system_id_map=None):
5998 """
5999 Maps id numbers from an original set of ids to a new set of ids.
6001 This function accepts id_map classes defining "from" and "to" ids
6002 Existing ids found in the "from" set are transformed to the corresponding
6003 id in the "to" set.
6005 Parameters
6006 ----------
6007 node_id_map : id_map, optional
6008 An id_map defining the mapping applied to node ids.
6009 The default is None, which results in the ids being unchanged
6010 traceline_id_map : id_map, optional
6011 An id_map defining the mapping applied to traceline ids.
6012 The default is None, which results in the ids being unchanged
6013 element_id_map : id_map, optional
6014 An id_map defining the mapping applied to element ids.
6015 The default is None, which results in the ids being unchanged
6016 coordinate_system_id_map : id_map, optional
6017 An id_map defining the mapping applied to coordinate system ids.
6018 The default is None, which results in the ids being unchanged
6020 Returns
6021 -------
6022 geom_out : Geometry
6023 A copy of the original geometry with id numbers modified
6025 """
6026 geom_out = self.copy()
6027 if node_id_map is not None:
6028 geom_out.node.id = node_id_map(self.node.id)
6029 for key, val in self.traceline.ndenumerate():
6030 geom_out.traceline.connectivity[key] = node_id_map(self.traceline.connectivity[key])
6031 for key, val in self.element.ndenumerate():
6032 geom_out.element.connectivity[key] = node_id_map(self.element.connectivity[key])
6033 if coordinate_system_id_map is not None:
6034 geom_out.coordinate_system.id = coordinate_system_id_map(self.coordinate_system.id)
6035 geom_out.node.disp_cs = coordinate_system_id_map(self.node.disp_cs)
6036 geom_out.node.def_cs = coordinate_system_id_map(self.node.def_cs)
6037 if traceline_id_map is not None:
6038 geom_out.traceline.id = traceline_id_map(self.traceline.id)
6039 if element_id_map is not None:
6040 geom_out.element.id = element_id_map(self.element.id)
6041 return geom_out
6043 def response_kinematic_transformation(self, response_coordinate,
6044 virtual_point_node_number,
6045 virtual_point_location=[0,0,0]):
6046 """
6047 Creates a kinematic response transformation from the rigid body body
6048 shapes for the associated geometry and response coordinates.
6050 Parameters
6051 ----------
6052 response_coordinate : CoordinateArray
6053 The "physical" response coordinates for the kinematic transformation.
6054 virtual_point_node_number : int
6055 The numeric ID for the virtual point in the kinematic transformation.
6056 virtual_point_location : list or ndarray
6057 The [X, Y, Z] coordinates that defined the location of the virtual
6058 point. The default is is [0, 0, 0].
6060 Returns
6061 -------
6062 transformation : Matrix
6063 The kinematic transformation as a matrix object. It is organized with
6064 the virtual point CoordinateArray on the rows and the physical
6065 response CoordinateArray on the columns.
6067 Notes
6068 -----
6069 The transformation automatically handles polarity differences in the geometry
6070 and response_coordinate.
6072 References
6073 ----------
6074 .. [1] M. Van der Seijs, D. van den Bosch, D. Rixen, and D. Klerk, "An improved
6075 methodology for the virtual point transformation of measured frequency
6076 response functions in dynamic substructuring," in Proceedings of the 4th
6077 International Conference on Computational Methods in Structural Dynamics
6078 and Earthquake Engineering, Kos Island, 2013, pp. 4334-4347,
6079 doi: 10.7712/120113.4816.C1539.
6080 """
6081 shapes = self.rigid_body_shapes(response_coordinate, cg = virtual_point_location)
6082 transformation_array = np.linalg.pinv(shapes.shape_matrix.T)
6084 virtual_point_coordinate = coordinate_array(node=virtual_point_node_number, direction=[1,2,3,4,5,6])
6086 return matrix(transformation_array, virtual_point_coordinate, response_coordinate)
6088 def force_kinematic_transformation(self, force_coordinate,
6089 virtual_point_node_number,
6090 virtual_point_location=[0,0,0]):
6091 """
6092 Creates a kinematic force transformation from the rigid body body
6093 shapes for the associated geometry and force coordinates.
6095 Parameters
6096 ----------
6097 force_coordinate : CoordinateArray
6098 The "physical" force coordinates for the kinematic transformation.
6099 virtual_point_node_number : int
6100 The numeric ID for the virtual point in the kinematic transformation.
6101 virtual_point_location : list or ndarray
6102 The [X, Y, Z] coordinates that defined the location of the virtual
6103 point. The default is is [0, 0, 0].
6105 Returns
6106 -------
6107 transformation : Matrix
6108 The kinematic transformation as a matrix object. It is organized with
6109 the virtual point CoordinateArray on the rows and the physical
6110 force CoordinateArray on the columns.
6112 Notes
6113 -----
6114 The transformation automatically handles polarity differences in the geometry
6115 and force_coordinate.
6117 References
6118 ----------
6119 .. [1] M. Van der Seijs, D. van den Bosch, D. Rixen, and D. Klerk, "An improved
6120 methodology for the virtual point transformation of measured frequency
6121 response functions in dynamic substructuring," in Proceedings of the 4th
6122 International Conference on Computational Methods in Structural Dynamics
6123 and Earthquake Engineering, Kos Island, 2013, pp. 4334-4347,
6124 doi: 10.7712/120113.4816.C1539.
6125 """
6126 shapes = self.rigid_body_shapes(force_coordinate, cg = virtual_point_location)
6127 virtual_point_coordinate = coordinate_array(node=virtual_point_node_number, direction=[1,2,3,4,5,6])
6129 return matrix(shapes.shape_matrix, virtual_point_coordinate, force_coordinate)
6131 def copy(self):
6132 """
6133 Return's a copy of the current Geometry
6135 Changes to the copy will not also be applied to the original geometry
6137 Returns
6138 -------
6139 Geometry
6140 A copy of the current Geometry
6142 """
6143 return Geometry(self.node.copy(), self.coordinate_system.copy(), self.traceline.copy(), self.element.copy())
6145 def save(self, filename):
6146 """Saves the geometry to a numpy .npz file
6148 The .npz file will have fields 'node', 'coordinate_system', 'element',
6149 and 'traceline', with each field storing the respective portion of the
6150 geometry
6152 Parameters
6153 ----------
6154 filename : str
6155 Filename to save the geometry to. If the filename doesn't end with
6156 .npz, it will be appended.
6158 Returns
6159 -------
6160 None.
6162 """
6163 np.savez(filename,
6164 node=self.node.view(np.ndarray),
6165 coordinate_system=self.coordinate_system.view(np.ndarray),
6166 element=self.element.view(np.ndarray),
6167 traceline=self.traceline.view(np.ndarray))
6169 def savemat(self, filename):
6170 """Saves the geometry to a .mat file
6172 The .mat file will have fields 'node', 'coordinate_system', 'element',
6173 and 'traceline', with each field storing the respective portion of the
6174 geometry
6176 Parameters
6177 ----------
6178 filename : str
6179 Filename to save the geometry to.
6181 Returns
6182 -------
6183 None.
6184 """
6186 mat_dict = {key: getattr(self, key).assemble_mat_dict() for key in
6187 ['node', 'coordinate_system', 'element', 'traceline']}
6188 scipy_savemat(filename,mat_dict)
6190 def mass_props_from_rigid_body_shapes(self, rigid_shapes, return_rectified_shapes = False,
6191 rectified_rotations_about_center_of_mass = False):
6192 """
6193 Computes mass properties (mass, inertia tensor, center of mass) given
6194 a set of mass-normalized rigid body shapes.
6196 Parameters
6197 ----------
6198 rigid_shapes : ShapeArray
6199 A set of mass normalized rigid body mode shapes. They need not be
6200 in any particular ordering, and may be linear combinations of
6201 translations and rotations. They must, however, span the space of
6202 rigid body modes.
6203 return_rectified_shapes : bool
6204 If True, a new set of shapes will be returned that will be rectified
6205 to the three pure translations, and the three pure rotations about
6206 the center of mass.
6207 rectified_rotations_about_center_of_mass : bool
6208 If True, the rectified shapes will be constructed such that they
6209 rotate about the test article's center of mass.
6211 Returns
6212 -------
6213 mass_properties : dict
6214 Contains:
6215 - 'A' : (6, 6) transformation matrix such that B ≈ Phi @ A
6216 - 'G' : (6, 6) rigid-body generalized mass matrix in basis B
6217 - 'mass' : scalar total mass
6218 - 'com_relative' : (3,) center of mass relative to rotation reference point
6219 - 'com_global' : (3,) center of mass in global coordinates (if center_of_rotation provided)
6220 - 'inertia_about_reference' : (3, 3) inertia tensor about the rotation reference point
6221 - 'inertia_about_com' : (3, 3) inertia tensor about center of mass
6222 - 'principal_inertias' : (3,) eigenvalues of inertia_about_com
6223 - 'principal_axes' : (3, 3) eigenvectors of inertia_about_com
6224 - 'reconstruction_error' : relative Frobenius norm of B - Phi @ A
6225 - 'G_tt', 'G_tr', 'G_rt', 'G_rr' : 3x3 blocks of G
6226 rectified_rigid_shapes : ShapeArray
6227 Scaled versions of the rigid body shapes consisting of three pure
6228 translation and three pure rotations.
6230 """
6231 if rigid_shapes.size != 6:
6232 raise ValueError('There must be six rigid shapes.')
6233 if not np.all(rigid_shapes.ravel()[0].coordinate == rigid_shapes.coordinate):
6234 raise ValueError('All rigid shapes must have the same coordinate.')
6235 coordinates = rigid_shapes.ravel()[0].coordinate
6236 rbm_modes = rigid_shapes.ravel().shape_matrix.T
6237 rb_shapes = self.rigid_body_shapes(coordinates)
6238 translation_shapes = rb_shapes[:3].shape_matrix.T
6239 rotation_shapes = rb_shapes[3:].shape_matrix.T
6240 scaled_shapes, mass_properties = rigid_body_shape_scaling(
6241 rbm_modes, translation_shapes, rotation_shapes,
6242 center_of_rotation=np.array([0,0,0]),
6243 return_mass_properties = True)
6244 if return_rectified_shapes:
6245 if rectified_rotations_about_center_of_mass:
6246 com = mass_properties['com_global']
6247 rb_shapes = self.rigid_body_shapes(coordinates, cg=com)
6248 translation_shapes = rb_shapes[:3].shape_matrix.T
6249 rotation_shapes = rb_shapes[3:].shape_matrix.T
6250 scaled_shapes, mass_properties_com = rigid_body_shape_scaling(
6251 rbm_modes, translation_shapes, rotation_shapes,
6252 center_of_rotation=com,
6253 return_mass_properties = True)
6254 rectified_rigid_shapes = rigid_shapes.copy()
6255 rectified_rigid_shapes.frequency = 0
6256 rectified_rigid_shapes.damping = 0
6257 rectified_rigid_shapes.comment1 = [
6258 'Translation X','Translation Y','Translation Z',
6259 'Rotation X','Rotation Y','Rotation Z']
6260 rectified_rigid_shapes.comment2 = ''
6261 rectified_rigid_shapes.comment3 = ''
6262 rectified_rigid_shapes.comment4 = ''
6263 rectified_rigid_shapes.comment5 = ''
6264 rectified_rigid_shapes.shape_matrix = scaled_shapes.T
6265 return mass_properties, rectified_rigid_shapes
6266 else:
6267 return mass_properties
6270 def rigid_body_shapes(self, coordinates, mass=1, inertia=np.eye(3), cg=np.zeros(3), principal_axes=False):
6271 """
6272 Creates a set of shapes corresponding to the rigid body motions
6274 Rigid body translation and rotation shapes are computed analytically
6275 from the Geoemtry
6277 Parameters
6278 ----------
6279 coordinates : CoordinateArray
6280 coordinates at which to compute deformations
6281 mass : float, optional
6282 The mass of the geometry used to scale the rigid body translation
6283 shapes. The default is 1.
6284 inertia : np.ndarray, optional
6285 A 3x3 array consisting of the mass moments of inertia, used to
6286 scale the rotation shapes. The default is np.eye(3).
6287 cg : np.ndarray, optional
6288 The center of gravity of the Geometry about which the rotations
6289 occur. The default is np.zeros(3).
6290 principal_axes : bool, optional
6291 If True, compute the principal axes of the test article and
6292 perform rotations about those axes. The default is False.
6294 Returns
6295 -------
6296 output_shape : ShapeArray
6297 A set of rigid body shapes for the current geometry
6299 """
6300 if principal_axes:
6301 rotation_scale, rotation_directions = np.linalg.eig(inertia)
6302 rotation_directions /= np.linalg.norm(rotation_directions, axis=0, keepdims=True)
6303 else:
6304 rotation_scale = np.diag(inertia)
6305 rotation_directions = np.eye(3)
6306 M = np.diag(np.concatenate((mass * np.ones(3), rotation_scale)))
6307 phi = np.sqrt(np.linalg.inv(M))
6308 # Construct the full coordinate array for all degrees of freedom for each node
6309 full_coordinates = coordinate_array(np.unique(coordinates.node)[:, np.newaxis],
6310 ['X+', 'Y+', 'Z+']).flatten()
6311 # Assume local coordinates for now, transform later
6312 translation_shape_matrix = full_coordinates.local_direction().T/np.sqrt(mass)
6313 rotation_shape_matrix = np.zeros(translation_shape_matrix.shape)
6314 for i in range(3):
6315 for j in range(len(full_coordinates)):
6316 direction = full_coordinates[j].local_direction()
6317 coord = self.global_node_coordinate(full_coordinates[j].node)
6318 rotation_shape_matrix[i, j] = phi[i + 3, i + 3] * \
6319 np.dot(direction, np.cross(rotation_directions[:, i], coord - cg))
6320 # Concatenate shapes together
6321 shape_matrix = np.concatenate((translation_shape_matrix, rotation_shape_matrix), axis=0)
6322 full_shapes = shape_array(full_coordinates, shape_matrix)
6323 # Now we need to transform to the local coordinate systems
6324 # First set up a geometry with just the global coordinate system in it
6325 global_geometry = Geometry(node=self.node.copy(),
6326 coordinate_system=coordinate_system_array((1,)))
6327 global_geometry.node.disp_cs = 1
6328 global_geometry.node.def_cs = 1
6329 transformed_shape = full_shapes.transform_coordinate_system(global_geometry, self)
6330 # Now turn it into our requested coordinates
6331 reduced_shape_matrix = transformed_shape[coordinates]
6332 output_shape = shape_array(coordinates, reduced_shape_matrix)
6333 return output_shape
6335 def coordinate_transformation_matrix(self, to_geometry : 'Geometry',
6336 nodes : np.ndarray = None,
6337 rotations : bool = False):
6338 """
6339 Creates a transformation matrix that transforms one geometry to a new
6340 geometry with different coordinate systems defined.
6342 Parameters
6343 ----------
6344 to_geometry : Geometry
6345 Geometry object with coordinate systems that the transformation
6346 matrix will transform into
6347 nodes : np.ndarray, optional
6348 An array of node id numbers to include in the transformation.
6349 The default is to include all nodes in the geometry.
6350 rotations : bool, optional
6351 If True, create degrees of freedom for rotations as well as
6352 translations.
6355 Returns
6356 -------
6357 Matrix
6358 A Matrix object that transforms data or shapes currently
6359 represented by this Geometry to the Geometry specified in the
6360 `to_geometry` argument.
6361 """
6362 if nodes is None:
6363 nodes = np.intersect1d(self.node.id,to_geometry.node.id)
6364 else:
6365 nodes = np.ravel(nodes)
6366 if rotations:
6367 dofs = coordinate_array(nodes[:,np.newaxis,np.newaxis],[[1, 2, 3],
6368 [4, 5, 6]])
6369 else:
6370 dofs = coordinate_array(nodes[:,np.newaxis],[1, 2, 3])
6371 transform_from_original = self.global_deflection(dofs)
6372 transform_to_new = to_geometry.global_deflection(dofs)
6373 # We need to transform this matrix into block-diagonal form.
6374 transform_from_original = transform_from_original.reshape(-1,3,3)
6375 transform_to_new = transform_to_new.reshape(-1,3,3)
6376 transform_from_original = block_diag(*transform_from_original)
6377 transform_to_new = block_diag(*transform_to_new)
6378 transform = transform_to_new @ np.linalg.pinv(transform_from_original)
6379 transformation_matrix = matrix(transform,dofs.flatten(),dofs.flatten())
6380 return transformation_matrix
6382 @classmethod
6383 def load(cls, filename):
6384 """
6385 Loads a geometry from a numpy .npz or .unv file
6387 The .npz file must have fields 'node', 'coordinate_system', 'element',
6388 and 'traceline', with each field storing the respective portion of the
6389 geometry
6391 The .unv file will need to have the proper datasets specified to define
6392 a geometry
6394 Parameters
6395 ----------
6396 filename : str
6397 Filename to load the geometry from.
6399 Raises
6400 ------
6401 AttributeError
6402 Raised if the calling class does not have a `from_unv` method
6403 defined
6405 Returns
6406 -------
6407 Geometry
6408 Geometry constructed from the data in the loaded file
6410 """
6411 if filename[-4:].lower() in ['.unv', '.uff']:
6412 try:
6413 from ..fileio.sdynpy_uff import readunv
6414 unv_dict = readunv(filename)
6415 return cls.from_unv(unv_dict)
6416 except AttributeError:
6417 raise AttributeError('Class {:} has no from_unv attribute defined'.format(cls))
6418 else:
6419 with np.load(filename, allow_pickle=True) as data:
6420 return Geometry(node=data['node'].view(NodeArray),
6421 coordinate_system=data['coordinate_system'].view(
6422 CoordinateSystemArray),
6423 traceline=data['traceline'].view(TracelineArray),
6424 element=data['element'].view(ElementArray))
6426 def write_to_unv(self, filename, write_nodes=True, write_coordinate_systems=True,
6427 write_tracelines=True, write_elements=True,
6428 dataset_2412_kwargs={},
6429 dataset_2420_kwargs={}):
6430 """
6431 Write the geometry to a unversal file format file
6433 Parameters
6434 ----------
6435 filename : str
6436 Filename to which the geometry will be written. If None, a
6437 unv data dictionary will be returned instead, similar to that
6438 obtained from the readunv function in sdynpy
6439 write_nodes : bool, optional
6440 If True, write the geometry's nodes to dataset 2411 in the output
6441 file. The default is True.
6442 write_coordinate_systems : True, optional
6443 If True, write the geometry's coordinate systems to dataset 2420
6444 in the output file. The default is True.
6445 write_tracelines : bool, optional
6446 If True, write the geometry's tracelines to dataset 82 in the
6447 output file. The default is True.
6448 write_elements : TYPE, optional
6449 If True, write the geometry's elements to dataset 2412 in the
6450 output file. The default is True.
6451 dataset_2412_kwargs : dict, optional
6452 Allows users to specify additional element parameters not stored
6453 by the Geometry. The default is {}.
6454 dataset_2420_kwargs : dict, optional
6455 Allows users to specify additional coordinate system parameters not
6456 stored by the Geometry. The default is {}.
6458 Returns
6459 -------
6460 unv_dict : dict
6461 Dictionary containing unv information, similar to that obtained from
6462 readunv. Only returned if filename is None.
6464 """
6465 from ..fileio.sdynpy_uff import dataset_82, dataset_2411, dataset_2412, dataset_2420
6466 # Write Coordinate systems
6467 if write_coordinate_systems:
6468 if self.coordinate_system.size > 0:
6469 cs_unv = dataset_2420.Sdynpy_UFF_Dataset_2420(
6470 cs_labels=self.coordinate_system.id.flatten(),
6471 cs_types=self.coordinate_system.cs_type.flatten(),
6472 cs_colors=self.coordinate_system.color.flatten(),
6473 cs_names=[cs.name for index, cs in self.coordinate_system.ndenumerate()],
6474 cs_matrices=self.coordinate_system.matrix.reshape(-1, 4, 3),
6475 **dataset_2420_kwargs)
6476 else:
6477 cs_unv = None
6478 if write_nodes:
6479 if self.node.size > 0:
6480 node_unv = dataset_2411.Sdynpy_UFF_Dataset_2411(
6481 self.node.id.flatten(),
6482 self.node.def_cs.flatten(),
6483 self.node.disp_cs.flatten(),
6484 self.node.color.flatten(),
6485 self.node.coordinate.reshape(-1, 3))
6486 else:
6487 node_unv = None
6488 if write_elements:
6489 if self.element.size > 0:
6490 if 'physical_property_table_numbers' not in dataset_2412_kwargs:
6491 dataset_2412_kwargs['physical_property_table_numbers'] = [1] * self.element.size
6492 if 'material_property_table_numbers' not in dataset_2412_kwargs:
6493 dataset_2412_kwargs['material_property_table_numbers'] = [1] * self.element.size
6494 if 'beam_orientations' not in dataset_2412_kwargs:
6495 dataset_2412_kwargs['beam_aft_cross_section_numbers'] = [
6496 None] * self.element.size
6497 if 'beam_aft_cross_section_numbers' not in dataset_2412_kwargs:
6498 dataset_2412_kwargs['beam_orientations'] = [None] * self.element.size
6499 if 'beam_fore_cross_section_numbers' not in dataset_2412_kwargs:
6500 dataset_2412_kwargs['beam_fore_cross_section_numbers'] = [
6501 None] * self.element.size
6502 elem_unv = dataset_2412.Sdynpy_UFF_Dataset_2412(
6503 element_labels=self.element.id.flatten(),
6504 fe_descriptor_ids=self.element.type.flatten(),
6505 colors=self.element.color.flatten(),
6506 connectivities=[elem.connectivity for index,
6507 elem in self.element.ndenumerate()],
6508 **dataset_2412_kwargs)
6509 else:
6510 elem_unv = None
6511 if write_tracelines:
6512 tl_unv = []
6513 for index, tl in self.traceline.ndenumerate():
6514 tl_unv.append(
6515 dataset_82.Sdynpy_UFF_Dataset_82(
6516 tl.id, tl.color, tl.description, tl.connectivity,
6517 ))
6518 if filename is None:
6519 unv_dict = {}
6520 unv_dict[2420] = [cs_unv]
6521 unv_dict[2411] = [node_unv]
6522 unv_dict[2412] = [elem_unv]
6523 unv_dict[82] = tl_unv
6524 return unv_dict
6525 else:
6526 with open(filename, 'w') as f:
6527 if cs_unv is not None:
6528 f.write(' -1\n')
6529 f.write(' 2420\n')
6530 f.write(cs_unv.write_string())
6531 f.write(' -1\n')
6532 if node_unv is not None:
6533 f.write(' -1\n')
6534 f.write(' 2411\n')
6535 f.write(node_unv.write_string())
6536 f.write(' -1\n')
6537 if elem_unv is not None:
6538 f.write(' -1\n')
6539 f.write(' 2412\n')
6540 f.write(elem_unv.write_string())
6541 f.write(' -1\n')
6542 for dataset in tl_unv:
6543 f.write(' -1\n')
6544 f.write(' 82\n')
6545 f.write(dataset.write_string())
6546 f.write(' -1\n')
6549from_excel_template = Geometry.from_excel_template
6550write_excel_template = Geometry.write_excel_template
6551from_imat_struct = Geometry.from_imat_struct
6552from_exodus = Geometry.from_exodus
6553from_unv = Geometry.from_unv
6554from_uff = Geometry.from_uff
6555load = Geometry.load
6558class id_map:
6559 """Class defining mapping between two sets of id numbers"""
6561 def __init__(self, from_ids, to_ids):
6562 """
6563 Initializes the id map
6565 Parameters
6566 ----------
6567 from_ids : np.ndarray
6568 Id numbers to map from
6569 to_ids : np.ndarray
6570 Id numbers to map to
6572 Returns
6573 -------
6574 None.
6576 """
6577 self.from_ids = from_ids
6578 self.to_ids = to_ids
6579 map_dict = {f: t for f, t in zip(from_ids, to_ids)}
6580 map_dict[0] = 0
6581 self.mapper = np.vectorize(map_dict.__getitem__, otypes=['uint64'])
6583 def __call__(self, val):
6584 """
6585 Map id numbers
6587 Parameters
6588 ----------
6589 val : np.ndarray
6590 A set of id numbers to map. These must all be contained in the
6591 id_map's from_ids array.
6593 Returns
6594 -------
6595 np.ndarray
6596 The set of id numbers from the id_map's to_ids array corresponding
6597 to those passed in to the function from the from_ids array.
6599 """
6600 return self.mapper(val)
6602 def inverse(self):
6603 """
6604 Returns an inverse map, swapping the from and to ids.
6606 If duplicate id entries exist, they will be overwritten.
6608 Returns
6609 -------
6610 id_map
6611 An id_map object with swapped from and to ids.
6613 """
6614 return id_map(self.to_ids,self.from_ids)