Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #pragma once
11 :
12 : #include "GeneralUserObject.h"
13 :
14 : // Local Includes
15 : #include "ElemIndexHelper.h"
16 : #include "ParallelRayStudy.h"
17 : #include "RayTracingAttributes.h"
18 : #include "TraceData.h"
19 : #include "TraceRayBndElement.h"
20 :
21 : // MOOSE Includes
22 : #include "TheWarehouse.h"
23 :
24 : // libMesh includes
25 : #include "libmesh/mesh.h"
26 : #include "libmesh/elem_side_builder.h"
27 :
28 : // Forward Declarations
29 : class RayBoundaryConditionBase;
30 : class RayKernelBase;
31 : class TraceRay;
32 : class RayTracingObject;
33 :
34 : /**
35 : * Base class for Ray tracing studies that will generate Rays and then propagate
36 : * all of them to termination.
37 : *
38 : * Subclasses _must_ override generateRays()
39 : */
40 : class RayTracingStudy : public GeneralUserObject
41 : {
42 : public:
43 : RayTracingStudy(const InputParameters & parameters);
44 :
45 : /**
46 : * Key that is used for restricting access to acquireRayInternal().
47 : *
48 : * In general, we restrict access to the construction of Rays to only protected
49 : * methods within RayTracingStudy. However, there are a few classes/functions
50 : * that have a legimate reason to construct Rays. Here, the friends are said
51 : * classes/functions.
52 : */
53 : class AcquireRayInternalKey
54 : {
55 : friend class Parallel::Packing<std::shared_ptr<Ray>>;
56 : friend void dataLoad(std::istream & stream, std::shared_ptr<Ray> & ray, void * context);
57 : AcquireRayInternalKey() {}
58 : AcquireRayInternalKey(const AcquireRayInternalKey &) {}
59 : };
60 :
61 : /**
62 : * Key that is used for restricting access to moveRayToBufferDuringTrace() and
63 : * acquireRayDuringTrace().
64 : */
65 : class AcquireMoveDuringTraceKey
66 : {
67 : friend class RayBoundaryConditionBase;
68 : friend class RayKernelBase;
69 : AcquireMoveDuringTraceKey() {}
70 : AcquireMoveDuringTraceKey(const AcquireMoveDuringTraceKey &) {}
71 : };
72 :
73 : static InputParameters validParams();
74 :
75 : virtual void initialSetup() override;
76 : virtual void residualSetup() override;
77 : virtual void jacobianSetup() override;
78 : virtual void meshChanged() override;
79 : virtual void timestepSetup() override;
80 8208 : virtual void initialize() override {}
81 8108 : virtual void finalize() override {}
82 :
83 : /**
84 : * Executes the study (generates and propagates Rays)
85 : */
86 : virtual void execute() override;
87 :
88 : /**
89 : * Setup for on subdomain change or subdomain AND ray change during ray tracing
90 : * @param subdomain The subdomain changed to
91 : * @param tid Thread id
92 : * @param ray_id ID of the ray initiating the change
93 : */
94 : virtual void
95 : segmentSubdomainSetup(const SubdomainID subdomain, const THREAD_ID tid, const RayID ray_id);
96 :
97 : /**
98 : * Reinitialize objects for a Ray segment for ray tracing
99 : * @param elem The elem the segment is in
100 : * @param start Start point of the segment
101 : * @param end End point of the segment
102 : * @param length The length of the start -> end segment
103 : * @param tid Thread id
104 : */
105 : virtual void reinitSegment(
106 : const Elem * elem, const Point & start, const Point & end, const Real length, THREAD_ID tid);
107 :
108 : /**
109 : * Called at the end of a Ray segment
110 : * @param tid Thread id
111 : * @param ray The ray
112 : */
113 : virtual void postOnSegment(const THREAD_ID tid, const std::shared_ptr<Ray> & ray);
114 :
115 : /**
116 : * Called at the beginning of a trace for a ray
117 : * @param tid Thread id
118 : * @param ray The ray
119 : */
120 5380196 : virtual void preTrace(const THREAD_ID /* tid */, const std::shared_ptr<Ray> & /* ray */) {}
121 :
122 : /**
123 : * Method for executing the study so that it can be called out of the standard UO execute()
124 : */
125 : void executeStudy();
126 :
127 : /**
128 : * Total number of processor crossings for Rays that finished on this processor
129 : */
130 : unsigned long long int endingProcessorCrossings() const { return _ending_processor_crossings; }
131 : /**
132 : * Max number of total processor crossings for Rays that finished on this processor
133 : */
134 : unsigned int endingMaxProcessorCrossings() const { return _ending_max_processor_crossings; }
135 : /**
136 : * Total number of processor crossings
137 : */
138 14 : unsigned long long int totalProcessorCrossings() const { return _total_processor_crossings; }
139 : /**
140 : * Max number of processor crossings for all Rays
141 : */
142 14 : unsigned int maxProcessorCrossings() const { return _max_processor_crossings; }
143 :
144 : /**
145 : * Total number of Ray/element intersections for Rays that finished on this processor
146 : */
147 : unsigned long long int endingIntersections() const { return _ending_intersections; }
148 : /**
149 : * Max number of intersections for Rays that finished on this processor
150 : */
151 : unsigned int endingMaxIntersections() const { return _ending_max_intersections; }
152 : /**
153 : * Total number of Ray/element intersections
154 : */
155 : unsigned long long int totalIntersections() const { return _total_intersections; }
156 : /**
157 : * Max number of intersections for a Ray
158 : */
159 : unsigned int maxIntersections() const { return _max_intersections; }
160 : /**
161 : * Max number of trajectory changes for a Ray
162 : */
163 : unsigned int maxTrajectoryChanges() const { return _max_trajectory_changes; }
164 :
165 : /**
166 : * Total amount of distance traveled by the rays that end on this processor
167 : */
168 : Real endingDistance() const { return _ending_distance; }
169 : /**
170 : * Total distance traveled by all Rays
171 : */
172 2266 : Real totalDistance() const { return _total_distance; }
173 :
174 : unsigned long long int localTraceRayResult(const int result) const
175 : {
176 196 : return _local_trace_ray_results[result];
177 : }
178 :
179 : const ParallelRayStudy & parallelRayStudy() const { return *_parallel_ray_study; }
180 :
181 : /**
182 : * Max distance any Ray can travel
183 : */
184 35293907 : Real rayMaxDistance() const { return _ray_max_distance; }
185 :
186 : /**
187 : * Duration for execute() in seconds
188 : */
189 : Real executionTime() { return std::chrono::duration<Real>(_execution_time).count(); }
190 : /**
191 : * Duration for execute() in nanoseconds
192 : */
193 : Real executionTimeNano()
194 : {
195 : return std::chrono::duration<Real, std::nano>(_execution_time).count();
196 : }
197 : /**
198 : * Duration for creation of all Rays in seconds
199 : */
200 : Real generationTime() const { return std::chrono::duration<Real>(_generation_time).count(); }
201 : /**
202 : * Duration for creation of all Rays in seconds
203 : */
204 : Real propagationTime() const { return std::chrono::duration<Real>(_propagation_time).count(); }
205 :
206 : /**
207 : * Whether or not to tolerate failure
208 : */
209 9 : bool tolerateFailure() const { return _tolerate_failure; }
210 :
211 : /**
212 : * Whether or not to bank Rays on completion
213 : */
214 810 : bool bankRaysOnCompletion() const { return _bank_rays_on_completion; }
215 :
216 : /**
217 : * Whether or not to use Ray dependent subdomain setup
218 : */
219 5370064 : bool rayDependentSubdomainSetup() const { return _ray_dependent_subdomain_setup; }
220 :
221 : /**
222 : * Register a value to be filled in the data on a Ray with a given name.
223 : * @param name The name to register
224 : * \return The Ray data index for the registered value
225 : *
226 : * Note that this does not actually allocate the storage. It is simply a registration system to
227 : * keep track of the values stored in the Ray data.
228 : */
229 : RayDataIndex registerRayData(const std::string & name);
230 : /**
231 : * Register values to be filled in the data on a Ray with a given name.
232 : * @param names The names to register
233 : * \return The Ray data indices for the registered value
234 : *
235 : * Note that this does not actually allocate the storage. It is simply a registration system to
236 : * keep track of the values stored in the Ray data.
237 : */
238 : std::vector<RayDataIndex> registerRayData(const std::vector<std::string> & names);
239 : /**
240 : * Gets the index associated with a registered value in the Ray data
241 : * @param name The value name to get the index of
242 : * @param graceful Whether or not to exit gracefully if none is found (return
243 : * INVALID_RAY_DATA_INDEX)
244 : * \return The index for the value
245 : */
246 : RayDataIndex getRayDataIndex(const std::string & name, const bool graceful = false) const;
247 : /**
248 : * Gets the indices associated with registered values in the Ray data
249 : * @param names The value names to get the indices of
250 : * @param graceful Whether or not to exit gracefully if none is found (index is filled with
251 : * INVALID_RAY_DATA_INDEX)
252 : * \return The indices for the values
253 : */
254 : std::vector<RayDataIndex> getRayDataIndices(const std::vector<std::string> & names,
255 : const bool graceful = false) const;
256 : /**
257 : * Gets the name associated with a registered value in the Ray data
258 : * @param index The index to get the name of
259 : * \return The name associated with index
260 : */
261 : const std::string & getRayDataName(const RayDataIndex index) const;
262 : /**
263 : * Gets the names associated with registered values in the Ray data
264 : * @param indices The indices to get the names of
265 : * \return The associated names
266 : */
267 : std::vector<std::string> getRayDataNames(const std::vector<RayDataIndex> & indices) const;
268 :
269 : /**
270 : * The registered size of values in the Ray data.
271 : */
272 : std::size_t rayDataSize() const { return _ray_data_names.size(); }
273 : /**
274 : * Whether or not any Ray data are registered
275 : */
276 : bool hasRayData() const { return _ray_data_names.size(); }
277 : /**
278 : * The Ray data names
279 : */
280 190 : const std::vector<std::string> & rayDataNames() const { return _ray_data_names; }
281 :
282 : /**
283 : * Register a value to be filled in the aux data on a Ray with a given name.
284 : * @param name The name to register
285 : * \returns The Ray aux data index for the registered value
286 : *
287 : * Note that this does not actually allocate the storage. It is simply a registration system to
288 : * keep track of the values stored in the Ray aux data.
289 : */
290 : RayDataIndex registerRayAuxData(const std::string & name);
291 : /**
292 : * Register values to be filled in the aux data on a Ray with a given name.
293 : * @param names The names to register
294 : * \returns The Ray aux data indices for the registered values
295 : *
296 : * Note that this does not actually allocate the storage. It is simply a registration system to
297 : * keep track of the values stored in the Ray aux data.
298 : */
299 : std::vector<RayDataIndex> registerRayAuxData(const std::vector<std::string> & names);
300 : /**
301 : * Gets the index associated with a registered value in the Ray aux data
302 : * @param name The value name to get the index of
303 : * @param graceful Whether or not to exit gracefully if none is found (return
304 : * INVALID_RAY_DATA_INDEX)
305 : * \returns The index for the value
306 : */
307 : RayDataIndex getRayAuxDataIndex(const std::string & name, const bool graceful = false) const;
308 : /**
309 : * Gets the indices associated with registered values in the Ray aux data
310 : * @param names The value names to get the indices of
311 : * @param graceful Whether or not to exit gracefully if none is found (index is filled with
312 : * INVALID_RAY_DATA_INDEX)
313 : * \return The indices for the values
314 : */
315 : std::vector<RayDataIndex> getRayAuxDataIndices(const std::vector<std::string> & names,
316 : const bool graceful = false) const;
317 : /**
318 : * Gets the name associated with a registered value in the Ray aux data
319 : * @param index The index to get the name of
320 : * \return The name associated with index
321 : */
322 : const std::string & getRayAuxDataName(const RayDataIndex index) const;
323 : /**
324 : * Gets the names associated with registered values in the Ray aux data
325 : * @param indices The indices to get the names of
326 : * \return The associated names
327 : */
328 : std::vector<std::string> getRayAuxDataNames(const std::vector<RayDataIndex> & indices) const;
329 : /**
330 : * The registered size of values in the Ray aux data.
331 : */
332 : std::size_t rayAuxDataSize() const { return _ray_aux_data_names.size(); }
333 : /**
334 : * Whether or not any Ray aux data are registered
335 : */
336 : bool hasRayAuxData() const { return _ray_aux_data_names.size(); }
337 : /**
338 : * The Ray aux data names
339 : */
340 82 : const std::vector<std::string> & rayAuxDataNames() const { return _ray_aux_data_names; }
341 :
342 : /**
343 : * Whether or not there are currently any active RayKernel objects
344 : */
345 : bool hasRayKernels(const THREAD_ID tid);
346 : /**
347 : * Fills the active RayKernels associated with this study and a block into result
348 : */
349 : void getRayKernels(std::vector<RayKernelBase *> & result, SubdomainID id, THREAD_ID tid);
350 : /**
351 : * Fills the active RayKernels associated with this study into result
352 : */
353 : template <typename T>
354 28717 : void getRayKernels(std::vector<T *> & result, THREAD_ID tid)
355 : {
356 28717 : _fe_problem.theWarehouse()
357 : .query()
358 57434 : .condition<AttribRayTracingStudy>(this)
359 28717 : .condition<AttribSystem>("RayKernel")
360 28717 : .condition<AttribThread>(tid)
361 : .queryInto(result);
362 28717 : }
363 : /**
364 : * Fills the active RayKernels associeted with this study, block, and potentially Ray into result
365 : */
366 : void
367 : getRayKernels(std::vector<RayKernelBase *> & result, SubdomainID id, THREAD_ID tid, RayID ray_id);
368 : /**
369 : * Fills the active RayBCs associated with this study and a boundary into result
370 : */
371 : void getRayBCs(std::vector<RayBoundaryConditionBase *> & result, BoundaryID id, THREAD_ID tid);
372 : /**
373 : * Fills the active RayBCs associated with this study and boundaries result
374 : */
375 : template <typename T>
376 86375 : void getRayBCs(std::vector<T *> & result, const std::vector<BoundaryID> & ids, THREAD_ID tid)
377 : {
378 86375 : _fe_problem.theWarehouse()
379 : .query()
380 172750 : .condition<AttribRayTracingStudy>(this)
381 86375 : .condition<AttribSystem>("RayBoundaryCondition")
382 86375 : .condition<AttribBoundaries>(ids)
383 86375 : .condition<AttribThread>(tid)
384 : .queryInto(result);
385 86375 : }
386 : /**
387 : * Fills the active RayBCs associated with this study into result
388 : */
389 : template <typename T>
390 3705 : void getRayBCs(std::vector<T *> & result, THREAD_ID tid)
391 : {
392 3705 : _fe_problem.theWarehouse()
393 : .query()
394 7410 : .condition<AttribRayTracingStudy>(this)
395 3705 : .condition<AttribSystem>("RayBoundaryCondition")
396 3705 : .condition<AttribThread>(tid)
397 : .queryInto(result);
398 3705 : }
399 : /**
400 : * Fills the active RayBCs associated with thie study, boundary elements, and potentially Ray into
401 : * result.
402 : *
403 : * This is purposely virtual because it allows derived studies to optimize the retrieval of RayBCs
404 : * during the trace in TraceRay.
405 : */
406 : virtual void getRayBCs(std::vector<RayBoundaryConditionBase *> & result,
407 : const std::vector<TraceRayBndElement> & bnd_elems,
408 : THREAD_ID tid,
409 : RayID ray_id);
410 :
411 : /**
412 : * Gets the current RayKernels for a thread, which are set in segmentSubdomainSetup()
413 : */
414 : const std::vector<RayKernelBase *> & currentRayKernels(THREAD_ID tid) const
415 : {
416 111920756 : return _threaded_current_ray_kernels[tid];
417 : }
418 :
419 : /**
420 : * Get the nodal bounding box for the domain
421 : */
422 4118378 : const BoundingBox & boundingBox() const { return _b_box; }
423 : /**
424 : * Get the loose nodal bounding box for the domain
425 : *
426 : * BoundingBox::contains_point does not have a tolerance, so we need to stretch the
427 : * box a little bit for some contains_point checks.
428 : */
429 13528053 : const BoundingBox & looseBoundingBox() const { return _loose_b_box; }
430 : /**
431 : * Get the inflated maximum length across the domain
432 : */
433 177110859 : Real domainMaxLength() const { return _domain_max_length; }
434 : /**
435 : * Get the current total volume of the domain
436 : */
437 12485 : Real totalVolume() const { return _total_volume; }
438 : /**
439 : * Whether or not the domain is rectangular (if it is prefectly encompassed by its bounding box)
440 : */
441 : bool isRectangularDomain() const;
442 :
443 : /**
444 : * Whether or not the local mesh has internal sidesets that have RayBCs on them
445 : *
446 : * NOTE: if useInternalSidesets() == false, this will be false even if the mesh does
447 : * have internal sidesets
448 : */
449 : bool hasInternalSidesets() const { return _internal_sidesets.size(); }
450 : /**
451 : * Get the internal sidesets (that have RayBC(s)) for each side for a given element
452 : *
453 : * This will be empty if the elem does not have any internal sidesets that have RayBC(s)
454 : */
455 : const std::vector<std::vector<BoundaryID>> & getInternalSidesets(const Elem * elem) const;
456 : /**
457 : * Gets the internal sidesets (that have RayBCs) within the local domain
458 : */
459 : const std::set<BoundaryID> & getInternalSidesets() const { return _internal_sidesets; }
460 : /**
461 : * Whether or not the side \s on elem \p elem is non-planar
462 : *
463 : * This is cached because checking whether or not a face is planar is costly
464 : */
465 14380501 : bool sideIsNonPlanar(const Elem * elem, const unsigned short s) const
466 : {
467 14380501 : return _has_non_planar_sides && _non_planar_sides[_elem_index_helper.getIndex(elem)][s] == 1;
468 : }
469 : /**
470 : * Whether or not the mesh has active elements of the same level
471 : *
472 : * Use this over sameLevelActiveElems(), which is for internally setting
473 : * _has_same_level_active_elems
474 : */
475 63749280 : bool hasSameLevelActiveElems() const { return _has_same_level_active_elems; }
476 :
477 : /**
478 : * INTERNAL methods for acquiring a Ray during a trace in RayKernels and RayBCs.
479 : *
480 : * You should not use these APIs directly. If you wish to acquire a Ray during generation
481 : * during generateRays()), use the protected RayTracingStudy::acquire{}Ray() methods.
482 : * If you wish to acquire a Ray during propagation in RayKernels and RayBC, use the
483 : * protected RayKernelBase::acquireRay() and RayBoundaryConditionBase::acquireRay(),
484 : * respectively.
485 : */
486 : ///@{
487 : std::shared_ptr<Ray> acquireRayDuringTrace(const THREAD_ID tid,
488 : const AcquireMoveDuringTraceKey &);
489 3925497 : std::shared_ptr<Ray> acquireRayInternal(const RayID id,
490 : const std::size_t data_size,
491 : const std::size_t aux_data_size,
492 : const bool reset,
493 : const AcquireRayInternalKey &)
494 : {
495 3925497 : return _parallel_ray_study->acquireParallelData(
496 7850994 : 0, this, id, data_size, aux_data_size, reset, Ray::ConstructRayKey());
497 : }
498 : ///@}
499 :
500 : /**
501 : * INTERNAL method for moving a Ray into the buffer during tracing.
502 : *
503 : * You should not and cannot use this API. It is protected by the
504 : * AcquireMoveDuringTraceKey, which is only constructable by
505 : * RayKernelBase and RayBoundaryConditionBase.
506 : */
507 : void moveRayToBufferDuringTrace(std::shared_ptr<Ray> & ray,
508 : const THREAD_ID tid,
509 : const AcquireMoveDuringTraceKey &);
510 : /**
511 : * Access to the libMesh MeshBase.
512 : *
513 : * This is needed for unpack routines for a Ray, which has a context of this study.
514 : */
515 5202717 : MeshBase & meshBase() const { return _mesh; }
516 :
517 : /**
518 : * @returns A reference to the MooseMesh associated with this study.
519 : */
520 3614 : MooseMesh & mesh() { return _mesh; }
521 :
522 : /**
523 : * Get the outward normal for a given element side.
524 : */
525 : virtual const Point &
526 : getSideNormal(const Elem * elem, const unsigned short side, const THREAD_ID tid);
527 : /**
528 : * Gets the outward normals for a given element.
529 : * Returns a pointer to the normal for the zeroth side.
530 : */
531 2 : virtual const Point * getElemNormals(const Elem * /* elem */, const THREAD_ID /* tid */)
532 : {
533 2 : mooseError("Unimplemented element normal caching in ", type(), "::getElemNormals()");
534 : }
535 :
536 : /**
537 : * Gets the data value for a banked ray with a given ID
538 : *
539 : * This will return the value replicated across all processors
540 : */
541 : RayData getBankedRayData(const RayID ray_id, const RayDataIndex index) const;
542 : /**
543 : * Gets the data value for a banked ray with a given ID
544 : *
545 : * This will return the value replicated across all processors
546 : */
547 : RayData getBankedRayAuxData(const RayID ray_id, const RayDataIndex index) const;
548 :
549 : /**
550 : * Gets the ID of a registered ray
551 : * @param name The name of said ray
552 : * @param graceful Whether or not to exit gracefully if none is found (with invalid_id)
553 : */
554 : RayID registeredRayID(const std::string & name, const bool graceful = false) const;
555 : /**
556 : * Gets the name of a registered ray
557 : * @param ray_id The ID of said ray
558 : */
559 : const std::string & registeredRayName(const RayID ray_id) const;
560 :
561 : /**
562 : * Whether or not ray registration is being used
563 : */
564 4501 : bool useRayRegistration() const { return _use_ray_registration; }
565 :
566 : /**
567 : * Whether or not to store the Ray data on the cached Ray traces
568 : */
569 72358 : bool dataOnCacheTraces() const { return _data_on_cache_traces; }
570 : /**
571 : * Whether or not to store the Ray aux data on the cached Ray traces
572 : */
573 72268 : bool auxDataOnCacheTraces() const { return _aux_data_on_cache_traces; }
574 : /**
575 : * Whether or not to cache individual element segments when _cache_traces = true
576 : */
577 62181 : bool segmentsOnCacheTraces() const { return _segments_on_cache_traces; }
578 :
579 : /**
580 : * Virtual that allows for selection in if a Ray should be cached or not (only used when
581 : * _cache_traces).
582 : */
583 5380196 : virtual bool shouldCacheTrace(const std::shared_ptr<Ray> & /* ray */) const
584 : {
585 5380196 : return _always_cache_traces;
586 : }
587 :
588 : /**
589 : * Initialize a Ray in the threaded cached trace map to be filled with segments
590 : */
591 : TraceData & initThreadedCachedTrace(const std::shared_ptr<Ray> & ray, THREAD_ID tid);
592 : /**
593 : * Get the cached trace data structure
594 : */
595 : const std::vector<TraceData> & getCachedTraces() const { return _cached_traces; }
596 :
597 : /**
598 : * Get the cached hmax for all elements in a subdomain
599 : *
600 : * Used for scaling tolerances in ray tracing.
601 : */
602 : Real subdomainHmax(const SubdomainID subdomain_id) const;
603 :
604 : /**
605 : * Entry point for acting on a ray when it is completed (shouldContinue() == false)
606 : */
607 : virtual void onCompleteRay(const std::shared_ptr<Ray> & ray);
608 :
609 : /**
610 : * Verifies that the Rays in the given range have unique Ray IDs.
611 : *
612 : * @param begin The beginning of the range of Rays to check
613 : * @param end The end of the range of Rays to check
614 : * @param global If true, this will complete the verification across all processors
615 : * @param error_suffix Entry point for additional information in the error message
616 : */
617 : void verifyUniqueRayIDs(const std::vector<std::shared_ptr<Ray>>::const_iterator begin,
618 : const std::vector<std::shared_ptr<Ray>>::const_iterator end,
619 : const bool global,
620 : const std::string & error_suffix) const;
621 :
622 : /**
623 : * Verifies that the Rays in the given range are unique. That is, that there are not multiple
624 : * shared_ptrs that point to the same Ray
625 : *
626 : * @param begin The beginning of the range of Rays to check
627 : * @param end The end of the range of Rays to check
628 : * @param error_suffix Entry point for additional information in the error message
629 : */
630 : void verifyUniqueRays(const std::vector<std::shared_ptr<Ray>>::const_iterator begin,
631 : const std::vector<std::shared_ptr<Ray>>::const_iterator end,
632 : const std::string & error_suffix);
633 :
634 : /**
635 : * Whether or not the study is propagating (tracing Rays)
636 : */
637 : bool currentlyPropagating() const { return _parallel_ray_study->currentlyExecuting(); }
638 : /**
639 : * Whether or not the study is generating
640 : */
641 : bool currentlyGenerating() const { return _parallel_ray_study->currentlyPreExecuting(); }
642 :
643 : /**
644 : * Gets the threaded TraceRay object for \p tid.
645 : */
646 : ///@{
647 44 : TraceRay & traceRay(const THREAD_ID tid) { return *_threaded_trace_ray[tid]; }
648 8443 : const TraceRay & traceRay(const THREAD_ID tid) const { return *_threaded_trace_ray[tid]; }
649 : ///@}
650 :
651 : /**
652 : * Whether or not to verify if Rays have valid information before being traced
653 : */
654 18926369 : bool verifyRays() const { return _verify_rays; }
655 : /**
656 : * Whether or not trace verification is enabled in devel/dbg modes
657 : */
658 : #ifndef NDEBUG
659 : bool verifyTraceIntersections() const { return _verify_trace_intersections; }
660 : #endif
661 :
662 : /**
663 : * Whether or not \p side is incoming on element \p elem in direction \p direction.
664 : */
665 : bool sideIsIncoming(const Elem * const elem,
666 : const unsigned short side,
667 : const Point & direction,
668 : const THREAD_ID tid);
669 :
670 : /**
671 : * Whether or not to produce a warning when interacting with a non-planar mesh.
672 : */
673 : bool warnNonPlanar() const { return _warn_non_planar; }
674 :
675 : /**
676 : * The underlying parallel study: used for the context for calling the packed range routines.
677 : */
678 : ParallelStudy<std::shared_ptr<Ray>, Ray> * parallelStudy() { return _parallel_ray_study.get(); }
679 :
680 : /**
681 : * Get an element's side pointer without excessive memory allocation
682 : *
683 : * @param elem The element to build a side for
684 : * @param s The side to build
685 : * @param tid The thread id
686 : * @return A pointer to the side element
687 : */
688 : const libMesh::Elem &
689 : elemSide(const libMesh::Elem & elem, const unsigned int s, const THREAD_ID tid = 0)
690 : {
691 3092126 : return _threaded_elem_side_builders[tid](elem, s);
692 : }
693 :
694 : protected:
695 : /**
696 : * Subclasses should override this to determine how to generate Rays.
697 : * This will be called within execute() and makes up the "generation phase"
698 : * of the algorithm.
699 : */
700 : virtual void generateRays() = 0;
701 :
702 : /**
703 : * Entry point before study execution
704 : */
705 8023 : virtual void preExecuteStudy() {}
706 : /**
707 : * Entry point after study execution
708 : */
709 7567 : virtual void postExecuteStudy() {}
710 :
711 : /**
712 : * Helper function for computing the total domain volume
713 : */
714 : Real computeTotalVolume();
715 :
716 : /**
717 : * Gets the writeable current RayKernels for a thread
718 : *
719 : * Allows for other ray studies to fill the current ray kernels in a custom manner
720 : */
721 : std::vector<RayKernelBase *> & currentRayKernelsWrite(THREAD_ID tid)
722 : {
723 : return _threaded_current_ray_kernels[tid];
724 : }
725 :
726 : /**
727 : * Reserve \p size entires in the Ray buffer.
728 : *
729 : * This can only be used within generateRays() and should be used when possible
730 : * to reserve entires in the buffer before adding Rays via moveRay(s)ToBuffer.
731 : */
732 : void reserveRayBuffer(const std::size_t size);
733 :
734 : /**
735 : * Determine whether or not the mesh currently has active elements that are
736 : * all the same level
737 : */
738 : bool sameLevelActiveElems() const;
739 :
740 : /**
741 : * Builds quadrature points for a given segment using the _segment_qrule
742 : * @param start Start point of the segment
743 : * @param end End point of the segment
744 : * @param length The lengh of the start -> end segment
745 : * @param points Points to fill into (should be sized ahead of time)
746 : * @param weights Weights to fill into (should be sized ahead of time)
747 : */
748 : virtual void buildSegmentQuadrature(const Point & start,
749 : const Point & end,
750 : const Real length,
751 : std::vector<Point> & points,
752 : std::vector<Real> & weights) const;
753 :
754 : /**
755 : * Get the Ray bank. This is the bank of Rays that have completed on this processor
756 : * after an execution of the study.
757 : *
758 : * This is only available when the private parameter _bank_rays_on_completion
759 : * is set to true.
760 : */
761 : const std::vector<std::shared_ptr<Ray>> & rayBank() const;
762 :
763 : /**
764 : * Gets the Ray with the ID \p ray_id from the Ray bank.
765 : *
766 : * If the Ray with \p ray_id is not found across all processors, this will error.
767 : *
768 : * This will ONLY return a valid Ray (not a null shared_ptr) on the processor that
769 : * has the Ray.
770 : */
771 : std::shared_ptr<Ray> getBankedRay(const RayID ray_id) const;
772 :
773 : /**
774 : * Resets the generation of unique RayIDs via generateUniqueRayID() to the beginning
775 : * of the range.
776 : */
777 : void resetUniqueRayIDs();
778 : /**
779 : * Resets the generation of unique replicated RayIDs accessed via generateReplicatedRayID().
780 : */
781 : void resetReplicatedRayIDs();
782 :
783 : /**
784 : * Gets all of the currently active RayTracingObjects
785 : */
786 : std::vector<RayTracingObject *> getRayTracingObjects();
787 :
788 : /**
789 : * Generates a unique RayID to be used for a Ray.
790 : *
791 : * This is used internally when acquiring new Rays.
792 : */
793 : virtual RayID generateUniqueRayID(const THREAD_ID tid);
794 : /**
795 : * Generates a Ray ID that is replicated across all processors.
796 : */
797 : RayID generateReplicatedRayID();
798 :
799 : /**
800 : * User APIs for constructing Rays within the RayTracingStudy.
801 : *
802 : * Rays can ONLY be constructed by users within the RayTracingStudy via the following methods.
803 : */
804 : ///@{
805 : /**
806 : * Acquire a Ray from the pool of Rays within generateRays().
807 : *
808 : * A unique ID is generated and assigned to the acquired Ray. The data and aux
809 : * data sizes are set according to the sizes required by the RayTracingStudy.
810 : */
811 : std::shared_ptr<Ray> acquireRay();
812 : /**
813 : * Acquire a Ray from the pool of Rays within generateRays(), without resizing
814 : * the data (sizes the data to zero).
815 : *
816 : * A unique ID is generated and assigned to the acquired Ray.
817 : */
818 : std::shared_ptr<Ray> acquireUnsizedRay();
819 : /**
820 : * Acquire a Ray from the pool of Rays within generateRays() in a replicated fashion.
821 : *
822 : * That is, this method must be called on all processors at the same time and the ID of the
823 : * resulting Ray is the same across all processors.
824 : *
825 : * The data and aux data sizes are set according to the sizes required by the RayTracingStudy.
826 : */
827 : std::shared_ptr<Ray> acquireReplicatedRay();
828 : /**
829 : * Acquires a Ray that that is copied from another Ray within generateRays().
830 : *
831 : * All of the information is copied except for the counters (intersections, processor crossings,
832 : * etc), which are reset.
833 : */
834 : std::shared_ptr<Ray> acquireCopiedRay(const Ray & ray);
835 : /**
836 : * Acquires a Ray with a given name within generateRays(). Used when ray registration is enabled,
837 : * that is, the private paramater '_use_ray_registration' == true.
838 : *
839 : * This method must be called on all processors at the same time with the same name.
840 : * This method can only be called on thread 0, which is why there is no thread argument.
841 : */
842 : std::shared_ptr<Ray> acquireRegisteredRay(const std::string & name);
843 : ///@}
844 :
845 : /**
846 : * Moves a ray to the buffer to be traced during generateRays().
847 : *
848 : * This moves the Ray into the buffer (with std::move), therefore \p ray will be nullptr after
849 : * this call.
850 : */
851 : void moveRayToBuffer(std::shared_ptr<Ray> & ray);
852 : /**
853 : * Moves rays to the buffer to be traced during generateRays().
854 : *
855 : * This moves the rays into the buffer (with std::move), therefore all valid rays
856 : * in \p rays will be nullptr after this call.
857 : */
858 : void moveRaysToBuffer(std::vector<std::shared_ptr<Ray>> & rays);
859 :
860 : /// The Mesh
861 : MooseMesh & _mesh;
862 : /// The Communicator
863 : const Parallel::Communicator & _comm;
864 : /// The rank of this processor (this actually takes time to lookup - so just do it once)
865 : const processor_id_type _pid;
866 :
867 : /// Whether or not to perform coverage checks on RayKernels
868 : const bool _ray_kernel_coverage_check;
869 : /// Whether not to warn if non-planar faces are found
870 : const bool _warn_non_planar;
871 : /// Whether or not to use Ray registration
872 : const bool _use_ray_registration;
873 : /// Whether or not to use the internal sidesets in ray tracing
874 : const bool _use_internal_sidesets;
875 : /// Whether or not to tolerate a Ray Tracing failure
876 : const bool _tolerate_failure;
877 : /// Whether or not to bank rays on completion
878 : const bool _bank_rays_on_completion;
879 : /// Whether or not subdomain setup is dependent on the Ray
880 : const bool _ray_dependent_subdomain_setup;
881 :
882 : /// Whether or not to cache traces on every trace execution
883 : bool _always_cache_traces;
884 : /// Whether or not to store the Ray data on the cache traces
885 : const bool _data_on_cache_traces;
886 : /// Whether or not to store the Ray aux data on the cache traces
887 : const bool _aux_data_on_cache_traces;
888 : /// Whether or not to cache individual element segments when caching
889 : const bool _segments_on_cache_traces;
890 : /// Max distance a Ray can travel before being killed (can change)
891 : const Real _ray_max_distance;
892 : /// Whether or not to verify if Rays have valid information before being traced
893 : const bool _verify_rays;
894 : /// Whether or not to verify the trace intersections in devel and dbg modes
895 : #ifndef NDEBUG
896 : const bool _verify_trace_intersections;
897 : #endif
898 :
899 : private:
900 : /**
901 : * Perform coverage checks (coverage of RayMaterials and RayKernels, if enabled)
902 : */
903 : void coverageChecks();
904 :
905 : /**
906 : * Perform checks to see if the listed dependencies in the RayTracingObjects exist
907 : */
908 : void dependencyChecks();
909 : /**
910 : * Verifies that the dependencies exist for a set of RayTracingObjects
911 : */
912 : void verifyDependenciesExist(const std::vector<RayTracingObject *> & rtos);
913 :
914 : /**
915 : * Check for if all of the element types in the mesh are supported by ray tracing
916 : */
917 : void traceableMeshChecks();
918 :
919 : /**
920 : * Does the setup for internal sidesets. This includes:
921 : * - Setting if the local mesh has internal sidesets that have RayBCs and on which
922 : * boundary (stored in _internal_sidesets)
923 : * - Sets up the internal sideset mapping that maps element sides to the BoundaryIDs
924 : * on said side that have internal sidesets that have RayBCs
925 : * (stored in _internal_sidesets_map)
926 : *
927 : * If _use_internal_sidesets == false, this will still check to make sure that the
928 : * user does not have internal sidesets with RayBCs on them, and will report an error
929 : * if so.
930 : */
931 : void internalSidesetSetup();
932 :
933 : /**
934 : * Sets up the caching of whether or not each element side is non-planar, which is
935 : * stored in _non_planar_sides.
936 : *
937 : * This is done because this check is made within TraceRay, and determining whether
938 : * or not a side is non planar is quite expensive and we don't want to call it
939 : * mid-trace.
940 : */
941 : void nonPlanarSideSetup();
942 :
943 : /**
944 : * Sets up the _elem_index_helper, which is used for obtaining a contiguous index
945 : * for all elements that this processor knows about.
946 : *
947 : * Said index is used for quickly accessing things like internal sidesets and checking
948 : * if an element side is non-planar during tracing.
949 : */
950 : void localElemIndexSetup();
951 :
952 : /**
953 : * Sets up the maps from Ray to associated RayTracingObjects if _use_ray_registration.
954 : *
955 : * If ray registration is disabled, this makes sure no RayTracingObjects provided rays.
956 : */
957 : void registeredRaySetup();
958 :
959 : /**
960 : * Zero the AuxVariables that the registered AuxRayKernels contribute to.
961 : */
962 : void zeroAuxVariables();
963 :
964 : /**
965 : * Caches the hmax for all elements in each subdomain
966 : */
967 : void subdomainHMaxSetup();
968 :
969 : /**
970 : * Internal method for registering Ray data or Ray aux data with a name.
971 : */
972 : RayDataIndex registerRayDataInternal(const std::string & name, const bool aux);
973 : /**
974 : * Internal method for registering Ray data or Ray aux data with names.
975 : */
976 : std::vector<RayDataIndex> registerRayDataInternal(const std::vector<std::string> & names,
977 : const bool aux);
978 : /**
979 : * Internal method for getting the index of Ray data or Ray aux data.
980 : */
981 : RayDataIndex
982 : getRayDataIndexInternal(const std::string & name, const bool aux, const bool graceful) const;
983 : /**
984 : * Internal method for getting the indicies of Ray data or Ray aux data.
985 : */
986 : std::vector<RayDataIndex> getRayDataIndicesInternal(const std::vector<std::string> & names,
987 : const bool aux,
988 : const bool graceful) const;
989 :
990 : /**
991 : * Internal method for getting the name of Ray data or Ray aux data.
992 : */
993 : const std::string & getRayDataNameInternal(const RayDataIndex index, const bool aux) const;
994 :
995 : /**
996 : * Internal method for getting the value (replicated across all processors) in a Ray's data
997 : * or aux data from the Ray banks.
998 : */
999 : RayData
1000 : getBankedRayDataInternal(const RayID ray_id, const RayDataIndex index, const bool aux) const;
1001 :
1002 : /**
1003 : * Registers a Ray with a given name.
1004 : *
1005 : * This is for internal use only and is an internal map for replicated Ray IDs -> names for
1006 : * easy access for studies that have a small number of rays to generate.
1007 : *
1008 : * @param name The name to register
1009 : * @return The allocated ID for the registered ray
1010 : */
1011 : RayID registerRay(const std::string & name);
1012 :
1013 : /// Threaded helpers for building element sides without extraneous allocation
1014 : std::vector<libMesh::ElemSideBuilder> _threaded_elem_side_builders;
1015 :
1016 : /// Timing
1017 : //@{
1018 : std::chrono::steady_clock::time_point _execution_start_time;
1019 : std::chrono::steady_clock::duration _execution_time;
1020 : std::chrono::steady_clock::duration _generation_time;
1021 : std::chrono::steady_clock::duration _propagation_time;
1022 : //@}
1023 :
1024 : /// The map from Ray data names to index
1025 : std::unordered_map<std::string, RayDataIndex> _ray_data_map;
1026 : /// The map from Ray aux data names to index
1027 : std::unordered_map<std::string, RayDataIndex> _ray_aux_data_map;
1028 :
1029 : /// The names for each Ray data entry
1030 : std::vector<std::string> _ray_data_names;
1031 : /// The names for each Ray aux data entry
1032 : std::vector<std::string> _ray_aux_data_names;
1033 :
1034 : /// Map from registered Ray name to ID
1035 : std::unordered_map<std::string, RayID> & _registered_ray_map;
1036 : /// Map from registered Ray ID to name
1037 : std::vector<std::string> & _reverse_registered_ray_map;
1038 :
1039 : /// Storage for the cached traces
1040 : std::vector<TraceData> _cached_traces;
1041 : /// The threaded storage for cached traces
1042 : std::vector<std::vector<TraceData>> _threaded_cached_traces;
1043 :
1044 : /// Number of currently cached objects for Jacobian/residual for each thread
1045 : std::vector<std::size_t> _num_cached;
1046 :
1047 : /// The BoundaryIDs on the local mesh that have internal RayBCs
1048 : std::set<BoundaryID> _internal_sidesets;
1049 : /// Internal sideset data, if internal sidesets exist (indexed with getLocalElemIndex())
1050 : std::vector<std::vector<std::vector<BoundaryID>>> _internal_sidesets_map;
1051 : /// Whether or not the local mesh has elements with non-planar sides
1052 : bool _has_non_planar_sides;
1053 : /// Non planar side data, which is for quick checking if an elem side is non-planar
1054 : /// We use unsigned short here to avoid a std::vector<bool>; 0 = false, otherwise true
1055 : std::vector<std::vector<unsigned short>> _non_planar_sides;
1056 : /// Whether or not the mesh has active elements of the same level
1057 : bool _has_same_level_active_elems;
1058 :
1059 : /// Nodal bounding box for the domain
1060 : libMesh::BoundingBox _b_box;
1061 : /// Loose nodal bounding box for the domain
1062 : libMesh::BoundingBox _loose_b_box;
1063 : /// An inflated max distance for the domain
1064 : Real _domain_max_length;
1065 : /// The total volume of the domain
1066 : Real _total_volume;
1067 :
1068 : /// Threaded cached subdomain query for RayKernelBase objects pertaining to this study
1069 : std::vector<TheWarehouse::QueryCache<AttribSubdomains>> _threaded_cache_ray_kernel;
1070 : /// Threaded cached boundary query for RayBC objects pertaining to this study
1071 : std::vector<TheWarehouse::QueryCache<AttribBoundaries>> _threaded_cache_ray_bc;
1072 : /// Threaded storage for all of the RayTracingObjects associated with a single Ray
1073 : std::vector<std::vector<std::set<const RayTracingObject *>>> _threaded_ray_object_registration;
1074 : /// The current RayKernel objects for each thread
1075 : std::vector<std::vector<RayKernelBase *>> _threaded_current_ray_kernels;
1076 : /// The TraceRay objects for each thread (they do the physical tracing)
1077 : std::vector<std::shared_ptr<TraceRay>> _threaded_trace_ray;
1078 : /// Face FE used for computing face normals for each thread
1079 : std::vector<std::unique_ptr<libMesh::FEBase>> _threaded_fe_face;
1080 : /// Face quadrature used for computing face normals for each thread
1081 : std::vector<std::unique_ptr<libMesh::QBase>> _threaded_q_face;
1082 : /// Threaded cache for side normals that have been computed already during tracing
1083 : std::vector<std::unordered_map<std::pair<const Elem *, unsigned short>, Point>>
1084 : _threaded_cached_normals;
1085 : /// Cumulative Ray bank - stored only when _bank_rays_on_completion
1086 : std::vector<std::shared_ptr<Ray>> _ray_bank;
1087 : /// Storage for the next available unique RayID, obtained via generateUniqueRayID()
1088 : std::vector<RayID> _threaded_next_ray_id;
1089 : /// Storage for the next available replicated RayID, obtained via generateReplicatedRayID()
1090 : RayID _replicated_next_ray_id;
1091 :
1092 : /// The study that used is to actually execute (trace) the Rays
1093 : const std::unique_ptr<ParallelRayStudy> _parallel_ray_study;
1094 :
1095 : /// Quadrature rule for laying points across a 1D ray segment
1096 : std::unique_ptr<libMesh::QBase> _segment_qrule;
1097 :
1098 : /// Total number of processor crossings for Rays that finished on this processor
1099 : unsigned long long int _ending_processor_crossings;
1100 : /// Max number of total processor crossings for Rays that finished on this processor
1101 : unsigned int _ending_max_processor_crossings;
1102 : /// Total number of processor crossings
1103 : unsigned long long int _total_processor_crossings;
1104 : /// Max number of processor crossings for all Rays
1105 : unsigned int _max_processor_crossings;
1106 :
1107 : /// Total number of Ray/element intersections for Rays that finished on this processor
1108 : unsigned long long int _ending_intersections;
1109 : /// Max number of intersections for Rays that finished on this processor
1110 : unsigned int _ending_max_intersections;
1111 : /// Max number of trajectory changes for Rays that finished on this processor
1112 : unsigned int _ending_max_trajectory_changes;
1113 : /// Total number of Ray/element intersections
1114 : unsigned long long int _total_intersections;
1115 : /// Max number of intersections for a single Ray
1116 : unsigned int _max_intersections;
1117 : /// Max number of trajectory changes for a single Ray
1118 : unsigned int _max_trajectory_changes;
1119 :
1120 : /// Total distance traveled by Rays that end on this processor
1121 : Real _ending_distance;
1122 : /// Total distance traveled by all Rays
1123 : Real _total_distance;
1124 :
1125 : /// Cumulative results on this processor from the threaded TraceRay objects
1126 : std::vector<unsigned long long int> _local_trace_ray_results;
1127 :
1128 : /// The cached hmax for all elements in a subdomain
1129 : std::unordered_map<SubdomainID, Real> _subdomain_hmax;
1130 :
1131 : /// Whether or not we've called initial setup - used to stop from late registration
1132 : bool _called_initial_setup;
1133 :
1134 : /// Helper for defining a local contiguous index for each element
1135 : ElemIndexHelper _elem_index_helper;
1136 :
1137 : /// Spin mutex object for locks
1138 : mutable Threads::spin_mutex _spin_mutex;
1139 : };
|
