Go to the documentation of this file.
57 const DynamicList<List<point>>& faces,
62 OBJstream
os(fileName +
".obj");
64 if (Pstream::parRun())
67 List<DynamicList<List<point>>> allProcFaces(Pstream::nProcs());
68 allProcFaces[Pstream::myProcNo()] = faces;
69 Pstream::gatherList(allProcFaces);
71 if (Pstream::master())
73 Info<<
"Writing file: " << fileName <<
endl;
75 for (
const DynamicList<List<point>>& procFaces : allProcFaces)
77 for (
const List<point>& facePts : procFaces)
79 os.write(face(
identity(facePts.size())), facePts);
86 Info<<
"Writing file: " << fileName <<
endl;
88 for (
const List<point>& facePts : faces)
90 os.write(face(
identity(facePts.size())), facePts);
98 DynamicList<List<point>>& facePts,
104 cutCellIso cutCell(
mesh,
f);
105 cutFaceIso cutFace(
mesh,
f);
109 cutCell.calcSubCell(cellI, 0.0);
111 alpha1[cellI] =
max(
min(cutCell.VolumeOfFluid(), 1), 0);
115 facePts.append(cutCell.facePoints());
122 if (
mesh.boundary()[patchi].size() > 0)
124 const label start =
mesh.boundary()[patchi].patch().start();
128 forAll(alphap, patchFacei)
130 const label facei = patchFacei + start;
131 cutFace.calcSubFace(facei, 0.0);
133 mag(cutFace.subFaceArea())/magSfp[patchFacei];
140 int main(
int argc,
char *argv[])
144 "Uses cutCellIso to create a volume fraction field from an "
154 const word
dictName(
"setAlphaFieldDict");
157 IOdictionary setAlphaFieldDict(
dictIO);
159 Info<<
"Reading " << setAlphaFieldDict.name() <<
endl;
161 const word fieldName = setAlphaFieldDict.get<word>(
"field");
162 const bool invert = setAlphaFieldDict.getOrDefault(
"invert",
false);
163 const bool writeOBJ = setAlphaFieldDict.getOrDefault(
"writeOBJ",
false);
165 Info<<
"Reading field " << fieldName <<
nl <<
endl;
181 setAlphaFieldDict.get<word>(
"type"),
192 DynamicList<List<point>> facePts;
198 isoFacesToFile(facePts, fieldName +
"0");
201 ISstream::defaultPrecision(18);
214 <<
", 1-max(alpha1): " << 1 -
gMax(
alpha)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static constexpr const zero Zero
Global zero (0)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const word dictName("faMeshDefinition")
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
const volScalarField & alpha1
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
labelList invert(const label len, const labelUList &map)
Create an inverse one-to-one mapping.
#define forAll(list, i)
Loop across all elements in list.
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
void func(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
OBJstream os(runTime.globalPath()/outputName)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
constexpr scalar pi(M_PI)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
Type gMin(const FieldField< Field, Type > &f)
Type gMax(const FieldField< Field, Type > &f)