1 // Copyright(C) 1999-2020 National Technology & Engineering Solutions 2 // of Sandia, LLC (NTESS). Under the terms of Contract DE-NA0003525 with 3 // NTESS, the U.S. Government retains certain rights in this software. 4 // 5 // See packages/seacas/LICENSE for details 6 7 #ifndef IOSS_Iogs_GeneratedMesh_h 8 #define IOSS_Iogs_GeneratedMesh_h 9 10 #include "vtk_ioss_mangle.h" 11 12 #include <Ioss_CodeTypes.h> 13 #include <Ioss_EntityType.h> // for EntityType 14 #include <array> 15 #include <cstddef> // for size_t 16 #include <cstdint> // for int64_t 17 #include <map> // for map, etc 18 #include <string> // for string 19 #include <utility> // for pair 20 #include <vector> // for vector 21 22 namespace Iogs { 23 using MapVector = std::vector<int64_t>; 24 25 class GeneratedMesh 26 { 27 public: 28 enum ShellLocation { MX = 0, PX = 1, MY = 2, PY = 3, MZ = 4, PZ = 5 }; 29 30 /** 31 Generate a cube mesh of size 'num_x' by 'num_y' by 'num_z' elements. 32 By default, the mesh is gen_struc on a single processor. If 'proc_count' is 33 greater than 1, then the mesh will be distributed over 'proc_count' processors 34 and this process will get the portion of the mesh for 'my_proc'. 35 The mesh will be decomposed along the 'Z' axis so 'num_z' must be greater than 36 or equal to 'proc_count' and for even distribution of the hexes 'num_z' mod 'proc_count' 37 should be zero. 38 39 The mesh can optionally include sidesets along each 40 face of the cube mesh. These are specified via the 41 'add_sidesets' function. 42 43 If the 'parameters' string constructor is used, the string 44 is parsed to determine the intervals in each direction and, 45 optionally, additional information. The form of the string 46 is "IxJxK" where I, J, and K are the number of intervals 47 in the X, Y, and Z directions respectively and the "x" are 48 literal 'x' characters. For example, the constructor 49 GeneratedMesh("10x12x14") will create the same mesh as 50 GeneratedMesh(10,12,14) 51 52 Additional valid options are: 53 - help -- no argument, shows valid options 54 - show -- no argument, prints out a summary of the 55 GeneratedMesh() parameters. The output will look similar 56 to: 57 \code 58 "10x12x8|bbox:-10,-10,-10,10,10,10|sideset:XYZ|show" 59 60 Mesh Parameters: 61 Intervals: 10 by 12 by 8 62 X = 2 * (0..10) + -10 Range: -10 <= X <= 10 63 Y = 1.66667 * (0..12) + -10 Range: -10 <= Y <= 10 64 Z = 2.5 * (0..8) + -10 Range: -10 <= Z <= 10 65 Node Count (total) = 1287 66 Element Count (total) = 1152 67 Block Count = 3 68 SideSet Count = 3 69 \endcode 70 71 - sideset -- argument = xXyYzZ which specifies whether there is 72 a sideset at that location. 'x' is minimum x face, 'X' is 73 maximum x face, similarly for y and z. Note that the argument 74 string is a single multicharacter string. You can add multiple 75 sidesets to a face, for example, sideset:xxx would add three 76 sidesets on the minimum x face. An error is output if a non 77 xXyYzZ character is found, but execution continues. 78 79 - zdecomp -- argument = n0, n1, n2, ..., n#proc-1 which are the number 80 of intervals in the z direction for each processor in a pallel run. 81 If this option is specified, then the total number of intervals in the 82 z direction is the sum of the n0, n1, ... An interval count must be 83 specified for each processor. If this option is not specified, then 84 the number of intervals on each processor in the z direction is 85 numZ/numProc with the extras added to the lower numbered processors. 86 87 - scale -- argument = xs, ys, zs which are the scale factors in the x, 88 y, and z directions. All three must be specified if this option is 89 present. 90 91 - offset -- argument = xoff, yoff, zoff which are the offsets in the 92 x, y, and z directions. All three must be specified if this option 93 is present. 94 95 - bbox -- argument = xmin, ymin, zmin, xmax, ymax, zmax 96 which specify the lower left and upper right corners of 97 the bounding box for the gen_struc mesh. This will 98 calculate the scale and offset which will fit the mesh in 99 the specified box. All calculations are based on the currently 100 active interval settings. If scale or offset or zdecomp 101 specified later in the option list, you may not get the 102 desired bounding box. 103 104 - rotate -- argument = axis,angle,axis,angle,... 105 where axis is 'x', 'y', or 'z' and angle is the rotation angle in 106 degrees. Multiple rotations are cumulative. The composite rotation 107 matrix is applied at the time the coordinates are retrieved after 108 scaling and offset are applied. 109 110 The unrotated coordinate of a node at grid location i,j,k is: 111 \code 112 x = x_scale * i + x_off, 113 y = z_scale * j + y_off, 114 z = z_scale * k + z_off, 115 \endcode 116 117 The extent of the unrotated mesh will be: 118 \code 119 x_off <= x <= x_scale * numX + x_off 120 y_off <= y <= y_scale * numY + y_off 121 z_off <= z <= z_scale * numZ + z_off 122 \endcode 123 124 If an unrecognized option is specified, an error message will be 125 output and execution will continue. 126 127 An example of valid input is: 128 \code 129 "10x20x40|scale:1,0.5,0.25|offset:-5,-5,-5" 130 \endcode 131 132 This would create a mesh with 10 intervals in x, 20 in y, 40 in z 133 The mesh would be centered on 0,0,0 with a range of 10 in each 134 direction. 135 136 NOTE: All options are processed in the order they appear in 137 the parameters string (except rotate which is applied at the 138 time the coordinates are gen_struc/retrieved) 139 */ 140 explicit GeneratedMesh(const std::string ¶meters, int proc_count = 1, int my_proc = 0); 141 GeneratedMesh(int64_t num_x, int64_t num_y, int64_t num_z, int proc_count = 1, int my_proc = 0); 142 GeneratedMesh(); 143 virtual ~GeneratedMesh(); 144 145 /** 146 * Add a sideset along the specified face of the hex mesh. 147 * The sidesets will maintain the order of definition. The 148 * first sideset defined will be sideset 1. 149 * The loc options are: 150 * - MX = add sideset on the face with minimum X 151 * - PX = add sideset on the face with maximum X 152 * - MY = add sideset on the face with minimum Y 153 * - PY = add sideset on the face with maximum Y 154 * - MZ = add sideset on the face with minimum Z 155 * - PZ = add sideset on the face with maximum Z 156 * 157 */ 158 int64_t add_sideset(ShellLocation loc); 159 160 /** 161 * Specify the coordinate scaling and offset in all three 162 * spatial dimensions. 163 * 164 * node location of node at (i,j,k) is 165 * \code 166 * X = scale X * i + offset X 167 * Y = scale Y * i + offset Y 168 * Z = scale Z * i + offset Z 169 * \endcode 170 * 171 * WARNING: Should be called before retrieving node 172 * coordinates. 173 */ 174 void set_scale(double scl_x, double scl_y, double scl_z); 175 void set_offset(double off_x, double off_y, double off_z); 176 void set_bbox(double xmin, double ymin, double zmin, double xmax, double ymax, double zmax); 177 178 /** 179 * Set rotation. Multiple calls are cumulative. 180 * Rotate 'angle_degrees' degrees about the axis 'axis' 181 * Center of rotation is about the origin and operates 182 * on the scaled/offset coordinates of the mesh. 183 */ 184 void set_rotation(const std::string &axis, double angle_degrees); 185 186 /** 187 * Return number of nodes in the entire model. 188 */ 189 virtual int64_t node_count() const; 190 191 /** 192 * Return number of nodes on this processor. 193 */ 194 virtual int64_t node_count_proc() const; 195 196 /** 197 * Return number of structured blocks in the entire model. 198 */ 199 virtual int64_t structured_block_count() const; 200 201 /** 202 * Return number of sidesets in the entire model. 203 */ 204 virtual int64_t sideset_count() const; 205 206 /** 207 * Return number of sideset 'sides' on sideset 'id' 208 */ 209 int64_t sideset_side_count(int64_t id) const; 210 211 /** 212 * Return number of sideset 'sides' on sideset 'id' on the current 213 * processor. 214 */ 215 virtual int64_t sideset_side_count_proc(int64_t id) const; 216 block_range(int64_t)217 Ioss::IJK_t block_range(int64_t /* id */) const 218 { 219 return Ioss::IJK_t{{(int)numX, (int)numY, (int)numZ}}; 220 } 221 Ioss::IJK_t block_range_proc(int64_t id) const; 222 Ioss::IJK_t block_offset_proc(int64_t id) const; 223 224 /** 225 * Return number of elements in all structured blocks in the model. 226 */ 227 virtual int64_t element_count() const; 228 229 /** 230 * Return number of elements in all structured blocks on this processor. 231 */ 232 virtual int64_t element_count_proc() const; 233 timestep_count()234 int64_t timestep_count() const { return timestepCount; } 235 /** 236 * Return number of elements in the structured block with id 237 * 'block_number'. The 'block_number' ranges from '1' to 238 * 'block_count()'. 239 */ 240 virtual int64_t element_count(int64_t block_number) const; 241 242 /** 243 * Return number of elements on this processor in the structured 244 * block with id 'block_number'. The 'block_number' ranges from 245 * '1' to 'block_count()'. 246 */ 247 virtual int64_t element_count_proc(int64_t block_number) const; 248 249 /** 250 * Returns pair containing "topology type string" and "number of 251 * nodes / element". The topology type string will be "hex8" for 252 * the hex element block 253 */ 254 virtual std::pair<std::string, int> topology_type(int64_t block_number) const; 255 256 void build_node_map(Ioss::Int64Vector &map, std::vector<int> &proc, int64_t slab, 257 size_t slabOffset, size_t adjacentProc, size_t index); 258 virtual int64_t communication_node_count_proc() const; 259 virtual void node_communication_map(MapVector &map, std::vector<int> &proc); 260 virtual void owning_processor(int *owner, int64_t num_node); 261 262 /** 263 * Fill the passed in 'map' argument with the node map 264 * "map[local_position] = global_id" for the nodes on this 265 * processor. 266 */ 267 virtual void node_map(MapVector &map) const; 268 virtual void node_map(Ioss::IntVector &map) const; 269 270 /** 271 * Fill the passed in 'map' argument with the element map 272 * "map[local_position] = global_id" for the elements on this 273 * processor in block "block_number". 274 */ 275 virtual void element_map(int64_t block_number, MapVector &map) const; 276 virtual void element_map(int64_t block_number, Ioss::IntVector &map) const; 277 278 /** 279 * Fill the passed in 'map' argument with the element map 280 * "map[local_position] = global_id" for all elements on this 281 * processor 282 */ 283 virtual void element_map(MapVector &map) const; 284 virtual void element_map(Ioss::IntVector &map) const; 285 286 /** 287 * Fill the passed in 'map' argument with the element map pair 288 * "map[local_position] = element global_id" and 289 * "map[local_position+1] = element local face id (0-based)" for 290 * all elements on the current processor having a face on the 291 * surface defined by ShellLocation. 292 */ 293 void element_surface_map(ShellLocation loc, MapVector &map) const; 294 295 /** 296 * Return the connectivity for the elements on this processor in 297 * the block with id 'block_number'. If the elements in this block 298 * have 'npe' nodes per element, then the first 'npe' entries in 299 * the 'conn' vector will be the nodal connectivity for the first 300 * element; the next 'npe' entries are the nodal connectivity for 301 * the second element. The 'connect' vector will be resized to the 302 * size required to contain the nodal connectivity for the 303 * specified block; all information in 'connect' will be overwritten. 304 */ 305 void connectivity(int64_t block_number, Ioss::Int64Vector &connect) const; 306 void connectivity(int64_t block_number, Ioss::IntVector &connect) const; 307 void connectivity(int64_t block_number, int64_t *connect) const; 308 virtual void connectivity(int64_t block_number, int *connect) const; 309 310 /** 311 * Return the coordinates for all nodes on this processor. The 312 * first 3 entries in the 'coord' vector are the x, y, and z 313 * coordinates of the first node, etc. The 'coord' vector will be 314 * resized to the size required to contain the nodal coordinates; 315 * all information in 'coord' will be overwritten. 316 */ 317 virtual void coordinates(std::vector<double> &coord) const; 318 virtual void coordinates(double *coord) const; 319 320 /** 321 * Return the coordinates for all nodes on this processor in 322 * separate vectors. The vectors will be resized to the size 323 * required to contain the nodal coordinates; all information in 324 * the vectors will be overwritten. 325 */ 326 virtual void coordinates(std::vector<double> &x, std::vector<double> &y, 327 std::vector<double> &z) const; 328 329 /** 330 * Return the coordinates for componenet 'comp' (1=x, 2=y, 3=z) 331 * for all nodes on this processor. The 332 * vector will be resized to the size required to contain the 333 * nodal coordinates; all information in the vector will be 334 * overwritten. 335 * It is an error to request the coordinates via this function 336 * if a rotation is defined. 337 */ 338 virtual void coordinates(int component, std::vector<double> &xyz) const; 339 340 /** 341 * Return the coordinates for componenet 'comp' (1=x, 2=y, 3=z, 0=all) 342 * for all nodes in zone `zone` on this processor. The 343 * vector will be resized to the size required to contain the 344 * nodal coordinates; all information in the vector will be 345 * overwritten. 346 * It is an error to request the coordinates via this function 347 * if a rotation is defined. 348 */ 349 void coordinates(int component, int zone, double *coord) const; 350 351 /** 352 * Return the list of the face/ordinal pairs 353 * "elem_sides[local_position] = element global_id" and 354 * "elem_sides[local_position+1] = element local face id (0-based)" 355 * for the faces in sideset 'id' on this 356 * processor. The 'elem_sides' vector will be resized to the size 357 * required to contain the list. The element ids are global ids, 358 * the side ordinal is 0-based. 359 */ 360 virtual void sideset_elem_sides(int64_t id, Ioss::Int64Vector &elem_sides) const; 361 362 virtual std::vector<std::string> sideset_touching_blocks(int64_t set_id) const; 363 get_num_x()364 int64_t get_num_x() const { return numX; } get_num_y()365 int64_t get_num_y() const { return numY; } get_num_z()366 int64_t get_num_z() const { return numZ; } 367 get_variable_count(Ioss::EntityType type)368 size_t get_variable_count(Ioss::EntityType type) const 369 { 370 return variableCount.find(type) != variableCount.end() ? variableCount.find(type)->second : 0; 371 } 372 373 private: 374 template <typename INT> void raw_element_map(int64_t block_number, std::vector<INT> &map) const; 375 template <typename INT> void raw_element_map(std::vector<INT> &map) const; 376 template <typename INT> void raw_connectivity(int64_t block_number, INT *connect) const; 377 378 GeneratedMesh(const GeneratedMesh &); 379 GeneratedMesh &operator=(const GeneratedMesh &); 380 381 void set_variable_count(const std::string &type, size_t count); 382 void parse_options(const std::vector<std::string> &groups); 383 void show_parameters() const; 384 void initialize(); 385 386 std::vector<ShellLocation> sidesets{}; 387 std::array<std::array<double, 3>, 3> rotmat; 388 size_t numX{0}, numY{0}, numZ{0}; 389 size_t myNumZ{0}, myStartZ{0}; 390 391 size_t processorCount{1}; 392 size_t myProcessor{0}; 393 394 size_t timestepCount{0}; 395 std::map<Ioss::EntityType, size_t> variableCount{}; 396 397 double offX{0}, offY{0}, offZ{0}; /** Offsets in X, Y, and Z directions */ 398 double sclX{1}, sclY{1}, sclZ{1}; /** Scale in X, Y, and Z directions 399 * location of node at (i,j,k) 400 * position is (sclX*i+offX, 401 * sclY*i+offY, sclZ*i+offZ) */ 402 bool doRotation{false}; 403 }; 404 } // namespace Iogs 405 #endif 406