IDWarp
IDWarp is a tool for performing mesh manipulation of 3-dimensional
multi-block, overset and unstructured meshes.
Building
IDWarp depends on the follow libraries: - CGNS Library - PETSc - MPI
See the MDO Lab installation guide here for the supported versions and installation instructions.
All the core computations in IDWarp are coded in Fortran. It is therefore necessary to build this library before using IDWarp.
To see a list of architectures that IDWarp has been known to compile on run:
make
from the root directory.
Follow the instructions and copy the closest architecture file and attempt a build using
make
If everything was successful, the following lines will be printed to the screen (near the end):
Testing if module idwarp can be imported...
Module idwarp was successfully imported.
If you don’t see this, it will be necessary modify the configure options in the config file.
It will most likely be necessary to modify the CGNS_INCLUDE_FLAGS
and the CGNS_LINKER_FLAGS variables. After changes to the
configuration file, run make clean before attempting a new build.
Lastly, to build the Python interface, go to the root directory and type:
pip install .
Some features require additional dependencies.
Using IDWarp with OpenFOAM meshes in DAFoam requires pyOFM.
This dependency can be checked with
pip install .[dafoam]
Using MultiUSMesh requires cgnsUtilities, which can be checked with
pip install .[multi]
Tutorial
This is a short tutorial explaining how to get started with IDWarp. We will start by presenting a simple stand-alone warping script and explaining the various parts.
from mpi4py import MPI
from idwarp import USMesh
options = {
'gridFile':'../../input_files/o_mesh.cgns',
'fileType':'cgns',
'specifiedSurfaces':None,
'symmetrySurfaces':None,
'symmetryPlanes':[],
'aExp': 3.0,
'bExp': 5.0,
'LdefFact':1.0,
'alpha':0.25,
'errTol':0.0001,
'evalMode':'fast',
'useRotations':True,
'zeroCornerRotations':True,
'cornerAngle':30.0,
'bucketSize':8,
}
# Create the mesh object
mesh = USMesh(options=options, comm=MPI.COMM_WORLD)
# Extract all coordinates
coords0 = mesh.getSurfaceCoordinates()
# Modify the coordinates as required
newCoords = coords0.copy()
for i in range(len(coords0)):
newCoords[i, :] *= 1.1
# Reset the newly computed surface coordinates
mesh.setSurfaceCoordinates(newCoords)
# Actually run the mesh warping
mesh.warpMesh()
# Write the new grid file.
mesh.writeGrid('warped.cgns')
The first two lines of code imports mpi4py and the mesh warping module
from mpi4py import MPI
from idwarp import USMesh
The next chunk lists the options for warping. The options are explained in Options.
options = {
'gridFile':'../../input_files/o_mesh.cgns',
'fileType':'cgns',
'specifiedSurfaces':None,
'symmetrySurfaces':None,
'symmetryPlanes':[],
'aExp': 3.0,
'bExp': 5.0,
'LdefFact':1.0,
'alpha':0.25,
'errTol':0.0001,
'evalMode':'fast',
'useRotations':True,
'zeroCornerRotations':True,
'cornerAngle':30.0,
'bucketSize':8,
}
Next we create the actual mesh object itself
# Create the mesh object
mesh = USMesh(options=options, comm=MPI.COMM_WORLD)
Note that we have explicitly passed in the MPI intracommunicator on which we want to create the object. If the ‘comm’ keyword argument is not given, it will default to MPI.COMM_WORLD. Therefor this example, mpi4py is not strictly required to be imported in the run script.
Next we request the surface coordinates from the mesh object. These will correspond to coordinates in the ‘specifiedSurfaces’ option.
# Extract all coordinates
coords0 = mesh.getSurfaceCoordinates()
coords0 is a numpy array of size (N,3). It is now up to the user to manipulate these coordinates however they wish for this example we simply loop over all coordinates and uniformly scale by a factor of 1.1:
new_coords = coords0.copy()
for i in range(len(coords0)):
new_coords[i, :] *= 1.1
Once the new set of coordinates have been determined, return them to the mesh warping object with the following command.
# Reset the newly computed surface coordiantes
mesh.setSurfaceCoordinates(new_coords)
Note that the shape of ‘new_coords’ must be identical to the coords0 array that was originally provided by the warping. Next we run the actual mesh warp using
# Actually run the mesh warping
mesh.warpMesh()
And finally to produce an updated grid file we can write the grid:
# Write the new grid file.
mesh.writeGrid('warped.cgns')
The warped grid file ‘warped.cgns’ will contain all the boundary condition/connectivity/auxiliary information as the original cgns file. Only the coordinates are updated to their new positions.
Options
Name |
Type |
Default value |
Description |
|---|---|---|---|
|
str or NoneType |
None |
This is the grid file to use. It must always be specified. It may be a structured or unstructured grid file or it may be an OpenFOAM directory containing a mesh specification. |
|
str |
|
Specify the type of grid file.
|
|
list or NoneType |
None |
This option is used to specify which surfaces are used to build the surface tree where deformations are to be specified. The default is None which will automatically use all wall-type surfaces in the grid file. For CGNS files this corresponds to the following boundary condtiions: |
|
list or NoneType |
None |
This option is used to specify which surfaces are used to determine symmetry planes. If None, IDwarp will automatically use the |
|
list or NoneType |
None |
This is sort of a “last-resort” option. It is used to overwrite and explicitly define symmetry conditions IDWarp is to use. It is the only method for specifying symmetry for unstructured CGNS files and OpenFOAM files. For the |
|
float |
3.0 |
This is the first exponent in the inverse distance weighting function. The default is usually the most robust, but changing this value slightly might help with some edge cases. |
|
float |
5.0 |
This is the second exponent in the inverse distance weighting function. The default is usually the most robust, but changing this value slightly might help with some edge cases. |
|
float |
1.0 |
This parameter is used to determine how “far” into the field a surface deformation propagates before it is attenuated. For small shape modifications, such as small changes to the shape of a airfoil, the default value of 1.0 is likely to be sufficient. However, for much larger changes in the surface such as wing planform changes, much larger values tend to be more robust. For these cases, values in the range 50-100 are common. |
|
float |
0.25 |
This value determines how the two different exponent terms are blended. It determines how much of the higher exponent |
|
float |
0.0005 |
This is the relative tolerance used for the fast sum approximation. A larger tolerance is faster, but may result in small mesh imperfections away from the surface. If mesh edge lengths grow uniformly away from the body, small “errors” in the node position are not an issue. However, if the mesh has small edge lengths a great distance from the body, these imperfections may cause issues and it may be required to lower the tolerance by an order of magnitude or two at the cost of more computational time. |
|
str |
|
How to compute the sums.
|
|
float |
1e-06 |
Possibly deprecated option. IDWarp currently uses the BC info in the CGNS meshes to figure out which nodes are on the symmetry plane. This option seems to be a leftover from older implementations. It was used to specify the distance from the symmetry plane for a node to tag it as a symmetry plane node. |
|
bool |
True |
Flag specifying if rotations are to be interpolated in addition to displacements. For small mesh changes it may not be necessary to interpolate rotations. However, if the surface is undergoing large changes in orientation, using rotations will help preserve boundary layer orthogonality which is generally desirable. |
|
bool |
True |
Flag specifying if rotations at corners (defined by |
|
float |
30.0 |
The minimum deviation in degrees between surface normals surrounding a node for it to be considered a corner point. |
|
str or NoneType |
None |
A restart file can be used to load the node coefficients to IDWarp to speed up initialization. There are 3 execution paths: If no restart file is provided, the code performs a full initialization. If the restart file is provided but the file does not exist or does not match with some internal checks, the code performs a full initialization and writes the load balancing and the coefficients to this restart file for future use. Finally, if the restart file is provided and it is good, the code loads the initialized coefficients and does not need to do a full initialization. |
|
int |
8 |
The size of the “buckets” at the last level of the KD-tree. A large bucket size reduces the number of levels in the tree and the overall tree size but may require more computation since a less fine granularity of leaves are available. Experiments have indicated there is little difference in run time for bucket sizes 1, 2, 4 and 8. |
API
USMesh
- class idwarp.UnstructuredMesh.USMesh(*args, **kwargs)[source]
This is the main Unstructured Mesh. This mesh object is designed to interact with an structured or unstructured CFD solver though a variety of interface functions.
Create the USMesh object.
- Parameters:
- optionsdictionary
A dictionary containing the options for the mesh movement strategies. THIS IS NOT OPTIONAL.
- commMPI_INTRA_COMM
MPI communication (as obtained from mpi4py) on which to create the USMesh object. If not provided, MPI_COMM_WORLD is used.
- debugbool
Flag specifying if the MExt import is automatically deleted. This needs to be true ONLY when a symbolic debugger is used.
- getCommonGrid()[source]
Return the grid in the original ordering. This is required for the OpenFOAM tecplot writer since the connectivity is only known in this ordering.
- getSolverGrid()[source]
Return the current grid in the order specified by setExternalMeshIndices(). This is the main routine for returning the defomed mesh to the external CFD solver.
- Returns:
- solverGrid, numpy array, real: The resulting grid.
The output is returned in flatted 1D coordinate format. The len of the array is 3*len(indices) as set by setExternalMeshIndices()
- getSurfaceCoordinates()[source]
Returns all defined surface coordinates on this processor
- Returns:
- coordsnumpy array size (N,3)
Specified surface coordinates residing on this processor. This may be empty array, size (0,3)
- getWarpGrid()[source]
Return the current grid. This funtion is typically unused. See getSolverGrid for the more useful interface functionality.
- Returns:
- warpGrid, numpy array, real: The resulting grid.
The output is returned in flatted 1D coordinate format.
- getdXs()[source]
Return the current values in dXs. This is the result from a mesh-warp derivative computation.
- Returns:
- dXsnumpy array
The specific components of dXs. size(N,3). This the same size as the array obtained with getSurfaceCoordiantes(). N may be zero if this processor does not have any surface coordinates.
- setExternalMeshIndices(ind)[source]
Set the indicies defining the transformation of an external solver grid to the original CGNS grid. This is required to use USMesh functions that involve the word “Solver” and warpDeriv. The indices must be zero-based.
- Parameters:
- indnumpy integer array
The list of indicies this processor needs from the common mesh file
- setSurfaceCoordinates(coordinates)[source]
Sets all surface coordinates on this processor
- Parameters:
- coordinatesnumpy array, size(N, 3)
The coordinate to set. This MUST be exactly the same size as the array obtained from getSurfaceCoordinates()
- setSurfaceDefinition(pts, conn, faceSizes, cgnsBlockID=None)[source]
This is the master function that determines the definition of the surface to be used for the mesh movement. This surface may be supplied from an external solver (such as SUmb) or it may be generated by IDWarp internally.
- Parameters:
- ptsarray, size (M, 3)
Nodes on this processor
- connint array, size (sum(faceSizes))
Connectivity of the nodes on this processor
- faceSizesint Array size (N)
Treat the conn array as a flat list with the faceSizes giving connectivity offset for each element.
- cgnsBlockIDdummy argument.
This argument is not used at all. It is here just to have the same API as the IDWarpMulti class.
- setSurfaceDefinitionFromFile(surfFile)[source]
Set the surface definition for the warping from a multiblock PLOT3D surface file
- Parameters:
- surfFilefilename of multiblock PLOT3D surface file.
- setSurfaceFromFile(surfFile)[source]
Update the internal surface surface coordinates using an external PLOT3D surface file. This can be used in an analogous way to setSurfaceDefinitionFromFile. The ‘sense’ of the file must be same as the file used with setSurfaceDefinitionFromFile. That means, the same number of blocks, with the same sizes, in the same order.
- Parameters:
- surfFile: filename of multiblock PLOT3D surface file’
- verifyWarpDeriv(dXv=None, solverVec=True, dofStart=0, dofEnd=10, h=1e-06, randomSeed=314)[source]
Run an internal verification of the solid warping derivatives
- warpDeriv(dXv, solverVec=True)[source]
Compute the warping derivative (dXv/dXs^T)*Vec (where vec is the dXv argument to this function.
This is the main routine to compute the mesh warping derivative.
- Parameters:
- dXvnumpy array
Vector of size external solver_grid. This is typically obtained from the external solver’s dRdx^T * psi calculation.
- solverVeclogical
Flag to indicate that the dXv vector is in the solver ordering and must be converted to the warp ordering first. This is the usual approach and thus defaults to True.
- Returns:
- None. The resulting calculation is available from the getdXs()
- function.
- warpDerivFwd(dXs, solverVec=True)[source]
Compute the forward mode warping derivative
This routine is not used for “regular” optimization; it is used for matrix-free type optimization. dXs is assumed to be the the peturbation on all the surface nodes.
- Parameters:
- dXsarray, size Nsx3
This is the forward mode peturbation seed. Same size as the surface mesh from getSurfaceCoordinates().
- solverVecbool
Whether or not to convert to the solver ordering.
- Returns:
- dXvarray
The peturbation on the volume meshes. It may be in warp ordering or solver ordering depending on the solverVec flag.
- warpMesh()[source]
Run the applicable mesh warping strategy.
This will update the volume coordinates to match surface coordinates set with setSurfaceCoordinates()
- writeGrid(fileName=None)[source]
Write the current grid to the correct format
- Parameters:
- fileNamestr or None
Filename for grid. Should end in .cgns for CGNS files. For PLOT3D whatever you want. It is not optional for CGNS/PLOT3D. It is not required for OpenFOAM meshes. This call will update the ‘points’ file.
MultiUSMesh
- class idwarp.MultiUnstructuredMesh.MultiUSMesh(CGNSFile, optionsDict, comm=None, dtype='d', debug=False)[source]
This mesh object is designed to support independent deformation of multiple overset component meshes.
Create the MultiUSMesh object.
INPUTS:
CGNSFile: string -> file name of the CGNS file. This CGNS file should be generated with cgns_utils combine, so that the domain names have the appropriate convention. That is, domains will have the same name as their original files. Domains that share the same name will be grouped to make an IDWarp instance.
optionsDict: dictionary of dictionaries -> Dictionary containing dictionaries that will be used to initialize multiple IDWarp instances. The keys are domain names and the values are dictionaries of standard IDWarp options that will be applied to this domain. The domains of the full CGNS file that do not have a corresponding entry in optionsDict will not be warped. For instance, if the CGNS file has the domains wing.00000, wing.00001, and wing.00002 associated with a wing mesh that we want to warp, then optionsDict should have an entry for ‘wing’.
Ney Secco 2017-02
- getSolverGrid()[source]
Return the current grid in the order specified by setExternalMeshIndices(). This is the main routine for returning the deformed mesh to the external CFD solver.
- Returns:
- solverGrid, numpy array, real: The resulting grid.
The output is returned in flatted 1D coordinate format. The len of the array is 3*len(indices) as set by setExternalMeshIndices()
- Ney Secco 2017-02
- getSurfaceCoordinates()[source]
Returns all defined surface coordinates on this processor, with the Solver ordering
- Returns:
- ptsnumpy array size (N,3)
Specified surface coordinates residing on this processor. This may be empty array, size (0,3)
- Ney Secco 2017-02
- getWarpGrid()[source]
Return the current grids. This function is typically unused. See getSolverGrid for the more useful interface functionality.
This only returns the nearfield meshes.
- Returns:
- volNodesList, list of 1D numpy arrays, real: These are the local volume nodes (in a flat 1D array)
- of each instance. That is, volNodesList[i] has the volume nodes stored in the local proc for
- the i-th IDWarp instance.
- numCoorTotal: The total number of coordinates, across all procs and IDWarp instances.
- Ney Secco 2017-02
- getdXs()[source]
Return the current values in dXs. This is the result from a mesh-warp derivative computation. Note that the same steps used in this function are done at warpDeriv, so in theory you could save time by just saving the output of warpDeriv.
Here we accumulate all seeds coming from each IDWarp instance
- Returns:
- dXsnumpy array
The specific components of dXs. size(N,3). This the same size as the array obtained with getSurfaceCoordiantes(). N may be zero if this processor does not have any surface coordinates. These can be, for instance, reverse AD seeds.
- Ney Secco 2017-03
- setExternalMeshIndices(ind)[source]
Set the indices defining the transformation of an external solver grid to the original CGNS grid. This is required to use USMesh functions that involve the word “Solver” and warpDeriv. The indices must be zero-based.
- Parameters:
- indnumpy integer array[3*n] where n is the number of nodes stored in the proc.
The list of indicies this processor needs from the common mesh file. ind[3*i:3*i+3] represents the position of each coordinate of the i-th ADflow node in the CGNS file. ind has 0-based indices.
- Ney Secco 2017-02
- setSurfaceCoordinates(pts)[source]
Sets all surface coordinates on this processor, with pts given in solver ordering
- Parameters:
- ptsnumpy array, size(N, 3)
The coordinate to set. This MUST be exactly the same size as the array obtained from getSurfaceCoordinates()
- Ney Secco 2017-02
- setSurfaceDefinition(pts, conn=None, faceSizes=None, cgnsBlockIDs=None, distTol=1e-06)[source]
This is the master function that determines the definition of the surface to be used for the mesh movement. This surface may be supplied from an external solver (such as ADflow) or it may be generated by IDWarp internally.
- Parameters:
- ptsarray, size (M, 3)
Nodes on this processor
- connint array, size (sum(faceSizes))
Connectivity of the nodes on this processor
- faceSizesint Array size (N)
Treat the conn array as a flat list with the faceSizes giving connectivity offset for each element.
- cgnsBlockIDsint Array size (N)
Block ID, in CGNS ordering, that contains each element.
- distTol: Distance tolerance to flag that a given surface node does not
belong to the current IDWarp surface definition in the current proc.
- Ney Secco 2017-02
- warpDeriv(dXv, solverVec=True)[source]
Compute the warping derivative (dXv/dXs^T)*Vec (where vec is the dXv argument to this function.
This is the main routine to compute the mesh warping derivative.
- Parameters:
- dXvnumpy array
Vector of size external solver_grid. This is typically obtained from the external solver’s dRdx^T * psi calculation.
- solverVeclogical
Flag to indicate that the dXv vector is in the solver ordering and must be converted to the warp ordering first. This is the usual approach and thus defaults to True.
- Returns:
- dXsnumpy array
The specific components of dXs. size(N,3). This the same size as the array obtained with getSurfaceCoordiantes(). N may be zero if this processor does not have any surface coordinates. These can be, for instance, reverse AD seeds.
- Ney Secco 2017-03
- warpDerivFwd(dXs)[source]
Compute the forward mode warping derivative
This routine is not used for “regular” optimization; it is used for matrix-free type optimization. dXs is assumed to be the the peturbation on all the surface nodes.
- Parameters:
- dXsarray, size Nsx3
This is the forward mode peturbation seed. Same size as the surface mesh from getSurfaceCoordinates().
- solverVecbool
Whether or not to convert to the solver ordering.
- Returns:
- dXvarray
The peturbation on the volume meshes. It may be in warp ordering or solver ordering depending on the solverVec flag.
- Ney Secco 2017-03
Limitations
Cannot be used where you have 2 symmetry planes e.g. in a 2D case. Use pywarp instead.