Appendix¶
Mint Application Code Example¶
Below is the complete Mint Application Code Example presented in
the Getting Started with Mint section. The code can be found in the Axom
source code under src/axom/mint/examples/user_guide/mint_getting_started.cpp
.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 | // sphinx_tutorial_walkthrough_includes_start
#include "axom/config.hpp" // compile time definitions
#include "axom/core/execution/execution_space.hpp" // for execution_space traits
#include "axom/mint.hpp" // for mint classes and functions
#include "axom/core/numerics/Matrix.hpp" // for numerics::Matrix
// sphinx_tutorial_walkthrough_includes_end
// namespace aliases
namespace mint = axom::mint;
namespace numerics = axom::numerics;
namespace xargs = mint::xargs;
using IndexType = axom::IndexType;
// compile-time switch for execution policy
#if defined(AXOM_USE_RAJA) && defined(AXOM_USE_CUDA)
constexpr int NUM_BLOCKS = 512;
using ExecPolicy = axom::CUDA_EXEC<NUM_BLOCKS>;
#elif defined(AXOM_USE_RAJA) && defined(AXOM_USE_OPENMP)
using ExecPolicy = axom::OMP_EXEC;
#else
using ExecPolicy = axom::SEQ_EXEC;
#endif
constexpr IndexType NUM_NODES_PER_CELL = 4;
constexpr double ONE_OVER_4 = 1. / static_cast<double>(NUM_NODES_PER_CELL);
/*!
* \brief Holds command-line arguments
*/
static struct
{
int res;
bool useUnstructured;
} Arguments;
//------------------------------------------------------------------------------
// FUNCTION PROTOTYPES
//------------------------------------------------------------------------------
void parse_args(int argc, char** argv);
mint::Mesh* getUniformMesh();
mint::Mesh* getUnstructuredMesh();
//------------------------------------------------------------------------------
// PROGRAM MAIN
//------------------------------------------------------------------------------
int main(int argc, char** argv)
{
parse_args(argc, argv);
// sphinx_tutorial_walkthrough_set_memory_start
// NOTE: use unified memory if we are using CUDA
const int allocID = axom::execution_space<ExecPolicy>::allocatorID();
axom::setDefaultAllocator(allocID);
// sphinx_tutorial_walkthrough_set_memory_end
// sphinx_tutorial_walkthrough_construct_mesh_start
mint::Mesh* mesh =
(Arguments.useUnstructured) ? getUnstructuredMesh() : getUniformMesh();
// sphinx_tutorial_walkthrough_construct_mesh_end
// sphinx_tutorial_walkthrough_add_fields_start
// add a cell-centered and a node-centered field
double* phi = mesh->createField<double>("phi", mint::NODE_CENTERED);
double* hc = mesh->createField<double>("hc", mint::CELL_CENTERED);
constexpr int NUM_COMPONENTS = 2;
double* xc =
mesh->createField<double>("xc", mint::CELL_CENTERED, NUM_COMPONENTS);
// sphinx_tutorial_walkthrough_add_fields_end
// sphinx_tutorial_walkthrough_compute_hf_start
// loop over the nodes and evaluate Himmelblaus Function
mint::for_all_nodes<ExecPolicy, xargs::xy>(
mesh,
AXOM_LAMBDA(IndexType nodeIdx, double x, double y) {
const double x_2 = x * x;
const double y_2 = y * y;
const double A = x_2 + y - 11.0;
const double B = x + y_2 - 7.0;
phi[nodeIdx] = A * A + B * B;
});
// sphinx_tutorial_walkthrough_compute_hf_end
// sphinx_tutorial_walkthrough_cell_centers_start
// loop over cells and compute cell centers
mint::for_all_cells<ExecPolicy, xargs::coords>(
mesh,
AXOM_LAMBDA(IndexType cellIdx,
const numerics::Matrix<double>& coords,
const IndexType* nodeIds) {
// NOTE: A column vector of the coords matrix corresponds to a nodes coords
// Sum the cell's nodal coordinates
double xsum = 0.0;
double ysum = 0.0;
double hsum = 0.0;
const IndexType numNodes = coords.getNumColumns();
for(IndexType inode = 0; inode < numNodes; ++inode)
{
const double* node = coords.getColumn(inode);
xsum += node[mint::X_COORDINATE];
ysum += node[mint::Y_COORDINATE];
hsum += phi[nodeIds[inode]];
} // END for all cell nodes
// compute cell centroid by averaging the nodal coordinate sums
const IndexType offset = cellIdx * NUM_COMPONENTS;
const double invnnodes = 1.f / static_cast<double>(numNodes);
xc[offset] = xsum * invnnodes;
xc[offset + 1] = ysum * invnnodes;
hc[cellIdx] = hsum * invnnodes;
});
// sphinx_tutorial_walkthrough_cell_centers_end
// sphinx_tutorial_walkthrough_vtk_start
// write the mesh in a VTK file for visualization
std::string vtkfile =
(Arguments.useUnstructured) ? "unstructured_mesh.vtk" : "uniform_mesh.vtk";
mint::write_vtk(mesh, vtkfile);
// sphinx_tutorial_walkthrough_vtk_end
delete mesh;
mesh = nullptr;
return 0;
}
//------------------------------------------------------------------------------
// FUNCTION PROTOTYPE IMPLEMENTATION
//------------------------------------------------------------------------------
void parse_args(int argc, char** argv)
{
Arguments.res = 25;
Arguments.useUnstructured = false;
for(int i = 1; i < argc; ++i)
{
if(strcmp(argv[i], "--unstructured") == 0)
{
Arguments.useUnstructured = true;
}
else if(strcmp(argv[i], "--resolution") == 0)
{
Arguments.res = std::atoi(argv[++i]);
}
} // END for all arguments
SLIC_ERROR_IF(
Arguments.res < 2,
"invalid mesh resolution! Please, pick a value greater than 2.");
}
//------------------------------------------------------------------------------
// sphinx_tutorial_walkthrough_construct_umesh_start
mint::Mesh* getUniformMesh()
{
// construct a N x N grid within a domain defined in [-5.0, 5.0]
const double lo[] = {-5.0, -5.0};
const double hi[] = {5.0, 5.0};
mint::Mesh* m = new mint::UniformMesh(lo, hi, Arguments.res, Arguments.res);
return (m);
}
// sphinx_tutorial_walkthrough_construct_umesh_end
//------------------------------------------------------------------------------
mint::Mesh* getUnstructuredMesh()
{
mint::Mesh* umesh = getUniformMesh();
const IndexType umesh_ncells = umesh->getNumberOfCells();
const IndexType umesh_nnodes = umesh->getNumberOfNodes();
const IndexType ncells = umesh_ncells * 4; // split each quad into 4 triangles
const IndexType nnodes = umesh_nnodes + umesh_ncells;
constexpr int DIMENSION = 2;
using MeshType = mint::UnstructuredMesh<mint::SINGLE_SHAPE>;
MeshType* m = new MeshType(DIMENSION, mint::TRIANGLE, nnodes, ncells);
m->resize(nnodes, ncells);
double* x = m->getCoordinateArray(mint::X_COORDINATE);
double* y = m->getCoordinateArray(mint::Y_COORDINATE);
IndexType* cells = m->getCellNodesArray();
// fill coordinates from uniform mesh
mint::for_all_nodes<ExecPolicy, xargs::xy>(
umesh,
AXOM_LAMBDA(IndexType nodeIdx, double nx, double ny) {
x[nodeIdx] = nx;
y[nodeIdx] = ny;
});
// loop over cells, compute cell centers and fill connectivity
mint::for_all_cells<ExecPolicy, xargs::coords>(
umesh,
AXOM_LAMBDA(IndexType cellIdx,
const numerics::Matrix<double>& coords,
const IndexType* nodeIds) {
// NOTE: A column vector of the coords matrix corresponds to a nodes coords
// Sum the cell's nodal coordinates
double xsum = 0.0;
double ysum = 0.0;
for(IndexType inode = 0; inode < NUM_NODES_PER_CELL; ++inode)
{
const double* node = coords.getColumn(inode);
xsum += node[mint::X_COORDINATE];
ysum += node[mint::Y_COORDINATE];
} // END for all cell nodes
// compute cell centroid by averaging the nodal coordinate sums
const IndexType nc = umesh_nnodes + cellIdx; /* centroid index */
x[nc] = xsum * ONE_OVER_4;
y[nc] = ysum * ONE_OVER_4;
// triangulate
const IndexType& n0 = nodeIds[0];
const IndexType& n1 = nodeIds[1];
const IndexType& n2 = nodeIds[2];
const IndexType& n3 = nodeIds[3];
const IndexType offset = cellIdx * 12;
cells[offset] = n0;
cells[offset + 1] = nc;
cells[offset + 2] = n3;
cells[offset + 3] = n0;
cells[offset + 4] = n1;
cells[offset + 5] = nc;
cells[offset + 6] = n1;
cells[offset + 7] = n2;
cells[offset + 8] = nc;
cells[offset + 9] = n2;
cells[offset + 10] = n3;
cells[offset + 11] = nc;
});
// delete uniform mesh
delete umesh;
umesh = nullptr;
return (m);
}
|
AXOM_LAMBDA Macro¶
The AXOM_LAMBDA
convenience macro expands to:
[=]
capture by value when the Axom Toolkit is compiled without CUDA.[=] __host__ __device__
when the Axom Toolkit is compiled with CUDA
Raw Sidre Data¶
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 | {
"name": "",
"groups":
{
"state":
{
"name": "state",
"groups":
{
"t1":
{
"name": "t1",
"views":
{
"block_id":
{
"name": "block_id",
"schema": "{\"dtype\":\"int32\", \"number_of_elements\": 1, \"offset\": 0, \"stride\": 4, \"element_bytes\": 4, \"endianness\": \"little\"}",
"value": "-1",
"state": "SCALAR",
"is_applied": 1
},
"partition_id":
{
"name": "partition_id",
"schema": "{\"dtype\":\"int32\", \"number_of_elements\": 1, \"offset\": 0, \"stride\": 4, \"element_bytes\": 4, \"endianness\": \"little\"}",
"value": "-1",
"state": "SCALAR",
"is_applied": 1
}
}
}
}
},
"coordsets":
{
"name": "coordsets",
"groups":
{
"c1":
{
"name": "c1",
"views":
{
"type":
{
"name": "type",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 9, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"explicit\"",
"state": "STRING",
"is_applied": 1
}
},
"groups":
{
"values":
{
"name": "values",
"views":
{
"x":
{
"name": "x",
"schema": "{\"dtype\":\"float64\", \"number_of_elements\": 6, \"offset\": 0, \"stride\": 8, \"element_bytes\": 8, \"endianness\": \"little\"}",
"value": "[0.0, 2.0, 1.0, 3.5, 2.5, 5.0]",
"state": "BUFFER",
"is_applied": 1
},
"y":
{
"name": "y",
"schema": "{\"dtype\":\"float64\", \"number_of_elements\": 6, \"offset\": 0, \"stride\": 8, \"element_bytes\": 8, \"endianness\": \"little\"}",
"value": "[0.0, 0.0, 1.0, 1.0, 2.0, 0.0]",
"state": "BUFFER",
"is_applied": 1
}
}
}
}
}
}
},
"topologies":
{
"name": "topologies",
"groups":
{
"t1":
{
"name": "t1",
"views":
{
"coordset":
{
"name": "coordset",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 3, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"c1\"",
"state": "STRING",
"is_applied": 1
},
"type":
{
"name": "type",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 13, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"unstructured\"",
"state": "STRING",
"is_applied": 1
}
},
"groups":
{
"elements":
{
"name": "elements",
"views":
{
"shape":
{
"name": "shape",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 4, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"tri\"",
"state": "STRING",
"is_applied": 1
},
"connectivity":
{
"name": "connectivity",
"schema": "{\"dtype\":\"int32\", \"number_of_elements\": 12, \"offset\": 0, \"stride\": 4, \"element_bytes\": 4, \"endianness\": \"little\"}",
"value": "[1, 3, 2, 2, 0, 1, 3, 4, 2, 1, 5, 3]",
"state": "BUFFER",
"is_applied": 1
},
"stride":
{
"name": "stride",
"schema": "{\"dtype\":\"int32\", \"number_of_elements\": 1, \"offset\": 0, \"stride\": 4, \"element_bytes\": 4, \"endianness\": \"little\"}",
"value": "3",
"state": "SCALAR",
"is_applied": 1
}
}
}
}
}
}
},
"fields":
{
"name": "fields",
"groups":
{
"den":
{
"name": "den",
"views":
{
"association":
{
"name": "association",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 8, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"element\"",
"state": "STRING",
"is_applied": 1
},
"volume_dependent":
{
"name": "volume_dependent",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 5, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"true\"",
"state": "STRING",
"is_applied": 1
},
"topology":
{
"name": "topology",
"schema": "{\"dtype\":\"char8_str\", \"number_of_elements\": 3, \"offset\": 0, \"stride\": 1, \"element_bytes\": 1, \"endianness\": \"little\"}",
"value": "\"t1\"",
"state": "STRING",
"is_applied": 1
},
"values":
{
"name": "values",
"schema": "{\"dtype\":\"float64\", \"number_of_elements\": 4, \"offset\": 0, \"stride\": 8, \"element_bytes\": 8, \"endianness\": \"little\"}",
"value": "[0.5, 1.2, 2.5, 0.9]",
"state": "BUFFER",
"is_applied": 1
}
}
}
}
}
}
}
|