1 /* ---------------------------------------------------------------------- 2 SPARTA - Stochastic PArallel Rarefied-gas Time-accurate Analyzer 3 http://sparta.sandia.gov 4 Steve Plimpton, sjplimp@sandia.gov 5 Michael Gallis, magalli@sandia.gov 6 Sandia National Laboratories 7 8 Copyright (2014) Sandia Corporation. Under the terms of Contract 9 DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains 10 certain rights in this software. This software is distributed under 11 the GNU General Public License. 12 13 See the README file in the top-level SPARTA directory. 14 ------------------------------------------------------------------------- */ 15 16 #ifndef SPARTA_SURF_H 17 #define SPARTA_SURF_H 18 19 #include "stdio.h" 20 #include "pointers.h" 21 #include "hash3.h" 22 #include "hashlittle.h" 23 24 namespace SPARTA_NS { 25 26 class Surf : protected Pointers { 27 public: 28 int exist; // 1 if any surfaces are defined, else 0 29 int implicit; // 1 = implicit surfs, 0 = explicit surfs 30 int distributed; // 1 = surfs spread across procs (exp or impl) 31 // 0 = each proc owns all 32 int surf_collision_check; // flag for whether init() check is required 33 // for assign of collision models to surfs 34 35 double bblo[3],bbhi[3]; // bounding box around surfs 36 int tally_comm; // style of comm for surf tallies 37 38 int nreact_one; // surface reactions in current step 39 bigint nreact_running; // running count of surface reactions 40 41 int ngroup; // # of defined groups 42 char **gnames; // name of each group 43 int *bitmask; // one-bit mask for each group 44 int *inversemask; // inverse mask for each group 45 46 // surf data structs 47 // explicit, all: each proc owns all surfs 48 // nlocal = Nsurf, nghost = 0 49 // tris = all surfs, nown = Nsurf/P, mytris = NULL 50 // explicit, distributed: each proc owns N/P surfs 51 // nlocal/nghost = surfs in owned/ghost grid cells 52 // tris = nloc+ngh surfs, nown = Nsurf/P, mytris = surfs I uniquely own 53 // implicit, distributed: each proc owns surfs in its owned grid cells 54 // nlocal = surfs in owned grid cells, nghost = surfs in ghost grid cells 55 // tris = nloc+ngh surfs, nown = nlocal, mytris = NULL 56 57 bigint nsurf; // total # of surf elements, lines or tris 58 59 struct Line { 60 surfint id; // unique ID for explicit surf 61 // cell ID for implicit surf 62 int type,mask; // type and mask of the surf 63 int isc,isr; // index of surface collision and reaction models 64 // -1 if unassigned 65 double p1[3],p2[3]; // end points of line segment 66 // rhand rule: Z x (p2-p1) = outward normal 67 double norm[3]; // outward normal to line segment 68 int transparent; // 1 if surf is transparent 69 }; 70 71 struct Tri { 72 surfint id; // unique ID for explicit surf 73 // cell ID for implicit surf 74 int type,mask; // type and mask of the surf 75 int isc,isr; // index of surface collision and reaction models 76 // -1 if unassigned 77 double p1[3],p2[3],p3[3]; // corner points of triangle 78 // rhand rule: (p2-p1) x (p3-p1) = outward normal 79 double norm[3]; // outward normal to triangle 80 int transparent; // 1 if surf is transparent 81 }; 82 83 Line *lines; // list of lines for surface collisions 84 Tri *tris; // list of tris for surface collisions 85 int nlocal; // # of lines or tris 86 // explicit, all: nlocal = nsurf 87 // explicit, distributed: 88 // surfs which overlap my owned grid cells 89 // implicit: surfs within my owned grid cells 90 int nghost; // # of ghost surfs I store for collisions 91 // explicit, all: nghost = 0 92 // explicit, distributed: 93 // surfs which overlap my ghost grid cells 94 // implicit: surfs within my ghost grid cells 95 int nmax; // max length of lines/tris vecs 96 97 Line *mylines; // list of lines assigned uniquely to me 98 // only for explicit, distributed 99 Tri *mytris; // list of tris assigned uniquely to me 100 // only for explicit, distributed 101 int nown; // # of lines or tris I own uniquely 102 int maxown; // max length of owned lines/tris vecs 103 104 Line *tmplines; // list of temporary lines, filled by ReadSurf 105 Tri *tmptris; // list of temporary tris, filled by ReadSurf 106 int ntmp; // # of temporary surfs 107 int nmaxtmp; // max size of tmplines/tmptris 108 109 int nsc,nsr; // # of surface collision and reaction models 110 class SurfCollide **sc; // list of surface collision models 111 class SurfReact **sr; // list of surface reaction models 112 113 int pushflag; // set to 1 to push surf pts near grid cell faces 114 double pushlo,pushhi; // lo/hi ranges to push on 115 double pushvalue; // new position to push to 116 117 // extra custom vectors/arrays for per-surf data 118 // ncustom > 0 if there are any extra arrays 119 // custom attributes are created by various commands 120 // these variables are public, others below are private 121 122 int ncustom; // # of custom attributes, some may be deleted 123 int *etype; // type = INT/DOUBLE of each attribute 124 int *esize; // size = 0 for vector, N for array columns 125 int *ewhich; // index into eivec,eiarray,edvec,edarray for data 126 127 int **eivec; // pointer to each integer vector 128 int ***eiarray; // pointer to each integer array 129 double **edvec; // pointer to each double vector 130 double ***edarray; // pointer to each double array 131 132 #include "hash_options.h" 133 134 #ifdef SPARTA_MAP 135 typedef std::map<surfint,int> MySurfHash; 136 typedef std::map<OnePoint2d,int> MyHashPoint; 137 typedef std::map<OnePoint2d,int>::iterator MyPointIt; 138 typedef std::map<TwoPoint3d,int> MyHash2Point; 139 typedef std::map<TwoPoint3d,int>::iterator My2PointIt; 140 typedef std::map<cellint,int> MyCellHash; 141 #elif defined SPARTA_UNORDERED_MAP 142 typedef std::unordered_map<surfint,int> MySurfHash; 143 typedef std::unordered_map<OnePoint2d,int,OnePoint2dHash> MyHashPoint; 144 typedef std::unordered_map<OnePoint2d,int,OnePoint2dHash>::iterator MyPointIt; 145 typedef std::unordered_map<TwoPoint3d,int,TwoPoint3dHash> MyHash2Point; 146 typedef std::unordered_map<TwoPoint3d,int,TwoPoint3dHash>::iterator My2PointIt; 147 typedef std::unordered_map<cellint,int> MyCellHash; 148 #else 149 typedef std::tr1::unordered_map<surfint,int> MySurfHash; 150 typedef std::tr1::unordered_map<OnePoint2d,int,OnePoint2dHash> MyHashPoint; 151 typedef std::tr1::unordered_map<OnePoint2d,int,OnePoint2dHash>:: 152 iterator MyPointIt; 153 typedef std::tr1::unordered_map<TwoPoint3d,int,TwoPoint3dHash> MyHash2Point; 154 typedef std::tr1::unordered_map<TwoPoint3d,int,TwoPoint3dHash>:: 155 iterator My2PointIt; 156 typedef std::tr1::unordered_map<cellint,int> MyCellHash; 157 #endif 158 159 MySurfHash *hash; // hash for nlocal surf IDs 160 int hashfilled; // 1 if hash is filled with surf IDs 161 162 Surf(class SPARTA *); 163 ~Surf(); 164 void global(char *); 165 void modify_params(int, char **); 166 void init(); 167 void clear(); 168 void remove_ghosts(); 169 void add_line(surfint, int, double *, double *); 170 void add_line_copy(int, Line *); 171 void add_line_own(surfint, int, double *, double *); 172 void add_line_temporary(surfint, int, double *, double *); 173 void add_tri(surfint, int, double *, double *, double *); 174 void add_tri_copy(int, Tri *); 175 void add_tri_own(surfint, int, double *, double *, double *); 176 void add_tri_own_clip(surfint, int, double *, double *, double *); 177 void add_tri_temporary(surfint, int, double *, double *, double *); 178 void rehash(); 179 int all_transparent(); 180 181 void setup_owned(); 182 void bbox_all(); 183 void bbox_one(void *, double *, double *); 184 185 void compute_line_normal(int); 186 void compute_tri_normal(int); 187 void quad_corner_point(int, double *, double *, double *); 188 void hex_corner_point(int, double *, double *, double *); 189 190 double line_size(int); 191 double line_size(Line *); 192 double line_size(double *, double *); 193 double axi_line_size(int); 194 double axi_line_size(Line *); 195 double tri_size(int, double &); 196 double tri_size(Tri *, double &); 197 double tri_size(double *, double *, double *, double &); 198 199 void check_watertight_2d(); 200 void check_watertight_3d(); 201 void check_point_inside(int); 202 void check_point_near_surf_2d(); 203 void check_point_near_surf_3d(); 204 205 void output_extent(int); 206 double shortest_line(int); 207 void smallest_tri(int, double &, double &); 208 209 void add_collide(int, char **); 210 int find_collide(const char *); 211 void add_react(int, char **); 212 int find_react(const char *); 213 214 void group(int, char **); 215 int add_group(const char *); 216 int find_group(const char *); 217 218 void compress_explicit(); 219 void compress_implicit(); 220 221 void collate_vector(int, surfint *, double *, int, double *); 222 void collate_vector_reduce(int, surfint *, double *, int, double *); 223 void collate_vector_rendezvous(int, surfint *, double *, int, double *); 224 225 void collate_array(int, int, surfint *, double **, double **); 226 void collate_array_reduce(int, int, surfint *, double **, double **); 227 void collate_array_rendezvous(int, int, surfint *, double **, double **); 228 void collate_vector_implicit(int, surfint *, double *, double *); 229 void collate_array_implicit(int, int, surfint *, double **, double **); 230 231 void redistribute_lines_clip(int, int); 232 void redistribute_lines_temporary(int); 233 void redistribute_tris_clip(int, int); 234 void redistribute_tris_temporary(int); 235 236 int find_custom(char *); 237 void error_custom(); 238 int add_custom(char *, int, int); 239 void allocate_custom(int, int); 240 void remove_custom(int); 241 242 void write_restart(FILE *); 243 void read_restart(FILE *); 244 virtual void grow(int); 245 virtual void grow_own(int); 246 virtual void grow_temporary(int); 247 bigint memory_usage(); 248 249 protected: 250 int me,nprocs; 251 int maxsc; // max # of models in sc 252 int maxsr; // max # of models in sr 253 254 // collate vector rendezvous data 255 256 struct InRvousVec { 257 surfint id; // surface ID 258 double value; // compute value 259 }; 260 261 double *out_rvous; 262 int ncol_rvous; 263 264 // watertight rendezvous data 265 266 struct InRvousPoint { 267 double x[2]; // 2d point coords 268 int which; // 1 for first endpoint, 2 for second endpoint 269 }; 270 271 struct InRvousEdge { 272 double x1[3],x2[3]; // 3d edge point coords 273 int which; // 1 for forward order, 2 for reverse order 274 }; 275 276 // union data struct for packing 32-bit and 64-bit ints into double bufs 277 // this avoids aliasing issues by having 2 pointers (double,int) 278 // to same buf memory 279 // constructor for 32-bit int prevents compiler 280 // from possibly calling the double constructor when passed an int 281 // copy to a double *buf: 282 // buf[m++] = ubuf(foo).d, where foo is a 32-bit or 64-bit int 283 // copy from a double *buf: 284 // foo = (int) ubuf(buf[m++]).i;, where (int) or (tagint) match foo 285 // the cast prevents compiler warnings about possible truncation 286 287 union ubuf { 288 double d; 289 int64_t i; ubuf(double arg)290 ubuf(double arg) : d(arg) {} ubuf(int64_t arg)291 ubuf(int64_t arg) : i(arg) {} ubuf(int arg)292 ubuf(int arg) : i(arg) {} 293 }; 294 295 // extra custom vectors/arrays for per-surf data 296 // these variables are private, others above are public 297 298 char **ename; // name of each attribute 299 300 int ncustom_ivec; // # of integer vector attributes 301 int ncustom_iarray; // # of integer array attributes 302 int *icustom_ivec; // index into ncustom for each integer vector 303 int *icustom_iarray; // index into ncustom for each integer array 304 int *eicol; // # of columns in each integer array (esize) 305 306 int ncustom_dvec; // # of double vector attributes 307 int ncustom_darray; // # of double array attributes 308 int *icustom_dvec; // index into ncustom for each double vector 309 int *icustom_darray; // index into ncustom for each double array 310 int *edcol; // # of columns in each double array (esize) 311 312 int *custom_restart_flag; // flag on each custom vec/array read from restart 313 // used to delete them if not redefined in 314 // restart script 315 316 // private methods 317 318 void point_line_compare(double *, double *, double *, double, int &, int &); 319 void point_tri_compare(double *, double *, double *, double *, double *, 320 double, int &, int &, int, int, int); 321 322 void collate_vector_allreduce(int, int *, double *, int, double *); 323 void collate_vector_irregular(int, int *, double *, int, double *); 324 void collate_array_allreduce(int, int, int *, double **, double **); 325 void collate_array_irregular(int, int, int *, double **, double **); 326 327 void check_watertight_2d_all(); 328 void check_watertight_2d_distributed(); 329 void check_watertight_3d_all(); 330 void check_watertight_3d_distributed(); 331 332 // callback functions for rendezvous communication 333 334 static int rendezvous_vector(int, char *, int &, int *&, char *&, void *); 335 static int rendezvous_array(int, char *, int &, int *&, char *&, void *); 336 static int rendezvous_implicit(int, char *, int &, int *&, char *&, void *); 337 static int rendezvous_watertight_2d(int, char *, 338 int &, int *&, char *&, void *); 339 static int rendezvous_watertight_3d(int, char *, 340 int &, int *&, char *&, void *); 341 static int rendezvous_lines(int, char *, 342 int &, int *&, char *&, void *); 343 static int rendezvous_tris(int, char *, 344 int &, int *&, char *&, void *); 345 }; 346 347 } 348 349 #endif 350 351 /* ERROR/WARNING messages: 352 353 E: Illegal ... command 354 355 Self-explanatory. Check the input script syntax and compare to the 356 documentation for the command. You can use -echo screen as a 357 command-line option when running SPARTA to see the offending line. 358 359 E: Could not find surf_modify surf-ID 360 361 Self-explanatory. 362 363 E: Could not find surf_modify sc-ID 364 365 Self-explanatory. 366 367 E: %d surface elements not assigned to a collision model 368 369 All surface elements must be assigned to a surface collision model via 370 the surf_modify command before a simulation is perforemd. 371 372 E: Reuse of surf_collide ID 373 374 A surface collision model ID cannot be used more than once. 375 376 E: Invalid surf_collide style 377 378 Self-explanatory. 379 380 */ 381