surfaceFieldValue.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) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2017-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 "surfaceFieldValue.H"
30#include "fvMesh.H"
31#include "emptyPolyPatch.H"
32#include "coupledPolyPatch.H"
33#include "sampledSurface.H"
35#include "PatchTools.H"
37
38// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39
40namespace Foam
41{
42namespace functionObjects
43{
44namespace fieldValues
45{
49}
50}
51}
52
53
54const Foam::Enum
55<
57>
59({
60 { regionTypes::stFaceZone, "faceZone" },
61 { regionTypes::stPatch, "patch" },
62 { regionTypes::stObject, "functionObjectSurface" },
63 { regionTypes::stSampled, "sampledSurface" },
64});
65
66
67const Foam::Enum
68<
70>
72({
73 // Normal operations
74 { operationType::opNone, "none" },
75 { operationType::opMin, "min" },
76 { operationType::opMax, "max" },
77 { operationType::opSum, "sum" },
78 { operationType::opSumMag, "sumMag" },
79 { operationType::opSumDirection, "sumDirection" },
80 { operationType::opSumDirectionBalance, "sumDirectionBalance" },
81 { operationType::opAverage, "average" },
82 { operationType::opAreaAverage, "areaAverage" },
83 { operationType::opAreaIntegrate, "areaIntegrate" },
84 { operationType::opCoV, "CoV" },
85 { operationType::opAreaNormalAverage, "areaNormalAverage" },
86 { operationType::opAreaNormalIntegrate, "areaNormalIntegrate" },
87 { operationType::opUniformity, "uniformity" },
88
89 // Using weighting
90 { operationType::opWeightedSum, "weightedSum" },
91 { operationType::opWeightedAverage, "weightedAverage" },
92 { operationType::opWeightedAreaAverage, "weightedAreaAverage" },
93 { operationType::opWeightedAreaIntegrate, "weightedAreaIntegrate" },
94 { operationType::opWeightedUniformity, "weightedUniformity" },
95
96 // Using absolute weighting
97 { operationType::opAbsWeightedSum, "absWeightedSum" },
98 { operationType::opAbsWeightedAverage, "absWeightedAverage" },
99 { operationType::opAbsWeightedAreaAverage, "absWeightedAreaAverage" },
100 { operationType::opAbsWeightedAreaIntegrate, "absWeightedAreaIntegrate" },
101 { operationType::opAbsWeightedUniformity, "absWeightedUniformity" },
102});
103
104const Foam::Enum
105<
107>
109({
110 { postOperationType::postOpNone, "none" },
111 { postOperationType::postOpMag, "mag" },
112 { postOperationType::postOpSqrt, "sqrt" },
113});
114
115
116// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
117
120{
121 if (stObject == regionType_)
122 {
124 }
125
126 return mesh_;
127}
128
129
130void Foam::functionObjects::fieldValues::surfaceFieldValue::setFaceZoneFaces()
131{
132 // Indices for all matches, already sorted
133 const labelList zoneIds
134 (
135 mesh_.faceZones().indices(selectionNames_)
136 );
137
138 // Total number of faces selected
139 label numFaces = 0;
140 for (const label zoneId : zoneIds)
141 {
142 numFaces += mesh_.faceZones()[zoneId].size();
143 }
144
145 if (zoneIds.empty())
146 {
148 << type() << ' ' << name() << ": "
149 << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
150 << " No matching face zone(s): "
151 << flatOutput(selectionNames_) << nl
152 << " Known face zones: "
153 << flatOutput(mesh_.faceZones().names()) << nl
154 << exit(FatalError);
155 }
156
157 // Could also check this
158 #if 0
159 if (!returnReduce(bool(numFaces), orOp<bool>()))
160 {
162 << type() << ' ' << name() << ": "
163 << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
164 << " The faceZone specification: "
165 << flatOutput(selectionNames_) << nl
166 << " resulted in 0 faces" << nl
167 << exit(FatalError);
168 }
169 #endif
170
171 faceId_.resize(numFaces);
172 facePatchId_.resize(numFaces);
173 faceFlip_.resize(numFaces);
174
175 numFaces = 0;
176
177 for (const label zoneId : zoneIds)
178 {
179 const faceZone& fZone = mesh_.faceZones()[zoneId];
180
181 forAll(fZone, i)
182 {
183 const label meshFacei = fZone[i];
184 const bool isFlip = fZone.flipMap()[i];
185
186 // Internal faces
187 label faceId = meshFacei;
188 label facePatchId = -1;
189
190 // Boundary faces
191 if (!mesh_.isInternalFace(meshFacei))
192 {
193 facePatchId = mesh_.boundaryMesh().whichPatch(meshFacei);
194 const polyPatch& pp = mesh_.boundaryMesh()[facePatchId];
195 const auto* cpp = isA<coupledPolyPatch>(pp);
196
197 if (cpp)
198 {
199 faceId = (cpp->owner() ? pp.whichFace(meshFacei) : -1);
200 }
201 else if (!isA<emptyPolyPatch>(pp))
202 {
203 faceId = pp.whichFace(meshFacei);
204 }
205 else
206 {
207 faceId = -1;
208 facePatchId = -1;
209 }
210 }
211
212 if (faceId >= 0)
213 {
214 faceId_[numFaces] = faceId;
215 facePatchId_[numFaces] = facePatchId;
216 faceFlip_[numFaces] = isFlip;
217
218 ++numFaces;
219 }
220 }
221 }
222
223 // Shrink to size used
224 faceId_.resize(numFaces);
225 facePatchId_.resize(numFaces);
226 faceFlip_.resize(numFaces);
227 nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
228}
229
230
231void Foam::functionObjects::fieldValues::surfaceFieldValue::setPatchFaces()
232{
233 // Patch indices for all matches
235
236 // Total number of faces selected
237 label numFaces = 0;
238
239 labelList selected
240 (
241 mesh_.boundaryMesh().patchSet
242 (
243 selectionNames_,
244 false // warnNotFound - we do that ourselves
245 ).sortedToc()
246 );
247
248 DynamicList<label> bad;
249 for (const label patchi : selected)
250 {
251 const polyPatch& pp = mesh_.boundaryMesh()[patchi];
252
253 if (isA<emptyPolyPatch>(pp))
254 {
255 bad.append(patchi);
256 }
257 else
258 {
259 numFaces += pp.size();
260 }
261 }
262
263 if (bad.size())
264 {
265 label nGood = (selected.size() - bad.size());
266
267 auto& os = (nGood > 0 ? WarningInFunction : FatalErrorInFunction);
268
269 os << "Cannot sample an empty patch" << nl;
270
271 for (const label patchi : bad)
272 {
273 os << " "
274 << mesh_.boundaryMesh()[patchi].name() << nl;
275 }
276
277 if (nGood)
278 {
279 os << "No non-empty patches selected" << endl
280 << exit(FatalError);
281 }
282 else
283 {
284 os << "Selected " << nGood << " non-empty patches" << nl;
285 }
286
287 patchIds.resize(nGood);
288 nGood = 0;
289 for (const label patchi : selected)
290 {
291 if (!bad.found(patchi))
292 {
293 patchIds[nGood] = patchi;
294 ++nGood;
295 }
296 }
297 }
298 else
299 {
300 patchIds = std::move(selected);
301 }
302
303 if (patchIds.empty())
304 {
306 << type() << ' ' << name() << ": "
307 << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
308 << " No matching patch name(s): "
309 << flatOutput(selectionNames_) << nl
310 << " Known patch names:" << nl
311 << mesh_.boundaryMesh().names() << nl
312 << exit(FatalError);
313 }
314
315 // Could also check this
316 #if 0
317 if (!returnReduce(bool(numFaces), orOp<bool>()))
318 {
320 << type() << ' ' << name() << ": "
321 << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
322 << " The patch specification: "
323 << flatOutput(selectionNames_) << nl
324 << " resulted in 0 faces" << nl
325 << exit(FatalError);
326 }
327 #endif
328
329 faceId_.resize(numFaces);
330 facePatchId_.resize(numFaces);
331 faceFlip_.resize(numFaces);
332 nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
333
334 numFaces = 0;
335 for (const label patchi : patchIds)
336 {
337 const polyPatch& pp = mesh_.boundaryMesh()[patchi];
338 const label len = pp.size();
339
340 SubList<label>(faceId_, len, numFaces) = identity(len);
341 SubList<label>(facePatchId_, len, numFaces) = patchi;
342 SubList<bool>(faceFlip_, len, numFaces) = false;
343
344 numFaces += len;
345 }
346}
347
348
349void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry
350(
351 faceList& faces,
353) const
354{
355 labelList whichFaces(faceId_);
356
357 // Remap patch-face ids to mesh face ids
358 forAll(whichFaces, i)
359 {
360 const label patchi = facePatchId_[i];
361 if (patchi != -1)
362 {
363 whichFaces[i] += mesh_.boundaryMesh()[patchi].start();
364 }
365 }
366
368 (
369 IndirectList<face>(mesh_.faces(), std::move(whichFaces)),
370 mesh_.points()
371 );
372
373
374 if (Pstream::parRun())
375 {
376 // Topological merge
377 labelList pointToGlobal;
378 labelList uniqueMeshPointLabels;
379 autoPtr<globalIndex> globalPoints;
380 autoPtr<globalIndex> globalFaces;
382 (
383 mesh_,
384 pp.localFaces(),
385 pp.meshPoints(),
386 pp.meshPointMap(),
387
388 pointToGlobal,
389 uniqueMeshPointLabels,
390 globalPoints,
391 globalFaces,
392
393 faces,
394 points
395 );
396 }
397 else
398 {
399 faces = pp.localFaces();
400 points = pp.localPoints();
401 }
402}
403
404
405void Foam::functionObjects::fieldValues::surfaceFieldValue::
406combineSurfaceGeometry
407(
408 faceList& faces,
410) const
411{
412 if (stObject == regionType_)
413 {
414 const polySurface& s = dynamicCast<const polySurface>(obr());
415
416 if (Pstream::parRun())
417 {
418 // Dimension as fraction of surface
419 const scalar mergeDim = 1e-10*boundBox(s.points(), true).mag();
420
421 labelList pointsMap;
422
424 (
425 mergeDim,
426 primitivePatch(SubList<face>(s.faces()), s.points()),
427 points,
428 faces,
429 pointsMap
430 );
431 }
432 else
433 {
434 faces = s.faces();
435 points = s.points();
436 }
437 }
438 else if (sampledPtr_)
439 {
440 const sampledSurface& s = *sampledPtr_;
441
442 if (Pstream::parRun())
443 {
444 // Dimension as fraction of mesh bounding box
445 const scalar mergeDim = 1e-10*mesh_.bounds().mag();
446
447 labelList pointsMap;
448
450 (
451 mergeDim,
452 primitivePatch(SubList<face>(s.faces()), s.points()),
453 points,
454 faces,
455 pointsMap
456 );
457 }
458 else
459 {
460 faces = s.faces();
461 points = s.points();
462 }
463 }
464}
465
466
467Foam::scalar
468Foam::functionObjects::fieldValues::surfaceFieldValue::totalArea() const
469{
470 scalar totalArea = 0;
471
472 if (stObject == regionType_)
473 {
474 const polySurface& s = dynamicCast<const polySurface>(obr());
475
476 totalArea = gSum(s.magSf());
477 }
478 else if (sampledPtr_)
479 {
480 totalArea = gSum(sampledPtr_->magSf());
481 }
482 else
483 {
484 totalArea = gSum(filterField(mesh_.magSf()));
485 }
486
487 return totalArea;
488}
489
490
491// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
492
494const noexcept
495{
496 // Few operations do not require the Sf field
497 switch (operation_)
498 {
499 case opNone:
500 case opMin:
501 case opMax:
502 case opSum:
503 case opSumMag:
504 case opAverage:
505 return false;
506
507 default:
508 return true;
509 }
510}
511
512
514{
515 if (sampledPtr_)
516 {
517 sampledPtr_->update();
518 }
519
520 if (!needsUpdate_)
521 {
522 return false;
523 }
524
525 switch (regionType_)
526 {
527 case stFaceZone:
528 {
529 setFaceZoneFaces();
530 break;
531 }
532 case stPatch:
533 {
534 setPatchFaces();
535 break;
536 }
537 case stObject:
538 {
539 const polySurface& s = dynamicCast<const polySurface>(obr());
540 nFaces_ = returnReduce(s.size(), sumOp<label>());
541 break;
542 }
543 case stSampled:
544 {
545 nFaces_ = returnReduce(sampledPtr_->faces().size(), sumOp<label>());
546 break;
547 }
548
549 // Compiler warning if we forgot an enumeration
550 }
551
552 if (nFaces_ == 0)
553 {
555 << type() << ' ' << name() << ": "
556 << regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
557 << " Region has no faces" << exit(FatalError);
558 }
559
560 totalArea_ = totalArea();
561
562 Log << " total faces = " << nFaces_ << nl
563 << " total area = " << totalArea_ << nl
564 << endl;
565
566 writeFileHeader(file());
567
568 needsUpdate_ = false;
569 return true;
570}
571
572
574(
575 Ostream& os
576)
577{
578 if (canWriteHeader() && (operation_ != opNone))
579 {
580 writeCommented(os, "Region type : ");
581 os << regionTypeNames_[regionType_] << ' ' << regionName_ << nl;
582
583 writeHeaderValue(os, "Faces", nFaces_);
584 writeHeaderValue(os, "Area", totalArea_);
585 writeHeaderValue(os, "Scale factor", scaleFactor_);
586
587 if (weightFieldNames_.size())
588 {
589 writeHeaderValue
590 (
591 os,
592 "Weight field",
593 flatOutput(weightFieldNames_, FlatOutput::BareComma{})
594 );
595 }
596
597 writeCommented(os, "Time");
598 if (writeArea_)
599 {
600 os << tab << "Area";
601 }
602
603 // TBD: add in postOperation information?
604
605 for (const word& fieldName : fields_)
606 {
607 os << tab << operationTypeNames_[operation_]
608 << '(' << fieldName << ')';
609 }
610
611 os << endl;
612 }
613
614 writtenHeader_ = true;
615}
616
617
618template<>
619Foam::scalar
621(
622 const Field<scalar>& values,
623 const vectorField& Sf,
624 const scalarField& weightField
625) const
626{
627 switch (operation_)
628 {
629 case opSumDirection:
630 {
631 const vector n(dict_.get<vector>("direction"));
632 return gSum(pos0(values*(Sf & n))*mag(values));
633 }
634 case opSumDirectionBalance:
635 {
636 const vector n(dict_.get<vector>("direction"));
637 const scalarField nv(values*(Sf & n));
638
639 return gSum(pos0(nv)*mag(values) - neg(nv)*mag(values));
640 }
641
642 case opUniformity:
643 case opWeightedUniformity:
644 case opAbsWeightedUniformity:
645 {
646 const scalar areaTotal = gSum(mag(Sf));
647 tmp<scalarField> areaVal(values * mag(Sf));
648
649 scalar mean, numer;
650
651 if (is_weightedOp() && canWeight(weightField))
652 {
653 // Weighted quantity = (Weight * phi * dA)
654
655 tmp<scalarField> weight
656 (
657 weightingFactor(weightField, is_magOp())
658 );
659
660 // Mean weighted value (area-averaged)
661 mean = gSum(weight()*areaVal()) / areaTotal;
662
663 // Abs. deviation from weighted mean value
664 numer = gSum(mag(weight*areaVal - (mean * mag(Sf))));
665 }
666 else
667 {
668 // Unweighted quantity = (1 * phi * dA)
669
670 // Mean value (area-averaged)
671 mean = gSum(areaVal()) / areaTotal;
672
673 // Abs. deviation from mean value
674 numer = gSum(mag(areaVal - (mean * mag(Sf))));
675 }
676
677 // Uniformity index
678 const scalar ui = 1 - numer/(2*mag(mean*areaTotal) + ROOTVSMALL);
679
680 return min(max(ui, 0), 1);
681 }
682
683 default:
684 {
685 // Fall through to other operations
686 return processSameTypeValues(values, Sf, weightField);
687 }
688 }
689}
690
691
692template<>
695(
696 const Field<vector>& values,
697 const vectorField& Sf,
698 const scalarField& weightField
699) const
700{
701 switch (operation_)
702 {
703 case opSumDirection:
704 {
705 const vector n(dict_.get<vector>("direction").normalise());
706
707 const scalarField nv(n & values);
708 return gSum(pos0(nv)*n*(nv));
709 }
710 case opSumDirectionBalance:
711 {
712 const vector n(dict_.get<vector>("direction").normalise());
713
714 const scalarField nv(n & values);
715 return gSum(pos0(nv)*n*(nv));
716 }
717 case opAreaNormalAverage:
718 {
719 const scalar val = gSum(values & Sf)/gSum(mag(Sf));
720 return vector(val, 0, 0);
721 }
722 case opAreaNormalIntegrate:
723 {
724 const scalar val = gSum(values & Sf);
725 return vector(val, 0, 0);
726 }
727
728 case opUniformity:
729 case opWeightedUniformity:
730 case opAbsWeightedUniformity:
731 {
732 const scalar areaTotal = gSum(mag(Sf));
733 tmp<scalarField> areaVal(values & Sf);
734
735 scalar mean, numer;
736
737 if (is_weightedOp() && canWeight(weightField))
738 {
739 // Weighted quantity = (Weight * phi . dA)
740
741 tmp<scalarField> weight
742 (
743 weightingFactor(weightField, is_magOp())
744 );
745
746 // Mean weighted value (area-averaged)
747 mean = gSum(weight()*areaVal()) / areaTotal;
748
749 // Abs. deviation from weighted mean value
750 numer = gSum(mag(weight*areaVal - (mean * mag(Sf))));
751 }
752 else
753 {
754 // Unweighted quantity = (1 * phi . dA)
755
756 // Mean value (area-averaged)
757 mean = gSum(areaVal()) / areaTotal;
758
759 // Abs. deviation from mean value
760 numer = gSum(mag(areaVal - (mean * mag(Sf))));
761 }
762
763 // Uniformity index
764 const scalar ui = 1 - numer/(2*mag(mean*areaTotal) + ROOTVSMALL);
765
766 return vector(min(max(ui, 0), 1), 0, 0);
767 }
768
769 default:
770 {
771 // Fall through to other operations
772 return processSameTypeValues(values, Sf, weightField);
773 }
774 }
775}
776
777
778// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
779
780template<>
783(
784 const Field<scalar>& weightField,
785 const bool useMag
786)
787{
788 if (useMag)
789 {
790 return mag(weightField);
791 }
792
793 // pass through
794 return weightField;
795}
796
797
798template<>
801(
802 const Field<scalar>& weightField,
803 const vectorField& Sf,
804 const bool useMag
805)
806{
807 // scalar * unit-normal
808
809 // Can skip this check - already used canWeight()
815
816 if (useMag)
817 {
818 return mag(weightField);
819 }
820
821 // pass through
822 return weightField;
823}
824
825
826template<>
829(
830 const Field<scalar>& weightField,
831 const vectorField& Sf,
832 const bool useMag
833)
834{
835 // scalar * Area
836
837 // Can skip this check - already used canWeight()
843
844 if (useMag)
845 {
846 return mag(weightField * mag(Sf));
847 }
848
849 return (weightField * mag(Sf));
850}
851
852
853template<>
856(
857 const Field<vector>& weightField,
858 const vectorField& Sf,
859 const bool useMag
860)
861{
862 // vector (dot) unit-normal
863
864 // Can skip this check - already used canWeight()
870
871 const label len = weightField.size();
872
873 auto tresult = tmp<scalarField>::New(weightField.size());
874 auto& result = tresult.ref();
875
876 for (label facei=0; facei < len; ++facei)
877 {
878 const vector unitNormal(normalised(Sf[facei]));
879 result[facei] = (weightField[facei] & unitNormal);
880 }
881
882 if (useMag)
883 {
884 for (scalar& val : result)
885 {
886 val = mag(val);
887 }
888 }
889
890 return tresult;
891}
892
893
894template<>
897(
898 const Field<vector>& weightField,
899 const vectorField& Sf,
900 const bool useMag
901)
902{
903 // vector (dot) Area
904
905 // Can skip this check - already used canWeight()
911
912 if (useMag)
913 {
914 return mag(weightField & Sf);
915 }
916
917 return (weightField & Sf);
918}
919
920
921// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
922
924(
925 const word& name,
926 const Time& runTime,
927 const dictionary& dict
928)
929:
930 fieldValue(name, runTime, dict, typeName),
931 regionType_(regionTypeNames_.get("regionType", dict)),
932 operation_(operationTypeNames_.get("operation", dict)),
933 postOperation_
934 (
935 postOperationTypeNames_.getOrDefault
936 (
937 "postOperation",
938 dict,
939 postOperationType::postOpNone,
940 true // Failsafe behaviour
941 )
942 ),
943 needsUpdate_(true),
944 writeArea_(false),
945 selectionNames_(),
946 weightFieldNames_(),
947 totalArea_(0),
948 nFaces_(0),
949 faceId_(),
950 facePatchId_(),
951 faceFlip_()
952{
953 read(dict);
954}
955
956
958(
959 const word& name,
960 const objectRegistry& obr,
961 const dictionary& dict
962)
963:
964 fieldValue(name, obr, dict, typeName),
965 regionType_(regionTypeNames_.get("regionType", dict)),
966 operation_(operationTypeNames_.get("operation", dict)),
967 postOperation_
968 (
969 postOperationTypeNames_.getOrDefault
970 (
971 "postOperation",
972 dict,
973 postOperationType::postOpNone,
974 true // Failsafe behaviour
975 )
976 ),
977 needsUpdate_(true),
978 writeArea_(false),
979 selectionNames_(),
980 weightFieldNames_(),
981 totalArea_(0),
982 nFaces_(0),
983 faceId_(),
984 facePatchId_(),
985 faceFlip_()
986{
987 read(dict);
988}
989
990
991// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
992
993// Needs completed sampledSurface, surfaceWriter
995{}
996
997
998// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
999
1001(
1002 const dictionary& dict
1003)
1004{
1006
1007 needsUpdate_ = true;
1008 writeArea_ = dict.getOrDefault("writeArea", false);
1009 weightFieldNames_.clear();
1010 // future?
1011 // sampleFaceScheme_ = dict.getOrDefault<word>("sampleScheme", "cell");
1012
1013 totalArea_ = 0;
1014 nFaces_ = 0;
1015 faceId_.clear();
1016 facePatchId_.clear();
1017 faceFlip_.clear();
1018 sampledPtr_.reset(nullptr);
1019 surfaceWriterPtr_.reset(nullptr);
1020
1021 // Can have "name" (word) and/or "names" (wordRes)
1022 //
1023 // If "names" exists AND contains a literal (non-regex) that can be used
1024 // as a suitable value for "name", the "name" entry becomes optional.
1025
1026 regionName_.clear();
1027 selectionNames_.clear();
1028
1029 {
1030 dict.readIfPresent("names", selectionNames_);
1031
1032 for (const auto& item : selectionNames_)
1033 {
1034 if (item.isLiteral())
1035 {
1036 regionName_ = item;
1037 break;
1038 }
1039 }
1040
1041 // Mandatory if we didn't pick up a value from selectionNames_
1042 dict.readEntry
1043 (
1044 "name",
1045 regionName_,
1047 regionName_.empty()
1048 );
1049
1050 // Ensure there is always content for selectionNames_
1051 if (selectionNames_.empty())
1052 {
1053 selectionNames_.resize(1);
1054 selectionNames_.first() = regionName_;
1055 }
1056 }
1057
1058
1059 // Create sampled surface, but leave 'expired' (ie, no update) since it
1060 // may depend on fields or data that do not yet exist
1061 if (stSampled == regionType_)
1062 {
1063 sampledPtr_ = sampledSurface::New
1064 (
1065 name(),
1066 mesh_,
1067 dict.subDict("sampledSurfaceDict")
1068 );
1069
1070 // Internal consistency. Want face values, never point values!
1071 sampledPtr_->isPointData(false);
1072 }
1073
1074 Info<< type() << ' ' << name() << ':' << nl
1075 << " operation = ";
1076
1077 if (postOperation_ != postOpNone)
1078 {
1079 Info<< postOperationTypeNames_[postOperation_] << '('
1080 << operationTypeNames_[operation_] << ')' << nl;
1081 }
1082 else
1083 {
1084 Info<< operationTypeNames_[operation_] << nl;
1085 }
1086
1087 if (is_weightedOp())
1088 {
1089 // Can have "weightFields" or "weightField"
1090
1091 bool missing = true;
1092 if (dict.readIfPresent("weightFields", weightFieldNames_))
1093 {
1094 missing = false;
1095 }
1096 else
1097 {
1098 weightFieldNames_.resize(1);
1099
1100 if (dict.readIfPresent("weightField", weightFieldNames_.first()))
1101 {
1102 missing = false;
1103 if ("none" == weightFieldNames_.first())
1104 {
1105 // "none" == no weighting
1106 weightFieldNames_.clear();
1107 }
1108 }
1109 }
1110
1111 if (missing)
1112 {
1113 // Suggest possible alternative unweighted operation?
1115 << "The '" << operationTypeNames_[operation_]
1116 << "' operation is missing a weightField." << nl
1117 << "Either provide the weightField, "
1118 << "use weightField 'none' to suppress weighting," << nl
1119 << "or use a different operation."
1120 << exit(FatalIOError);
1121 }
1122
1123 Info<< " weight field = ";
1124 if (weightFieldNames_.empty())
1125 {
1126 Info<< "none" << nl;
1127 }
1128 else
1129 {
1130 Info<< flatOutput(weightFieldNames_) << nl;
1131 }
1132 }
1133
1134 if (stSampled == regionType_ && sampledPtr_)
1135 {
1136 Info<< " sampled surface: ";
1137 sampledPtr_->print(Info, 0);
1138 Info<< nl;
1139 }
1140
1141 if (writeFields_)
1142 {
1143 const word formatName(dict.get<word>("surfaceFormat"));
1144
1145 surfaceWriterPtr_.reset
1146 (
1148 (
1149 formatName,
1150 dict.subOrEmptyDict("formatOptions").subOrEmptyDict(formatName)
1151 )
1152 );
1153
1154 // Propagate field counts (per surface)
1155 surfaceWriterPtr_->nFields(fields_.size());
1156
1157 if (debug)
1158 {
1159 surfaceWriterPtr_->verbose(true);
1160 }
1161
1162 if (surfaceWriterPtr_->enabled())
1163 {
1164 Info<< " surfaceFormat = " << formatName << nl;
1165 }
1166 else
1167 {
1168 surfaceWriterPtr_->clear();
1169 }
1170 }
1171
1172 Info<< nl << endl;
1173
1174 return true;
1175}
1176
1177
1179{
1180 if (needsUpdate_ || operation_ != opNone)
1181 {
1183 }
1184
1185 update();
1186
1187 if (operation_ != opNone)
1188 {
1189 writeCurrentTime(file());
1190 }
1191
1192 if (writeArea_)
1193 {
1194 totalArea_ = totalArea();
1195 Log << " total area = " << totalArea_ << endl;
1196
1197 if (operation_ != opNone && Pstream::master())
1198 {
1199 file() << tab << totalArea_;
1200 }
1201 }
1202
1203 // Many operations use the Sf field
1204 vectorField Sf;
1205 if (usesSf())
1206 {
1207 if (stObject == regionType_)
1208 {
1209 const polySurface& s = dynamicCast<const polySurface>(obr());
1210 Sf = s.Sf();
1211 }
1212 else if (sampledPtr_)
1213 {
1214 Sf = sampledPtr_->Sf();
1215 }
1216 else
1217 {
1218 Sf = filterField(mesh_.Sf());
1219 }
1220 }
1221
1222 // Faces and points for surface format (if specified)
1223 faceList faces;
1225
1226 if (surfaceWriterPtr_)
1227 {
1228 if (withTopologicalMerge())
1229 {
1230 combineMeshGeometry(faces, points);
1231 }
1232 else
1233 {
1234 combineSurfaceGeometry(faces, points);
1235 }
1236 }
1237
1238
1239 // Check availability and type of weight field
1240 // Only support a few weight types:
1241 // scalar: 0-N fields
1242 // vector: 0-1 fields
1243
1244 // Default is a zero-size scalar weight field (ie, weight = 1)
1245 scalarField scalarWeights;
1246 vectorField vectorWeights;
1247
1248 for (const word& weightName : weightFieldNames_)
1249 {
1250 if (validField<scalar>(weightName))
1251 {
1252 tmp<scalarField> tfld = getFieldValues<scalar>(weightName, true);
1253
1254 if (scalarWeights.empty())
1255 {
1256 scalarWeights = tfld;
1257 }
1258 else
1259 {
1260 scalarWeights *= tfld;
1261 }
1262 }
1263 else if (validField<vector>(weightName))
1264 {
1265 tmp<vectorField> tfld = getFieldValues<vector>(weightName, true);
1266
1267 if (vectorWeights.empty())
1268 {
1269 vectorWeights = tfld;
1270 }
1271 else
1272 {
1274 << "weightField " << weightName
1275 << " - only one vector weight field allowed. " << nl
1276 << "weights: " << flatOutput(weightFieldNames_) << nl
1277 << abort(FatalError);
1278 }
1279 }
1280 else if (weightName != "none")
1281 {
1282 // Silently ignore "none", flag everything else as an error
1283
1284 // TBD: treat missing "rho" like incompressible with rho = 1
1285 // and/or provided rhoRef value
1286
1288 << "weightField " << weightName
1289 << " not found or an unsupported type" << nl
1290 << abort(FatalError);
1291 }
1292 }
1293
1294
1295 // Process the fields
1296 if (returnReduce(!vectorWeights.empty(), orOp<bool>()))
1297 {
1298 if (scalarWeights.size())
1299 {
1300 vectorWeights *= scalarWeights;
1301 }
1302
1303 writeAll(Sf, vectorWeights, points, faces);
1304 }
1305 else
1306 {
1307 writeAll(Sf, scalarWeights, points, faces);
1308 }
1309
1310
1311 if (operation_ != opNone)
1312 {
1313 file() << endl;
1314 Log << endl;
1315 }
1316
1317
1318 return true;
1319}
1320
1321
1323(
1324 const mapPolyMesh& mpm
1325)
1326{
1327 needsUpdate_ = true;
1328}
1329
1330
1332(
1333 const polyMesh& mesh
1334)
1335{
1336 needsUpdate_ = true;
1337}
1338
1339
1340// ************************************************************************* //
#define Log
Definition: PDRblock.C:35
label n
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
Generic templated field type.
Definition: Field.H:82
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
static void gatherAndMerge(const scalar mergeDist, const PrimitivePatch< FaceList, PointField > &pp, Field< typename PrimitivePatch< FaceList, PointField >::point_type > &mergedPoints, List< typename PrimitivePatch< FaceList, PointField >::face_type > &mergedFaces, labelList &pointMergeMap, const bool useLocal=false)
Gather points and faces onto master and merge into single patch.
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
bool empty() const noexcept
True if the UList is empty (ie, size() is zero)
Definition: UListI.H:427
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
Vector< Cmpt > & normalise(const scalar tol=ROOTVSMALL)
Inplace normalise the vector by its magnitude.
Definition: VectorI.H:123
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Abstract base-class for Time/database function objects.
Watches for presence of the named trigger file in the case directory and signals a simulation stop (o...
Definition: abort.H:128
Intermediate class for handling field value-based function objects.
Definition: fieldValue.H:123
word regionName_
Name of region (patch, zone, etc.)
Definition: fieldValue.H:132
const dictionary & dict() const noexcept
Return the reference to the construction dictionary.
Definition: fieldValueI.H:31
virtual bool write()
Write.
Definition: fieldValue.C:115
A face regionType variant of the fieldValues function object.
static const Enum< regionTypes > regionTypeNames_
Region type names.
@ stObject
Calculate with function object surface.
static tmp< scalarField > weightingFactor(const Field< WeightType > &weightField, const bool useMag)
Weighting factor.
bool usesSf() const noexcept
True if the operation needs a surface Sf.
const objectRegistry & obr() const
The volume mesh or surface registry being used.
virtual bool read(const dictionary &dict)
Read from dictionary.
bool update()
Update the surface and surface information as required.
static const Enum< postOperationType > postOperationTypeNames_
Operation type names.
virtual void writeFileHeader(Ostream &os)
Output file header information.
static tmp< scalarField > areaWeightingFactor(const Field< WeightType > &weightField, const vectorField &Sf, const bool useMag)
Weighting factor, weight field with area factor.
static const Enum< operationType > operationTypeNames_
Operation type names.
Type processValues(const Field< Type > &values, const vectorField &Sf, const Field< WeightType > &weightField) const
Apply the 'operation' to the values. Wrapper around.
const fvMesh & mesh_
Reference to the fvMesh.
Computes the magnitude of an input field.
Definition: mag.H:153
@ LITERAL
String literal.
Definition: keyType.H:81
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
void movePoints()
Update for new mesh geometry.
void updateMesh()
Update for new mesh topology.
Registry of regIOobjects.
const Type & lookupObject(const word &name, const bool recursive=false) const
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
A surface mesh consisting of general polygon faces and capable of holding fields.
Definition: polySurface.H:71
splitCell * master() const
Definition: splitCell.H:113
A class for managing temporary objects.
Definition: tmp.H:65
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
mesh update()
labelList patchIds
dynamicFvMesh & mesh
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
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
label faceId(-1)
#define WarningInFunction
Report a warning using Foam::Warning.
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Type gSum(const FieldField< Field, Type > &f)
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Definition: labelList.C:38
dimensionedScalar pos0(const dimensionedScalar &ds)
List< label > labelList
A List of labels.
Definition: List.H:66
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:680
messageStream Info
Information stream (stdout output on master, null elsewhere)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
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
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:215
IOerror FatalIOError
dimensionedScalar neg(const dimensionedScalar &ds)
const direction noexcept
Definition: Scalar.H:223
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
Vector< scalar > vector
Definition: vector.H:61
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field.
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
constexpr char tab
The tab '\t' character(0x09)
Definition: Ostream.H:52
dictionary dict
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Surround with '\0' and '\0' separate with ','.
Definition: FlatOutput.H:84