Non-Conformal Coupled (NCC): Conservative coupling of non-conforming patches

This major development provides coupling of patches which are
non-conformal, i.e. where the faces of one patch do not match the faces
of the other. The coupling is fully conservative and second order
accurate in space, unlike the Arbitrary Mesh Interface (AMI) and
associated ACMI and Repeat AMI methods which NCC replaces.

Description:

A non-conformal couple is a connection between a pair of boundary
patches formed by projecting one patch onto the other in a way that
fills the space between them. The intersection between the projected
surface and patch forms new faces that are incorporated into the finite
volume mesh. These new faces are created identically on both sides of
the couple, and therefore become equivalent to internal faces within the
mesh. The affected cells remain closed, meaning that the area vectors
sum to zero for all the faces of each cell. Consequently, the main
benefits of the finite volume method, i.e. conservation and accuracy,
are not undermined by the coupling.

A couple connects parts of mesh that are otherwise disconnected and can
be used in the following ways:

+ to simulate rotating geometries, e.g. a propeller or stirrer, in which
  a part of the mesh rotates with the geometry and connects to a
  surrounding mesh which is not moving;
+ to connect meshes that are generated separately, which do not conform
  at their boundaries;
+ to connect patches which only partially overlap, in which the
  non-overlapped section forms another boundary, e.g. a wall;
+ to simulate a case with a geometry which is periodically repeating by
  creating multiple couples with different transformations between
  patches.

The capability for simulating partial overlaps replaces the ACMI
functionality, currently provided by the 'cyclicACMI' patch type, and
which is unreliable unless the couple is perfectly flat. The capability
for simulating periodically repeating geometry replaces the Repeat AMI
functionality currently provided by the 'cyclicRepeatAMI' patch type.

Usage:

The process of meshing for NCC is very similar to existing processes for
meshing for AMI. Typically, a mesh is generated with an identifiable set
of internal faces which coincide with the surface through which the mesh
will be coupled. These faces are then duplicated by running the
'createBaffles' utility to create two boundary patches. The points are
then split using 'splitBaffles' in order to permit independent motion of
the patches.

In AMI, these patches are assigned the 'cyclicAMI' patch type, which
couples them using AMI interpolation methods.

With NCC, the patches remain non-coupled, e.g. a 'wall' type. Coupling
is instead achieved by running the new 'createNonConformalCouples'
utility, which creates additional coupled patches of type
'nonConformalCyclic'. These appear in the 'constant/polyMesh/boundary'
file with zero faces; they are populated with faces in the finite volume
mesh during the connection process in NCC.

For a single couple, such as that which separates the rotating and
stationary sections of a mesh, the utility can be called using the
non-coupled patch names as arguments, e.g.

    createNonConformalCouples -overwrite rotatingZoneInner rotatingZoneOuter

where 'rotatingZoneInner' and 'rotatingZoneOuter' are the names of the
patches.

For multiple couples, and/or couples with transformations,
'createNonConformalCouples' should be run without arguments. Settings
will then be read from a configuration file named
'system/createNonConformalCouplesDict'. See
'$FOAM_ETC/caseDicts/annotated/createNonConformalCouplesDict' for
examples.

Boundary conditions must be specified for the non-coupled patches. For a
couple where the patches fully overlap, boundary conditions
corresponding to a slip wall are typically applied to fields, i.e
'movingWallSlipVelocity' (or 'slip' if the mesh is stationary) for
velocity U, 'zeroGradient' or 'fixedFluxPressure' for pressure p, and
'zeroGradient' for other fields.  For a couple with
partially-overlapping patches, boundary conditions are applied which
physically represent the non-overlapped region, e.g. a no-slip wall.

Boundary conditions also need to be specified for the
'nonConformalCyclic' patches created by 'createNonConformalCouples'. It
is generally recommended that this is done by including the
'$FOAM_ETC/caseDicts/setConstraintTypes' file in the 'boundaryField'
section of each of the field files, e.g.

    boundaryField
    {
        #includeEtc "caseDicts/setConstraintTypes"

        inlet
        {
             ...
        }

        ...
    }

For moving mesh cases, it may be necessary to correct the mesh fluxes
that are changed as a result of the connection procedure. If the
connected patches do not conform perfectly to the mesh motion, then
failure to correct the fluxes can result in noise in the pressure
solution.

Correction for the mesh fluxes is enabled by the 'correctMeshPhi' switch
in the 'PIMPLE' (or equivalent) section of 'system/fvSolution'. When it
is enabled, solver settings are required for 'MeshPhi'. The solution
just needs to distribute the error enough to dissipate the noise. A
smooth solver with a loose tolerance is typically sufficient, e.g. the
settings in 'system/fvSolution' shown below:

    solvers
    {
        MeshPhi
        {
            solver          smoothSolver;
            smoother        symGaussSeidel;
            tolerance       1e-2;
            relTol          0;
        }
        ...
    }

    PIMPLE
    {
         correctMeshPhi      yes;
         ...
    }

The solution of 'MeshPhi' is an inexpensive computation since it is
applied only to a small subset of the mesh adjacent to the
couple. Conservation is maintained whether or not the mesh flux
correction is enabled, and regardless of the solution tolerance for
'MeshPhi'.

Advantages of NCC:

+ NCC maintains conservation which is required for many numerical
  schemes and algorithms to operate effectively, in particular those
  designed to maintain boundedness of a solution.

+ Closed-volume systems no longer suffer from accumulation or loss of
  mass, poor convergence of the pressure equation, and/or concentration
  of error in the reference cell.

+ Partially overlapped simulations are now possible on surfaces that are
  not perfectly flat. The projection fills space so no overlaps or
  spaces are generated inside contiguously overlapping sections, even if
  those sections have sharp angles.

+ The finite volume faces created by NCC have geometrically accurate
  centres. This makes the method second order accurate in space.

+ The polyhedral mesh no longer requires duplicate boundary faces to be
  generated in order to run a partially overlapped simulation.

+ Lagrangian elements can now transfer across non-conformal couplings in
  parallel.

+ Once the intersection has been computed and applied to the finite
  volume mesh, it can use standard cyclic or processor cyclic finite
  volume boundary conditions, with no need for additional patch types or
  matrix interfaces.

+ Parallel communication is done using the standard
  processor-patch-field system. This is more efficient than alternative
  systems since it has been carefully optimised for use within the
  linear solvers.

+ Coupled patches are disconnected prior to mesh motion and topology
  change and reconnected afterwards. This simplifies the boundary
  condition specification for mesh motion fields.

Resolved Bug Reports:

+ https://bugs.openfoam.org/view.php?id=663
+ https://bugs.openfoam.org/view.php?id=883
+ https://bugs.openfoam.org/view.php?id=887
+ https://bugs.openfoam.org/view.php?id=1337
+ https://bugs.openfoam.org/view.php?id=1388
+ https://bugs.openfoam.org/view.php?id=1422
+ https://bugs.openfoam.org/view.php?id=1829
+ https://bugs.openfoam.org/view.php?id=1841
+ https://bugs.openfoam.org/view.php?id=2274
+ https://bugs.openfoam.org/view.php?id=2561
+ https://bugs.openfoam.org/view.php?id=3817

Deprecation:

NCC replaces the functionality provided by AMI, ACMI and Repeat AMI.
ACMI and Repeat AMI are insufficiently reliable to warrant further
maintenance so are removed in an accompanying commit to OpenFOAM-dev.
AMI is more widely used so will be retained alongside NCC for the next
version release of OpenFOAM and then subsequently removed from
OpenFOAM-dev.
This commit is contained in:
Will Bainbridge
2022-05-09 14:35:11 +01:00
parent 94679fa88d
commit 569fa31d09
254 changed files with 18751 additions and 4327 deletions

View File

@ -1,3 +1,6 @@
reconstructPar.C
fvFieldReconstructor.C
pointFieldReconstructor.C
reconstructLagrangianPositions.C
EXE = $(FOAM_APPBIN)/reconstructPar

View File

@ -1,16 +1,14 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/parallel/reconstruct/reconstruct/lnInclude \
-I$(LIB_SRC)/parallel/parallel/lnInclude \
-I$(LIB_SRC)/lagrangian/basic/lnInclude \
-I$(LIB_SRC)/regionModels/regionModel/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-lgenericPatchFields \
-llagrangian \
-ldynamicMesh \
-lmeshTools \
-lreconstruct \
-lregionModels
-lparallel \
-llagrangian \
-lregionModels \
-lgenericPatchFields

View File

@ -1,48 +0,0 @@
{
// Foam version 2.1 changes the addressing of faces in faceProcAddressing
// The following code checks and modifies the addressing for cases where
// the decomposition has been done with the foam2.0 and earlier tools, but
// the reconstruction is attempted with version 2.1 or later
label minFaceIndex = labelMax;
PtrList<labelIOList>& faceProcAddressing = procMeshes.faceProcAddressing();
forAll(faceProcAddressing, proci)
{
const labelList& curFaceAddr = faceProcAddressing[proci];
forAll(curFaceAddr, facei)
{
if (mag(curFaceAddr[facei]) < minFaceIndex)
{
minFaceIndex = mag(curFaceAddr[facei]);
}
}
}
if (minFaceIndex < 1)
{
WarningInFunction
<< "parallel decomposition addressing." << endl
<< "It looks like you are trying to reconstruct the case "
<< "decomposed with an earlier version of FOAM, which could\n"
<< "potentially cause compatibility problems. The code will "
<< "attempt to update the addressing automatically; in case of\n"
<< "failure, please repeat the decomposition of the case using "
<< "the current version fo decomposePar"
<< endl;
forAll(faceProcAddressing, proci)
{
labelList& curFaceAddr = faceProcAddressing[proci];
forAll(curFaceAddr, facei)
{
curFaceAddr[facei] += sign(curFaceAddr[facei]);
}
faceProcAddressing[proci].write();
}
}
}

View File

@ -0,0 +1,96 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "fvFieldReconstructor.H"
#include "processorCyclicFvPatch.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::label Foam::fvFieldReconstructor::completePatchID
(
const label proci,
const label procPatchi
) const
{
const fvPatch& procPatch = procMeshes_[proci].boundary()[procPatchi];
if (procPatchi < completeMesh_.boundary().size())
{
return procPatchi;
}
else if (isA<processorCyclicFvPatch>(procPatch))
{
return refCast<const processorCyclicFvPatch>(procPatch).referPatchID();
}
else
{
return -1;
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::fvFieldReconstructor::fvFieldReconstructor
(
const fvMesh& completeMesh,
const PtrList<fvMesh>& procMeshes,
const labelListList& faceProcAddressing,
const labelListList& cellProcAddressing,
const PtrList<surfaceLabelField::Boundary>& faceProcAddressingBf
)
:
completeMesh_(completeMesh),
procMeshes_(procMeshes),
faceProcAddressing_(faceProcAddressing),
cellProcAddressing_(cellProcAddressing),
faceProcAddressingBf_(faceProcAddressingBf),
nReconstructed_(0)
{
forAll(procMeshes_, proci)
{
const fvMesh& procMesh = procMeshes_[proci];
if
(
faceProcAddressing[proci].size() != procMesh.nFaces()
|| cellProcAddressing[proci].size() != procMesh.nCells()
)
{
FatalErrorInFunction
<< "Size of maps does not correspond to size of mesh"
<< " for processor " << proci << endl
<< "faceProcAddressing : " << faceProcAddressing[proci].size()
<< " nFaces : " << procMesh.nFaces() << endl
<< "cellProcAddressing : " << cellProcAddressing[proci].size()
<< " nCell : " << procMesh.nCells() << endl
<< " nFaces : " << procMesh.boundary().size()
<< exit(FatalError);
}
}
}
// ************************************************************************* //

View File

@ -0,0 +1,229 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::fvFieldReconstructor
Description
Finite volume reconstructor for volume and surface fields.
SourceFiles
fvFieldReconstructor.C
fvFieldReconstructorReconstructFields.C
\*---------------------------------------------------------------------------*/
#ifndef fvFieldReconstructor_H
#define fvFieldReconstructor_H
#include "PtrList.H"
#include "surfaceFields.H"
#include "IOobjectList.H"
#include "fvPatchFieldMapper.H"
#include "setSizeFieldMapper.H"
#include "labelIOList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class fvFieldReconstructor Declaration
\*---------------------------------------------------------------------------*/
class fvFieldReconstructor
{
// Private Data
//- Reconstructed mesh reference
const fvMesh& completeMesh_;
//- List of processor meshes
const PtrList<fvMesh>& procMeshes_;
//- List of processor face addressing lists
const labelListList& faceProcAddressing_;
//- List of processor cell addressing lists
const labelListList& cellProcAddressing_;
//- Boundary field of face addressing
const PtrList<surfaceLabelField::Boundary>& faceProcAddressingBf_;
//- Number of fields reconstructed
label nReconstructed_;
// Private Member Functions
//- Convert a processor patch to the corresponding complete patch index
label completePatchID(const label proci, const label procPatchi) const;
//- ...
template<class Type>
static void rmapFaceToFace
(
Field<Type>& toField,
const Field<Type>& fromField,
const labelUList& addressing,
const bool isFlux
);
public:
//- Mapper for sizing only - does not do any actual mapping.
class fvPatchFieldReconstructor
:
public fvPatchFieldMapper,
public setSizeFieldMapper
{
public:
// Constructors
//- Construct given size
fvPatchFieldReconstructor(const label size)
:
setSizeFieldMapper(size)
{}
};
// Constructors
//- Construct from components
fvFieldReconstructor
(
const fvMesh& mesh,
const PtrList<fvMesh>& procMeshes,
const labelListList& faceProcAddressing,
const labelListList& cellProcAddressing,
const PtrList<surfaceLabelField::Boundary>& faceProcAddressingBf
);
//- Disallow default bitwise copy construction
fvFieldReconstructor(const fvFieldReconstructor&) = delete;
// Member Functions
//- Return number of fields reconstructed
label nReconstructed() const
{
return nReconstructed_;
}
//- Reconstruct volume internal field
template<class Type>
tmp<DimensionedField<Type, volMesh>>
reconstructFvVolumeInternalField
(
const IOobject& fieldIoObject,
const PtrList<DimensionedField<Type, volMesh>>& procFields
) const;
//- Read and reconstruct volume internal field
template<class Type>
tmp<DimensionedField<Type, volMesh>>
reconstructFvVolumeInternalField(const IOobject& fieldIoObject) const;
//- Reconstruct volume field
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
reconstructFvVolumeField
(
const IOobject& fieldIoObject,
const PtrList<GeometricField<Type, fvPatchField, volMesh>>&
) const;
//- Read and reconstruct volume field
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
reconstructFvVolumeField(const IOobject& fieldIoObject) const;
//- Reconstruct surface field
template<class Type>
tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
reconstructFvSurfaceField
(
const IOobject& fieldIoObject,
const PtrList<GeometricField<Type, fvsPatchField, surfaceMesh>>&
) const;
//- Read and reconstruct surface field
template<class Type>
tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
reconstructFvSurfaceField(const IOobject& fieldIoObject) const;
//- Read, reconstruct and write all/selected volume internal fields
template<class Type>
void reconstructFvVolumeInternalFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
);
//- Read, reconstruct and write all/selected volume fields
template<class Type>
void reconstructFvVolumeFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
);
//- Read, reconstruct and write all/selected surface fields
template<class Type>
void reconstructFvSurfaceFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
);
// Member Operators
//- Disallow default bitwise assignment
void operator=(const fvFieldReconstructor&) = delete;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "fvFieldReconstructorReconstructFields.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,556 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "fvFieldReconstructor.H"
#include "Time.H"
#include "PtrList.H"
#include "fvPatchFields.H"
#include "emptyFvPatch.H"
#include "emptyFvPatchField.H"
#include "emptyFvsPatchField.H"
#include "processorCyclicFvPatch.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
void Foam::fvFieldReconstructor::rmapFaceToFace
(
Field<Type>& toField,
const Field<Type>& fromField,
const labelUList& addressing,
const bool isFlux
)
{
forAll(addressing, i)
{
toField[mag(addressing[i]) - 1] =
(isFlux && addressing[i] < 0 ? -1 : +1)*fromField[i];
}
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::DimensionedField<Type, Foam::volMesh>>
Foam::fvFieldReconstructor::reconstructFvVolumeInternalField
(
const IOobject& fieldIoObject,
const PtrList<DimensionedField<Type, volMesh>>& procFields
) const
{
// Create the internalField
Field<Type> internalField(completeMesh_.nCells());
forAll(procMeshes_, proci)
{
const DimensionedField<Type, volMesh>& procField = procFields[proci];
// Set the cell values in the reconstructed field
internalField.rmap
(
procField.field(),
cellProcAddressing_[proci]
);
}
return tmp<DimensionedField<Type, volMesh>>
(
new DimensionedField<Type, volMesh>
(
fieldIoObject,
completeMesh_,
procFields[0].dimensions(),
internalField
)
);
}
template<class Type>
Foam::tmp<Foam::DimensionedField<Type, Foam::volMesh>>
Foam::fvFieldReconstructor::reconstructFvVolumeInternalField
(
const IOobject& fieldIoObject
) const
{
PtrList<DimensionedField<Type, volMesh>>
procFields(procMeshes_.size());
forAll(procMeshes_, proci)
{
procFields.set
(
proci,
new DimensionedField<Type, volMesh>
(
IOobject
(
fieldIoObject.name(),
procMeshes_[proci].time().timeName(),
procMeshes_[proci],
IOobject::MUST_READ,
IOobject::NO_WRITE,
false
),
procMeshes_[proci]
)
);
}
return reconstructFvVolumeInternalField
(
IOobject
(
fieldIoObject.name(),
completeMesh_.time().timeName(),
completeMesh_,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
procFields
);
}
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>>
Foam::fvFieldReconstructor::reconstructFvVolumeField
(
const IOobject& fieldIoObject,
const PtrList<GeometricField<Type, fvPatchField, volMesh>>& procFields
) const
{
// Create the internalField
Field<Type> internalField(completeMesh_.nCells());
// Create the patch fields
PtrList<fvPatchField<Type>> patchFields(completeMesh_.boundary().size());
forAll(procFields, proci)
{
const GeometricField<Type, fvPatchField, volMesh>& procField =
procFields[proci];
// Set the cell values in the reconstructed field
internalField.rmap
(
procField.primitiveField(),
cellProcAddressing_[proci]
);
// Set the boundary patch values in the reconstructed field
forAll(procField.boundaryField(), procPatchi)
{
const fvPatch& procPatch =
procMeshes_[proci].boundary()[procPatchi];
const label completePatchi = completePatchID(proci, procPatchi);
if (completePatchi == procPatchi)
{
if (!patchFields(completePatchi))
{
patchFields.set
(
completePatchi,
fvPatchField<Type>::New
(
procField.boundaryField()[procPatchi],
completeMesh_.boundary()[completePatchi],
DimensionedField<Type, volMesh>::null(),
fvPatchFieldReconstructor
(
completeMesh_.boundary()[completePatchi].size()
)
)
);
}
patchFields[completePatchi].rmap
(
procField.boundaryField()[procPatchi],
faceProcAddressingBf_[proci][procPatchi] - 1
);
}
else if (isA<processorCyclicFvPatch>(procPatch))
{
if (!patchFields(completePatchi))
{
patchFields.set
(
completePatchi,
fvPatchField<Type>::New
(
completeMesh_.boundary()[completePatchi].type(),
completeMesh_.boundary()[completePatchi],
DimensionedField<Type, volMesh>::null()
)
);
}
patchFields[completePatchi].rmap
(
procField.boundaryField()[procPatchi],
faceProcAddressingBf_[proci][procPatchi] - 1
);
}
}
}
// Construct and return the field
return tmp<GeometricField<Type, fvPatchField, volMesh>>
(
new GeometricField<Type, fvPatchField, volMesh>
(
fieldIoObject,
completeMesh_,
procFields[0].dimensions(),
internalField,
patchFields
)
);
}
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>>
Foam::fvFieldReconstructor::reconstructFvVolumeField
(
const IOobject& fieldIoObject
) const
{
PtrList<GeometricField<Type, fvPatchField, volMesh>>
procFields(procMeshes_.size());
forAll(procMeshes_, proci)
{
procFields.set
(
proci,
new GeometricField<Type, fvPatchField, volMesh>
(
IOobject
(
fieldIoObject.name(),
procMeshes_[proci].time().timeName(),
procMeshes_[proci],
IOobject::MUST_READ,
IOobject::NO_WRITE,
false
),
procMeshes_[proci]
)
);
}
return reconstructFvVolumeField
(
IOobject
(
fieldIoObject.name(),
completeMesh_.time().timeName(),
completeMesh_,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
procFields
);
}
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
Foam::fvFieldReconstructor::reconstructFvSurfaceField
(
const IOobject& fieldIoObject,
const PtrList<GeometricField<Type, fvsPatchField, surfaceMesh>>& procFields
) const
{
// Create the internalField
Field<Type> internalField(completeMesh_.nInternalFaces());
// Create the patch fields
PtrList<fvsPatchField<Type>> patchFields(completeMesh_.boundary().size());
forAll(procMeshes_, proci)
{
const GeometricField<Type, fvsPatchField, surfaceMesh>& procField =
procFields[proci];
// Set the internal face values in the reconstructed field
rmapFaceToFace
(
internalField,
procField.primitiveField(),
SubList<label>
(
faceProcAddressing_[proci],
procMeshes_[proci].nInternalFaces()
),
isFlux(procFields[proci])
);
// Set the boundary patch values in the reconstructed field
forAll(procField.boundaryField(), procPatchi)
{
const fvPatch& procPatch =
procMeshes_[proci].boundary()[procPatchi];
const label completePatchi = completePatchID(proci, procPatchi);
if (completePatchi == procPatchi)
{
if (!patchFields(completePatchi))
{
patchFields.set
(
completePatchi,
fvsPatchField<Type>::New
(
procField.boundaryField()[procPatchi],
completeMesh_.boundary()[completePatchi],
DimensionedField<Type, surfaceMesh>::null(),
fvPatchFieldReconstructor
(
completeMesh_.boundary()[completePatchi].size()
)
)
);
}
patchFields[completePatchi].rmap
(
procField.boundaryField()[procPatchi],
faceProcAddressingBf_[proci][procPatchi] - 1
);
}
else if (isA<processorCyclicFvPatch>(procPatch))
{
if (!patchFields(completePatchi))
{
patchFields.set
(
completePatchi,
fvsPatchField<Type>::New
(
completeMesh_.boundary()[completePatchi].type(),
completeMesh_.boundary()[completePatchi],
DimensionedField<Type, surfaceMesh>::null()
)
);
}
patchFields[completePatchi].rmap
(
procField.boundaryField()[procPatchi],
faceProcAddressingBf_[proci][procPatchi] - 1
);
}
else if (isA<processorFvPatch>(procPatch))
{
rmapFaceToFace
(
internalField,
procField.boundaryField()[procPatchi],
faceProcAddressingBf_[proci][procPatchi],
isFlux(procFields[proci])
);
}
}
}
// Construct and return the field
return tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
(
new GeometricField<Type, fvsPatchField, surfaceMesh>
(
fieldIoObject,
completeMesh_,
procFields[0].dimensions(),
internalField,
patchFields
)
);
}
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
Foam::fvFieldReconstructor::reconstructFvSurfaceField
(
const IOobject& fieldIoObject
) const
{
PtrList<GeometricField<Type, fvsPatchField, surfaceMesh>>
procFields(procMeshes_.size());
forAll(procMeshes_, proci)
{
procFields.set
(
proci,
new GeometricField<Type, fvsPatchField, surfaceMesh>
(
IOobject
(
fieldIoObject.name(),
procMeshes_[proci].time().timeName(),
procMeshes_[proci],
IOobject::MUST_READ,
IOobject::NO_WRITE,
false
),
procMeshes_[proci]
)
);
}
return reconstructFvSurfaceField
(
IOobject
(
fieldIoObject.name(),
completeMesh_.time().timeName(),
completeMesh_,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
procFields
);
}
template<class Type>
void Foam::fvFieldReconstructor::reconstructFvVolumeInternalFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
)
{
const word& fieldClassName = DimensionedField<Type, volMesh>::typeName;
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
selectedFields.empty()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructFvVolumeInternalField<Type>(*fieldIter())().write();
nReconstructed_++;
}
}
Info<< endl;
}
}
template<class Type>
void Foam::fvFieldReconstructor::reconstructFvVolumeFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
)
{
const word& fieldClassName =
GeometricField<Type, fvPatchField, volMesh>::typeName;
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
selectedFields.empty()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructFvVolumeField<Type>(*fieldIter())().write();
nReconstructed_++;
}
}
Info<< endl;
}
}
template<class Type>
void Foam::fvFieldReconstructor::reconstructFvSurfaceFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
)
{
const word& fieldClassName =
GeometricField<Type, fvsPatchField, surfaceMesh>::typeName;
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
selectedFields.empty()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructFvSurfaceField<Type>(*fieldIter())().write();
nReconstructed_++;
}
}
Info<< endl;
}
}
// ************************************************************************* //

View File

@ -0,0 +1,93 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "pointFieldReconstructor.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::pointFieldReconstructor::pointFieldReconstructor
(
const pointMesh& completeMesh,
const PtrList<pointMesh>& procMeshes,
const labelListList& pointProcAddressing
)
:
completeMesh_(completeMesh),
procMeshes_(procMeshes),
pointProcAddressing_(pointProcAddressing),
patchPointAddressing_(procMeshes.size()),
nReconstructed_(0)
{
// Inverse-addressing of the patch point labels.
labelList pointMap(completeMesh_.size(), -1);
// Create the pointPatch addressing
forAll(procMeshes_, proci)
{
const pointMesh& procMesh = procMeshes_[proci];
patchPointAddressing_[proci].setSize(procMesh.boundary().size());
forAll(procMesh.boundary(), patchi)
{
if (patchi < completeMesh_.boundary().size())
{
labelList& procPatchAddr = patchPointAddressing_[proci][patchi];
procPatchAddr.setSize(procMesh.boundary()[patchi].size(), -1);
const labelList& patchPointLabels =
completeMesh_.boundary()[patchi].meshPoints();
// Create the inverse-addressing of the patch point labels.
forAll(patchPointLabels, pointi)
{
pointMap[patchPointLabels[pointi]] = pointi;
}
const labelList& procPatchPoints =
procMesh.boundary()[patchi].meshPoints();
forAll(procPatchPoints, pointi)
{
procPatchAddr[pointi] =
pointMap
[
pointProcAddressing_[proci][procPatchPoints[pointi]]
];
}
if (procPatchAddr.size() && min(procPatchAddr) < 0)
{
FatalErrorInFunction
<< "Incomplete patch point addressing"
<< abort(FatalError);
}
}
}
}
}
// ************************************************************************* //

View File

@ -0,0 +1,150 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::pointFieldReconstructor
Description
Point field reconstructor.
SourceFiles
pointFieldReconstructor.C
\*---------------------------------------------------------------------------*/
#ifndef pointFieldReconstructor_H
#define pointFieldReconstructor_H
#include "pointMesh.H"
#include "pointFields.H"
#include "pointPatchFieldMapper.H"
#include "setSizeFieldMapper.H"
#include "IOobjectList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class pointFieldReconstructor Declaration
\*---------------------------------------------------------------------------*/
class pointFieldReconstructor
{
// Private Data
//- Reconstructed mesh reference
const pointMesh& completeMesh_;
//- List of processor meshes
const PtrList<pointMesh>& procMeshes_;
//- List of processor point addressing lists
const labelListList& pointProcAddressing_;
//- Point patch addressing
labelListListList patchPointAddressing_;
//- Number of fields reconstructed
label nReconstructed_;
public:
//- Mapper for sizing only - does not do any actual mapping.
class pointPatchFieldReconstructor
:
public pointPatchFieldMapper,
public setSizeFieldMapper
{
public:
// Constructors
//- Construct given size
pointPatchFieldReconstructor(const label size)
:
setSizeFieldMapper(size)
{}
};
// Constructors
//- Construct from components
pointFieldReconstructor
(
const pointMesh& mesh,
const PtrList<pointMesh>& procMeshes,
const labelListList& pointProcAddressing
);
//- Disallow default bitwise copy construction
pointFieldReconstructor(const pointFieldReconstructor&) = delete;
// Member Functions
//- Return number of fields reconstructed
label nReconstructed() const
{
return nReconstructed_;
}
//- Reconstruct field
template<class Type>
tmp<GeometricField<Type, pointPatchField, pointMesh>>
reconstructField(const IOobject& fieldIoObject);
//- Reconstruct and write all fields
template<class Type>
void reconstructFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
);
// Member Operators
//- Disallow default bitwise assignment
void operator=(const pointFieldReconstructor&) = delete;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "pointFieldReconstructorReconstructFields.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,180 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "pointFieldReconstructor.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::pointPatchField, Foam::pointMesh>>
Foam::pointFieldReconstructor::reconstructField(const IOobject& fieldIoObject)
{
// Read the field for all the processors
PtrList<GeometricField<Type, pointPatchField, pointMesh>> procFields
(
procMeshes_.size()
);
forAll(procMeshes_, proci)
{
procFields.set
(
proci,
new GeometricField<Type, pointPatchField, pointMesh>
(
IOobject
(
fieldIoObject.name(),
procMeshes_[proci]().time().timeName(),
procMeshes_[proci](),
IOobject::MUST_READ,
IOobject::NO_WRITE
),
procMeshes_[proci]
)
);
}
// Create the internalField
Field<Type> internalField(completeMesh_.size());
// Create the patch fields
PtrList<pointPatchField<Type>> patchFields(completeMesh_.boundary().size());
forAll(procMeshes_, proci)
{
const GeometricField<Type, pointPatchField, pointMesh>&
procField = procFields[proci];
// Get processor-to-global addressing for use in rmap
const labelList& procToGlobalAddr = pointProcAddressing_[proci];
// Set the cell values in the reconstructed field
internalField.rmap
(
procField.primitiveField(),
procToGlobalAddr
);
// Set the boundary patch values in the reconstructed field
forAll(procField.boundaryField(), patchi)
{
// Get patch index of the original patch
const label curBPatch =
patchi < completeMesh_.boundary().size() ? patchi : -1;
// check if the boundary patch is not a processor patch
if (curBPatch != -1)
{
if (!patchFields(curBPatch))
{
patchFields.set(
curBPatch,
pointPatchField<Type>::New
(
procField.boundaryField()[patchi],
completeMesh_.boundary()[curBPatch],
DimensionedField<Type, pointMesh>::null(),
pointPatchFieldReconstructor
(
completeMesh_.boundary()[curBPatch].size()
)
)
);
}
patchFields[curBPatch].rmap
(
procField.boundaryField()[patchi],
patchPointAddressing_[proci][patchi]
);
}
}
}
// Construct and write the field
// setting the internalField and patchFields
return tmp<GeometricField<Type, pointPatchField, pointMesh>>
(
new GeometricField<Type, pointPatchField, pointMesh>
(
IOobject
(
fieldIoObject.name(),
completeMesh_().time().timeName(),
completeMesh_(),
IOobject::NO_READ,
IOobject::NO_WRITE
),
completeMesh_,
procFields[0].dimensions(),
internalField,
patchFields
)
);
}
// Reconstruct and write all point fields
template<class Type>
void Foam::pointFieldReconstructor::reconstructFields
(
const IOobjectList& objects,
const HashSet<word>& selectedFields
)
{
word fieldClassName
(
GeometricField<Type, pointPatchField, pointMesh>::typeName
);
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
!selectedFields.size()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructField<Type>(*fieldIter())().write();
}
}
Info<< endl;
}
}
// ************************************************************************* //

View File

@ -0,0 +1,117 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
InClass
Foam::reconstructLagrangian
Description
SourceFiles
reconstructLagrangianPositions.C
reconstructLagrangianFields.C
\*---------------------------------------------------------------------------*/
#ifndef reconstructLagrangian_H
#define reconstructLagrangian_H
#include "cloud.H"
#include "polyMesh.H"
#include "IOobjectList.H"
#include "CompactIOField.H"
#include "fvMesh.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
void reconstructLagrangianPositions
(
const polyMesh& mesh,
const word& cloudName,
const PtrList<fvMesh>& meshes,
const labelListList& faceProcAddressing,
const labelListList& cellProcAddressing
);
template<class Type>
tmp<IOField<Type>> reconstructLagrangianField
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const word& fieldName
);
template<class Type>
tmp<CompactIOField<Field<Type>, Type>> reconstructLagrangianFieldField
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const word& fieldName
);
template<class Type>
void reconstructLagrangianFields
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const IOobjectList& objects,
const HashSet<word>& selectedFields
);
template<class Type>
void reconstructLagrangianFieldFields
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const IOobjectList& objects,
const HashSet<word>& selectedFields
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "reconstructLagrangianFields.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,273 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "IOField.H"
#include "CompactIOField.H"
#include "Time.H"
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
template<class Type>
Foam::tmp<Foam::IOField<Type>> Foam::reconstructLagrangianField
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const word& fieldName
)
{
// Construct empty field on mesh
tmp<IOField<Type>> tfield
(
new IOField<Type>
(
IOobject
(
fieldName,
mesh.time().timeName(),
cloud::prefix/cloudName,
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
Field<Type>(0)
)
);
Field<Type>& field = tfield.ref();
forAll(meshes, i)
{
// Check object on local mesh
typeIOobject<IOField<Type>> localIOobject
(
fieldName,
meshes[i].time().timeName(),
cloud::prefix/cloudName,
meshes[i],
IOobject::MUST_READ,
IOobject::NO_WRITE
);
if (localIOobject.headerOk())
{
IOField<Type> fieldi(localIOobject);
label offset = field.size();
field.setSize(offset + fieldi.size());
forAll(fieldi, j)
{
field[offset + j] = fieldi[j];
}
}
}
return tfield;
}
template<class Type>
Foam::tmp<Foam::CompactIOField<Foam::Field<Type>, Type>>
Foam::reconstructLagrangianFieldField
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const word& fieldName
)
{
// Construct empty field on mesh
tmp<CompactIOField<Field<Type>, Type >> tfield
(
new CompactIOField<Field<Type>, Type>
(
IOobject
(
fieldName,
mesh.time().timeName(),
cloud::prefix/cloudName,
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
Field<Field<Type>>(0)
)
);
Field<Field<Type>>& field = tfield.ref();
forAll(meshes, i)
{
// No type checking is done to handle CompactIOField and IOField
IOobject localIOobject
(
fieldName,
meshes[i].time().timeName(),
cloud::prefix/cloudName,
meshes[i],
IOobject::MUST_READ,
IOobject::NO_WRITE
);
if (localIOobject.headerOk())
{
CompactIOField<Field<Type>, Type> fieldi(localIOobject);
label offset = field.size();
field.setSize(offset + fieldi.size());
forAll(fieldi, j)
{
field[offset + j] = fieldi[j];
}
}
}
return tfield;
}
template<class Type>
void Foam::reconstructLagrangianFields
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const IOobjectList& objects,
const HashSet<word>& selectedFields
)
{
const word fieldClassName(IOField<Type>::typeName);
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing lagrangian "
<< fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
selectedFields.empty()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructLagrangianField<Type>
(
cloudName,
mesh,
meshes,
fieldIter()->name()
)().write();
}
}
Info<< endl;
}
}
template<class Type>
void Foam::reconstructLagrangianFieldFields
(
const word& cloudName,
const polyMesh& mesh,
const PtrList<fvMesh>& meshes,
const IOobjectList& objects,
const HashSet<word>& selectedFields
)
{
{
const word fieldClassName(CompactIOField<Field<Type>, Type>::typeName);
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing lagrangian "
<< fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
selectedFields.empty()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructLagrangianFieldField<Type>
(
cloudName,
mesh,
meshes,
fieldIter()->name()
)().write();
}
}
Info<< endl;
}
}
{
const word fieldClassName(IOField<Field<Type>>::typeName);
IOobjectList fields = objects.lookupClass(fieldClassName);
if (fields.size())
{
Info<< " Reconstructing lagrangian "
<< fieldClassName << "s\n" << endl;
forAllConstIter(IOobjectList, fields, fieldIter)
{
if
(
selectedFields.empty()
|| selectedFields.found(fieldIter()->name())
)
{
Info<< " " << fieldIter()->name() << endl;
reconstructLagrangianFieldField<Type>
(
cloudName,
mesh,
meshes,
fieldIter()->name()
)().write();
}
}
Info<< endl;
}
}
}
// ************************************************************************* //

View File

@ -0,0 +1,83 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "reconstructLagrangian.H"
#include "labelIOList.H"
#include "passiveParticleCloud.H"
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
void Foam::reconstructLagrangianPositions
(
const polyMesh& mesh,
const word& cloudName,
const PtrList<fvMesh>& meshes,
const labelListList& faceProcAddressing,
const labelListList& cellProcAddressing
)
{
passiveParticleCloud lagrangianPositions
(
mesh,
cloudName,
IDLList<passiveParticle>()
);
forAll(meshes, i)
{
const labelList& cellMap = cellProcAddressing[i];
const labelList& faceMap = faceProcAddressing[i];
Cloud<passiveParticle> lpi(meshes[i], cloudName, false);
forAllConstIter(Cloud<passiveParticle>, lpi, iter)
{
const passiveParticle& ppi = iter();
const label mappedCell = cellMap[ppi.cell()];
// Inverting sign if necessary and subtracting 1 from
// faceProcAddressing
label mappedTetFace = mag(faceMap[ppi.tetFace()]) - 1;
lagrangianPositions.append
(
new passiveParticle
(
mesh,
ppi.coordinates(),
mappedCell,
mappedTetFace,
ppi.procTetPt(mesh, mappedCell, mappedTetFace)
)
);
}
}
IOPosition<Cloud<passiveParticle>>(lagrangianPositions).write();
}
// ************************************************************************* //

View File

@ -32,21 +32,15 @@ Description
#include "argList.H"
#include "timeSelector.H"
#include "fvCFD.H"
#include "IOobjectList.H"
#include "processorMeshes.H"
#include "processorRunTimes.H"
#include "domainDecomposition.H"
#include "regionProperties.H"
#include "fvFieldReconstructor.H"
#include "pointFieldReconstructor.H"
#include "reconstructLagrangian.H"
#include "cellSet.H"
#include "faceSet.H"
#include "pointSet.H"
#include "hexRef8Data.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
@ -120,7 +114,6 @@ int main(int argc, char *argv[])
);
#include "setRootCase.H"
#include "createTime.H"
HashSet<word> selectedFields;
if (args.optionFound("fields"))
@ -144,7 +137,6 @@ int main(int argc, char *argv[])
<< nl << endl;
}
const bool noReconstructSets = args.optionFound("noSets");
if (noReconstructSets)
@ -153,7 +145,6 @@ int main(int argc, char *argv[])
<< nl << endl;
}
HashSet<word> selectedLagrangianFields;
if (args.optionFound("lagrangianFields"))
{
@ -168,12 +159,20 @@ int main(int argc, char *argv[])
args.optionLookup("lagrangianFields")() >> selectedLagrangianFields;
}
const wordList regionNames(selectRegionNames(args, runTime));
// Set time from database
Info<< "Create time\n" << endl;
processorRunTimes runTimes(Foam::Time::controlDictName, args);
// Allow override of time
const instantList times = runTimes.selectProc(args);
// Get region names
const wordList regionNames =
selectRegionNames(args, runTimes.procTimes()[0]);
// Determine the processor count
const label nProcs =
fileHandler().nProcs(args.path(), regionDir(regionNames[0]));
if (!nProcs)
{
FatalErrorInFunction
@ -184,38 +183,12 @@ int main(int argc, char *argv[])
// Warn fileHandler of number of processors
const_cast<fileOperation&>(fileHandler()).setNProcs(nProcs);
// Create the processor databases
PtrList<Time> databases(nProcs);
forAll(databases, proci)
{
databases.set
(
proci,
new Time
(
Time::controlDictName,
args.rootPath(),
args.caseName()/fileName(word("processor") + name(proci))
)
);
}
// Use the times list from the master processor
// and select a subset based on the command-line options
instantList timeDirs = timeSelector::select
(
databases[0].times(),
args
);
// Note that we do not set the runTime time so it is still the
// one set through the controlDict. The -time option
// only affects the selected set of times from processor0.
// - can be illogical
// + any point motion handled through mesh.readUpdate
if (timeDirs.empty())
if (times.empty())
{
WarningInFunction << "No times selected" << endl;
exit(1);
@ -226,130 +199,76 @@ int main(int argc, char *argv[])
instantList masterTimeDirs;
if (newTimes)
{
masterTimeDirs = runTime.times();
masterTimeDirs = runTimes.completeTime().times();
}
HashSet<word> masterTimeDirSet(2*masterTimeDirs.size());
forAll(masterTimeDirs, i)
{
masterTimeDirSet.insert(masterTimeDirs[i].name());
}
if
(
newTimes
&& regionNames.size() == 1
&& regionNames[0] == fvMesh::defaultRegion
&& haveAllTimes(masterTimeDirSet, timeDirs)
&& haveAllTimes(masterTimeDirSet, times)
)
{
Info<< "All times already reconstructed.\n\nEnd\n" << endl;
return 0;
}
// Set all times on processor meshes equal to reconstructed mesh
forAll(databases, proci)
{
databases[proci].setTime(runTime);
}
// Reconstruct all regions
forAll(regionNames, regioni)
{
const word& regionName = regionNames[regioni];
const word& regionDir = Foam::regionDir(regionName);
Info<< "\n\nReconstructing fields for mesh " << regionName << nl
<< endl;
fvMesh mesh
(
IOobject
(
regionName,
runTime.timeName(),
runTime,
IOobject::MUST_READ
),
false
);
// Read all meshes and addressing to reconstructed mesh
processorMeshes procMeshes(databases, regionName);
// Check face addressing for meshes that have been decomposed
// with a very old foam version
#include "checkFaceAddressingComp.H"
// Create meshes
Info<< "\n\nReconstructing fields for mesh " << regionName
<< nl << endl;
domainDecomposition meshes(runTimes, regionName);
meshes.readComplete();
meshes.readProcs();
meshes.readAddressing();
meshes.readUpdate();
// Loop over all times
forAll(timeDirs, timei)
forAll(times, timei)
{
if (newTimes && masterTimeDirSet.found(timeDirs[timei].name()))
if (newTimes && masterTimeDirSet.found(times[timei].name()))
{
Info<< "Skipping time " << timeDirs[timei].name()
Info<< "Skipping time " << times[timei].name()
<< endl << endl;
continue;
}
// Set the time
runTimes.setTime(times[timei], timei);
// Set time for global database
runTime.setTime(timeDirs[timei], timei);
Info<< "Time = " << runTimes.completeTime().userTimeName()
<< nl << endl;
Info<< "Time = " << runTime.userTimeName() << endl << endl;
// Set time for all databases
forAll(databases, proci)
// Update the meshes
const fvMesh::readUpdateState state = meshes.readUpdate();
if (state == fvMesh::POINTS_MOVED)
{
databases[proci].setTime(timeDirs[timei], timei);
meshes.writeComplete(false);
}
// Check if any new meshes need to be read.
fvMesh::readUpdateState meshStat = mesh.readUpdate();
fvMesh::readUpdateState procStat = procMeshes.readUpdate();
if
(
meshStat != procStat
// Allowed status difference. Mesh has moved in parallel, but
// the change has not yet been reconstructed.
&& !(
meshStat == fvMesh::UNCHANGED
&& procStat == fvMesh::POINTS_MOVED
)
// Allowed status difference. Mesh has changed topology in both
// cases, and this has lead to new processor interfaces
// (probably). These interfaces are compatible with a
// reconstructed mesh in which patches have not changed.
&& !(
meshStat == fvMesh::TOPO_CHANGE
&& procStat == fvMesh::TOPO_PATCH_CHANGE
)
state == fvMesh::TOPO_CHANGE
|| state == fvMesh::TOPO_PATCH_CHANGE
)
{
WarningInFunction
<< "readUpdate for the reconstructed mesh:"
<< meshStat << nl
<< "readUpdate for the processor meshes :"
<< procStat << nl
<< "These should be equal or your addressing"
<< " might be incorrect."
<< " Please check your time directories for any "
<< "mesh directories." << endl;
}
// Reconstruct and write the points for moving mesh cases
if (procStat == fvMesh::POINTS_MOVED)
{
procMeshes.reconstructPoints(mesh);
meshes.writeComplete(!noReconstructSets);
}
// Get list of objects from processor0 database
IOobjectList objects
(
procMeshes.meshes()[0],
databases[0].timeName()
meshes.procMeshes()[0],
runTimes.procTimes()[0].timeName()
);
if (!noFields)
@ -359,10 +278,11 @@ int main(int argc, char *argv[])
fvFieldReconstructor fvReconstructor
(
mesh,
procMeshes.meshes(),
procMeshes.faceProcAddressing(),
procMeshes.cellProcAddressing()
meshes.completeMesh(),
meshes.procMeshes(),
meshes.procFaceAddressing(),
meshes.procCellAddressing(),
meshes.procFaceAddressingBf()
);
fvReconstructor.reconstructFvVolumeInternalFields<scalar>
@ -454,23 +374,23 @@ int main(int argc, char *argv[])
{
Info<< "Reconstructing point fields" << nl << endl;
const pointMesh& pMesh = pointMesh::New(mesh);
PtrList<pointMesh> pMeshes(procMeshes.meshes().size());
forAll(pMeshes, proci)
const pointMesh& completePMesh =
pointMesh::New(meshes.completeMesh());
PtrList<pointMesh> procPMeshes(nProcs);
forAll(procPMeshes, proci)
{
pMeshes.set
procPMeshes.set
(
proci,
new pointMesh(procMeshes.meshes()[proci])
new pointMesh(meshes.procMeshes()[proci])
);
}
pointFieldReconstructor pointReconstructor
(
pMesh,
pMeshes,
procMeshes.pointProcAddressing()
completePMesh,
procPMeshes,
meshes.procPointAddressing()
);
pointReconstructor.reconstructFields<scalar>
@ -517,15 +437,15 @@ int main(int argc, char *argv[])
{
HashTable<IOobjectList> cloudObjects;
forAll(databases, proci)
forAll(runTimes.procTimes(), proci)
{
fileName lagrangianDir
(
fileHandler().filePath
(
databases[proci].timePath()
/ regionDir
/ cloud::prefix
runTimes.procTimes()[proci].timePath()
/regionDir
/cloud::prefix
)
);
@ -551,8 +471,8 @@ int main(int argc, char *argv[])
// Do local scan for valid cloud objects
IOobjectList sprayObjs
(
procMeshes.meshes()[proci],
databases[proci].timeName(),
meshes.procMeshes()[proci],
runTimes.procTimes()[proci].timeName(),
cloud::prefix/cloudDirs[i]
);
@ -567,16 +487,13 @@ int main(int argc, char *argv[])
}
}
if (cloudObjects.size())
{
// Pass2: reconstruct the cloud
forAllConstIter(HashTable<IOobjectList>, cloudObjects, iter)
{
const word cloudName = string::validate<word>
(
iter.key()
);
const word cloudName =
string::validate<word>(iter.key());
// Objects (on arbitrary processor)
const IOobjectList& sprayObjs = iter();
@ -586,105 +503,105 @@ int main(int argc, char *argv[])
reconstructLagrangianPositions
(
mesh,
meshes.completeMesh(),
cloudName,
procMeshes.meshes(),
procMeshes.faceProcAddressing(),
procMeshes.cellProcAddressing()
meshes.procMeshes(),
meshes.procFaceAddressing(),
meshes.procCellAddressing()
);
reconstructLagrangianFields<label>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFieldFields<label>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFields<scalar>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFieldFields<scalar>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFields<vector>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFieldFields<vector>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFields<sphericalTensor>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFieldFields<sphericalTensor>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFields<symmTensor>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFieldFields<symmTensor>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFields<tensor>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
reconstructLagrangianFieldFields<tensor>
(
cloudName,
mesh,
procMeshes.meshes(),
meshes.completeMesh(),
meshes.procMeshes(),
sprayObjs,
selectedLagrangianFields
);
@ -696,278 +613,6 @@ int main(int argc, char *argv[])
}
}
if (!noReconstructSets)
{
// Scan to find all sets
HashTable<label> cSetNames;
HashTable<label> fSetNames;
HashTable<label> pSetNames;
forAll(procMeshes.meshes(), proci)
{
const fvMesh& procMesh = procMeshes.meshes()[proci];
// Note: look at sets in current time only or between
// mesh and current time?. For now current time. This will
// miss out on sets in intermediate times that have not
// been reconstructed.
IOobjectList objects
(
procMesh,
databases[0].timeName(), // procMesh.facesInstance()
polyMesh::meshSubDir/"sets"
);
IOobjectList cSets(objects.lookupClass(cellSet::typeName));
forAllConstIter(IOobjectList, cSets, iter)
{
cSetNames.insert(iter.key(), cSetNames.size());
}
IOobjectList fSets(objects.lookupClass(faceSet::typeName));
forAllConstIter(IOobjectList, fSets, iter)
{
fSetNames.insert(iter.key(), fSetNames.size());
}
IOobjectList pSets(objects.lookupClass(pointSet::typeName));
forAllConstIter(IOobjectList, pSets, iter)
{
pSetNames.insert(iter.key(), pSetNames.size());
}
}
if (cSetNames.size() || fSetNames.size() || pSetNames.size())
{
// Construct all sets
PtrList<cellSet> cellSets(cSetNames.size());
PtrList<faceSet> faceSets(fSetNames.size());
PtrList<pointSet> pointSets(pSetNames.size());
Info<< "Reconstructing sets:" << endl;
if (cSetNames.size())
{
Info<< " cellSets "
<< cSetNames.sortedToc() << endl;
}
if (fSetNames.size())
{
Info<< " faceSets "
<< fSetNames.sortedToc() << endl;
}
if (pSetNames.size())
{
Info<< " pointSets "
<< pSetNames.sortedToc() << endl;
}
// Load sets
forAll(procMeshes.meshes(), proci)
{
const fvMesh& procMesh = procMeshes.meshes()[proci];
IOobjectList objects
(
procMesh,
databases[0].timeName(),
polyMesh::meshSubDir/"sets"
);
// cellSets
const labelList& cellMap =
procMeshes.cellProcAddressing()[proci];
IOobjectList cSets
(
objects.lookupClass(cellSet::typeName)
);
forAllConstIter(IOobjectList, cSets, iter)
{
// Load cellSet
const cellSet procSet(*iter());
label setI = cSetNames[iter.key()];
if (!cellSets.set(setI))
{
cellSets.set
(
setI,
new cellSet
(
mesh,
iter.key(),
procSet.size()
)
);
}
cellSet& cSet = cellSets[setI];
cSet.instance() = runTime.timeName();
forAllConstIter(cellSet, procSet, iter)
{
cSet.insert(cellMap[iter.key()]);
}
}
// faceSets
const labelList& faceMap =
procMeshes.faceProcAddressing()[proci];
IOobjectList fSets
(
objects.lookupClass(faceSet::typeName)
);
forAllConstIter(IOobjectList, fSets, iter)
{
// Load faceSet
const faceSet procSet(*iter());
label setI = fSetNames[iter.key()];
if (!faceSets.set(setI))
{
faceSets.set
(
setI,
new faceSet
(
mesh,
iter.key(),
procSet.size()
)
);
}
faceSet& fSet = faceSets[setI];
fSet.instance() = runTime.timeName();
forAllConstIter(faceSet, procSet, iter)
{
fSet.insert(mag(faceMap[iter.key()])-1);
}
}
// pointSets
const labelList& pointMap =
procMeshes.pointProcAddressing()[proci];
IOobjectList pSets
(
objects.lookupClass(pointSet::typeName)
);
forAllConstIter(IOobjectList, pSets, iter)
{
// Load pointSet
const pointSet propSet(*iter());
label setI = pSetNames[iter.key()];
if (!pointSets.set(setI))
{
pointSets.set
(
setI,
new pointSet
(
mesh,
iter.key(),
propSet.size()
)
);
}
pointSet& pSet = pointSets[setI];
pSet.instance() = runTime.timeName();
forAllConstIter(pointSet, propSet, iter)
{
pSet.insert(pointMap[iter.key()]);
}
}
}
// Write sets
forAll(cellSets, i)
{
cellSets[i].write();
}
forAll(faceSets, i)
{
faceSets[i].write();
}
forAll(pointSets, i)
{
pointSets[i].write();
}
}
}
// Reconstruct refinement data
{
PtrList<hexRef8Data> procData(procMeshes.meshes().size());
forAll(procMeshes.meshes(), procI)
{
const fvMesh& procMesh = procMeshes.meshes()[procI];
procData.set
(
procI,
new hexRef8Data
(
IOobject
(
"dummy",
procMesh.time().timeName(),
polyMesh::meshSubDir,
procMesh,
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE,
false
)
)
);
}
// Combine individual parts
const PtrList<labelIOList>& cellAddr =
procMeshes.cellProcAddressing();
UPtrList<const labelList> cellMaps(cellAddr.size());
forAll(cellAddr, i)
{
cellMaps.set(i, &cellAddr[i]);
}
const PtrList<labelIOList>& pointAddr =
procMeshes.pointProcAddressing();
UPtrList<const labelList> pointMaps(pointAddr.size());
forAll(pointAddr, i)
{
pointMaps.set(i, &pointAddr[i]);
}
UPtrList<const hexRef8Data> procRefs(procData.size());
forAll(procData, i)
{
procRefs.set(i, &procData[i]);
}
hexRef8Data
(
IOobject
(
"dummy",
mesh.time().timeName(),
polyMesh::meshSubDir,
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
cellMaps,
pointMaps,
procRefs
).write();
}
// If there is a "uniform" directory in the time region
// directory copy from the master processor
{
@ -975,13 +620,17 @@ int main(int argc, char *argv[])
(
fileHandler().filePath
(
databases[0].timePath()/regionDir/"uniform"
runTimes.procTimes()[0].timePath()/regionDir/"uniform"
)
);
if (!uniformDir0.empty() && fileHandler().isDir(uniformDir0))
{
fileHandler().cp(uniformDir0, runTime.timePath()/regionDir);
fileHandler().cp
(
uniformDir0,
runTimes.completeTime().timePath()/regionDir
);
}
}
@ -993,13 +642,17 @@ int main(int argc, char *argv[])
(
fileHandler().filePath
(
databases[0].timePath()/"uniform"
runTimes.procTimes()[0].timePath()/"uniform"
)
);
if (!uniformDir0.empty() && fileHandler().isDir(uniformDir0))
{
fileHandler().cp(uniformDir0, runTime.timePath());
fileHandler().cp
(
uniformDir0,
runTimes.completeTime().timePath()
);
}
}
}