Merge branch 'master' of github.com:OpenFOAM/OpenFOAM-2.3.x

This commit is contained in:
mattijs 2014-07-07 11:49:20 +01:00
commit 05976b0a22
9 changed files with 296 additions and 69 deletions

View file

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -47,6 +47,7 @@ Description
#include "fvMesh.H"
#include "MeshedSurfaces.H"
#include "globalIndex.H"
#include "cellSet.H"
#include "extrudedMesh.H"
#include "extrudeModel.H"
@ -192,6 +193,23 @@ void updateFaceLabels(const mapPolyMesh& map, labelList& faceLabels)
}
void updateCellSet(const mapPolyMesh& map, labelHashSet& cellLabels)
{
const labelList& reverseMap = map.reverseCellMap();
labelHashSet newCellLabels(2*cellLabels.size());
forAll(cellLabels, i)
{
label oldCellI = cellLabels[i];
if (reverseMap[oldCellI] >= 0)
{
newCellLabels.insert(reverseMap[oldCellI]);
}
}
cellLabels.transfer(newCellLabels);
}
int main(int argc, char *argv[])
@ -270,6 +288,9 @@ int main(int argc, char *argv[])
word backPatchName;
labelList backPatchFaces;
// Optional added cells (get written to cellSet)
labelHashSet addedCellsSet;
if (mode == PATCH || mode == MESH)
{
if (flipNormals && mode == MESH)
@ -670,11 +691,33 @@ int main(int argc, char *argv[])
map().reverseFaceMap(),
backPatchFaces
);
// Store added cells
if (mode == MESH)
{
const labelListList addedCells
(
layerExtrude.addedCells
(
meshFromMesh,
layerExtrude.layerFaces()
)
);
forAll(addedCells, faceI)
{
const labelList& aCells = addedCells[faceI];
forAll(aCells, i)
{
addedCellsSet.insert(aCells[i]);
}
}
}
}
else
{
// Read from surface
fileName surfName(dict.lookup("surface"));
surfName.expand();
Info<< "Extruding surfaceMesh read from file " << surfName << nl
<< endl;
@ -810,6 +853,7 @@ int main(int argc, char *argv[])
// Update stored data
updateFaceLabels(map(), frontPatchFaces);
updateFaceLabels(map(), backPatchFaces);
updateCellSet(map(), addedCellsSet);
// Move mesh (if inflation used)
if (map().hasMotionPoints())
@ -913,6 +957,9 @@ int main(int argc, char *argv[])
// Update fields
mesh.updateMesh(map);
// Update local data
updateCellSet(map(), addedCellsSet);
// Move mesh (if inflation used)
if (map().hasMotionPoints())
{
@ -929,6 +976,21 @@ int main(int argc, char *argv[])
<< exit(FatalError);
}
// Need writing cellSet
label nAdded = returnReduce(addedCellsSet.size(), sumOp<label>());
if (nAdded > 0)
{
cellSet addedCells(mesh, "addedCells", addedCellsSet);
Info<< "Writing added cells to cellSet " << addedCells.name()
<< nl << endl;
if (!addedCells.write())
{
FatalErrorIn(args.executable()) << "Failed writing cellSet"
<< addedCells.name()
<< exit(FatalError);
}
}
Info<< "End\n" << endl;
return 0;

View file

@ -46,15 +46,6 @@ wedgeFvPatch::wedgeFvPatch(const polyPatch& patch, const fvBoundaryMesh& bm)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::tmp<Foam::vectorField> Foam::wedgeFvPatch::delta() const
{
const vectorField nHat(nf());
return nHat*(nHat & (Cf() - Cn()));
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

View file

@ -83,9 +83,6 @@ public:
{
return wedgePolyPatch_.cellT();
}
//- Return cell-centre to face normal vector
virtual tmp<vectorField> delta() const;
};

View file

@ -140,7 +140,9 @@ const Foam::scalarField& Foam::fvPatch::magSf() const
Foam::tmp<Foam::vectorField> Foam::fvPatch::delta() const
{
return Cf() - Cn();
// Use patch-normal delta for all non-coupled BCs
const vectorField nHat(nf());
return nHat*(nHat & (Cf() - Cn()));
}

View file

@ -30,7 +30,7 @@ License
template<class SourcePatch, class TargetPatch>
void Foam::directAMI<SourcePatch, TargetPatch>::appendToDirectSeeds
(
boolList& mapFlag,
labelList& mapFlag,
labelList& srcTgtSeed,
DynamicList<label>& srcSeeds,
DynamicList<label>& nonOverlapFaces,
@ -50,19 +50,32 @@ void Foam::directAMI<SourcePatch, TargetPatch>::appendToDirectSeeds
{
label srcI = srcNbr[i];
if (mapFlag[srcI] && srcTgtSeed[srcI] == -1)
if ((mapFlag[srcI] == 0) && (srcTgtSeed[srcI] == -1))
{
// first attempt: match by comparing face centres
const face& srcF = this->srcPatch_[srcI];
const vector srcN = srcF.normal(srcPoints);
const point& srcC = srcCf[srcI];
scalar tol = GREAT;
forAll(srcF, fpI)
{
const point& p = srcPoints[srcF[fpI]];
scalar d2 = magSqr(p - srcC);
if (d2 < tol)
{
tol = d2;
}
}
tol = max(SMALL, 0.0001*sqrt(tol));
bool found = false;
forAll(tgtNbr, j)
{
label tgtI = tgtNbr[j];
const face& tgtF = this->tgtPatch_[tgtI];
pointHit ray = tgtF.ray(srcCf[srcI], srcN, tgtPoints);
const point tgtC = tgtF.centre(tgtPoints);
if (ray.hit())
if (mag(srcC - tgtC) < tol)
{
// new match - append to lists
found = true;
@ -74,11 +87,57 @@ void Foam::directAMI<SourcePatch, TargetPatch>::appendToDirectSeeds
}
}
// second attempt: match by shooting a ray into the tgt face
if (!found)
{
// no match available for source face srcI
mapFlag[srcI] = false;
const vector srcN = srcF.normal(srcPoints);
forAll(tgtNbr, j)
{
label tgtI = tgtNbr[j];
const face& tgtF = this->tgtPatch_[tgtI];
pointHit ray = tgtF.ray(srcCf[srcI], srcN, tgtPoints);
if (ray.hit())
{
// new match - append to lists
found = true;
srcTgtSeed[srcI] = tgtI;
srcSeeds.append(srcI);
break;
}
}
}
// no match available for source face srcI
if (!found)
{
mapFlag[srcI] = -1;
nonOverlapFaces.append(srcI);
if (debug)
{
Pout<< "source face not found: id=" << srcI
<< " centre=" << srcCf[srcI]
<< " normal=" << srcF.normal(srcPoints)
<< " points=" << srcF.points(srcPoints)
<< endl;
Pout<< "target neighbours:" << nl;
forAll(tgtNbr, j)
{
label tgtI = tgtNbr[j];
const face& tgtF = this->tgtPatch_[tgtI];
Pout<< "face id: " << tgtI
<< " centre=" << tgtF.centre(tgtPoints)
<< " normal=" << tgtF.normal(tgtPoints)
<< " points=" << tgtF.points(tgtPoints)
<< endl;
}
}
}
}
}
@ -96,6 +155,36 @@ void Foam::directAMI<SourcePatch, TargetPatch>::appendToDirectSeeds
}
template<class SourcePatch, class TargetPatch>
void Foam::directAMI<SourcePatch, TargetPatch>::restartAdvancingFront
(
labelList& mapFlag,
DynamicList<label>& nonOverlapFaces,
label& srcFaceI,
label& tgtFaceI
) const
{
forAll(mapFlag, faceI)
{
if (mapFlag[faceI] == 0)
{
tgtFaceI = this->findTargetFace(faceI);
if (tgtFaceI < 0)
{
mapFlag[faceI] = -1;
nonOverlapFaces.append(faceI);
}
else
{
srcFaceI = faceI;
break;
}
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
@ -176,15 +265,19 @@ void Foam::directAMI<SourcePatch, TargetPatch>::calculate
srcTgtSeed[srcFaceI] = tgtFaceI;
// list to keep track of whether src face can be mapped
boolList mapFlag(srcAddr.size(), true);
// 1 = mapped, 0 = untested, -1 = cannot map
labelList mapFlag(srcAddr.size(), 0);
label nTested = 0;
DynamicList<label> nonOverlapFaces;
do
{
srcAddr[srcFaceI].append(tgtFaceI);
tgtAddr[tgtFaceI].append(srcFaceI);
mapFlag[srcFaceI] = false;
mapFlag[srcFaceI] = 1;
nTested++;
// Do advancing front starting from srcFaceI, tgtFaceI
appendToDirectSeeds
@ -196,6 +289,12 @@ void Foam::directAMI<SourcePatch, TargetPatch>::calculate
srcFaceI,
tgtFaceI
);
if (srcFaceI < 0 && nTested < this->srcPatch_.size())
{
restartAdvancingFront(mapFlag, nonOverlapFaces, srcFaceI, tgtFaceI);
}
} while (srcFaceI >= 0);
if (nonOverlapFaces.size() != 0)
@ -211,14 +310,16 @@ void Foam::directAMI<SourcePatch, TargetPatch>::calculate
forAll(srcAddr, i)
{
scalar magSf = this->srcMagSf_[i];
srcAddress[i].transfer(srcAddr[i]);
// srcWeights[i] = scalarList(srcAddr[i].size(), magSf);
srcWeights[i] = scalarList(1, magSf);
srcAddress[i].transfer(srcAddr[i]);
}
forAll(tgtAddr, i)
{
scalar magSf = this->tgtMagSf_[i];
tgtAddress[i].transfer(tgtAddr[i]);
// tgtWeights[i] = scalarList(tgtAddr[i].size(), magSf);
tgtWeights[i] = scalarList(1, magSf);
tgtAddress[i].transfer(tgtAddr[i]);
}
}

View file

@ -67,7 +67,7 @@ private:
//- Append to list of src face seed indices
void appendToDirectSeeds
(
boolList& mapFlag,
labelList& mapFlag,
labelList& srcTgtSeed,
DynamicList<label>& srcSeeds,
DynamicList<label>& nonOverlapFaces,
@ -75,6 +75,16 @@ private:
label& tgtFaceI
) const;
//- Restart the advancing front - typically happens for
// disconnected regions
void restartAdvancingFront
(
labelList& mapFlag,
DynamicList<label>& nonOverlapFaces,
label& srcFaceI,
label& tgtFaceI
) const;
// Evaluation

View file

@ -99,11 +99,10 @@ void Foam::fieldValues::faceSource::setFaceZoneFaces()
const faceZone& fZone = mesh().faceZones()[zoneId];
faceId_.setSize(fZone.size());
facePatchId_.setSize(fZone.size());
faceSign_.setSize(fZone.size());
DynamicList<label> faceIds(fZone.size());
DynamicList<label> facePatchIds(fZone.size());
DynamicList<label> faceSigns(fZone.size());
label count = 0;
forAll(fZone, i)
{
label faceI = fZone[i];
@ -145,27 +144,26 @@ void Foam::fieldValues::faceSource::setFaceZoneFaces()
{
if (fZone.flipMap()[i])
{
faceSign_[count] = -1;
faceSigns.append(-1);
}
else
{
faceSign_[count] = 1;
faceSigns.append(1);
}
faceId_[count] = faceId;
facePatchId_[count] = facePatchId;
count++;
faceIds.append(faceId);
facePatchIds.append(facePatchId);
}
}
faceId_.setSize(count);
facePatchId_.setSize(count);
faceSign_.setSize(count);
faceId_.transfer(faceIds);
facePatchId_.transfer(facePatchIds);
faceSign_.transfer(faceSigns);
nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
if (debug)
{
Pout<< "Original face zone size = " << fZone.size()
<< ", new size = " << count << endl;
<< ", new size = " << faceId_.size() << endl;
}
}
@ -445,12 +443,43 @@ void Foam::fieldValues::faceSource::initialise(const dictionary& dict)
if (dict.readIfPresent("weightField", weightFieldName_))
{
Info<< " weight field = " << weightFieldName_;
Info<< " weight field = " << weightFieldName_ << nl;
}
if (dict.found("orientedWeightField"))
{
if (weightFieldName_ == "none")
{
dict.lookup("orientedWeightField") >> weightFieldName_;
Info<< " weight field = " << weightFieldName_ << nl;
orientWeightField_ = true;
}
else
{
FatalIOErrorIn
(
"void Foam::fieldValues::faceSource::initialise"
"("
"const dictionary&"
")",
dict
)
<< "Either weightField or orientedWeightField can be supplied, "
<< "but not both"
<< exit(FatalIOError);
}
}
List<word> orientedFields;
if (dict.readIfPresent("orientedFields", orientedFields))
{
orientedFieldsStart_ = fields_.size();
fields_.append(orientedFields);
}
if (dict.readIfPresent("scaleFactor", scaleFactor_))
{
Info<< " scale factor = " << scaleFactor_;
Info<< " scale factor = " << scaleFactor_ << nl;
}
Info<< nl << endl;
@ -581,6 +610,8 @@ Foam::fieldValues::faceSource::faceSource
source_(sourceTypeNames_.read(dict.lookup("source"))),
operation_(operationTypeNames_.read(dict.lookup("operation"))),
weightFieldName_("none"),
orientWeightField_(false),
orientedFieldsStart_(labelMax),
scaleFactor_(1.0),
nFaces_(0),
faceId_(),
@ -628,24 +659,42 @@ void Foam::fieldValues::faceSource::write()
totalArea = gSum(filterField(mesh().magSf(), false));
}
if (Pstream::master())
{
file() << obr_.time().value() << tab << totalArea;
}
// construct weight field
scalarField weightField(faceId_.size(), 1.0);
if (weightFieldName_ != "none")
{
weightField =
getFieldValues<scalar>
(
weightFieldName_,
true,
orientWeightField_
);
}
// Combine onto master
combineFields(weightField);
// process the fields
forAll(fields_, i)
{
const word& fieldName = fields_[i];
bool processed = false;
bool ok = false;
processed = processed || writeValues<scalar>(fieldName);
processed = processed || writeValues<vector>(fieldName);
processed = processed || writeValues<sphericalTensor>(fieldName);
processed = processed || writeValues<symmTensor>(fieldName);
processed = processed || writeValues<tensor>(fieldName);
bool orient = i >= orientedFieldsStart_;
ok = ok || writeValues<scalar>(fieldName, weightField, orient);
ok = ok || writeValues<vector>(fieldName, weightField, orient);
ok = ok
|| writeValues<sphericalTensor>(fieldName, weightField, orient);
ok = ok || writeValues<symmTensor>(fieldName, weightField, orient);
ok = ok || writeValues<tensor>(fieldName, weightField, orient);
if (!processed)
if (!ok)
{
WarningIn("void Foam::fieldValues::faceSource::write()")
<< "Requested field " << fieldName

View file

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -71,8 +71,10 @@ Description
sourceName | name of face source if required | no |
operation | operation to perform | yes |
weightField | name of field to apply weighting | no |
orientedWeightField | name of oriented field to apply weighting | no |
scaleFactor | scale factor | no | 1
fields | list of fields to operate on | yes |
orientedFields | list of oriented fields to operate on | no |
\endtable
\linebreak
@ -109,8 +111,8 @@ Note
- faces on empty patches get ignored
- if the field is a volField the \c faceZone can only consist of boundary
faces
- all fields are oriented according to the \c faceZone (so you might
e.g. see negative pressure)
- the `oriented' entries relate to mesh-oriented fields, such as the
flux, phi. These fields will be oriented according to the face normals.
- using \c sampledSurfaces:
- not available for surface fields
- if interpolate=true they use \c interpolationCellPoint
@ -128,6 +130,7 @@ SeeAlso
SourceFiles
faceSource.C
faceSourceTemplates.C
\*---------------------------------------------------------------------------*/
@ -242,6 +245,12 @@ protected:
//- Weight field name - optional
word weightFieldName_;
//- Flag to indicate if flipMap should be applied to the weight field
bool orientWeightField_;
//- Start index of fields that require application of flipMap
label orientedFieldsStart_;
//- Scale factor - optional
scalar scaleFactor_;
@ -282,7 +291,8 @@ protected:
tmp<Field<Type> > getFieldValues
(
const word& fieldName,
const bool mustGet = false
const bool mustGet = false,
const bool applyOrientation = false
) const;
//- Apply the 'operation' to the values. Operation has to
@ -356,7 +366,12 @@ public:
//- Templated helper function to output field values
template<class Type>
bool writeValues(const word& fieldName);
bool writeValues
(
const word& fieldName,
const scalarField& weightField,
const bool orient
);
//- Filter a surface field according to faceIds
template<class Type>

View file

@ -54,7 +54,8 @@ template<class Type>
Foam::tmp<Foam::Field<Type> > Foam::fieldValues::faceSource::getFieldValues
(
const word& fieldName,
const bool mustGet
const bool mustGet,
const bool applyOrientation
) const
{
typedef GeometricField<Type, fvsPatchField, surfaceMesh> sf;
@ -62,7 +63,7 @@ Foam::tmp<Foam::Field<Type> > Foam::fieldValues::faceSource::getFieldValues
if (source_ != stSampledSurface && obr_.foundObject<sf>(fieldName))
{
return filterField(obr_.lookupObject<sf>(fieldName), true);
return filterField(obr_.lookupObject<sf>(fieldName), applyOrientation);
}
else if (obr_.foundObject<vf>(fieldName))
{
@ -103,7 +104,7 @@ Foam::tmp<Foam::Field<Type> > Foam::fieldValues::faceSource::getFieldValues
}
else
{
return filterField(fld, true);
return filterField(fld, applyOrientation);
}
}
@ -115,6 +116,7 @@ Foam::tmp<Foam::Field<Type> > Foam::fieldValues::faceSource::getFieldValues
"Foam::fieldValues::faceSource::getFieldValues"
"("
"const word&, "
"const bool, "
"const bool"
") const"
) << "Field " << fieldName << " not found in database"
@ -267,22 +269,20 @@ Type Foam::fieldValues::faceSource::processValues
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
bool Foam::fieldValues::faceSource::writeValues(const word& fieldName)
bool Foam::fieldValues::faceSource::writeValues
(
const word& fieldName,
const scalarField& weightField,
const bool orient
)
{
const bool ok = validField<Type>(fieldName);
if (ok)
{
Field<Type> values(getFieldValues<Type>(fieldName));
scalarField weightField(values.size(), 1.0);
if (weightFieldName_ != "none")
{
weightField = getFieldValues<scalar>(weightFieldName_, true);
}
Field<Type> values(getFieldValues<Type>(fieldName, true, orient));
vectorField Sf;
if (surfacePtr_.valid())
{
// Get oriented Sf
@ -291,13 +291,12 @@ bool Foam::fieldValues::faceSource::writeValues(const word& fieldName)
else
{
// Get oriented Sf
Sf = filterField(mesh().Sf(), false);
Sf = filterField(mesh().Sf(), true);
}
// Combine onto master
combineFields(values);
combineFields(Sf);
combineFields(weightField);
// Write raw values on surface if specified
if (surfaceWriterPtr_.valid())
@ -378,7 +377,8 @@ Foam::tmp<Foam::Field<Type> > Foam::fieldValues::faceSource::filterField
(
"fieldValues::faceSource::filterField"
"("
"const GeometricField<Type, fvPatchField, volMesh>&"
"const GeometricField<Type, fvPatchField, volMesh>&, "
"const bool"
") const"
) << type() << " " << name_ << ": "
<< sourceTypeNames_[source_] << "(" << sourceName_ << "):"