Geometry

In structural dynamics, geometry generally refers to the spatial layout of the measurements or analysis points. For example, in a modal test, we will generally have sensors distributed throughout the test article. The physical locations and orientations of those sensors represent the test geometry, and must be documented in order to map measurements taken during the test to physical locations on the test article. Similarly, in a finite element analysis, locations of nodes in the model will represent the analysis geometry.

This document will describe how geometry is defined and used in SDynPy. To start, a generic geometry will be created to demonstrate the process and describe the objects used. Subsequently, several example geometries of varying complexity will be created.

Let’s import SDynPy and start looking at geometry.

SDynPy Geometry Objects¶

In SDynPy, the primary data type used to represent test or analysis geometry is the Geometry class defined in the sdynpy.geometry module. A typical SDynPy Geometry has four components:

• Nodes: These are locations defined within the geometry and are represented by the NodeArray class. Associated with each node is a unique identification number, a coordinate in 3D space, and references to the identification numbers of the coordinate systems used to define the nodes position and orientation.

• Coordinate Systems: These are coordinate systems that are used to define node positions and orientations and are represented by the CoordinateSystemArray class. A simple geometry may have only one global coordinate used to define all node positions and orientations. A more complex geometry may have local coordinate systems defined relative to specific coordinates or measurements.

• Tracelines: These are lines drawn between nodes to enhance visualization of the geometry and are represented by the TracelineArray class. Without visualization aids, the geometry would simply be a point cloud, which can be difficult to interpret, as it is not clear which points are in front of or behind other points.

• Elements: These are lines (e.g., bars), faces (e.g., triangles, quadrilaterals), or volumes (e.g., hexahedra, tetrahedra) drawn between nodes to enhance visualization of the geometry and are represented by the ElementArray class. Without visualization aids, the geometry would simply be a point cloud, which can be difficult to interpret, as it is not clear which points are in front of or behind other points.

The general process to create a Geometry object is to first construct the desired components (NodeArray, CoordinateSystemArray, TracelineArray, and ElementArray) and then to pass them to the Geometry constructor.

Creating Geometry from Scratch¶

When performing a test, personnel instrumenting the test article should ideally be making measurements of the positions and orientations of each sensor. Similarly, for a finite element analysis positions of the nodes in the model can usually be extracted from the mesh file. These positions and orientations may be defined in many different ways depending on the test or analysis type. A cylindrical test article may most conveniently represent node positions in a cylindrical coordinate system, whereas test articles of other shapes may be better represented in cartesian or cylindrical coordinate systems. Orientation data may be represented in a multitude of formats: unit vectors of sensor directions, rotation matrices, quaternions, euler angles, etc. This information may be entered directly into a Python script to build the geometry, or it may come from an external file that must be read into Python. While SDynPy has several readers of common file formats which will be covered later in this section, the parsing of arbitrary file formats to extract desired information is outside the scope of this document. Therefore, for examples in this section, it will be assumed that the user can extract geometry data into the Python workspace in order to build these objects. This document will then show how that data can be packaged into a geometry object to represent the data.

Defining the Coordinate Systems in the Geometry¶

The first step to define a geometry is might intuitively be to define the nodes in the model; however, node positions are all defined with respect to some coordinate system, so the true first step should generally be to define the coordinate systems that the geometry will use. Simple geometries may only have one global coordinate system. More complex geometries may have an arbitrary number of coordinate systems defined.

Every node is assigned both a definition coordinate system and a displacement coordinate system. The definition coordinate system is the coordinate system in which the node positions are defined. The displacement coordinate system is the coordinate system in which the node deforms. The rationale for having two coordinate systems for each node is to allow the node to displace in a separate coordinate system from which it is defined. For example, a conical test article may have its nodes defined most conveniently in a cylindrical coordinate system with its radius from the cylindrical axis defined by the cone angle. However, if sensors are placed normal to the surface of the cone, they will not be aligned with the coordinate system used to define their location.

In the above photo, the test article is the green cone, and the blue square represents a sensor glued to the surface of the cone. The sensor will move in its own local coordinate system due to it being mounted obliquely to the part’s coordinate system, so it must have it’s own coordinate system defined. However, it may be significant additional effort to compute the position of the sensor in this local coordinate system, so it is allowable to define the sensor’s position in the global part coordinate system.

For simple geometries containing only a global coordinate system, the definition and displacement coordinate systems may be the same coordinate system.

For this initial demonstration, we will construct two coordinate systems, the first is a global cartesian coordinate system, the second will be a cylindrical coordinate system. Coordinate systems are defined in SDynPy using the CoordinateSystemArray object. Like most SdynpyArray objects, it has a helper function to allow easier definition, which is the same name as the object but in snake_case capitalization. In this case, the helper function is coordinate_system_array. We will define the coordinate systems unique identification numbers, as well as their type. SDynPy uses the integer 0 to represent a cartesian coordinate system, 1 to represent a cylindrical or polar coordinate system, and 2 to represent a spherical coordinate system. We also assign a name to each coordinate system, to make it more obvious what they represent. We could also assign a matrix representing the coordinate systems’ transformations from the global coordinate system; however, for simplicity at this point we will not do that. Later examples in this section will creating coordinate systems with custom matrices.

Let’s briefly explore the CoordinateSystemArray object before continuing. We can see by simply typing in the variable name into the terminal, we get a representation of the CoordinateSystemArray object.

To see the data fields of a CoordinateSystemArray, we can use the fields property or examine its dtype for more information.

To access one of the fields, we can simply use the field name like any other attribute of the object.

These are the data that are stored within each coordinate system within a CoordinateSystemArray. Each coordinate system must have a unique positive integer identifier with which to reference it, and this is stored in the id field. Additionally, the name field allows a 40 character string to be defined that provides a name for the coordinate system. A color can also be assigned to a coordinate system; SDynPy does not actually apply colors to its coordinate system, but this can be defined to maintain consistency with Universal File Format Dataset 2420. Colors are defined using an integer from 0 to 15; see sdynpy.colors for the mapping of integers to colors. As stated previously, the coordinate system type is encoded as an integer in the cs_type field. Finally, the transformation matrix of the coordinate system is stored in the matrix field. Note in the dtype output, we see that the size of the matrix field is actually (4,3), meaning for each coordinate system in the CoordinateSystemArray is a 4 x 3 matrix consisting of a 3 x 3 rotation matrix above a 1 x 3 translation vector. See Coordinate System Transformations for more information on defining coordinate system matrices.

Note that when a field is defined containing a shape, the shape of this field will be appended to the shape of the object itself. For example, we have a shape (2,) CoordinateSystemArray:

If we look at the shape of the matrix field, we will see that the field is then (2,4,3):

meaning that there is a (4,3) matrix for each of the (2,) coordinate systems.

Defining Nodes in the Geometry¶

With the coordinate systems defined, we will now define nodes in the geometry. Nodes are represented in SDynPy by the NodeArray object, which again is more conveniently constructed using the node_array function. Let’s construct a NodeArray for the nodes in the Cartesiean coordinate system. Depending on the complexity of the model, it may make sense to construct nodes in chunks and then combine them at the end, and we will show that here. One thing to keep in mind is that each node must have a unique identification number, so if combining nodes at a later point, ensure these identification numbers do not overlap.

We will start by defining the nodes for the cartesian section of the Geometry. We will simply make a 3 x 3 x 3 grid of points.

We see that our positions are a num_nodes by 3 array. Let’s also create the node identification numbers to go along with these points. To help us avoid numbering conflicts, we will start the cartesian nodes at 100.

Now let’s create the NodeArray:

Now let’s do the same for our cylindrical coordinate system. We will start the node numbers for these nodes with 200.

Recall that we can typically utilize NumPy functions with SdynpyArray objects. In this case, we will use np.concatenate to combine the two sets of nodes into one.

Before we go further, let’s also investigate the NodeArray object. Typing in the variable in the terminal will present a formatted node table for convenience.

To see the data fields of a NodeArray, we can use the fields property or examine its dtype for more information.

Similar to the CoordinateSystemArray, we see that the id field, which contains the node ID number, is a 8-byte (64-bit) unsigned integer. The geometry.node.disp_cs and geometry.node.def_cs arrays, which contain references to the coordinate systems in which the node is defined and in which the node displaces, respectively, are also this data type. The color array, while still an unsigned integer, is only 2 bytes, or 16 bits, and again represents the color assigned to each node in the geometry. Finally, the coordinate field, which contains the 3D position of the node as defined in the def_cs coordinate system, consists of 8-byte (64-bit) floating-point data, and also has a shape of (3,), which signifies there are three values of the coordinate for each entry in the NodeArray. These extra dimensions of the field arrays are appended at the end of dimension of the object itself. For example, if we compare the shape of the NodeArray array and its coordinate field, we will see that the shapes are identical except for the appending of the length-3 extra dimension on the latter array.

Again, this means that there are three coordinates for each of the 51 nodes in our NodeArray.

Creating a Geometry¶

With a NodeArray and a CoordinateSystemArray defined, we now have the minimum dataset required to create a Geometry object. We could also create an ElementArray and a TracelineArray prior to constructing the Geometry; however, since these latter types are used only for visualization, it may be useful to first construct and plot the Geometry to understand which nodes should be connected using elements or lines.

We construct a Geometry object using its constructor directly. Note that unlike NodeArray and CoordinateSystemArray, we generally do not work with collections of geometries like we work with collections of nodes or collections of coordinate systems. Therefore, the Geometry object is not a subclass of SdynpyArray, and we use its constructor directly rather than some helper function. Note that this does not preclude storing, for example, multiple Geometry objects in a list if the application requires it. However Geometry objects are not optimized for such operations like NodeArray or CoordinateSystemArray objects.

If we type in the variable name, we get a representation of the Geometry object.

Note that since the Geometry is not a SdynpyArray or numpy.ndarray, there are no fields or dtype. There are instead only regular attributes.

We can access the main properties of the geometry using the node, coordinate_system, element, or traceline attributes.

We see that because we have constructed our geometry with only nodes and coordinate systems, the element and traceline fields are currently empty.

Since Geometry objects represent the spatial distribution of measurement or analysis points in a test or analysis, probably the most common thing to do with a geometry is to plot it.

This will bring up an interactive window in which the geometry can be rotated, panned, zoomed, or otherwise interacted with. There are numerous plotting options to supply to the Geometry method (see e.g. label_nodes in the previous commant). However, as the geometry currently stands, it is rather hard to interpret due to the lack of elements or lines connecting the nodes. It is simply a point-cloud currently, so many of the 3D cues that enable our brain to interpret a 2D image as a 3D scene are missing. We should presently add in those cues to allow us to easier interpret the geometry.

Adding Tracelines¶

Tracelines are lines that connect nodes in the Geometry to aid in visualization. There are two approaches to creating tracelines in the model. The first approach is to call the Geometry method to add tracelines one-by-one. The second approach is to construct a Geometry.plot explicitly and then assign it to Geometry.traceline (alternatively, one could pass it to the traceline argument of the Geometry.add_traceline constructor when initially creating the TracelineArray object). Given that we often would like to be looking at our Geometry object when figuring out which nodes to connect with a traceline, the former approach is usually the most straightforward, though there are case where the latter may be preferred, for instance when there is a pattern in the node numbers that can be exploited.

In this case, we will add tracelines to the cartesian portion of the model. SDynPy uses the convention that a 0 in the traceline connectivity array means to “pick up the pen” and will mean that there will be a gap in the line at that position.

Note that when we call the Geometry.add_traceline method, we can explicitly assign an identification number with the id argument. If we don’t, SDynPy will automatically assign a unique identifier to the traceline that is one larger than the last largest identification number to avoid conflicts.

Now if we plot the geometry, we can see that the added lines help with the interpretation of the cartesian portion of the geometry.

Let’s look at what the traceline attribute to see what we’ve added to the Geometry.

Geometry.traceline gives a TracelineArray object that contains the traceline definitions. We can again look at its fields or dtype.

Similarly to the other geometry objects, TracelineArray objects have id and color. The description field stores a name or description of each item in the TracelineArray as a string with less than 40 characters. Finally, the traceline connectivity field is stored as an object array, where each entry in the array is a NumPy ndarray with length equal to the number of nodes in the traceline. This construction is necessary as each traceline might have a different number of nodes, so a single array of fixed size is not possible.

Note that due to how object arrays are used in NumPy, investigating the shape of the connectivity field will not immediately tell the user how many nodes are in each connectivity array, but will rather just return the shape of the TracelineArray itself (note the dtype definition previously, where the connectivity field has no additional shape associated with it).

However, if we actually index into a single connectivity array, we can then see how big it is.

Alternatively to the Geometry.add_traceline method, we could have constructed a TracelineArray explicitly and then used it in the construction of the Geometry or assigned it to the Geometry.traceline attribute. Like other SdynpyArray objects, there is a traceline_array function to help us construct the object.

We could then assign them to our geometry using the traceline attribute.

Or if we knew the tracelines ahead of time, we could use them in the construction of the geometry directly.

Adding Elements¶

The last portion of the geometry is elements. As opposed to a traceline, which can only connect between two nodes, elements can connect from three nodes in a triangle up to 20+ nodes in a hexahedral element. Elements in SDynPy follow closely the definition in Universal File Format Dataset 2412.

Because elements of different types may have the same number of nodes (e.g. a 4-node element could be a quadrilateral or a tetrahedral element), the type of each element must be explicitly defined. The element types used by SDynPy can be found in the sdpy.geometry._element_types attribute, and are encoded by an integer.

Note that SDynPy’s elements are only used for visualization, therefore SDynPy makes no distinction between, e.g., elements 51, 61, 71, 81, and 91, which are all variants of a triangular element. SDynPy will simply draw a triangle connecting the three nodes for any of these element types.

SDynPy stores element data in the ElementArray object. Similar to adding tracelines, elements can be added via a Geometry.add_element method, or alternatively by explicitly constructing an ElementArray object and assigning it to the Geometry.element attribute or by using it in the Geometry constructor. Again, given that we often would like to be looking at our Geometry object when figuring out which nodes to connect with an element, the former approach is usually the most straightforward, though there are case where the latter may be preferred, for instance when there is a pattern in the node numbers that can be exploited.

For this example, we will add triangular elements on the top row, and quadrilateral elements on the bottom row.

With elements, more interesting plotting options are available to us.

Let’s look at the Geometry.element attribute to see what we have added to the Geometry object.

Geometry.element gives an ElementArray object. We can look at its fields and dtype attributes to see the underlying data.

Similarly to the TracelineArray, there is an id field containing a unique identification number, a color field that specifies the color applied to the element, and a connectivity field stored in an object array, where each entry in the array is a NumPy ndarray with length equal to the number of nodes in the traceline. This construction is necessary as each element might have a different number of nodes, so a single array of fixed size is not possible. The element_type field is also present to store the type of element encoded as an integer.

Note that due to how object arrays are used in NumPy, investigating the shape of the connectivity field will not immediately tell the user how many nodes are in each connectivity array, but will rather just return the shape of the ElementArray itself (note the dtype definition previously, where the connectivity field has no additional shape associated with it).

However, if we index into one of the element’s connectivity, we can see the number of nodes in that element.

Alternatively to the add_element method, we could have constructed a ElementArray explicitly and then used it in the construction of the Geometry or assigned it to the Geometry.element attribute. Like other SdynpyArray objects, there is a element_array function to help us construct the object.

We could then assign the ElementArray explicitly to the geometry by assigning to the element attribute.

Alternatively, we could have used an explicitly constructed ElementArray in the Geometry constructor.

Creating Geometry from an Excel File Template¶

As an alternative to constructing a Geometry object using only code, SDynPy also has the capability to read in geometry data from a formatted Excel workbook using the from_excel_template class method. To generate an empty Excel workbook in the correct format, one can use the write_excel_template static method.

The template file contains a worksheet for the four main parts of the Geometry object; Coordinate Systems, Nodes, Elements, and Trace Lines. Many of the cells in the template contain data validation to ensure only appropriate content is entered.

Defining the Coordinate Systems¶

When defining the coordinate systems for a Geometry, the fields are much the same as when defining the items by code. The most significant different is that rather than defining the coordinate system transformation matrix by a 4 x 3 matrix, it is instead defined by a translation vector and a set of three rotations where the user can specify angle and axes. This could result in a more intuitive definition of rotated coordinate systems.

Defining the Nodes¶

Nodes are defined on the second sheet of the workbook. Users can manually enter identification number, color, location, and coordinate system information.

Defining the Elements¶

Elements are defined on the third sheet of the workbook. Identification number, type, color, and connectivity are all defined. In the workbook, the type is limited to beam, triangle, quadrilateral, or tetrahedron.

Defining the Tracelines¶

Tracelines can be defined on the last sheed of the workbook. Identification number, description, color, and connectivity are defined.

Loading the Geometry from the Completed Template¶

Once the template is filled out, the user can simply load the completed template into a Geometry object using the Geometry.from_excel_template class method.

Reading and Writing Universal Files¶

A common storage mechanism for structural dynamics geometries is the Universal File Format. This file format consists of several dataset definitions. The datasets pertaining to Geometry are Datasets 82 (Tracelines), 2411 (Nodes), 2412 (Elements), and 2420 (Coordinate Systems). SDynPy’s geometry objects are modeled after these datasets and largely contain the same fields. Therefore it should come as no surprise that SDynPy’s geometry is largely compatible with the universal file format.

A geometry can be written to the universal file format by using its write_to_unv method.

This will write a file containing the following text (with the nodes in dataset 2411 and elements in dataset 2412 truncated for brevity).

    -1
  2420
         1                                                                      
                                                                                
         1         0         1                                                  
Cartesian CS                                                                    
   1.0000000000000000e+00   0.0000000000000000e+00   0.0000000000000000e+00     
   0.0000000000000000e+00   1.0000000000000000e+00   0.0000000000000000e+00     
   0.0000000000000000e+00   0.0000000000000000e+00   1.0000000000000000e+00     
   0.0000000000000000e+00   0.0000000000000000e+00   0.0000000000000000e+00     
         2         1         1                                                  
Spherical CS                                                                    
   1.0000000000000000e+00   0.0000000000000000e+00   0.0000000000000000e+00     
   0.0000000000000000e+00   1.0000000000000000e+00   0.0000000000000000e+00     
   0.0000000000000000e+00   0.0000000000000000e+00   1.0000000000000000e+00     
   0.0000000000000000e+00   0.0000000000000000e+00   0.0000000000000000e+00     
    -1
    -1
  2411
       100         1         1         1                                        
  -1.0000000000000000e+00  -1.0000000000000000e+00  -3.0000000000000000e+00     
       101         1         1         1                                        
  -1.0000000000000000e+00  -1.0000000000000000e+00  -2.0000000000000000e+00     

. . .

       221         2         2         7                                        
   1.0000000000000000e+00   3.3500000000000000e+02   1.0000000000000000e+00     
       222         2         2         7                                        
   1.0000000000000000e+00   3.3500000000000000e+02   2.0000000000000000e+00     
       223         2         2         7                                        
   1.0000000000000000e+00   3.3500000000000000e+02   3.0000000000000000e+00     
    -1
    -1
  2412
         1        61         1         1         7         3                    
       222       223       201                                                  
         2        61         1         1         7         3                    
       223       201       202                                                  
         3        64         1         1         7         4                    
       221       222       201       200                                        
         4        61         1         1         7         3                    
       201       202       204                                                  

. . .

        23        61         1         1         7         3                    
       220       222       223                                                  
        24        64         1         1         7         4                    
       218       219       222       221                                        
    -1
    -1
    82
         1        29         1                                                  
Lines in the XY Plane                                                           
       100       103       106       115       124       121       118       109
       100         0       101       104       107       116       125       122
       119       110       101         0       102       105       108       117
       126       123       120       111       102                              
    -1
    -1
    82
         2        29         2                                                  
Lines in the YZ Plane                                                           
       120       123       126       125       124       121       118       119
       120         0       111       114       117       116       115       112
       109       110       111         0       102       105       108       107
       106       103       100       101       102                              
    -1
    -1
    82
         3        29         3                                                  
Lines in the XZ Plane                                                           
       126       117       108       107       106       115       124       125
       126         0       123       114       105       104       103       112
       121       122       123         0       120       111       102       101
       100       109       118       119       120                              
    -1

Similarly, a geometry can be read from a universal file. This can be done in multiple ways. If the universal file only contains geometry information, it can be easier to simply call the Geometry.load class method on the file to construct the Geometry object. The Geometry.load method will detect the .unv or .uff file extension and automatically pass the file to the SDynPy Universal File Format reader readunv. This will read the entire universal file and only keep datasets associated with the geometry. Note that the class method Geometry.load is also aliased to the module-level geometry.load function, so either is acceptable.

A secondary approach will call the readunv function explicitly. This approach is useful where the universal file may contain both geometry information and other types of information. Because Geometry.load will read the entire file, calling a separate function to read the remainder of the data will result in the entire file being read twice. Instead, the entire function can be read one time and the datasets parsed from the file can be passed to various functions to create SDynPy objects. A dictionary of dataset numbers and their contents is the output of readunv.

This dictionary can be passed to the Geometry.from_unv class method to construct the geometry. Again, this class method is aliased to a module function geometry.from_unv for convenience.

Reading and Writing to NumPy Files¶

SDynPy does not have a native storage format for its data types. However, being built mostly on NumPy arrays, it is almost trivial to use NumPy format for storage. A Geometry object can be saved using its Geometry.save method. It will be written to a NumPy .npz file containing the nodes, elements, tracelines, and coordinate systems as data members.

To load a Geometry object from this file, one can simply call the class method Geometry.load or its alias geometry.load, which will recognize the file extension and pass it to the NumPy loader.

Computing and Comparing Global Node Positions¶

Geometry in SDynPy is stored in the local coordinate system. For example, in the geometry we constructed, the coordinate for the cylindrical portion of the geometry is stored as ($r$, $\theta$, $z$) triples. However, it can be useful to quickly be able to compute the positions of these nodes in the global coordinate system, for example to compute the distance between two nodes. The global positions of specific nodes or all nodes can be computed using the Geometry method.

An optional node_ids argument can be passed to compute the positions of just specific nodes.

It may also be necessary to select nodes by global positions, meaning to find the node closest to a given position in space. For this, we can use the Geometry.node_by_global_position

And to get the closest nodes positions:

And the distances between the points:

Coordinate System Transformations¶

The coordinate systems presented up until this point have all been aligned with the global cartesian coordinate system. SDynPy also offers the ability to create coordinate systems which are shifted or rotated with respect to the global coordinate system.

As stated previously, a CoordinateSystemArray has a field called matrix which is size 4 x 3, meaning each coordinate system defined in the Geometry contains a 4 x 3 transformation matrix. This defines the coordinate system’s rotation and translation with respect to the global coordinate system. The first three rows of the matrix field are defined as a 3 x 3 rotation matrix $\mathbf{R}$ and the last row is defined as a 1 x 3 translation vector $\mathbf{T}$.

In SDynPy, a global coordinate $\mathbf{x}$ can be transformed to a local coordinate $\mathbf{x}'$ via

Similarly, a local coordinate can be transformed to a global coordinate via

With this definition, the rows of the rotation matrix end up being the local $x$, $y$, and $z$ axes of the coordinate system represented in the global coordinate system.

Note that any cylindrical or spherical transformation be applied after this transformation matrix. For example, a cylindrical coordinate system with rotation matrix

will have its cylindrical axis pointing in the $-y$ direction.

For the deflections, which are relative, the position component $\mathbf{T}$ will cancel out. For example, to compute a global deflection $\Delta\mathbf{x}$ from a local deflection $\Delta\mathbf{x}'$, the equation is

Coordinate system transformations can be a relatively abstract concept, therefore SDynPy makes it easy to visualize coordinate systems to help users to understand if they have set theirs up correctly. The plot_coordinate method will plot a coordinate system triad at each node showing the displacement directions of that node as defined by its disp_cs coordinate system. This function will be described more in depth in Coordinates

In the plot, it is clear that the displacement directions of the cylindrical portion of the geometry follow the curvature of the cylinder, as that is their defined displacement coordinate system. The cartesian portion of the geometry follow the global coordinate system.

Reading and Writing to Exodus Files¶

Exodus is a file format used at Sandia National Laboratories and elsewhere for finite element analysis mesh definition and results storage. Because modal test and other structural dynamics datasets are often used to calibrate or validate models, it can be useful to bring the Exodus data and test data into SDynPy for comparison, or otherwise export the SDynPy data to an Exodus file for comparison using some other toolset.

Two types of exodus files exist in SDynPy. The first is the Exodus class, which is the more traditional Exodus interface. In this interface, data remains on disk until it is requested by calling a method of the Exodus class. The second is the ExodusInMemory class, which loads all data from disk and stores it in memory. In this class, data can be accessed via attributes.

In SDynPy, the Geometry is analogous to the node positions and element connectivity in the Exodus file. In the Exodus file, elements are stored in blocks which contain elements with the same element type, material, and other properties. Since SDynPy does not have the concept of element blocks, an element block is defined for each element type in the ElementArray. Similarly, a block is defined for each traceline in the TracelineArray.

To write a SDynPy geometry to an Exodus file, we will use the ExodusInMemory.from_sdynpy class method and pass it a Geometry object. This will create an ExodusInMemory object from the data in the Geometry object. Note that Exodus files generally do not store local coordinate system information, so all data is transformed to the global coordinate system prior to export. The ExodusInMemory.from_sdynpy class method can also accept results in the form of a ShapeArray or NDDataArray object via the displacement_data optional argument, which will be written to the Exodus file as nodal variables; however, in this example we will only supply the geometry data.

Once we have an ExodusInMemory object, we can write it to a file using ExodusInMemory.write_to_file. An optional clobber argument will overwrite an already existing file. If the file name exists and clobber=True is not specified, an error will occur.

One the data is stored to an exodus file, it can be read in any number of Exodus readers. For example, the open-source software Paraview can read exodus files.

To read an Exodus file, we first use the Exodus class to open the file, then we can call Exodus.read_into_memory to generate an ExodusInMemory object if we desire. However, with the main goal to bring the Exodus data into a SDynPy geometry, we will generally use the Geometry.from_exodus class method to produce a SDynPy Geometry from either the Exodus or ExodusInMemory objects.

It is important to note that SDynPy Geometry and Exodus files do not have the same feature set, and therefore we should not expect the Geometry we loaded from the Exodus file to be equivalent to the Geometry we originally saved to the Exodus file. For example, we see that the Geometry object that we loaded from the Exodus file only has a single coordinate system (the global coordinate system), and all nodes are defined in that coordinate system. Similarly, Exodus does not have the concept of tracelines, so the lines are stored as beam elements. Therefore, the Geometry loaded from the Exodus file has zero tracelines and several elements of type 21 (Linear Beam). Finally, the color information stored in the original Geometry is lost.

Geometry Reduction¶

Many times when working with geometries, it is useful to reduce the geometry to just a subset of the original geometry. For example, we may want to select only the cylindrical portion of the geometry we have created. However, such a reduction is more complicated than simply throwing away the nodes we do not want, because those nodes may be referenced by elements or tracelines. The safest way to remove a portion of a geometry is to use the Geometry.reduce method. This will automatically reduce to a Geometry object that contains only the specified nodes and only the elements and tracelines that contain those specified nodes.

Here we will select nodes. We recall that in our geometry, the cylinder nodes have $z$-coordinate > 0 and the cube nodes have $z$-coordinate < 0, so we can use that as a selection criterion. We are after the id field of each node that matches our criterion. Beware when working with NodeArray objects that the coordinate field is a local coordinate. This means that we can incorrectly select for nodes if we use position criteria that reference global positions while comparing to the local coordinate field. Unless you are sure that your geometry only has a single, global coordinate system, it is usually safer to use the global_node_coordinate method to ensure you are accessing the correct position.

Once we have the nodes we want, we can pass the list of identification numbers to the Geometry.reduce method. This will not modify the Geometry in-place, but will rather return a copy of the Geometry object.

We can plot these to ensure the operations have been performed successfully.

Combining Geometries¶

Similarly to reducing Geometry objects, there can be times where we want to combine geometry objects into a single object. SDynPy has two ways to do this. The first way is to simply add the two Geometry objects together. In this case, the NodeArray, CoordinateSystemArray, TracelineArray, and ElementArray objects from the respective Geometry objects will simply be concatenated together.

Note that SDynPy will check to ensure there aren’t conflicts between the Geometry objects. For example, if two nodes are labeled 100 but they have different data associated with them, SDynPy will throw an error.

Note that if two items share the same identification number but they have identical underlying data, then the concatenation will succeed. The identical data between the two geometries will not be duplicated; only one copy will remain. This is useful when, for example, concatenating two geometries which both have a global coordinate system with identification number 1.

The second way to combine two geometries is to use the Geometry.overlay_geometries method. This method is typically used when two geometries from different sources are to be compared. Instead of combining the identification numbers of the items in the Geometry objects, it instead computes an offset value that gets applied to the identification numbers to ensure there are no conflicts. These offsets can be returned by the method as well by passing optional arguments.

Comparing the combined_geometry from the concatenation and the overlaid_geometry from Geometry.overlay_geometries, we can see how the offset is applied. Recalling the original node numbers:

We see that when we concatenate, those numbers are maintained.

However, when we overlay the geometries, an offset is applied.

SDynPy computed the offset of 1000, and will add 1*1000 to the first geometry, 2*1000 to the second geometry, etc. In this way, it can be ensured that no conflicts occur.

In the overlaid geometry, we see that the cube nodes start with 1000 and the cylinder nodes start with 2000.

Geometry Mapping¶

Many times in Structural Dynamics applications, we will want to compare two datasets. Perhaps a test is being compared to its equivalent simulation, or two simulations are compared to evaluate modeling assumptions. Regardless, these two datasets will very often not have compatible geometry. A test may have a few hundred sensors, whereas a finite element model may have a few million degrees of freedom. Regardless, we often wish to compare degrees of freedom, and for that reason, SDynPy defines an id_map object. The id_map object is essentially a vectorized dictionary that accepts an identification number and returns a new identification number that the first identification number is mapped to.

To show an example, let’s randomly generate some large numbers that we can pretend are node identification numbers from a finite element analysis.

To set up this map, we then supply the nodes we are mapping from (e.g. our geometry’s nodes) and the nodes we are mapping to (e.g. these fake finite element node numbers).

To explicitly map the nodes, we can then call the Geometry.map_ids method, and pass it this id_map object.

We can see that the node identification numbers no longer are the 100s and 200s they were previously, but are now our 5-digit finite element node numbers. Similarly, you can see that in the connectivity arrays for the TracelineArray and ElementArray portions of the geometry, the node numbers have also changed:

The Geometry.map_ids method is most often used to change node identification numbers for comparing between different datasets. However, the Geometry.map_ids method also can take optional arguments traceline_id_map, element_id_map, and coordinate_system_id_map to map the traceline, element, or coordinate system identification numbers.

Cone Geometry Examples¶

This section will contain multiple examples of creating a cone geometry using different approaches to defining the position and displacement directions of the nodes. Conical geometries are interesting in that they will generally force the use of local coordinate systems as sensors mounted to the cone’s surface will not be aligned with any cylindrical or cartesian coordinate system.

Positions Defined by a Cylindrical Coordinate System¶

When creating a cone, perhaps the most straightfoward way to define the node positions is with a cylindrical coordinate system.

First we will define the cone geometry.

We now have constructed our positions of our nodes and the identification numbers that go along with them. Note that the node_ids is a 8 x 12 array, and the node_positions is a 8 x 12 x 3 array, because there are three coordinates for each of the 8 x 12 nodes.

We normally work with 1D arrays of data when dealing with geometry, so it may seem to be a bit strange to do, but we can actually compute a NodeArray directly from these arrays. We won’t worry about the color or coordinate system definitions yet.

The NodeArray that we have made is now an 8 x 12 NodeArray.

And again, because the dtype of the coordinate field has shape (3,), it means that the shape of the coordinate field will be 8 x 12 x 3, which is why we were able to perform the assignment in the first place.

When we type the variable cone_nodes into the terminal, it will again print a table of the nodes; however we see the index is now a 2D index due to our NodeArray being two-dimensional.

Note that we did not assign a color or any coordinate systems to the nodes when they were defined, so they are automatically assigned to a default value. We must now define a CoordinateSystemArray containing a coordinate system with identification number 1. Based on how the node positions were defined, this should be a cylindrical coordinate system.

We will finally draw elements on the model. Because of the algorithmic way we defined our nodes, it should be easy to figure out which nodes should connect to which other nodes to produce a mesh of elements. Indeed, our nodes are already arranged that way in the NodeArray which we have constructed. We simply need to connect each node to its neighbors in the array. There is one additional complexity, and that is that the cone “wraps around” in the circumferential direction. This is rectified by concatenating the first column of identification numbers after the last column of identification numbers, meaning the first column will appear twice in the array.

Now we can simply loop through the array dimensions and construct an element for each entry.

We can then create the ElementArray manually.

We will forgo any tracelines in this geometry as the elements adequately convey the cone geometry. We can therefore just create the Geometry object.

We have successfully constructed a cone Geometry. However, if our Geometry represents a test in the sense that the sensors are mounted directly to the cone’s surface, we still have more work to do. Currently, the node displacement coordinate systems are defined as the cylindrical coordinate system.

If the sensors are truly mounted to the conical surface, they would have a direction normal to the cone’s surface, as well as the other two directions tangential to the cone’s surface.

Positions Defined by a Cartesian Coordinate System¶

A second approach to defining node positions is to use a cartesian coordinate system to define the node positions. This will follow a very similar path to the previous example. While initially defining the node positions in cartesian coordinates may require slightly more up-front work than simply using a cylindrical coordinate system, there could be an advantage in having the node locations already defined in a global coordinate system.

We can then again define a Geometry object.

Now when we construct the coordinate system, we specify a cartesian coordinate system cs_type=0 rather than a cylindrical coordinate system cs_type=1

Accepting that the nodes are in the same places, with the same node identification numbers, we can use the same ElementArray that we used previously to finish out the Geometry object.

We can verify that the geometries are equivalent by overlaying them.

But again, if we look at the displacement coordinate systems, they are not aligned with the cone. In this case, they are aligned with the global coordinate system.

This Geometry object may be similar to that extracted from a finite element analysis of a cone; finite element nodes are typically defined in a global coordinate system, and node displacements are also often represented in a global coordinate system. However, it is unlikely that a test would be set up with this Geometry, as it would require the machining of blocks to mount each sensor on so the sensor could be oriented in the global direction. Compared to a little extra math, the machining of blocks seems to be a great deal more work, so we will demonstrate how we can define local coordinate systems in our Geometry objects.

Defining Local Coordinate Systems with Rotation Matrices¶

We have seen in the previous two attempts to create Geometry objects that because the cone’s surface is not aligned with a global cartesian or cylindrical coordinate system, we must create local coordinate systems to enable our sensors to be mounted directly to the surface. This will involve the computation of the rotation matrix $\mathbf{R}$ that will define the local coordinate system that the sensors will use to displace. Recall that for displacements, the position of the origin of the coordinate system $\mathbf{T}$ will cancel out, so we can ignore that for now.

When working with rotations, users are directed to the rotation module in SDynPy. This contains functions for computing rotation matrices from axis/angle, quaternions, or Rodrigues parameters. We will use the R function to generate a rotation matrix given an angle and an axis. We can multiply rotation matrices together to allow for compound rotations.

In this case, when we go to compute the node positions, we will also compute the rotation matrix for the coordinate system. We will compute the rotation matrix, which is 3 x 3, and then we will append a row of zeros to it to make it the complete 4 x 3 transformation matrix.

If we recall, we have an 8 x 12 set of nodes. Thinking about what the shape of our transformation array should be, we might expect them to have the same shape of 8 x 12 nodes, each having a 4 x 3 transformation matrix. Checking the shape of the node_rotation_matrices variable shows us that is indeed what we have.

Don’t panic when you see arrays with higher dimensionality! NumPy can be incredibly efficient when working with multidimensional arrays rather than a for loop, so take the time to understand what is going on here. The way we have set up the node_rotation_matrix variable is precisely the correct size for the eventual creation of a 8 x 12 CoordinateSystemArray with its matrix field having a 4 x 3 transformation matrix.

Let’s go ahead and do that now. To simplify bookkeeping, we will assign the coordinate system identification numbers to be the same numbers as the node identification numbers. This way, the coordinate system associated with node 101 will also be called node 101, meaning the disp_cs field of the NodeArray object can be the same array as its id field.

Currently our CoordinateSystemArray contains all of the local coordinate systems, but not our global coordinate system, which we previously label 1. We have referenced this coordinate system as the definition coordinate system in the NodeArray object, so it must be defined as a coordinate system in our geometry. We should then concatenate our previous global coordinate system (defined in the variable cone_cs) to our local coordinate systems. Note that our current coordinate system array has a shape of 8 x 12, so we cannot simply add a single coordinate system to this array in its current shape. We will therefore first flatten the coordinate system arrays and then concatenate them together.

Now if we type the CoordinateSystemArray variable name in the terminal, we will see that we now have our global coordinate system 1 defined along with all of the local coordinate systems.

Let’s go ahead and create the geometry now.

Now if we plot the local coordinate systems, we will see that the coordinate system triads are aligned with the surface of the cone.

Summary¶

In summary, SDynPy’s Geometry object and constituent objects NodeArray, CoordinateSystemArray, TracelineArray, and ElementArray allow us to represent the spatial information relating to a test or analysis. NodeArray objects contain the nodes, which are defined points in space. The coordinate field of the NodeArray defines where each node is positioned. Each node has two coordinate systems associated with it. The first is the definition coordinate system, which determines the coordinate system in which the node’s coordinate field is interpreted, and the displacement coordinate system, which determines the coordinate system in which the node displaces. These are stored as references to identification numbers in the CoordinateSystemArray; the displacement coordinate system is stored in the disp_cs field, and the definition coordinate system is stored in the def_cs field.

The CoordinateSystemArray objects contain coordinate system definitions. These can be cartesian, cylindrical, or spherical coordinate systems, and they can have a local transformation applied.

Because a simple point cloud is missing 3D cues like shadows and occlusion, it can be difficult to intepret. The ElementArray and TracelineArray objects, therefore, exist to define visualization aids that connect nodes in a logical way to aid the user in interpreting the Geometry. Both of these objects contain a connectivity field which define how the nodes are connected.

The real advantage of using SDynPy Geometry objects is the ease in visualizing the geometry. A simple Geometry.plot method call will bring up an interactive 3D representation of the Geometry. Nodes defined in local coordinate systems are automatically converted to and plotted in a single global coordinate system for intuitive viewing. A call to Geometry.plot_coordinate will plot displacement coordinate systems as triads at each node. Each part of the Geometry object has a color field that defines how it is rendered, and the plotting functions are full-featured in that the user can customize the camera view, transparency, shading, and many other aspects of the visualization.

SDynPy Geometry objects can be easily read from or written to the universal file format, Exodus files, or NumPy files, making them easy to share with customers or collaborators.