(Looking for the newer parallel C++ version? It's here.)
Overset meshes are a method for representing complex geometry in computational fluid dynamics and other types of PDE-based simulations. An overset mesh is a composite mesh made up of a set of overlapping component grids (structured or unstructured), in which the individual grids exchange data with each other via interpolation.
Overset meshes have several nice properties, including:
- Simplicity: Component grids are (often) structured, so operating on them is similar to operating on Cartesian grids
- Flexibility: Overset meshes' modular nature make them convenient for moving geometry problems and adaptive mesh refinement, as well as simulations in which different physical phenomena or numerical methods are relevant in different regions
- Efficiency: Structured grids enable the use of efficient and accurate numerical methods.
However, some processing is required to convert a set of overlapping grids into a working overset mesh. Overkit is designed to automate the tasks involved in this processing, such as:
- Hole cutting of boundary geometry and overlapped coarse grids
- Determination of donors and receivers for inter-grid communication
- Generation of interpolation weights.
Overkit is also capable of generating interpolation data for mesh remapping (more on this below).
Note: Overkit is currently designed for structured grids with node-centered data. Support for cell-centered data and unstructured grids may be added in future versions.
Overkit requires Fortran/C compilers and CMake. It currently has no other third-party library/tool dependencies.
Overkit uses CMake, so the build process is similar to that of other CMake-based build systems. An example build procedure is:
cd ~/overkit-src
mkdir build && cd build
cmake .. -DCMAKE_INSTALL_PREFIX=/usr/local/overkit
make
make installSee CMake's documentation for more details.
Use the flag -DEXAMPLES=ON when configuring to build some example programs that use
Overkit. Source code for the examples can be found in examples/; the built programs are placed
in <cmake-build-dir>/examples/. Run <example-program-name> --help for information on how to
use each one.
Overkit has optional OpenMP support that can be enabled by configuring with -DOPENMP=ON.
What follows is a brief(ish) Overkit user's guide.
Note: It may be helpful to look at some of the examples' sources while reading this.
In an overset mesh, each component grid has fringes along its non-boundary edges. A fringe is made up of one or more layers of receiver points, which are points that receive data from other grids. Each receiver has a corresponding donor on the sending grid which consists of the set of points that makes up the interpolation stencil and each point's associated interpolation weight. Overkit calls the combination of a donor and a receiver a connection. Most of the time, a pair of grids in an overset mesh will communicate in both directions, so an interface between them will have two adjacent fringes (one on each grid).
In some cases when assembling an overset mesh, suitable donors may not be found for certain receiver points (e.g. if there is insufficient overlap between grids); these receiver points are referred to as orphan points.
Overset grids often have regions in which solutions are not computed called holes. These result from several processes during assembly (details in the assembly section below).
In a typical explicit solver the solution procedure for an overset mesh may be something like:
Set u to some initial value for all grid points
For istep in 1..nsteps
Solve for u at non-hole, non-fringe points
Interpolate to update u at fringe points
Note: When designing overset meshes for solvers that follow the above procedure, it is important for the grids to overlap enough to allow donors to contain only non-fringe points (otherwise data will be used in the interpolation stencils that has not yet been updated to the current step).
Overkit has three main data structures: domains, grids, and connectivities.
The domain is the primary Overkit data structure. It contains a set of grids and their associated connectivities.
The grid represents a single component grid. It holds the grid's metadata (size, periodicity, etc.), its coordinates, and a field representing the state of each grid point (interior point, boundary point, hole, etc.).
The connectivity holds inter-grid connectivity data for a given pair of grids. For each connection (donor/receiver pair), the following data is stored:
- Donor extents: the range of points on the donor grid from which data will be interpolated, stored as a lower and upper point
- Donor interpolation coefficients: the interpolation weights for each point in the interpolation stencil, stored in tensor product format (i.e., an array of one-dimensional stencils, with the first index indicating the point in the stencil and the second index indicating the dimension)
- Donor coordinates: the mapped coordinates of the receiver point inside the donor stencil (the coordinate range depends on the interpolation scheme; e.g., for linear interpolation the coordinates range from 0 to 1, with (0,0) representing the lower left vertex and (1,1) representing the upper right vertex). Note: These are typically only needed when generating interpolation coefficients separately from Overkit
- Receiver point: the point on the receiver grid to which the interpolated data will be sent.
Cart (ovk_cart): Type representing a grid's structure, namely its size and periodicity.
Field (ovk_field_<int,large_int,real,logical>): Type for grid data. Data is stored in the
values member.
Array (ovk_array_<int,large_int,real,logical>): Type for linear arrays. Data is stored in the
values member.
Domains are created with the following command:
call ovkCreateDomain(Domain, <numdims>, <numgrids>, StatusLogFile=<statusunit>, &
ErrorLogFile=<errorunit>)<numdims> is the dimension of the problem (2 or 3).
<numgrids> is the number of grids in the domain.
<statusunit> and <errorunit> are files to which to send the status and error output. They can
be set to OUTPUT_UNIT and ERROR_UNIT from iso_fortran_env to output to the command line, or
the arguments can be omitted to suppress output.
Note: There is also a command to clean up once finished:
call ovkDestroyDomain(Domain)Grids are created with the following command:
call ovkCreateGrid(Domain, <gridid>, <numpoints>, Periodic=<periodic>, &
PeriodicStorage=<periodicstorage>, PeriodicLength=<periodiclength>, &
GeometryType=<geometrytype>)<gridid> is the index of the grid in the domain.
<numpoints> is the number of points in each direction.
<periodic> specifies whether the grid is periodic or not (.true. or .false.) in each direction.
(optional; default=.false.).
<periodicstorage> can either be OVK_PERIODIC_STORAGE_UNIQUE or
OVK_PERIODIC_STORAGE_DUPLICATED. The former indicates that the grid contains only one instance
of each point along any of its periodic dimensions; the latter indicates that the first point is
duplicated at the end in any periodic dimension. (optional; default=OVK_PERIODIC_STORAGE_UNIQUE).
<periodiclength> specifies the offset that must be applied to the grid's coordinates
when crossing over from one period to the next. The length can be set to 0 for O-shaped grids.
Note: Overkit currently assumes that only the corresponding dimension's coordinate changes when
crossing periods, e.g., if crossing from point (i,N,k) to (i,1,k), the X and Z coordinates will
remain the same, and the Y coordinate will be offset by PeriodicLength(2). (optional; default=0.0).
<geometrytype> can be used to provide hints to Overkit about whether the grid's coordinates have
any special properties that can be used to improve performance. Possible values are:
OVK_GEOMETRY_UNIFORM- an axis-aligned uniformly-spaced gridOVK_GEOMETRY_ORIENTED_UNIFORM- a rotated uniform gridOVK_GEOMETRY_RECTILINEAR- a non-uniform grid in which all of the cells are rectangular and axis-alignedOVK_GEOMETRY_ORIENTED_RECTILINEAR- a rotated rectilinear gridOVK_GEOMETRY_CURVILINEAR- a general non-uniform structured grid.
(optional; default=OVK_GEOMETRY_CURVILINEAR).
Note: There is also a command to destroy a grid:
call ovkDestroyGrid(Domain, <gridid>)(but most of the time this isn't necessary, as the domain will clean up all of its contained grids
and connectivities when ovkDestroyDomain is called).
Grids' coordinates are set using the following sequence of commands:
! Request to edit the grid (returns pointer to ovk_grid)
call ovkEditGrid(Domain, <gridid>, Grid)
! Request to edit the grid's coordinates (returns pointers to ovk_field_real)
! 2D is similar but with Z omitted
call ovkEditGridCoords(Grid, 1, X)
call ovkEditGridCoords(Grid, 2, Y)
call ovkEditGridCoords(Grid, 3, Z)
! Set the coordinates
! 2D is similar but with Z omitted and k equal to 1
do k = 1, Nz
do j = 1, Ny
do i = 1, Nx
X%values(i,j,k) = <xvalue>
Y%values(i,j,k) = <yvalue>
Z%values(i,j,k) = <zvalue>
end do
end do
end do
! Tell the grid that you're done editing the coordinates
call ovkReleaseGridCoords(Grid, X)
call ovkReleaseGridCoords(Grid, Y)
call ovkReleaseGridCoords(Grid, Z)
! Tell the domain that you're done editing the grid
call ovkReleaseGrid(Domain, Grid)Grid coordinates are initially set to uniformly vary from 0 to N-1 in each direction (which is sometimes convenient as a starting point for generating grids).
Grids' states are set with the following sequence of commands:
! Request to edit the grid (returns pointer to ovk_grid)
call ovkEditGrid(Domain, <gridid>, Grid)
! Request to edit the grid's state (returns pointer to ovk_field_int)
call ovkEditGridState(Grid, State)
! Set the state
! 2D is similar but with k equal to 1
do k = 1, Nz
do j = 1, Ny
do i = 1, Nx
State%values(i,j,k) = <statevalue>
end do
end do
end do
! Tell the grid that you're done editing the state
call ovkReleaseGridState(Grid, State)
! Tell the domain that you're done editing the grid
call ovkReleaseGrid(Domain, Grid)<statevalue> can be one of the following:
OVK_INTERIOR_POINT- (the default) a normal pointOVK_DOMAIN_BOUNDARY_POINT- a point that lies along a boundary of the simulation domainOVK_INTERNAL_BOUNDARY_POINT- a point along a grid edge that isn't a domain boundary, but where a fringe is undesirable (perhaps due to the presence of some other type of inter-grid interface)OVK_EXTERIOR_POINT- a point that is outside of the simulation domain.
Note 1: Coordinates and state can be edited at the same time (in fact it's slightly more efficient to do this).
Note 2: Overkit expects valid grid coordinates even for points set to OVK_EXTERIOR_POINT.
Note 3: The above state values are actually combinations of individual state flags of the form
OVK_STATE_* (more on these later).
Once the domain and grids are set up, assembly can be run with the following commands:
AssemblyOptions = ovk_assembly_options_(<numdims>, <numgrids>)
call ovkAssemble(Domain, AssemblyOptions)However, for anything to actually happen we will need to first modify the assembly options. The
sections below describe the assembly process in Overkit and explain what each of the options in the
ovk_assembly_options structure do. Commands for modifying assembly options will have one of
the following forms:
! Per-grid settings
call ovkSetAssemblyOption<optionname>(AssemblyOptions, <gridid>, <optionvalue>)
! Per-grid-pair settings
call ovkSetAssemblyOption<optionname>(AssemblyOptions, <gridid1>, <gridid2>, <optionvalue>)A value of OVK_ALL_GRIDS can be passed to any of the <gridid*> parameters to modify the
option value for multiple grids at once.
Overkit first needs to determine how the component grids are situated with respect to each other. This consists of looking at every grid point and determining within which cell(s) on the other grids the point resides (if any). Naively, overlap detection has quadratic computational complexity in the number of grid points (because it involves testing each point against each cell among all pairs of grids), but this can be improved to linear complexity through the use of spatial partitioning data structures.
Overkit currently uses a hybrid approach consisting of a kd-tree at coarse levels (to prune out empty space and to separate regions of vastly different grid resolution) and a spatial hash grid at finer levels.
Relevant options:
Overlappable(m,n)- whether gridmis considered to potentially overlap with gridn(default=.false.)OverlapTolerance(m,n)- tolerance for points on gridnthat lie slightly outside of gridmbut should be considered to overlap (most often occurs along curved boundaries). Specified as a fraction of the boundary cells' sizes (e.g. a value of 0.1 means that a point can be up to 10% outside of a nearby cell and still be considered to overlap) (default=0.0)OverlapAccelDepthAdjust(m)- adjusts the desired minimum number of cells in each kd-tree node on gridm. Each increment by 1.0 halves the number (increasing the tree depth), and each decrement by 1.0 doubles the number (decreasing the tree depth) (default=0.0)OverlapAccelResolutionAdjust(m)- adjusts the size of hash grid bins on gridm. Each increment by 1.0 halves the bin size (resulting in fewer cells overlapping with each bin, i.e., better performance at the cost of more memory) and each decrement by 1.0 doubles the bin size (default=0.0).
The collection of domain boundary points on all of the grids defines the overall shape of the domain. Depending on the input grid configuration, some points may lie outside of this region and will need to be removed. (For example, in a two-grid airfoil configuration with a uniform background grid and a second grid that wraps around the airfoil geometry, points on the background grid will lie inside the non-domain region defined by inner edge of the second grid.)
Overkit handles this by (1) projecting a representation of other grids' boundaries onto the grid being cut, (2) detecting which side is outside of the domain, then (3) flooding this side and removing the detected points from the domain.
Relevant options:
InferBoundaries(n)- assume grid edges in non-overlapping regions are domain boundaries on gridn(default=.false.)CutBoundaryHoles(m,n)- whether the boundary points on gridmare allowed to be used in cutting gridn(default=.false.).
Next we need to determine which points should be receiving data. Overkit divides the points detected in this step into two categories, outer fringe points and occluded points. Outer fringe points are points near non-boundary grid edges. Occluded points are overlapped points for which one of the overlapping grids has been determined to be a better representation of that part of the domain than the present one.
Relevant options:
FringeSize(n)- how many layers of points away from non-boundary edges of gridnshould be counted as fringe (default=0)Occludes(m,n)- how to treat points on gridnthat are overlapped by gridm. Possible values areOVK_OCCLUDES_ALL, meaning all such points are occluded,OVK_OCCLUDES_NONE, meaning none are occluded, andOVK_OCCLUDES_COARSE, meaning only points with a lower resolution are occluded (default=OVK_OCCLUDES_NONE)EdgePadding(m,n)- in some cases it may be desirable to prevent points from receiving data from points that are within some distance of the donor grid's edge. This option specifies how far from the edge cells should be on the overlapping gridmin order to occlude points on the overlapped gridn(default=0)EdgeSmoothing(n)- used in conjunction withEdgePadding, helps to remove any sharp edges or disconnectedness in the set of occluded points on gridn. Higher values increase smoothing strength (default=0).
For overset meshes we typically don't need to keep all of the receiver points detected in the previous step. Instead we can cut out most of the occluded points and leave only the few layers nearest to the unoccluded part of the grid. Overkit calls these remaining points the inner fringe.
The EdgePadding option described in the previous section reduces the size of a grid's
occluded region, which has the side effect of increasing the amount of overlap remaining after
overlap minimization. This can be useful in certain situations such as that of the explicit solver
described above.
Relevant options:
MinimizeOverlap(m,n)- whether points on gridnthat are occluded by gridmand aren't part of the inner fringe should be removed (default=.false.)FringeSize(n)- (defined above)EdgePadding(m,n)- (defined above).
By this point we have our final set of receiver points and the only thing left to do is to locate the corresponding donors. This involves choosing which overlapping grid will donate to each receiver and then actually generating the interpolation stencils on the donor grids. In determining the best donor grid for a given receiver, both the local grid resolution and edge distance are used. If a suitable donor cannot be found, the point will be marked as an orphan and a warning message will be produced.
Relevant options:
ConnectionType(m,n)- how donors on gridmwill communicate with receivers on gridn. Possible values areOVK_CONNECTION_NONE, meaning communication is not allowed, andOVK_CONNECTION_NEAREST,OVK_CONNECTION_LINEAR, andOVK_CONNECTION_CUBICfor nearest-neighbor, linear, and cubic interpolation respectively (default=OVK_CONNECTION_NONE)EdgePadding(m,n)- (defined above; used again here as a threshold for deciding whether to prioritize edge distance or resolution).
In addition to overset meshes, Overkit can also generate interpolation data for multi-block or overset mesh remapping. For this application one would create a single domain containing both the source and destination grids, and set assembly options as follows (omitting options left at the default value):
Overlappable(<srcgrids>,<destgrids>)=.true.Occludes(<srcgrids>,<destgrids>)=OVK_OCCLUDES_ALLConnectionType(<srcgrids>,<destgrids>)=<desired interp scheme>.
Once assembly is completed, interpolation data can be accessed using the following commands:
! Get the connectivity data structure from the domain (returns pointer to ovk_connectivity)
call ovkGetConnectivity(Domain, <gridid1>, <gridid2>, Connectivity)
! Access the donor extents (returns pointer to integer array of shape (3,2,<numconnections>))
call ovkGetConnectivityDonorExtents(Connectivity, DonorExtents)
! Access the donor interpolation coefficients (returns pointer to real(kind=ovk_rk) array of
! shape (<stencilsize>,<numdims>,<numconnections>))
call ovkGetConnectivityDonorInterpCoefs(Connectivity, DonorInterpCoefs)
! Access the donor coordinates (returns pointer to real(kind=ovk_rk) array of shape
! (<numdims>,<numconnections>))
call ovkGetConnectivityDonorCoords(Connectivity, DonorCoords)
! Access the receiver points (returns pointer to integer array of shape (3,<numconnections>))
call ovkGetConnectivityReceiverPoints(Connectivity, ReceiverPoints)Occasionally after running the assembly, you may get a result that doesn't match what you expected.
In these instances it can be helpful to look at the intermediate results of the steps described
above. For this purpose Overkit records some of these results in the grid state and provides a
command ovkFilterGridState to access them:
! Get the grid from the domain
call ovkGetGrid(Domain, <gridid>, Grid)
! Returns ovk_field_logical
call ovkFilterGridState(Grid, <stateflags>, <filtermode>, Mask)<stateflags> consists of one or more state flags combined using Fortran's bitwise 'or' (ior)
<filtermode> tells the command how to treat these flags. OVK_ANY returns points matching
any of the specified flags and OVK_ALL returns points matching all of them. OVK_NONE
and OVK_NOT_ALL are the respective complements of OVK_ANY and OVK_ALL.
Currently-supported state flags are:
OVK_STATE_GRID- point is part of the domainOVK_STATE_INTERIOR- point is on the interiorOVK_STATE_BOUNDARY- point is a boundaryOVK_STATE_EXTERIOR- point is outside of the domainOVK_STATE_DOMAIN_BOUNDARY- point is a domain boundaryOVK_STATE_INTERNAL_BOUNDARY- point is an internal boundaryOVK_STATE_OVERLAPPED- point is overlapped by another gridOVK_STATE_INFERRED_DOMAIN_BOUNDARY- point is a domain boundary that was detected during boundary inferenceOVK_STATE_BOUNDARY_HOLE- point was removed during boundary hole cuttingOVK_STATE_OCCLUDED- point is occluded by another gridOVK_STATE_FRINGE- point is part of the fringeOVK_STATE_OUTER_FRINGE- point is part of the outer fringeOVK_STATE_INNER_FRINGE- point is part of the inner fringeOVK_STATE_OVERLAP_MINIMIZED- point was removed during overlap minimizationOVK_STATE_RECEIVER- point is a receiverOVK_STATE_ORPHAN- point was a receiver but no suitable donors were found.
Overkit also provides commands to write the grids and states to a file (PLOT3D format) that can be opened in visualization tools such as TecPlot and ParaView. See the examples' source code for details.
@misc{overkitfortran,
author = {Smith, M.},
title = {Overkit},
year = {2019},
publisher = {GitHub},
journal = {GitHub repository},
howpublished = {\url{https://github.com/majosm/overkit-fortran}},
doi = {<get from Zenodo link above>}
}Overkit is one of several software projects produced by the University of Illinois Center for Exascale Simulation of Plasma-Coupled Combustion (XPACC): http://xpacc.illinois.edu.
This material is based in part upon work supported by the Department of Energy, National Nuclear Security Administration, under Award Number DE-NA0002374.