streamLineBase.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015 OpenFOAM Foundation
9 Copyright (C) 2015-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "streamLineBase.H"
30#include "fvMesh.H"
31#include "ReadFields.H"
32#include "OFstream.H"
33#include "sampledSet.H"
34#include "globalIndex.H"
35#include "mapDistribute.H"
37#include "wallPolyPatch.H"
39#include "mapPolyMesh.H"
40
41// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42
43namespace Foam
44{
45namespace functionObjects
46{
48}
49}
50
51
52const Foam::Enum
53<
55>
57({
58 { trackDirType::FORWARD, "forward" },
59 { trackDirType::BACKWARD, "backward" },
60 { trackDirType::BIDIRECTIONAL, "bidirectional" }
61});
62
63
64// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
65
66const Foam::word&
68{
69 if (!sampledSetPtr_)
70 {
72 }
73
74 return sampledSetAxis_;
75}
76
77
80{
81 if (!sampledSetPtr_)
82 {
83 sampledSetPtr_ = sampledSet::New
84 (
85 "seedSampleSet",
86 mesh_,
88 dict_.subDict("seedSampleSet")
89 );
90
91 sampledSetAxis_ = sampledSetPtr_->axis();
92 }
93
94 return *sampledSetPtr_;
95}
96
97
100{
101 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
102
103 label nFaces = 0;
104
105 for (const polyPatch& pp : patches)
106 {
107 //if (!polyPatch::constraintType(pp.type()))
108 if (isA<wallPolyPatch>(pp))
109 {
110 nFaces += pp.size();
111 }
112 }
113
114 labelList addressing(nFaces);
115
116 nFaces = 0;
117
118 for (const polyPatch& pp : patches)
119 {
120 //if (!polyPatch::constraintType(pp.type()))
121 if (isA<wallPolyPatch>(pp))
122 {
123 forAll(pp, i)
124 {
125 addressing[nFaces++] = pp.start()+i;
126 }
127 }
128 }
129
131 (
133 (
134 mesh_.faces(),
135 addressing
136 ),
137 mesh_.points()
138 );
139}
140
141
143(
144 const label nSeeds,
145 label& UIndex,
150)
151{
152 // Read fields
153 label nScalar = 0;
154 label nVector = 0;
155
156 for (const word& fieldName : fields_)
157 {
158 if (foundObject<volScalarField>(fieldName))
159 {
160 ++nScalar;
161 }
162 else if (foundObject<volVectorField>(fieldName))
163 {
164 ++nVector;
165 }
166 else
167 {
169 << "Cannot find field " << fieldName << nl
170 << "Valid scalar fields: "
171 << flatOutput(mesh_.sortedNames<volScalarField>()) << nl
172 << "Valid vector fields: "
173 << flatOutput(mesh_.sortedNames<volVectorField>()) << nl
174 << exit(FatalError);
175 }
176 }
177 vsInterp.setSize(nScalar);
178 nScalar = 0;
179 vvInterp.setSize(nVector);
180 nVector = 0;
181
182 for (const word& fieldName : fields_)
183 {
184 if (foundObject<volScalarField>(fieldName))
185 {
186 const volScalarField& f = lookupObject<volScalarField>(fieldName);
187 vsInterp.set
188 (
189 nScalar++,
191 (
192 interpolationScheme_,
193 f
194 )
195 );
196 }
197 else if (foundObject<volVectorField>(fieldName))
198 {
199 const volVectorField& f = lookupObject<volVectorField>(fieldName);
200
201 if (f.name() == UName_)
202 {
203 UIndex = nVector;
204 }
205
206 vvInterp.set
207 (
208 nVector++,
210 (
211 interpolationScheme_,
212 f
213 )
214 );
215 }
216 }
217
218 // Store the names
219 scalarNames_.setSize(vsInterp.size());
220 forAll(vsInterp, i)
221 {
222 scalarNames_[i] = vsInterp[i].psi().name();
223 }
224 vectorNames_.setSize(vvInterp.size());
225 forAll(vvInterp, i)
226 {
227 vectorNames_[i] = vvInterp[i].psi().name();
228 }
229
230 // Check that we know the index of U in the interpolators.
231
232 if (UIndex == -1)
233 {
235 << "Cannot find field to move particles with : " << UName_ << nl
236 << "This field has to be present in the sampled fields " << fields_
237 << " and in the objectRegistry."
238 << exit(FatalError);
239 }
240
241 // Sampled data
242 // ~~~~~~~~~~~~
243
244 // Size to maximum expected sizes.
245 allTracks_.clear();
246 allTracks_.setCapacity(nSeeds);
247 allScalars_.setSize(vsInterp.size());
248 forAll(allScalars_, i)
249 {
250 allScalars_[i].clear();
251 allScalars_[i].setCapacity(nSeeds);
252 }
253 allVectors_.setSize(vvInterp.size());
254 forAll(allVectors_, i)
255 {
256 allVectors_[i].clear();
257 allVectors_[i].setCapacity(nSeeds);
258 }
259}
260
261
263(
264 const label tracki,
265
266 const scalar w,
267 const label lefti,
268 const label righti,
269
270 DynamicList<point>& newTrack,
271 DynamicList<scalarList>& newScalars,
272 DynamicList<vectorList>& newVectors
273) const
274{
275 const label sz = newTrack.size();
276
277 const List<point>& track = allTracks_[tracki];
278
279 newTrack.append((1.0-w)*track[lefti] + w*track[righti]);
280
281 // Scalars
282 {
283 newScalars.append(scalarList(allScalars_.size()));
284 scalarList& newVals = newScalars[sz];
285
286 forAll(allScalars_, scalari)
287 {
288 const scalarList& trackVals = allScalars_[scalari][tracki];
289 newVals[scalari] = (1.0-w)*trackVals[lefti] + w*trackVals[righti];
290 }
291 }
292
293 // Vectors
294 {
295 newVectors.append(vectorList(allVectors_.size()));
296 vectorList& newVals = newVectors[sz];
297
298 forAll(allVectors_, vectori)
299 {
300 const vectorList& trackVals = allVectors_[vectori][tracki];
301 newVals[vectori] = (1.0-w)*trackVals[lefti] + w*trackVals[righti];
302 }
303 }
304}
305
306
307// Can split a track into multiple tracks
309(
310 const treeBoundBox& bb,
311 const label tracki,
312 PtrList<DynamicList<point>>& newTracks,
313 PtrList<DynamicList<scalarList>>& newScalars,
315) const
316{
317 const List<point>& track = allTracks_[tracki];
318
319 if (track.size())
320 {
321 for
322 (
323 label segmenti = 1;
324 segmenti < track.size();
325 segmenti++
326 )
327 {
328 const point& startPt = track[segmenti-1];
329 const point& endPt = track[segmenti];
330
331 const vector d(endPt-startPt);
332 const scalar magD = mag(d);
333 if (magD > ROOTVSMALL)
334 {
335 if (bb.contains(startPt))
336 {
337 // Store 1.0*track[segmenti-1]+0*track[segmenti]
338 storePoint
339 (
340 tracki,
341
342 0.0,
343 segmenti-1,
344 segmenti,
345
346 newTracks.last(),
347 newScalars.last(),
348 newVectors.last()
349 );
350
351 if (!bb.contains(endPt))
352 {
353 point clipPt;
354 if (bb.intersects(endPt, startPt, clipPt))
355 {
356 // End of track. Store point and interpolated
357 // values
358 storePoint
359 (
360 tracki,
361
362 mag(clipPt-startPt)/magD,
363 segmenti-1,
364 segmenti,
365
366 newTracks.last(),
367 newScalars.last(),
368 newVectors.last()
369 );
370
371 newTracks.last().shrink();
372 newScalars.last().shrink();
373 newVectors.last().shrink();
374 }
375 }
376 }
377 else
378 {
379 // startPt outside box. New track. Get starting point
380
381 point clipPt;
382 if (bb.intersects(startPt, endPt, clipPt))
383 {
384 // New track
385 newTracks.append
386 (
387 new DynamicList<point>(track.size()/10)
388 );
389 newScalars.append
390 (
391 new DynamicList<scalarList>(track.size()/10)
392 );
393 newVectors.append
394 (
395 new DynamicList<vectorList>(track.size()/10)
396 );
397
398 // Store point and interpolated values
399 storePoint
400 (
401 tracki,
402
403 mag(clipPt-startPt)/magD,
404 segmenti-1,
405 segmenti,
406
407 newTracks.last(),
408 newScalars.last(),
409 newVectors.last()
410 );
411
412 if (!bb.contains(endPt))
413 {
414 bb.intersects
415 (
416 endPt,
417 point(clipPt),
418 clipPt
419 );
420
421 // Store point and interpolated values
422 storePoint
423 (
424 tracki,
425
426 mag(clipPt-startPt)/magD,
427 segmenti-1,
428 segmenti,
429
430 newTracks.last(),
431 newScalars.last(),
432 newVectors.last()
433 );
434
435 newTracks.last().shrink();
436 newScalars.last().shrink();
437 newVectors.last().shrink();
438 }
439 }
440 }
441 }
442 }
443
444 // Last point
445 if (bb.contains(track.last()))
446 {
447 storePoint
448 (
449 tracki,
450
451 1.0,
452 track.size()-2,
453 track.size()-1,
454
455 newTracks.last(),
456 newScalars.last(),
457 newVectors.last()
458 );
459 }
460 }
461}
462
463
465{
466 // Storage for new tracks. Per track, per sample the coordinate (newTracks)
467 // or values for all the sampled fields (newScalars, newVectors)
471
472 forAll(allTracks_, tracki)
473 {
474 const List<point>& track = allTracks_[tracki];
475
476 if (track.size())
477 {
478 // New track. Assume it consists of the whole track
479 newTracks.append(new DynamicList<point>(track.size()));
480 newScalars.append(new DynamicList<scalarList>(track.size()));
481 newVectors.append(new DynamicList<vectorList>(track.size()));
482
483 // Trim, split and append to newTracks
484 trimToBox(bb, tracki, newTracks, newScalars, newVectors);
485 }
486 }
487
488 // Transfer newTracks to allTracks_
489 allTracks_.setSize(newTracks.size());
490 forAll(allTracks_, tracki)
491 {
492 allTracks_[tracki].transfer(newTracks[tracki]);
493 }
494 // Replace track scalars
495 forAll(allScalars_, scalari)
496 {
497 DynamicList<scalarList>& fieldVals = allScalars_[scalari];
498 fieldVals.setSize(newTracks.size());
499
500 forAll(fieldVals, tracki)
501 {
502 scalarList& trackVals = allScalars_[scalari][tracki];
503 trackVals.setSize(newScalars[tracki].size());
504 forAll(trackVals, samplei)
505 {
506 trackVals[samplei] = newScalars[tracki][samplei][scalari];
507 }
508 }
509 }
510 // Replace track vectors
511 forAll(allVectors_, vectori)
512 {
513 DynamicList<vectorList>& fieldVals = allVectors_[vectori];
514 fieldVals.setSize(newTracks.size());
515 forAll(fieldVals, tracki)
516 {
517 vectorList& trackVals = allVectors_[vectori][tracki];
518 trackVals.setSize(newVectors[tracki].size());
519 forAll(trackVals, samplei)
520 {
521 trackVals[samplei] = newVectors[tracki][samplei][vectori];
522 }
523 }
524 }
525}
526
527
529{
530 if (Pstream::parRun())
531 {
532 // Append slave tracks to master ones
533 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
534
535 globalIndex globalTrackIDs(allTracks_.size());
536
537 // Construct a distribution map to pull all to the master.
540
541 if (Pstream::master())
542 {
543 // Master: receive all. My own first, then consecutive
544 // processors.
545 label tracki = 0;
546
547 forAll(recvMap, proci)
548 {
549 labelList& fromProc = recvMap[proci];
550 fromProc.setSize(globalTrackIDs.localSize(proci));
551 forAll(fromProc, i)
552 {
553 fromProc[i] = tracki++;
554 }
555 }
556 }
557
558 labelList& toMaster = sendMap[0];
559 toMaster.setSize(globalTrackIDs.localSize());
560 forAll(toMaster, i)
561 {
562 toMaster[i] = i;
563 }
564
565 const mapDistribute distMap
566 (
567 globalTrackIDs.totalSize(),
568 std::move(sendMap),
569 std::move(recvMap)
570 );
571
572
573 // Distribute the track positions. Note: use scheduled comms
574 // to prevent buffering.
575 allTracks_.shrink();
576 mapDistributeBase::distribute
577 (
579 distMap.schedule(),
580 distMap.constructSize(),
581 distMap.subMap(),
582 false,
583 distMap.constructMap(),
584 false,
585 allTracks_,
586 flipOp()
587 );
588 allTracks_.setCapacity(allTracks_.size());
589
590 // Distribute the scalars
591 forAll(allScalars_, scalari)
592 {
593 allScalars_[scalari].shrink();
594 mapDistributeBase::distribute
595 (
597 distMap.schedule(),
598 distMap.constructSize(),
599 distMap.subMap(),
600 false,
601 distMap.constructMap(),
602 false,
603 allScalars_[scalari],
604 flipOp()
605 );
606 allScalars_[scalari].setCapacity(allScalars_[scalari].size());
607 }
608 // Distribute the vectors
609 forAll(allVectors_, vectori)
610 {
611 allVectors_[vectori].shrink();
612 mapDistributeBase::distribute
613 (
615 distMap.schedule(),
616 distMap.constructSize(),
617 distMap.subMap(),
618 false,
619 distMap.constructMap(),
620 false,
621 allVectors_[vectori],
622 flipOp()
623 );
624 allVectors_[vectori].setCapacity(allVectors_[vectori].size());
625 }
626 }
627
628
629 // Note: filenames scattered below since used in global call
630 HashTable<fileName> outputFileNames;
631
632 if (Pstream::master())
633 {
634 if (!bounds_.empty())
635 {
636 // Clip to bounding box
637 trimToBox(treeBoundBox(bounds_));
638 }
639
640
641 label nTracks = 0;
642 label n = 0;
643 forAll(allTracks_, tracki)
644 {
645 if (allTracks_[tracki].size())
646 {
647 nTracks++;
648 n += allTracks_[tracki].size();
649 }
650 }
651
652 Log << " Tracks:" << nTracks << nl
653 << " Total samples:" << n
654 << endl;
655
656
657 // Massage into form suitable for writers
658 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
659
660 // Make output directory
661
662 fileName vtkPath
663 (
664 time_.globalPath()/functionObject::outputPrefix/"sets"
665 / name()/mesh_.regionName()
666 / mesh_.time().timeName()
667 );
668
669 mkDir(vtkPath);
670
671 // Convert track positions (and compact out empty tracks)
672
673 PtrList<coordSet> tracks(nTracks);
674 nTracks = 0;
675 labelList oldToNewTrack(allTracks_.size(), -1);
676
677 forAll(allTracks_, tracki)
678 {
679 if (allTracks_[tracki].size())
680 {
681 List<point>& points = allTracks_[tracki];
682 scalarList dist(points.size());
683 dist[0] = 0;
684 for (label pointi = 1; pointi < points.size(); ++pointi)
685 {
686 dist[pointi] =
687 dist[pointi-1] + mag(points[pointi] - points[pointi-1]);
688 }
689
690 tracks.set
691 (
692 nTracks,
693 new coordSet
694 (
695 "track" + Foam::name(nTracks),
696 sampledSetAxis(), // "xyz"
697 std::move(allTracks_[tracki]),
698 std::move(dist)
699 )
700 );
701 oldToNewTrack[tracki] = nTracks;
702 ++nTracks;
703 }
704 }
705
706
707 const bool canWrite =
708 (
709 !tracks.empty()
710 && trackWriterPtr_
711 && trackWriterPtr_->enabled()
712 && (!allScalars_.empty() || !allVectors_.empty())
713 );
714
715 if (canWrite)
716 {
717 auto& writer = trackWriterPtr_();
718
719 writer.nFields(allScalars_.size() + allVectors_.size());
720
721 writer.open
722 (
723 tracks,
724 (vtkPath / tracks[0].name())
725 );
726
727
728 // Temporary measure
729 if (!allScalars_.empty())
730 {
731 List<List<scalarField>> scalarValues(allScalars_.size());
732
733 forAll(allScalars_, scalari)
734 {
735 auto& allTrackVals = allScalars_[scalari];
736 scalarValues[scalari].resize(nTracks);
737
738 forAll(allTrackVals, tracki)
739 {
740 scalarList& vals = allTrackVals[tracki];
741 if (vals.size())
742 {
743 const label newTracki = oldToNewTrack[tracki];
744 scalarValues[scalari][newTracki].transfer(vals);
745 }
746 }
747 }
748
749 forAll(scalarNames_, i)
750 {
751 fileName outFile =
752 writer.write(scalarNames_[i], scalarValues[i]);
753
754 outputFileNames.insert
755 (
756 scalarNames_[i],
757 time_.relativePath(outFile, true)
758 );
759 }
760 }
761
762 if (!allVectors_.empty())
763 {
764 List<List<vectorField>> vectorValues(allVectors_.size());
765
766 forAll(allVectors_, vectori)
767 {
768 auto& allTrackVals = allVectors_[vectori];
769 vectorValues[vectori].setSize(nTracks);
770
771 forAll(allTrackVals, tracki)
772 {
773 vectorList& vals = allTrackVals[tracki];
774 if (vals.size())
775 {
776 const label newTracki = oldToNewTrack[tracki];
777 vectorValues[vectori][newTracki].transfer(vals);
778 }
779 }
780 }
781
782 forAll(vectorNames_, i)
783 {
784 fileName outFile =
785 writer.write(vectorNames_[i], vectorValues[i]);
786
787 outputFileNames.insert
788 (
789 vectorNames_[i],
790 time_.relativePath(outFile, true)
791 );
792 }
793 }
794
795 writer.close(true);
796 }
797
798 // Log << " Writing data to " << scalarVtkFile.path() << endl;
799 }
800
801
802 // File names generated on the master but setProperty needed everywher
803 Pstream::scatter(outputFileNames);
804
805 forAllConstIters(outputFileNames, iter)
806 {
807 const word& fieldName = iter.key();
808 const fileName& outputName = iter.val();
809
811 propsDict.add("file", outputName);
812 setProperty(fieldName, propsDict);
813 }
814
815 return true;
816}
817
818
820(
821 const word& newUName,
822 const wordList& newFieldNames
823)
824{
825 UName_ = newUName;
826 fields_ = newFieldNames;
827}
828
829
830// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
831
833(
834 const word& name,
835 const Time& runTime,
836 const dictionary& dict
837)
838:
839 functionObjects::fvMeshFunctionObject(name, runTime, dict),
840 dict_(dict),
841 fields_()
842{}
843
844
846(
847 const word& name,
848 const Time& runTime,
849 const dictionary& dict,
850 const wordList& fieldNames
851)
852:
853 functionObjects::fvMeshFunctionObject(name, runTime, dict),
854 dict_(dict),
855 fields_(fieldNames)
856{}
857
858
859// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
860
862{}
863
864
865// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
866
868{
869 if (&dict_ != &dict)
870 {
871 // Update local copy of dictionary:
872 dict_ = dict;
873 }
874
876
877 Info<< type() << " " << name() << ":" << nl;
878
879 UName_ = dict.getOrDefault<word>("U", "U");
880
881 if (fields_.empty())
882 {
883 dict.readEntry("fields", fields_);
884
885 if (!fields_.found(UName_))
886 {
888 << "Velocity field for tracking " << UName_
889 << " should be present in the list of fields " << fields_
890 << exit(FatalIOError);
891 }
892 }
893
894 Info<< " Employing velocity field " << UName_ << endl;
895
896 bool trackForward;
897 if (dict.readIfPresent("trackForward", trackForward))
898 {
899 trackDir_ =
900 (
901 trackForward
902 ? trackDirType::FORWARD
903 : trackDirType::BACKWARD
904 );
905
906 if (dict.found("direction"))
907 {
909 << "Cannot specify both 'trackForward' and 'direction'" << nl
910 << exit(FatalIOError);
911 }
912 }
913 else
914 {
915 trackDir_ = trackDirTypeNames.get("direction", dict);
916 }
917 dict.readEntry("lifeTime", lifeTime_);
918 if (lifeTime_ < 1)
919 {
921 << "Illegal value " << lifeTime_ << " for lifeTime"
922 << exit(FatalIOError);
923 }
924
925
926 trackLength_ = VGREAT;
927 if (dict.readIfPresent("trackLength", trackLength_))
928 {
929 Info<< type() << " : fixed track length specified : "
930 << trackLength_ << nl << endl;
931 }
932
933
934 bounds_ = boundBox::invertedBox;
935 if (dict.readIfPresent("bounds", bounds_) && !bounds_.empty())
936 {
937 Info<< " clipping all segments to " << bounds_ << nl << endl;
938 }
939
940
941 interpolationScheme_ = dict.getOrDefault
942 (
943 "interpolationScheme",
945 );
946
947 //Info<< " using interpolation " << interpolationScheme_ << endl;
948
949 cloudName_ = dict.getOrDefault<word>("cloud", type());
950
951 sampledSetPtr_.clear();
952 sampledSetAxis_.clear();
953
954 const word setFormat(dict.get<word>("setFormat"));
955
956 trackWriterPtr_ = coordSetWriter::New
957 (
958 setFormat,
959 dict.subOrEmptyDict("formatOptions").optionalSubDict(setFormat)
960 );
961
962 return true;
963}
964
965
967{
968 return true;
969}
970
971
973{
974 Log << type() << " " << name() << " write:" << nl;
975
976 // Do all injection and tracking
977 track();
978
979 writeToFile();
980
981 return true;
982}
983
984
986{
987 if (&mpm.mesh() == &mesh_)
988 {
989 read(dict_);
990 }
991}
992
993
995{
996 // Moving mesh affects the search tree
997 read(dict_);
998}
999
1000
1001// ************************************************************************* //
#define Log
Definition: PDRblock.C:35
Field reading functions for post-processing utilities.
label n
vtk::internalMeshWriter writer(topoMesh, topoCells, vtk::formatType::INLINE_ASCII, runTime.path()/"blockTopology")
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
void setSize(const label n)
Same as resize()
Definition: DynamicList.H:221
void append(const T &val)
Copy append an element to the end of this list.
Definition: DynamicListI.H:503
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
A HashTable similar to std::unordered_map.
Definition: HashTable.H:123
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Definition: HashTableI.H:180
A List with indirect addressing.
Definition: IndirectList.H:119
void transfer(List< T > &list)
Definition: List.C:447
void setSize(const label n)
Alias for resize()
Definition: List.H:218
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
label nProcs() const noexcept
Number of ranks associated with PstreamBuffers.
static void scatter(const List< commsStruct > &comms, T &value, const int tag, const label comm)
Broadcast data: Distribute without modification.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
const T * set(const label i) const
Definition: PtrList.H:138
void append(T *ptr)
Append an element to the end of the list.
Definition: PtrListI.H:113
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
T & last()
Return the last element of the list.
Definition: UListI.H:216
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
bool empty() const noexcept
True if the list is empty (ie, size() is zero)
Definition: UPtrListI.H:113
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
static const boundBox invertedBox
A large inverted boundBox: min/max == +/- ROOTVGREAT.
Definition: boundBox.H:86
Holds list of sampling positions.
Definition: coordSet.H:56
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A class for handling file names.
Definition: fileName.H:76
static word outputPrefix
Directory prefix.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Computes the magnitude of an input field.
Definition: mag.H:153
autoPtr< sampledSet > sampledSetPtr_
Seed set engine.
const sampledSet & sampledSetPoints() const
Demand driven construction of the sampledSet.
word sampledSetAxis_
Axis of the sampled points to output.
void initInterpolations(const label nSeeds, label &UIndex, PtrList< volScalarField > &vsFlds, PtrList< interpolation< scalar > > &vsInterp, PtrList< volVectorField > &vvFlds, PtrList< interpolation< vector > > &vvInterp)
Initialise fields, interpolators and track storage.
autoPtr< indirectPrimitivePatch > wallPatch() const
Construct patch out of all wall patch faces.
const word & sampledSetAxis() const
The axis of the sampledSet. Creates sampledSet if required.
void storePoint(const label tracki, const scalar w, const label lefti, const label righti, DynamicList< point > &newTrack, DynamicList< List< scalar > > &newScalars, DynamicList< List< vector > > &newVectors) const
Generate point and values by interpolating from existing values.
virtual void resetFieldNames(const word &newUName, const wordList &newFieldNames)
Reset the field names.
trackDirType
Enumeration defining the track direction.
static const Enum< trackDirType > trackDirTypeNames
Names for the trackDir.
virtual bool writeToFile()
Write tracks to file.
void trimToBox(const treeBoundBox &bb, const label tracki, PtrList< DynamicList< point > > &newTracks, PtrList< DynamicList< scalarList > > &newScalars, PtrList< DynamicList< vectorList > > &newVectors) const
Trim and possibly split a track.
virtual bool execute()
Execute the averaging.
virtual bool write()
Track and write.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
label localSize() const
My local size.
Definition: globalIndexI.H:207
label totalSize() const
Global sum of localSizes.
Definition: globalIndexI.H:125
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
Abstract base class for volume field interpolation.
Definition: interpolation.H:60
const labelListList & constructMap() const noexcept
From subsetted data to new reconstructed data.
static List< labelPair > schedule(const labelListList &subMap, const labelListList &constructMap, const int tag, const label comm=UPstream::worldComm)
Calculate a communication schedule. See above.
const labelListList & subMap() const noexcept
From subsetted data back to original data.
label constructSize() const noexcept
Constructed data size.
Class containing processor-to-processor mapping information.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
const polyMesh & mesh() const
Return polyMesh.
Definition: mapPolyMesh.H:363
void movePoints()
Update for new mesh geometry.
void updateMesh()
Update for new mesh topology.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
Holds list of sampling points which is filled at construction time. Various implementations of this b...
Definition: sampledSet.H:86
splitCell * master() const
Definition: splitCell.H:113
Standard boundBox with extra functionality for use in octree.
Definition: treeBoundBox.H:89
bool intersects(const point &overallStart, const vector &overallVec, const point &start, const point &end, point &pt, direction &ptBits) const
Intersects segment; set point to intersection position and face,.
Definition: treeBoundBox.C:162
bool contains(const vector &dir, const point &) const
Contains point (inside or on edge) and moving in direction.
Definition: treeBoundBox.C:320
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
const polyBoundaryMesh & patches
engineTime & runTime
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
word outputName("finiteArea-edges.obj")
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
const pointField & points
Namespace for OpenFOAM.
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
Definition: MSwindows.C:515
messageStream Info
Information stream (stdout output on master, null elsewhere)
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:64
vector point
Point is a vector.
Definition: point.H:43
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
List< vector > vectorList
A List of vectors.
Definition: vectorList.H:64
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:215
IOerror FatalIOError
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList f(nPoints)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278
Functor to negate primitives. Dummy for most other types.
Definition: flipOp.H:69
IOdictionary propsDict(dictIO)
word setFormat(propsDict.getOrDefault< word >("setFormat", "vtk"))