Compare commits

..

18 Commits

Author SHA1 Message Date
c0f56215ec ENH: AMI caching - addeduser-selectable stride 2025-10-23 20:06:28 +01:00
0ba3a91e0c WIP: removed initialisation to zero - not applicable to all use cases 2025-10-16 09:42:16 +01:00
af8bf54b3c WIP: small style changes 2025-10-14 10:39:48 +01:00
5012aa60c8 WIP: AMI cache - further updates - TO SQUASH 2025-10-13 17:36:13 +01:00
b9214f01c4 WIP: cyclicAMIFvPatchField - disable cacheNeighbourField when using AMI cache 2025-10-13 17:36:13 +01:00
2340bd0d74 TUT: mixerVesselAMI2D updated for MD24 cache weights and addressing 2025-10-13 17:36:13 +01:00
ae521bfca2 WIP: AMI caching - updates for parallel operation 2025-10-13 17:36:13 +01:00
a40ad0fa71 BUG: AMI local communicator allocated every step. Fixes #3372 2025-10-13 17:36:13 +01:00
6ba17c03b9 ENH: AMI - refactored caching 2025-10-13 17:36:13 +01:00
c34c720c6a ENH: FaceCellWave: enable through cyclicAMI 2025-10-13 17:36:13 +01:00
e7b9d57158 ENH: AMI - added caching of weights and addressing
Applicable to rotational cases:

- stores AMI weights and addressing on the first revolution
- cached evaluations performed on subsequent revolutions to reduce computational
  costs

Cached values are stored in angular bins, specified using the [optional]
`cacheSize` entry when defining the patch in the polyMesh/boundary file, e.g.

    AMI1
    {
        type            cyclicAMI;
        AMIMethod       faceAreaWeightAMI;
        neighbourPatch  AMI2;

        cacheSize       360; // New entry
        transform       rotational;
        rotationAxis    (0 0 1);
        rotationCentre  (0 0 0);
    }

Note that the transform must also be set to rotational; the additional
`rotationAxis` and `rotationCentre` entries are used to construct a local AMI
co-ordinate system to determine the rotation angle as the mesh moves.

360 bins are created in the example above, equating to a uniform bin width
of 1 degree.
2025-10-13 17:36:13 +01:00
09521d1304 ENH: add default "constant" name for fileOperations::findTimes()
- makes fileHandler calling resemble Time more closely

ENH: provide natural_sort less/greater functions

- more succinct than (compare < 0) or (compare > 0).
  The corresponding wrappers for UList renamed as list_less, etc
  but they were only used in unit tests

ENH: handle (-allRegions | -all-regions, ..) directly when retrieving args

- makes handling independent of aliasing order

STYLE: include <string_view> when including <string>

- the standard does not guarantee which headers (if any) actually
  include string_view
2025-10-13 17:36:13 +01:00
e02b4be7ca ENH: allow delayed construction of regionFaModel in volume bcs (#3419)
Ideally wish to delay construction of the finite-area mesh until
  updateCoeffs(), when it is actually needed.

  This helps avoid difficult to trace errors and avoids inadvertent
  parallel communication during construction from a dictionary.

  However, lazy evaluation fails at the later stage while attempting
  to load the corresponding initial fields, since the timeIndex has
  already advanced.

  Use the newly introduced regionModels::allowFaModels() switch
  to locally enable or disable lazy evaluation as needed.

  This allows selective disabling (eg, in post-processing utilities)
  where loading the volume field should not activate the associated
  regionFaModel and loading a finite-area mesh too (which may not
  exist or be defined at that point).
2025-10-13 17:36:13 +01:00
c7b5f1e3eb ENH: added clean up function remove0DirFields (RunFunctions)
- less typing than before and avoids relying on bash-specific behaviour
  (fixes #3448)

ENH: add -region support for cleanFaMesh and cleanPolyMesh

CONFIG: add bash completion help for -area-region

ENH: general improvements for regionProperties

- robustness and failsafe for foamListRegions, regionProperties
- additional global model switches for regionModels
2025-10-13 17:36:13 +01:00
ccb57c0499 ENH: increase some string_view coverage
- remove some stdFoam::span<char> handling, which was just a stop-gap
  measure
2025-10-13 17:36:12 +01:00
c83bdc1422 STYLE: surface/volume writing function objects
- further hide implementation (cxx ending)

STYLE: better distinction between code/templates for fileFormats/conversion

- for ensight and vtk, a large part of the code is header only.
  Make this easier to identify
2025-10-13 17:36:12 +01:00
22fd0b3e72 ENH: update globalIndex::mpiGather to use MPI intrinsic/user types
- add provisional support for selecting MPI_Gatherv
  when merging fields in the surface writers.
  Uses the 'gatherv' keyword [experimental]

BUG: gatherv/scatterv wrappers used the incorrect data type

- incorrectly checked against UPstream_basic_dataType trait instead of
  UPstream_dataType, which resulted in a count mismatch for
  user-defined types (eg, vector, tensor,...).

  This was not visibly active bug, since previous internal uses of
  gatherv were restricted to primitive types.

COMP: make UPstream dataType traits size(...) constexpr

- allows static asserts and/or compile-time selection
  of code based on size(1)
2025-10-13 17:36:12 +01:00
de91806a24 ENH: GAMGAgglomeration: increase repeatability. Fixes #3450 2025-10-13 17:36:12 +01:00
160 changed files with 4903 additions and 4626 deletions

View File

@ -49,7 +49,8 @@ int main(int argc, char *argv[])
Info<< "mesh 0: " << mesh.sortedNames() << nl;
auto& reg = const_cast<faMeshesRegistry&>(faMeshesRegistry::New(mesh));
faMeshesRegistry& reg =
const_cast<faMeshesRegistry&>(faMeshesRegistry::New(mesh));
// faMeshesRegistry faReg = faMeshesRegistry(mesh);

View File

@ -1,3 +1,3 @@
foamToEnsight-check.cxx
foamToEnsight-check.C
EXE = $(FOAM_USER_APPBIN)/foamToEnsight-check

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022-2025 OpenCFD Ltd.
Copyright (C) 2022-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -203,7 +203,6 @@ int main(int argc, char *argv[])
argList::addVerboseOption();
#include "addAllRegionOptions.H"
#include "addAllFaRegionOptions.H"
argList::addBoolOption
(
@ -252,14 +251,6 @@ int main(int argc, char *argv[])
// Handle -allRegions, -regions, -region
#include "getAllRegionOptions.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
if (!doFiniteArea)
{
areaRegionNames.clear(); // For consistency
}
// ------------------------------------------------------------------------
#include "createNamedMeshes.H"
@ -268,8 +259,8 @@ int main(int argc, char *argv[])
/// #include "createMeshAccounting.H"
PtrList<ensightMesh> ensightMeshes(regionNames.size());
List<PtrList<faMesh>> meshesFa(regionNames.size());
List<PtrList<ensightFaMesh>> ensightMeshesFa(regionNames.size());
PtrList<faMesh> meshesFa(regionNames.size());
PtrList<ensightFaMesh> ensightMeshesFa(regionNames.size());
forAll(regionNames, regioni)
{
@ -282,32 +273,21 @@ int main(int argc, char *argv[])
);
ensightMeshes[regioni].verbose(optVerbose);
if (!doFiniteArea)
if (doFiniteArea)
{
continue;
}
meshesFa[regioni].resize_null(areaRegionNames.size());
ensightMeshesFa[regioni].resize_null(areaRegionNames.size());
autoPtr<faMesh> faMeshPtr(faMesh::TryNew(mesh));
forAll(areaRegionNames, areai)
{
const word& areaName = areaRegionNames[areai];
autoPtr<faMesh> faMeshPtr(faMesh::TryNew(areaName, mesh));
autoPtr<faMesh> faMeshPtr(faMesh::TryNew(mesh));
if (faMeshPtr)
{
meshesFa[regioni].set(areai, std::move(faMeshPtr));
meshesFa.set(regioni, std::move(faMeshPtr));
ensightMeshesFa[regioni].set
ensightMeshesFa.set
(
areai,
new ensightFaMesh(meshesFa[regioni][areai])
regioni,
new ensightFaMesh(meshesFa[regioni])
);
ensightMeshesFa[regioni][areai].verbose(optVerbose);
ensightMeshesFa[regioni].verbose(optVerbose);
}
}
}
@ -315,7 +295,7 @@ int main(int argc, char *argv[])
// ------------------------------------------------------------------------
if (UPstream::master())
if (Pstream::master())
{
Info<< "Checking " << timeDirs.size() << " time steps" << nl;
}
@ -337,20 +317,17 @@ int main(int argc, char *argv[])
auto& ensMesh = ensightMeshes[regioni];
// Finite-area (can be missing)
auto& ensFaMeshes = ensightMeshesFa[regioni];
auto* ensFaMeshPtr = ensightMeshesFa.get(regioni);
if (moving)
{
ensMesh.expire();
ensMesh.correct();
forAll(areaRegionNames, areai)
if (ensFaMeshPtr)
{
if (auto* ensFaMeshPtr = ensFaMeshes.get(areai))
{
ensFaMeshPtr->expire();
ensFaMeshPtr->correct();
}
ensFaMeshPtr->expire();
ensFaMeshPtr->correct();
}
}
@ -363,15 +340,9 @@ int main(int argc, char *argv[])
printInfo(ensMesh, optVerbose);
// finite-area
forAll(areaRegionNames, areai)
if (ensFaMeshPtr)
{
auto* ensFaMeshPtr = ensFaMeshes.get(areai);
if (ensFaMeshPtr)
{
printInfo(*ensFaMeshPtr, optVerbose);
}
printInfo(*ensFaMeshPtr, optVerbose);
}
}
}

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2017 Wikki Ltd
Copyright (C) 2021-2025 OpenCFD Ltd.
Copyright (C) 2021-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -28,7 +28,7 @@ Application
makeFaMesh
Description
Check a finite-area mesh
Check a finiteArea mesh
Original Authors
Zeljko Tukovic, FAMENA
@ -39,7 +39,6 @@ Original Authors
#include "Time.H"
#include "argList.H"
#include "faMesh.H"
#include "faMeshTools.H"
#include "polyMesh.H"
#include "areaFaMesh.H"
#include "edgeFaMesh.H"
@ -48,7 +47,7 @@ Original Authors
#include "processorFaPatch.H"
#include "foamVtkIndPatchWriter.H"
#include "foamVtkLineWriter.H"
#include "regionProperties.H"
#include "faMeshTools.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -60,7 +59,7 @@ int main(int argc, char *argv[])
{
argList::addNote
(
"Check a finite-area mesh"
"Check a finiteArea mesh"
);
argList::addBoolOption
@ -78,15 +77,9 @@ int main(int argc, char *argv[])
);
#include "addRegionOption.H"
#include "addAllFaRegionOptions.H"
#include "addFaRegionOption.H"
#include "setRootCase.H"
#include "createTime.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
// ------------------------------------------------------------------------
#include "createNamedPolyMesh.H"
int geometryOrder(1);
@ -98,41 +91,16 @@ int main(int argc, char *argv[])
faMesh::geometryOrder(geometryOrder);
}
for (const word& areaName : areaRegionNames)
#include "createNamedFaMesh.H"
Info<< "Time = " << runTime.timeName() << nl << endl;
// Mesh information (verbose)
faMeshTools::printMeshChecks(aMesh);
if (args.found("write-vtk"))
{
Info<< "Create faMesh";
if (!polyMesh::regionName(areaName).empty())
{
Info<< " [" << areaName << "]";
}
Info<< " for time = " << runTime.timeName() << nl;
autoPtr<faMesh> faMeshPtr(faMesh::TryNew(areaName, mesh));
if (!faMeshPtr)
{
Info<< " ...failed to create area-mesh";
if (!polyMesh::regionName(areaName).empty())
{
Info<< " [" << areaName << "]";
}
Info<< endl;
continue;
}
else
{
Info<< endl;
}
const auto& aMesh = faMeshPtr();
// Mesh information (verbose)
faMeshTools::printMeshChecks(aMesh);
if (args.found("write-vtk"))
{
#include "faMeshWriteVTK.H"
}
#include "faMeshWriteVTK.H"
}
Info<< "\nEnd\n" << endl;

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2021-2025 OpenCFD Ltd.
Copyright (C) 2021-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
@ -15,16 +15,6 @@ Description
\*---------------------------------------------------------------------------*/
// The base name for output files
word vtkBaseFileName(aMesh.regionName());
if (polyMesh::regionName(vtkBaseFileName).empty())
{
vtkBaseFileName = "finiteArea";
}
Info<< nl;
Info<< "Write faMesh in vtk format:" << nl;
{
// finiteArea - faces
vtk::uindirectPatchWriter writer
@ -33,7 +23,7 @@ Info<< "Write faMesh in vtk format:" << nl;
// vtk::formatType::INLINE_ASCII,
fileName
(
aMesh.time().globalPath() / vtkBaseFileName
aMesh.time().globalPath() / "finiteArea"
)
);
@ -42,7 +32,7 @@ Info<< "Write faMesh in vtk format:" << nl;
globalIndex procAddr(aMesh.nFaces());
labelList cellIDs;
if (UPstream::master())
if (Pstream::master())
{
cellIDs.resize(procAddr.totalSize());
for (const labelRange& range : procAddr.ranges())
@ -63,7 +53,8 @@ Info<< "Write faMesh in vtk format:" << nl;
writer.beginPointData(1);
writer.write("normal", aMesh.pointAreaNormals());
Info<< " " << writer.output().name() << nl;
Info<< nl
<< "Wrote faMesh in vtk format: " << writer.output().name() << nl;
}
{
@ -75,7 +66,7 @@ Info<< "Write faMesh in vtk format:" << nl;
// vtk::formatType::INLINE_ASCII,
fileName
(
aMesh.time().globalPath() / (vtkBaseFileName + "-edges")
aMesh.time().globalPath() / "finiteArea-edges"
)
);
@ -97,7 +88,8 @@ Info<< "Write faMesh in vtk format:" << nl;
writer.beginPointData(1);
writer.write("normal", aMesh.pointAreaNormals());
Info<< " " << writer.output().name() << nl;
Info<< nl
<< "Wrote faMesh in vtk format: " << writer.output().name() << nl;
}
{
@ -116,7 +108,7 @@ Info<< "Write faMesh in vtk format:" << nl;
// vtk::formatType::INLINE_ASCII,
fileName
(
aMesh.time().globalPath() / (vtkBaseFileName + "-edgesCentres")
aMesh.time().globalPath() / "finiteArea-edgesCentres"
)
);
@ -142,7 +134,8 @@ Info<< "Write faMesh in vtk format:" << nl;
writer.write("normal", fld);
}
Info<< " " << writer.output().name() << nl;
Info<< nl
<< "Wrote faMesh in vtk format: " << writer.output().name() << nl;
}

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2023-2025 OpenCFD Ltd.
Copyright (C) 2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
@ -15,7 +15,7 @@ Description
Input
mesh (polyMesh)
meshDefDict (system/finite-area/faMeshDefinition)
meshDefDict (system/faMeshDefintion)
\*---------------------------------------------------------------------------*/
@ -25,16 +25,13 @@ Input
const polyBoundaryMesh& pbm = mesh.boundaryMesh();
labelList patchIDs;
wordRes polyPatchNames;
meshDefDict.readIfPresent("polyMeshPatches", polyPatchNames);
if
const labelList patchIDs
(
wordRes select;
meshDefDict.readIfPresent("polyMeshPatches", select)
)
{
patchIDs = pbm.indices(select, true); // useGroups
}
pbm.indices(polyPatchNames, true) // useGroups
);
label nFaceLabels = 0;
for (const label patchi : patchIDs)

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2021-2025 OpenCFD Ltd.
Copyright (C) 2021-2024 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
@ -13,18 +13,6 @@ License
Description
Search for the appropriate faMeshDefinition dictionary...
Expects to find a faMeshDefinition file in a location specified
by the command-line option, or one of the standard locations.
For default area region (region0):
(1) system/finite-area/faMeshDefinition
(2) system/faMeshDefinition [legacy location - v2312 and earlier]
For a named area region:
(1) system/finite-area/<area-name>/faMeshDefinition
(2) system/finite-area/faMeshDefinition.<area-name>
[symmetric with faOptions]
Required Classes
- Foam::polyMesh
- Foam::IOdictionary
@ -36,185 +24,81 @@ Required Variables
- runTime [Time]
Provided Variables
- faMeshDefinitionPtr [autoPtr<IOdictionary>]
If the dictionary cannot be found, exits with an error message or
reports a warning (dry-run mode)
- meshDefDict [IOdictionary]
- meshDictPtr [autoPtr<IOdictionary>]
\*---------------------------------------------------------------------------*/
const word dictName("faMeshDefinition");
autoPtr<IOdictionary> faMeshDefinitionPtr;
autoPtr<IOdictionary> meshDictPtr;
{
// Exit unless dry-run?
const bool exitOnFailure = (!args.dryRun());
fileName dictPath;
const word& regionDir = Foam::polyMesh::regionName(regionName);
const word& areaRegionDir = Foam::polyMesh::regionName(areaRegionName);
const word& regionDir = polyMesh::regionName(regionName);
const word& areaRegionDir = polyMesh::regionName(areaRegionName);
// Canonical location
fileName dictPath
(
runTime.system()/regionDir
/ faMesh::prefix()/areaRegionDir/dictName
);
// Dictionary specified on the command-line ...
const bool specified = args.readIfPresent("dict", dictPath);
if (specified)
if (args.readIfPresent("dict", dictPath))
{
// Dictionary specified on the command-line ...
if (isDir(dictPath))
{
dictPath /= dictName;
}
}
else if
(
// Dictionary under system/faMeshDefinition ?
// (v2312 and earlier)
areaRegionDir.empty()
&& exists
(
runTime.path()/runTime.caseSystem()
/ regionDir/faMesh::meshSubDir/dictName
)
)
{
// Dictionary present directly in system/ (v2312 and earlier)
dictPath = runTime.system()/regionDir/dictName;
}
else
{
// Use system/finite-area/ directory, with region qualifications
dictPath =
(
runTime.system()/regionDir
/ faMesh::prefix()/areaRegionDir/dictName
);
}
IOobject meshDictIO
(
dictPath,
runTime,
IOobjectOption::MUST_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::NO_REGISTER,
IOobject::MUST_READ,
IOobject::NO_WRITE,
IOobject::NO_REGISTER,
true // is globalObject
);
refPtr<IOobject> meshDictIOPtr;
bool foundIOobject = meshDictIO.typeHeaderOk<IOdictionary>(true);
// For reporting any alternative locations
dictPath.clear();
if (!foundIOobject && !specified)
if (!meshDictIO.typeHeaderOk<IOdictionary>(true))
{
if (!areaRegionDir.empty())
{
// Alternative location
// - system/finite-area/faMeshDefinition.<area-name>
dictPath =
(
runTime.system()/regionDir
/ faMesh::prefix()/IOobject::groupName(dictName, areaRegionDir)
);
}
else
{
// legacy location (v2312 and earlier)
// - system/faMeshDefinition
dictPath = runTime.system()/regionDir/dictName;
}
if (!dictPath.empty())
{
auto ioptr = autoPtr<IOobject>::New
(
dictPath,
runTime,
IOobjectOption::MUST_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::NO_REGISTER,
true // is globalObject
);
foundIOobject =
(
ioptr && ioptr->typeHeaderOk<IOdictionary>(true)
);
if (foundIOobject)
{
// Use if found
meshDictIOPtr = std::move(ioptr);
}
// Generally retain dictPath for error messages,
// but don't mention the really old legacy location
if (areaRegionDir.empty())
{
dictPath.clear();
}
}
FatalErrorInFunction
<< meshDictIO.objectPath() << nl
<< exit(FatalError);
}
// Alternative location or regular location?
if (!meshDictIOPtr)
{
meshDictIOPtr.cref(meshDictIO);
}
const auto& io = meshDictIOPtr();
Info<< "Creating faMesh from definition: "
<< meshDictIO.objectRelPath() << endl;
// Handle final success or failure
if (foundIOobject)
{
Info<< "----------------" << nl
<< "Using area-mesh definition";
if (!areaRegionDir.empty())
{
Info<< " [" << areaRegionName << "]";
}
Info<< nl
<< " " << io.objectRelPath() << nl << endl;
faMeshDefinitionPtr = autoPtr<IOdictionary>::New(io);
}
else
{
// Failure
if (exitOnFailure)
{
auto& err =
FatalErrorInFunction
<< "Did not find area-mesh definition";
if (!areaRegionDir.empty())
{
err<< " [" << areaRegionName << "]";
}
err << nl
<< " " << io.objectRelPath() << nl;
if (!dictPath.empty())
{
err << " " << dictPath << nl;
}
FatalError<< exit(FatalError);
}
else
{
// dry-run: warning
auto& err =
Warning << nl
<< "----------------" << nl
<< "Did not find area-mesh definition";
if (!areaRegionDir.empty())
{
err << " [" << areaRegionName << "]";
}
err << nl
<< " " << io.objectRelPath() << nl;
if (!dictPath.empty())
{
err << " " << dictPath << nl;
}
}
}
meshDictPtr = autoPtr<IOdictionary>::New(meshDictIO);
}
IOdictionary& meshDefDict = *meshDictPtr;
// ************************************************************************* //

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2017 Wikki Ltd
Copyright (C) 2021-2025 OpenCFD Ltd.
Copyright (C) 2021-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -40,6 +40,7 @@ Original Authors
#include "Time.H"
#include "argList.H"
#include "OSspecific.H"
#include "faMesh.H"
#include "faMeshTools.H"
#include "IOdictionary.H"
@ -52,7 +53,6 @@ Original Authors
#include "PtrListOps.H"
#include "foamVtkLineWriter.H"
#include "foamVtkIndPatchWriter.H"
#include "regionProperties.H"
#include "syncTools.H"
#include "OBJstream.H"
@ -104,14 +104,12 @@ int main(int argc, char *argv[])
);
#include "addRegionOption.H"
#include "addAllFaRegionOptions.H"
#include "addFaRegionOption.H"
#include "setRootCase.H"
#include "createTime.H"
#include "createNamedPolyMesh.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
#include "getFaRegionOption.H"
const bool doDecompose = !args.found("no-decompose");
const bool doDecompFields = !args.found("no-fields");
@ -125,83 +123,56 @@ int main(int argc, char *argv[])
Info<< "Skip decompose of finiteArea fields" << nl;
}
// ------------------------------------------------------------------------
// Reading faMeshDefinition dictionary
#include "findMeshDefinitionDict.H"
for (const word& areaRegionName : areaRegionNames)
// Inject/overwrite name for optional 'empty' patch
word patchName;
if (args.readIfPresent("empty-patch", patchName))
{
// Reading faMeshDefinition dictionary
#include "findMeshDefinitionDict.H"
meshDefDict.add("emptyPatch", patchName, true);
}
if (!faMeshDefinitionPtr)
{
if (args.dryRun())
{
break;
}
else
{
FatalErrorInFunction
<< "Did not find area-mesh definition"<< nl
<< exit(FatalError);
}
}
// Preliminary checks
#include "checkPatchTopology.H"
auto& meshDefDict = (*faMeshDefinitionPtr);
Info << "Create areaMesh";
if (!Foam::polyMesh::regionName(areaRegionName).empty())
{
Foam::Info << ' ' << areaRegionName;
}
Info << " for polyMesh at time = " << runTime.timeName() << nl;
// Create
faMesh aMesh(areaRegionName, mesh, meshDefDict);
// Inject/overwrite name for optional 'empty' patch
if (word name; args.readIfPresent("empty-patch", name))
{
meshDefDict.add("emptyPatch", name, true);
}
// Mesh information (less verbose)
faMeshTools::printMeshChecks(aMesh, 0);
// Preliminary checks
#include "checkPatchTopology.H"
if (args.found("write-edges-obj"))
{
#include "faMeshWriteEdgesOBJ.H"
}
Info<< "Create area-mesh";
if (!polyMesh::regionName(areaRegionName).empty())
{
Info<< " [" << areaRegionName << "]";
}
Info<< " with polyMesh at time = " << runTime.timeName() << nl;
if (args.found("write-vtk"))
{
#include "faMeshWriteVTK.H"
}
// Create
faMesh aMesh(areaRegionName, mesh, meshDefDict);
if (args.dryRun())
{
Info<< "\ndry-run: not writing mesh or decomposing fields\n" << nl;
}
else
{
// More precision (for points data)
IOstream::minPrecision(10);
// Mesh information (less verbose)
faMeshTools::printMeshChecks(aMesh, 0);
Info<< nl << "Write finite area mesh." << nl;
aMesh.write();
if (args.found("write-edges-obj"))
{
#include "faMeshWriteEdgesOBJ.H"
}
if (args.found("write-vtk"))
{
#include "faMeshWriteVTK.H"
}
if (args.dryRun())
{
Info<< "\ndry-run: not writing mesh or decomposing fields\n" << nl;
}
else
{
// More precision (for points data)
IOstream::minPrecision(10);
Info<< nl << "Write finite area mesh";
if (const word& name = aMesh.regionName(); !name.empty())
{
Info<< " [" << name << "]";
}
Info<< nl;
aMesh.write();
Info<< endl;
#include "decomposeFaFields.H"
}
Info<< endl;
#include "decomposeFaFields.H"
}
Info<< "\nEnd\n" << endl;

View File

@ -334,7 +334,6 @@ int main(int argc, char *argv[])
);
#include "addAllRegionOptions.H"
#include "addAllFaRegionOptions.H"
argList::addDryRunOption
(
@ -443,6 +442,7 @@ int main(int argc, char *argv[])
bool decomposeFieldsOnly = args.found("fields");
bool forceOverwrite = args.found("force");
// Set time from database
#include "createTime.H"
@ -491,17 +491,9 @@ int main(int argc, char *argv[])
decompDictFile = runTime.globalPath()/decompDictFile;
}
// Handle -allRegions, -regions, -region
// Get region names
#include "getAllRegionOptions.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
if (!doFiniteArea)
{
areaRegionNames.clear(); // For consistency
}
const bool optRegions =
(regionNames.size() != 1 || regionNames[0] != polyMesh::defaultRegion);
@ -821,76 +813,56 @@ int main(int argc, char *argv[])
// Field objects at this time
IOobjectList objects;
// faMesh fields - can have multiple finite-area per volume
HashTable<IOobjectList> faObjects;
IOobjectList faObjects;
if (doDecompFields)
{
// List of volume mesh objects for this instance
objects = IOobjectList(mesh, runTime.timeName());
// List of area mesh objects (assuming single region)
faObjects = IOobjectList
(
mesh.time(),
runTime.timeName(),
faMesh::dbDir(mesh, word::null),
IOobjectOption::NO_REGISTER
);
// Ignore generated fields: (cellDist)
objects.remove("cellDist");
// Lists of finite-area fields
faObjects.reserve(areaRegionNames.size());
for (const word& areaName : areaRegionNames)
{
// The finite-area objects for this area region
IOobjectList objs
(
faMesh::Registry(mesh),
runTime.timeName(),
polyMesh::regionName(areaName),
IOobjectOption::NO_REGISTER
);
if (!objs.empty())
{
faObjects.emplace_set(areaName, std::move(objs));
}
}
}
// Finite area handling
// - all area regions use the same volume decomposition
HashPtrTable<faMeshDecomposition> faMeshDecompHashes;
autoPtr<faMeshDecomposition> faMeshDecompPtr;
if (doFiniteArea)
{
const word boundaryInst =
mesh.time().findInstance(mesh.meshDir(), "boundary");
for (const word& areaName : areaRegionNames)
{
IOobject io
(
"faBoundary",
boundaryInst,
faMesh::meshDir(mesh, areaName),
mesh.time(),
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE,
IOobject::NO_REGISTER
);
IOobject io
(
"faBoundary",
boundaryInst,
faMesh::meshDir(mesh, word::null),
mesh.time(),
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE,
IOobject::NO_REGISTER
);
if (io.typeHeaderOk<faBoundaryMesh>(true))
{
// Always based on the volume decomposition!
faMeshDecompHashes.set
if (io.typeHeaderOk<faBoundaryMesh>(true))
{
// Always based on the volume decomposition!
faMeshDecompPtr.reset
(
new faMeshDecomposition
(
areaName,
autoPtr<faMeshDecomposition>::New
(
areaName,
mesh,
mesh.nProcs(),
mesh.model()
)
);
}
mesh,
mesh.nProcs(),
mesh.model()
)
);
}
}
@ -1094,7 +1066,7 @@ int main(int argc, char *argv[])
processorDb.setTime(runTime);
// Read the mesh
// read the mesh
if (!procMeshList.set(proci))
{
procMeshList.set
@ -1301,15 +1273,12 @@ int main(int argc, char *argv[])
}
// Finite-area mesh and field decomposition
for (auto& iter : faMeshDecompHashes.sorted())
// Finite area mesh and field decomposition
if (faMeshDecompPtr)
{
const word& areaName = iter.key();
Info<< "\nFinite area mesh decomposition" << endl;
faMeshDecomposition& aMesh = *(iter.val());
Info<< "\nFinite area mesh decomposition: "
<< areaName << endl;
faMeshDecomposition& aMesh = faMeshDecompPtr();
aMesh.decomposeMesh();
aMesh.writeDecomposition();
@ -1320,13 +1289,9 @@ int main(int argc, char *argv[])
faFieldDecomposer::fieldsCache areaFieldCache;
if
(
const auto objs = faObjects.cfind(areaName);
doDecompFields && objs.good()
)
if (doDecompFields)
{
areaFieldCache.readAllFields(aMesh, objs.val());
areaFieldCache.readAllFields(aMesh, faObjects);
}
const label nAreaFields = areaFieldCache.size();
@ -1357,7 +1322,7 @@ int main(int argc, char *argv[])
processorDb.setTime(runTime);
// Read the volume mesh
// Read the mesh
fvMesh procFvMesh
(
IOobject
@ -1368,7 +1333,7 @@ int main(int argc, char *argv[])
)
);
faMesh procMesh(areaName, procFvMesh);
faMesh procMesh(procFvMesh);
// // Does not work. HJ, 15/Aug/2017
// const labelIOList& faceProcAddressing =
@ -1427,11 +1392,7 @@ int main(int argc, char *argv[])
boundaryProcAddressing
);
areaFieldCache.decomposeAllFields
(
fieldDecomposer,
args.verbose() // report
);
areaFieldCache.decomposeAllFields(fieldDecomposer);
}
}
}

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2015-2025 OpenCFD Ltd.
Copyright (C) 2015-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -91,7 +91,6 @@ int main(int argc, char *argv[])
argList::noParallel();
#include "addAllRegionOptions.H"
#include "addAllFaRegionOptions.H"
argList::addVerboseOption();
argList::addOption
@ -202,17 +201,9 @@ int main(int argc, char *argv[])
const bool newTimes = args.found("newTimes");
// Handle -allRegions, -regions, -region
// Get region names
#include "getAllRegionOptions.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
if (!doFiniteArea)
{
areaRegionNames.clear(); // For consistency
}
// Determine the processor count
label nProcs{0};
@ -400,30 +391,19 @@ int main(int argc, char *argv[])
IOobjectOption::NO_REGISTER
);
// The finite-area fields (multiple finite-area per volume)
HashTable<IOobjectList> faObjects;
IOobjectList faObjects;
if (doFiniteArea && doFields)
{
// Lists of finite-area fields - scan on processor0
faObjects.reserve(areaRegionNames.size());
for (const word& areaName : areaRegionNames)
{
// The finite-area objects for this area region
IOobjectList objs
(
procMeshes.meshes()[0],
databases[0].timeName(),
faMesh::dbDir(areaName), // local relative to mesh
IOobjectOption::NO_REGISTER
);
if (!objs.empty())
{
faObjects.emplace_set(areaName, std::move(objs));
}
}
// List of area mesh objects (assuming single region)
// - scan on processor0
faObjects = IOobjectList
(
procMeshes.meshes()[0],
databases[0].timeName(),
faMesh::dbDir(word::null), // local relative to mesh
IOobjectOption::NO_REGISTER
);
}
if (doFields)
@ -573,60 +553,42 @@ int main(int argc, char *argv[])
}
// Reconstruct any finite-area fields
if (doFiniteArea)
// If there are any FA fields, reconstruct them
if (!doFiniteArea)
{
bool hadFaFields = false;
for (const word& areaName : areaRegionNames)
{
const auto objs = faObjects.cfind(areaName);
if (!objs.good())
{
continue;
}
const auto& faObjs = objs.val();
if
(
!faObjs.count(fieldTypes::is_area)
&& !faObjs.count<edgeScalarField>()
)
{
continue;
}
hadFaFields = true;
Info<< "Reconstructing finite-area fields ["
<< polyMesh::regionName(areaName)
<< "]" << nl << endl;
const faMesh aMesh(areaName, mesh);
processorFaMeshes procFaMeshes
(
procMeshes.meshes(),
areaName
);
faFieldReconstructor reconstructor
(
aMesh,
procFaMeshes.meshes(),
procFaMeshes.edgeProcAddressing(),
procFaMeshes.faceProcAddressing(),
procFaMeshes.boundaryProcAddressing()
);
reconstructor.reconstructAllFields(faObjs);
}
if (!hadFaFields)
{
Info << "No finite-area fields" << nl << endl;
}
}
else if
(
faObjects.count<areaScalarField>()
|| faObjects.count<areaVectorField>()
|| faObjects.count<areaSphericalTensorField>()
|| faObjects.count<areaSymmTensorField>()
|| faObjects.count<areaTensorField>()
|| faObjects.count<edgeScalarField>()
)
{
Info << "Reconstructing FA fields" << nl << endl;
faMesh aMesh(mesh);
processorFaMeshes procFaMeshes(procMeshes.meshes());
faFieldReconstructor reconstructor
(
aMesh,
procFaMeshes.meshes(),
procFaMeshes.edgeProcAddressing(),
procFaMeshes.faceProcAddressing(),
procFaMeshes.boundaryProcAddressing()
);
reconstructor.reconstructAllFields(faObjects);
}
else
{
Info << "No FA fields" << nl << endl;
}
if (doReconstructSets)
{

View File

@ -119,10 +119,17 @@ autoPtr<faceCoupleInfo> determineCoupledFaces
const polyBoundaryMesh& masterPatches = masterMesh.boundaryMesh();
DynamicList<label> masterFaces(masterMesh.nBoundaryFaces());
DynamicList<label> masterFaces
(
masterMesh.nFaces()
- masterMesh.nInternalFaces()
);
for (const polyPatch& pp : masterPatches)
forAll(masterPatches, patchi)
{
const polyPatch& pp = masterPatches[patchi];
if (isA<processorPolyPatch>(pp))
{
for
@ -132,8 +139,11 @@ autoPtr<faceCoupleInfo> determineCoupledFaces
proci++
)
{
const string toProcString("to" + Foam::name(proci));
if (pp.name().ends_with(toProcString))
const string toProcString("to" + name(proci));
if (
pp.name().rfind(toProcString)
== (pp.name().size()-toProcString.size())
)
{
label meshFacei = pp.start();
forAll(pp, i)
@ -146,6 +156,7 @@ autoPtr<faceCoupleInfo> determineCoupledFaces
}
}
masterFaces.shrink();
// Pick up all patches on meshToAdd ending in "procBoundaryDDDtoYYY"
@ -153,10 +164,16 @@ autoPtr<faceCoupleInfo> determineCoupledFaces
const polyBoundaryMesh& addPatches = meshToAdd.boundaryMesh();
DynamicList<label> addFaces(meshToAdd.nBoundaryFaces());
DynamicList<label> addFaces
(
meshToAdd.nFaces()
- meshToAdd.nInternalFaces()
);
for (const polyPatch& pp : addPatches)
forAll(addPatches, patchi)
{
const polyPatch& pp = addPatches[patchi];
if (isA<processorPolyPatch>(pp))
{
bool isConnected = false;
@ -198,6 +215,7 @@ autoPtr<faceCoupleInfo> determineCoupledFaces
}
}
}
addFaces.shrink();
return autoPtr<faceCoupleInfo>::New
(
@ -484,14 +502,15 @@ void writeMaps
Info<< " pointProcAddressing" << endl;
ioAddr.rename("pointProcAddressing");
labelIOList::writeContents(ioAddr, pointProcAddressing);
IOList<label>::writeContents(ioAddr, pointProcAddressing);
// From processor face to reconstructed mesh face
// - add turning index to faceProcAddressing.
// See reconstructPar for meaning of turning index.
labelList faceProcAddr(faceProcAddressing);
Info<< " faceProcAddressing" << endl;
ioAddr.rename("faceProcAddressing");
labelIOList faceProcAddr(ioAddr, faceProcAddressing);
// Now add turning index to faceProcAddressing.
// See reconstructPar for meaning of turning index.
forAll(faceProcAddr, procFacei)
{
const label masterFacei = faceProcAddr[procFacei];
@ -526,21 +545,19 @@ void writeMaps
}
}
Info<< " faceProcAddressing" << endl;
ioAddr.rename("faceProcAddressing");
labelIOList::writeContents(ioAddr, faceProcAddr);
faceProcAddr.clear();
faceProcAddr.write();
// From processor cell to reconstructed mesh cell
Info<< " cellProcAddressing" << endl;
ioAddr.rename("cellProcAddressing");
labelIOList::writeContents(ioAddr, cellProcAddressing);
IOList<label>::writeContents(ioAddr, cellProcAddressing);
// From processor patch to reconstructed mesh patch
Info<< " boundaryProcAddressing" << endl;
ioAddr.rename("boundaryProcAddressing");
labelIOList::writeContents(ioAddr, boundProcAddressing);
IOList<label>::writeContents(ioAddr, boundProcAddressing);
Info<< endl;
}
@ -628,8 +645,7 @@ void sortFaEdgeMapping
{
// From faMeshReconstructor.C - edge shuffling on patches
Map<label> remapGlobal;
remapGlobal.reserve(onePatch.nEdges());
Map<label> remapGlobal(2*onePatch.nEdges());
for (label edgei = 0; edgei < onePatch.nInternalEdges(); ++edgei)
{
remapGlobal.insert(edgei, remapGlobal.size());
@ -763,11 +779,6 @@ int main(int argc, char *argv[])
);
#include "addAllRegionOptions.H"
#include "addAllFaRegionOptions.H"
// Prevent volume fields [with regionFaModels] from incidental
// triggering finite-area
regionModels::allowFaModels(false);
// Prevent volume BCs from triggering finite-area
regionModels::allowFaModels(false);
@ -830,17 +841,9 @@ int main(int argc, char *argv[])
}
}
// Handle -allRegions, -regions, -region
// Get region names
#include "getAllRegionOptions.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
if (!doFiniteArea)
{
areaRegionNames.clear(); // For consistency
}
// Determine the processor count
label nProcs{0};
@ -932,8 +935,7 @@ int main(int argc, char *argv[])
polyMesh::meshDir(regionName),
databases[0],
IOobject::NO_READ,
IOobject::NO_WRITE,
IOobject::NO_REGISTER
IOobject::NO_WRITE
);
// Problem: faceCompactIOList recognises both 'faceList' and
@ -1403,12 +1405,8 @@ int main(int argc, char *argv[])
// Read finite-area
// For each named area region that exists create a HashTable
// entry that will contain the PtrList for all processors
PtrList<faMesh> procFaMeshes(databases.size());
PtrList<polyMesh> procMeshes(databases.size());
HashTable<PtrList<faMesh>> procAreaRegionMeshes;
procAreaRegionMeshes.reserve(areaRegionNames.size());
forAll(databases, proci)
{
@ -1416,16 +1414,20 @@ int main(int argc, char *argv[])
<< "Read processor mesh: "
<< (databases[proci].caseName()/regionDir) << endl;
const polyMesh& procMesh = procMeshes.emplace
procMeshes.set
(
proci,
IOobject
new polyMesh
(
regionName,
databases[proci].timeName(),
databases[proci]
IOobject
(
regionName,
databases[proci].timeName(),
databases[proci]
)
)
);
const polyMesh& procMesh = procMeshes[proci];
writeMaps
(
@ -1448,67 +1450,38 @@ int main(int argc, char *argv[])
"boundary"
);
for (const word& areaName : areaRegionNames)
IOobject io
(
"faBoundary",
boundaryInst,
faMesh::meshDir(procMesh, word::null),
procMesh.time(),
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE,
IOobject::NO_REGISTER
);
if (io.typeHeaderOk<faBoundaryMesh>(true))
{
IOobject io
(
"faBoundary",
boundaryInst,
faMesh::meshDir(procMesh, areaName),
procMesh.time(),
IOobject::READ_IF_PRESENT,
IOobject::NO_WRITE,
IOobject::NO_REGISTER
);
if (io.typeHeaderOk<faBoundaryMesh>(true))
{
// Always based on the volume decomposition!
auto& procFaMeshes = procAreaRegionMeshes(areaName);
procFaMeshes.resize(databases.size());
procFaMeshes.set
(
proci,
new faMesh(areaName, procMesh)
);
}
// Always based on the volume decomposition!
procFaMeshes.set(proci, new faMesh(procMesh));
}
}
}
// A finite-area mapping exists if procFaMeshes was filled
// Re-read reconstructed polyMesh. Note: could probably be avoided
// by merging into loops above.
refPtr<polyMesh> masterPolyMeshPtr;
if (!procAreaRegionMeshes.empty())
// Finite-area mapping
doFiniteArea = false;
forAll(procFaMeshes, proci)
{
masterPolyMeshPtr.reset
(
new polyMesh
(
IOobject
(
regionName,
runTime.timeName(),
runTime
),
true
)
);
if (procFaMeshes.set(proci))
{
doFiniteArea = true;
}
}
// Process any finite-area meshes
for (const auto& iter : procAreaRegionMeshes.csorted())
if (doFiniteArea)
{
const auto& areaName = iter.key();
const auto& procFaMeshes = iter.val();
const polyMesh& masterMesh = masterPolyMeshPtr();
// Addressing from processor to reconstructed case
labelListList faFaceProcAddressing(nProcs);
labelListList faEdgeProcAddressing(nProcs);
@ -1526,10 +1499,32 @@ int main(int argc, char *argv[])
faBoundProcAddressing[proci] = identity(bm.size());
// Mark processor patches
faBoundProcAddressing[proci].slice(bm.nNonProcessor()) = -1;
for
(
label patchi = bm.nNonProcessor();
patchi < bm.size();
++patchi
)
{
faBoundProcAddressing[proci][patchi] = -1;
}
}
// Re-read reconstructed polyMesh. Note: could probably be avoided
// by merging into loops above.
const polyMesh masterMesh
(
IOobject
(
regionName,
runTime.timeName(),
runTime
),
true
);
// faceProcAddressing
// ~~~~~~~~~~~~~~~~~~
@ -1564,13 +1559,13 @@ int main(int argc, char *argv[])
// Construct without patches
faMesh masterFaMesh
(
areaName,
masterMesh,
std::move(masterFaceLabels),
io
);
const auto& masterPatch = masterFaMesh.patch();
const uindirectPrimitivePatch& masterPatch =
masterFaMesh.patch();
// pointProcAddressing
@ -1582,7 +1577,7 @@ int main(int argc, char *argv[])
const auto& procPatch = procFaMeshes[proci].patch();
const auto& mp = procPatch.meshPoints();
auto& pointAddr = faPointProcAddressing[proci];
labelList& pointAddr = faPointProcAddressing[proci];
pointAddr.resize_nocopy(mp.size());
forAll(mp, i)
@ -1631,7 +1626,8 @@ int main(int argc, char *argv[])
label nPatches = 0;
forAll(completePatches, patchi)
{
const auto& patchEdgeLabels = singlePatchEdgeLabels[patchi];
const labelList& patchEdgeLabels =
singlePatchEdgeLabels[patchi];
const faPatch& fap = procMesh0.boundary()[patchi];
@ -1662,11 +1658,11 @@ int main(int argc, char *argv[])
// Serial mesh - no parallel communication
const bool oldParRun = UPstream::parRun(false);
const bool oldParRun = Pstream::parRun(false);
masterFaMesh.addFaPatches(completePatches);
UPstream::parRun(oldParRun); // Restore parallel state
Pstream::parRun(oldParRun); // Restore parallel state
// Write mesh & individual addressing

View File

@ -33,266 +33,99 @@ License
#include "decomposedBlockData.H"
#include "IFstream.H"
// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
// Trimmed-down version of lookupAndCacheProcessorsPath
// with Foam::exists() check. No caching.
// Check for two conditions:
// - file has to exist
// - if collated the entry has to exist inside the file
// Note: bypass fileOperation::filePath(IOobject&) since has problems
// with going to a different number of processors
// (in collated format). Use file-based searching instead
namespace Foam
bool Foam::checkFileExistence(const fileName& fName)
{
// Trimmed-down version of lookupAndCacheProcessorsPath
// with Foam::exists() check. No caching.
// If indeed collated format:
// Collect block-number in individual filenames
// (might differ on different processors)
static bool checkFileExistenceCollated
(
const Foam::fileOperation& handler,
const Foam::fileName& fName
)
{
// using namespace Foam;
// Check for two conditions:
// - file has to exist
// - if collated the entry has to exist inside the file
// Note: bypass fileOperation::filePath(IOobject&) since has problems
// with going to a different number of processors
// (in collated format). Use file-based searching instead
const auto& handler = Foam::fileHandler();
typedef fileOperation::procRangeType procRangeType;
fileName path, pDir, local;
procRangeType group;
label numProcs;
const label proci =
fileOperation::splitProcessorPath
(fName, path, pDir, local, group, numProcs);
bool found = false;
if (proci != -1)
{
const label handlerComm = handler.comm();
// Read all directories to see any beginning with processor
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
const label globalProci = UPstream::myProcNo(UPstream::worldComm);
const label handlerProci = UPstream::myProcNo(handlerComm);
const label nHandlerProcs = UPstream::nProcs(handlerComm);
// Determine my local block number
label myBlockNumber = -1;
{
fileOperation::procRangeType group;
label proci = fileOperation::detectProcessorPath(fName, group);
if (proci == -1 && group.empty())
{
// 'processorsXXX' format so contains all ranks
// according to worldComm
myBlockNumber = globalProci;
}
else
{
// 'processorsXXX_n-m' format so check for relative rank
myBlockNumber = handlerProci;
}
}
// Since we are streaming anyhow, could also pack as tuple:
// Tuple2<fileName, label>
// Collect file names on master of local communicator
const fileNameList fNames
const fileNameList dirEntries
(
Pstream::listGatherValues
(
fName,
handlerComm,
UPstream::msgType()
)
handler.readDir(path, fileName::Type::DIRECTORY)
);
// Collect block numbers on master of local communicator
const labelList myBlockNumbers
(
Pstream::listGatherValues
(
myBlockNumber,
handlerComm,
UPstream::msgType()
)
);
// Extract info from processorN or processorsNN
// - highest processor number
// - directory+offset containing data for proci
// Determine for all whether the filename exists in the collated file.
boolList allFound;
if (UPstream::master(handlerComm))
// label nProcs = 0;
for (const fileName& dirN : dirEntries)
{
allFound.resize(nHandlerProcs, false);
// Analyse directory name
label rNum(-1);
const label readProci =
fileOperation::detectProcessorPath(dirN, group, &rNum);
// Store nBlocks and index of file that was used for nBlocks
label nBlocks = -1;
label blockRanki = -1;
forAll(fNames, ranki)
if (proci == readProci)
{
if
(
blockRanki == -1
|| (fNames[ranki] != fNames[blockRanki])
)
// Found "processorN"
if (Foam::exists(path/dirN/local))
{
blockRanki = ranki;
IFstream is(fNames[ranki]);
nBlocks = decomposedBlockData::getNumBlocks(is);
found = true;
break;
}
}
else if (rNum != -1)
{
// "processorsNN" or "processorsNN_start-end"
if (group.empty())
{
// "processorsNN"
if (proci < rNum && Foam::exists(path/dirN/local))
{
found = true;
break;
}
}
else if (group.contains(proci))
{
// "processorsNN_start-end"
// - save the local proc offset
allFound[ranki] = (myBlockNumbers[ranki] < nBlocks);
if (Foam::exists(path/dirN/local))
{
found = true;
break;
}
}
}
}
}
// Scatter using the handler communicator
found = Pstream::listScatterValues
(
allFound,
handlerComm,
UPstream::msgType()
);
if (!found)
{
found = Foam::exists(fName);
}
return found;
}
} // End namespace
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
Foam::bitSet Foam::haveProcessorFile
(
const word& name, // eg "faces"
const fileName& instance, // eg "constant"
const fileName& local, // eg, polyMesh
const Time& runTime,
const bool verbose
)
{
const auto& handler = Foam::fileHandler();
const fileName fName
(
handler.filePath(runTime.path()/instance/local/name)
);
bool found = handler.isFile(fName);
// Assume non-collated (as fallback value).
// If everyone claims to have the file, use master to verify if
// collated is involved.
bool isCollated = false;
if (returnReduceAnd(found, UPstream::worldComm))
{
// Test for collated format.
// Note: can test only world-master. Since even host-collated will have
// same file format type for all processors
if (UPstream::master(UPstream::worldComm))
{
const bool oldParRun = UPstream::parRun(false);
if (IFstream is(fName); is.good())
{
IOobject io(name, instance, local, runTime);
io.readHeader(is);
isCollated = decomposedBlockData::isCollatedType(io);
}
UPstream::parRun(oldParRun);
}
Pstream::broadcast(isCollated, UPstream::worldComm);
}
// For collated, check that the corresponding blocks exist
if (isCollated)
{
found = checkFileExistenceCollated(handler, fName);
}
// Globally consistent information about who has the file
bitSet haveFileOnProc = bitSet::allGather(found, UPstream::worldComm);
if (verbose)
{
Info<< "Per processor availability of \""
<< name << "\" file in " << instance/local << nl
<< " " << flatOutput(haveFileOnProc) << nl << endl;
}
return haveFileOnProc;
}
Foam::boolList Foam::haveMeshFile
(
const word& name, // eg "faces"
const fileName& instance, // eg "constant"
const fileName& local, // eg, polyMesh
const Time& runTime,
const bool verbose
)
{
const auto& handler = Foam::fileHandler();
const fileName fName
(
handler.filePath(runTime.path()/instance/local/name)
);
bool found = handler.isFile(fName);
// Assume non-collated (as fallback value).
// If everyone claims to have the file, use master to verify if
// collated is involved.
bool isCollated = false;
if (returnReduceAnd(found, UPstream::worldComm))
{
// Test for collated format.
// Note: can test only world-master. Since even host-collated will have
// same file format type for all processors
if (UPstream::master(UPstream::worldComm))
{
const bool oldParRun = UPstream::parRun(false);
if (IFstream is(fName); is.good())
{
IOobject io(name, instance, local, runTime);
io.readHeader(is);
isCollated = decomposedBlockData::isCollatedType(io);
}
UPstream::parRun(oldParRun);
}
Pstream::broadcast(isCollated, UPstream::worldComm);
}
// For collated, check that the corresponding blocks exist
if (isCollated)
{
found = checkFileExistenceCollated(handler, fName);
}
// Globally consistent information about who has a mesh
boolList haveFileOnProc
(
UPstream::allGatherValues<bool>(found, UPstream::worldComm)
);
if (verbose)
{
Info<< "Per processor availability of \""
<< name << "\" file in " << instance/local << nl
<< " " << flatOutput(haveFileOnProc) << nl << endl;
}
return haveFileOnProc;
}
Foam::boolList Foam::haveMeshFile
(
@ -302,49 +135,143 @@ Foam::boolList Foam::haveMeshFile
const bool verbose
)
{
#if 0
// Simple directory scanning - too fragile
bool found = checkFileExistence(runTime.path()/meshPath/meshFile);
#else
// Trimmed-down version of lookupAndCacheProcessorsPath
// with Foam::exists() check. No caching.
// Check for two conditions:
// - file has to exist
// - if collated the entry has to exist inside the file
// Note: bypass fileOperation::filePath(IOobject&) since has problems
// with going to a different number of processors
// (in collated format). Use file-based searching instead
const auto& handler = Foam::fileHandler();
typedef fileOperation::procRangeType procRangeType;
const fileName fName
(
handler.filePath(runTime.path()/meshPath/meshFile)
);
bool found = handler.isFile(fName);
// Assume non-collated (as fallback value).
// If everyone claims to have the file, use master to verify if
// collated is involved.
bool isCollated = false;
if (returnReduceAnd(found, UPstream::worldComm))
if (returnReduceAnd(found)) // worldComm
{
// Test for collated format.
// Bit tricky: avoid having all slaves open file since this involves
// reading it on master and broadcasting it. This fails if file > 2G.
// So instead only read on master
bool isCollated = false;
// Note: can test only world-master. Since even host-collated will have
// same file format type for all processors
if (UPstream::master(UPstream::worldComm))
{
const bool oldParRun = UPstream::parRun(false);
if (IFstream is(fName); is.good())
IFstream is(fName);
if (is.good())
{
IOobject io(meshFile, meshPath, runTime);
io.readHeader(is);
isCollated = decomposedBlockData::isCollatedType(io);
}
UPstream::parRun(oldParRun);
}
Pstream::broadcast(isCollated, UPstream::worldComm);
}
Pstream::broadcast(isCollated); //UPstream::worldComm
// For collated, check that the corresponding blocks exist
if (isCollated)
{
found = checkFileExistenceCollated(handler, fName);
}
// Collect block-number in individual filenames (might differ
// on different processors)
if (isCollated)
{
const label nProcs = UPstream::nProcs(fileHandler().comm());
const label myProcNo = UPstream::myProcNo(fileHandler().comm());
// Collect file names on master of local communicator
const fileNameList fNames
(
Pstream::listGatherValues
(
fName,
fileHandler().comm(),
UPstream::msgType()
)
);
// Collect local block number
label myBlockNumber = -1;
{
procRangeType group;
label proci = fileOperation::detectProcessorPath(fName, group);
if (proci == -1 && group.empty())
{
// 'processorsXXX' format so contains all ranks
// according to worldComm
myBlockNumber = UPstream::myProcNo(UPstream::worldComm);
}
else
{
// 'processorsXXX_n-m' format so check for the
// relative rank
myBlockNumber = myProcNo;
}
}
const labelList myBlockNumbers
(
Pstream::listGatherValues
(
myBlockNumber,
fileHandler().comm(),
UPstream::msgType()
)
);
// Determine for all whether the filename exists in the collated
// file.
boolList allFound(nProcs, false);
if (UPstream::master(fileHandler().comm()))
{
// Store nBlocks and index of file that was used for nBlocks
label nBlocks = -1;
label blockRanki = -1;
forAll(fNames, ranki)
{
if
(
blockRanki == -1
|| (fNames[ranki] != fNames[blockRanki])
)
{
blockRanki = ranki;
IFstream is(fNames[ranki]);
nBlocks = decomposedBlockData::getNumBlocks(is);
}
allFound[ranki] = (myBlockNumbers[ranki] < nBlocks);
}
}
found = Pstream::listScatterValues
(
allFound,
fileHandler().comm(),
UPstream::msgType()
);
}
}
#endif
// Globally consistent information about who has a mesh
boolList haveFileOnProc

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2012 OpenFOAM Foundation
Copyright (C) 2022-2025 OpenCFD Ltd.
Copyright (C) 2022-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -48,32 +48,15 @@ namespace Foam
// Forward Declarations
class faMesh;
//- Check for availability of specified file (on all ranks)
bitSet haveProcessorFile
(
const word& name, // eg, "faces"
const fileName& instance, // eg, "constant"
const fileName& local, // eg, "polyMesh"
const Time& runTime,
const bool verbose = true
);
//- Check for availability of given file
bool checkFileExistence(const fileName& fName);
//- Check for availability of specified file (on all ranks)
boolList haveMeshFile
(
const word& name, // eg, "faces"
const fileName& instance, // eg, "constant"
const fileName& local, // eg, "polyMesh"
const Time& runTime,
const bool verbose = true
);
//- Check for availability of specified mesh file
//- Check for availability of specified mesh file (default: "faces")
boolList haveMeshFile
(
const Time& runTime,
const fileName& meshPath,
const word& meshFile,
const word& meshFile = "faces",
const bool verbose = true
);
@ -92,7 +75,7 @@ void masterMeshInstance
fileName& pointsInstance
);
} // End namespace Foam
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

View File

@ -21,59 +21,44 @@ Requires
// Initially all possible objects that are available at the final time
List<wordHashSet> availableRegionObjectNames(meshes.size());
List<HashTable<wordHashSet>> availableFaRegionObjectNames(meshes.size());
List<wordHashSet> availableFaRegionObjectNames(meshes.size());
forAll(meshes, regioni)
{
const auto& mesh = meshes[regioni];
IOobjectList objects;
HashTable<IOobjectList> faObjects;
IOobjectList faObjects;
if (doConvertFields && !timeDirs.empty())
{
// const word& checkTimeDir = timeDirs.back().name();
// List of volume mesh objects for this instance
objects = IOobjectList(mesh, timeDirs.back().name());
// List of area mesh objects (assuming single region)
faObjects = IOobjectList
(
mesh.time(),
timeDirs.back().name(),
faMesh::dbDir(mesh, word::null),
IOobjectOption::NO_REGISTER
);
if (fieldSelector)
{
objects.filterObjects(fieldSelector);
faObjects.filterObjects(fieldSelector);
}
// Remove "*_0" restart fields
objects.prune_0();
faObjects.prune_0();
if (!doPointValues)
{
// Prune point fields if disabled
objects.filterClasses(Foam::fieldTypes::is_point, true);
}
// The finite-area regions and their fields for this volume region
// and instance
faObjects.reserve(areaRegionNames.size());
for (const word& areaName : areaRegionNames)
{
IOobjectList objs
(
faMesh::Registry(mesh),
timeDirs.back().name(),
polyMesh::regionName(areaName),
IOobjectOption::NO_REGISTER
);
if (fieldSelector)
{
objs.filterObjects(fieldSelector);
}
// Remove "*_0" restart fields
objs.prune_0();
faObjects(areaName) = std::move(objs);
}
}
// Volume fields
@ -93,29 +78,20 @@ forAll(meshes, regioni)
}
// Area fields
for (const word& areaName : areaRegionNames)
if (!faObjects.empty())
{
wordList objectNames;
wordList objectNames(faObjects.sortedNames());
if (const auto iter = faObjects.cfind(areaName); iter.good())
{
objectNames = iter.val().sortedNames();
// Check availability for all times...
checkData
(
faMesh::Registry(mesh),
timeDirs,
objectNames,
polyMesh::regionName(areaName)
);
}
availableFaRegionObjectNames[regioni].emplace_set
// Check availability for all times... (assuming single region)
checkData
(
areaName,
std::move(objectNames)
mesh.time(),
timeDirs,
objectNames,
faMesh::dbDir(mesh, word::null)
);
availableFaRegionObjectNames[regioni] = objectNames;
}
}

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2021-2025 OpenCFD Ltd.
Copyright (C) 2021-2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
@ -15,14 +15,12 @@ Description
\*---------------------------------------------------------------------------*/
// Cases and meshes per volume region
PtrList<ensightCase> ensightCases(regionNames.size());
PtrList<ensightMesh> ensightMeshes(regionNames.size());
// Per volume region can have multiple finite-area regions
List<PtrDynList<faMesh>> meshesFa(regionNames.size());
List<PtrDynList<ensightCase>> ensightCasesFa(regionNames.size());
List<PtrDynList<ensightFaMesh>> ensightMeshesFa(regionNames.size());
PtrList<faMesh> meshesFa(regionNames.size());
PtrList<ensightCase> ensightCasesFa(regionNames.size());
PtrList<ensightFaMesh> ensightMeshesFa(regionNames.size());
{
forAll(regionNames, regioni)
@ -47,7 +45,6 @@ List<PtrDynList<ensightFaMesh>> ensightMeshesFa(regionNames.size());
}
}
// Volume mesh
ensightMeshes.set
(
regioni,
@ -62,59 +59,32 @@ List<PtrDynList<ensightFaMesh>> ensightMeshesFa(regionNames.size());
new ensightCase(ensCasePath, ensCaseName, caseOpts)
);
if (!doFiniteArea)
if (doFiniteArea)
{
continue;
}
// Note: not every volume region is guaranteed to have
// any or all area region(s)
meshesFa[regioni].reserve_exact(areaRegionNames.size());
for (const word& areaName : areaRegionNames)
{
auto faMeshPtr = faMesh::TryNew(areaName, mesh);
autoPtr<faMesh> faMeshPtr(faMesh::TryNew(mesh));
if (faMeshPtr)
{
meshesFa[regioni].push_back(std::move(faMeshPtr));
}
}
// Setup the ensight components
const label nAreas = meshesFa[regioni].size();
ensightCasesFa[regioni].reserve_exact(nAreas);
ensightMeshesFa[regioni].reserve_exact(nAreas);
for (const faMesh& areaMesh : meshesFa[regioni])
{
// Use regionName() - automatically filters for defaultRegion
const word& areaName = areaMesh.regionName();
if (areaName.empty())
{
// single-region
ensightCasesFa[regioni].emplace_back
ensightCasesFa.set
(
ensCasePath/"finite-area",
"finite-area",
caseOpts
regioni,
new ensightCase
(
ensCasePath/"finite-area",
"finite-area",
caseOpts
)
);
}
else
{
// multi-region
ensightCasesFa[regioni].emplace_back
(
ensCasePath/"finite-area"/areaName,
areaName,
caseOpts
);
}
auto& ensFaMesh = ensightMeshesFa[regioni].emplace_back(areaMesh);
ensFaMesh.verbose(optVerbose);
meshesFa.set(regioni, std::move(faMeshPtr));
ensightMeshesFa.set
(
regioni,
new ensightFaMesh(meshesFa[regioni])
);
ensightMeshesFa[regioni].verbose(optVerbose);
}
}
}
}

View File

@ -124,7 +124,6 @@ Usage
#include "OFstream.H"
#include "Pstream.H"
#include "HashOps.H"
#include "PtrDynList.H"
#include "regionProperties.H"
#include "fvc.H"
@ -173,7 +172,6 @@ int main(int argc, char *argv[])
argList::addVerboseOption();
#include "addAllRegionOptions.H"
#include "addAllFaRegionOptions.H"
argList::addBoolOption
(
@ -441,14 +439,6 @@ int main(int argc, char *argv[])
// Handle -allRegions, -regions, -region
#include "getAllRegionOptions.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
if (!doFiniteArea)
{
areaRegionNames.clear(); // For consistency
}
// ------------------------------------------------------------------------
// Directory management
@ -528,32 +518,29 @@ int main(int argc, char *argv[])
polyMesh::readUpdateState meshState = mesh.readUpdate();
const bool moving = (meshState != polyMesh::UNCHANGED);
// Ensight (fvMesh)
// Ensight
auto& ensCase = ensightCases[regioni];
auto& ensMesh = ensightMeshes[regioni];
// Finite-area (can be missing)
auto* ensFaCasePtr = ensightCasesFa.get(regioni);
auto* ensFaMeshPtr = ensightMeshesFa.get(regioni);
ensCase.setTime(timeDirs[timei], timeIndex);
// Finite-area (optional)
// Accounting exists for each volume region but may be empty
auto& ensFaCases = ensightCasesFa[regioni];
auto& ensFaMeshes = ensightMeshesFa[regioni];
for (auto& ensFaCase : ensFaCases)
if (ensFaCasePtr)
{
ensFaCase.setTime(timeDirs[timei], timeIndex);
ensFaCasePtr->setTime(timeDirs[timei], timeIndex);
}
// Movement
if (moving)
{
ensMesh.expire();
ensMesh.correct();
for (auto& ensFaMesh : ensFaMeshes)
if (ensFaMeshPtr)
{
ensFaMesh.expire();
ensFaMesh.correct();
ensFaMeshPtr->expire();
ensFaMeshPtr->correct();
}
}
@ -568,19 +555,15 @@ int main(int argc, char *argv[])
}
// finite-area
forAll(ensFaMeshes, areai)
if (ensFaCasePtr && ensFaMeshPtr)
{
const auto& ensFaCase = ensFaCases[areai];
const auto& ensFaMesh = ensFaMeshes[areai];
autoPtr<ensightGeoFile> os =
ensFaCase.newGeometry(hasMovingMesh);
ensFaCasePtr->newGeometry(hasMovingMesh);
ensFaMesh.write(os.ref());
ensFaMeshPtr->write(os.ref());
}
}
// Objects at this time
IOobjectList objects(mesh, runTime.timeName());
@ -592,37 +575,23 @@ int main(int argc, char *argv[])
// Volume, internal, point fields
#include "convertVolumeFields.H"
// finite-area
forAll(ensFaMeshes, areai)
// The finite-area objects at this time
IOobjectList faObjects;
if (ensFaMeshPtr)
{
auto* ensFaCasePtr = ensFaCases.get(areai);
auto* ensFaMeshPtr = ensFaMeshes.get(areai);
faObjects =
IOobjectList(ensFaMeshPtr->mesh(), runTime.timeName());
// The finite-area region objects at this time
IOobjectList faObjects;
if (ensFaMeshPtr)
{
const word& areaName = ensFaMeshPtr->mesh().name();
faObjects = IOobjectList
(
faMesh::Registry(mesh),
runTime.timeName(),
polyMesh::regionName(areaName),
IOobjectOption::NO_REGISTER
);
faObjects.filterObjects
(
availableFaRegionObjectNames[regioni](areaName)
);
}
// The finiteArea fields
#include "convertAreaFields.H"
faObjects.filterObjects
(
availableFaRegionObjectNames[regioni]
);
}
// The finiteArea fields
#include "convertAreaFields.H"
// Lagrangian fields
#include "convertLagrangian.H"
}
@ -635,13 +604,14 @@ int main(int argc, char *argv[])
// Write cases
forAll(ensightCases, regioni)
{
// finite-volume
ensightCases[regioni].write();
}
// Finite-area (if any)
for (const auto& ensFaCase : ensightCasesFa[regioni])
forAll(ensightCasesFa, regioni)
{
if (ensightCasesFa.set(regioni))
{
ensFaCase.write();
ensightCasesFa[regioni].write();
}
}

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2018-2025 OpenCFD Ltd.
Copyright (C) 2018-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -26,6 +26,7 @@ License
\*---------------------------------------------------------------------------*/
#include "readFields.H"
#include "volFields.H"
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
@ -41,30 +42,26 @@ Foam::label Foam::checkData
wordHashSet goodFields;
IOobject io
(
"any-name", // placeholder
"constant", // placeholder
local,
obr,
IOobjectOption::NO_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::NO_REGISTER
);
for (const word& fieldName : objectNames)
{
// In case prune_0() not previously used...
if (fieldName.ends_with("_0")) continue;
// // If prune_0() not previously used...
// if (objectNames.ends_with("_0")) continue;
bool good = false;
for (const instant& inst : timeDirs)
{
io.resetHeader(fieldName);
io.instance() = inst.name();
good = io.typeHeaderOk<regIOobject>(false, false);
good =
IOobject
(
fieldName,
inst.name(),
local,
obr,
IOobject::NO_READ,
IOobject::NO_WRITE,
IOobject::NO_REGISTER
).typeHeaderOk<volScalarField>(false, false);
if (!good)
{

View File

@ -22,66 +22,24 @@ Description
//
// No subsetting!
if (doFiniteArea && !areaRegionNames.empty())
if (doFiniteArea)
{
using reportFields = foamToVtkReportFields;
// The (region) polyMesh being used. No subsetting possible
const auto& basePolyMesh = meshProxy.baseMesh();
autoPtr<faMesh> faMeshPtr;
for (const word& areaName : areaRegionNames)
const label nAreaFields = faObjects.count(Foam::fieldTypes::is_area);
if (nAreaFields || withMeshIds)
{
const bool isDefaultRegion(polyMesh::regionName(areaName).empty());
faMeshPtr = faMesh::TryNew(meshProxy.baseMesh());
}
// CAUTION
// If we want to have constant access to the HashTable:
//
// (SEGFAULT)
// const auto& faObjs = faObjects.lookup(areaName, IOobjectList());
//
// Use an empty fallback to avoid binding to a temporary:
//
// const IOobjectList emptyObjectList;
// const auto& faObjs = faObjects.lookup(areaName, emptyObjectList);
if (faMeshPtr && (nAreaFields || withMeshIds))
{
const faMesh& areaMesh = faMeshPtr();
// Since we do not need the area fields afterwards,
// just move them out from the HashTable
IOobjectList faObjs;
if (auto iter = faObjects.find(areaName); iter.good())
{
faObjs = std::move(iter.val());
}
const label nAreaFields = faObjs.count(Foam::fieldTypes::is_area);
autoPtr<faMesh> faMeshPtr;
if (nAreaFields || withMeshIds)
{
faMeshPtr = faMesh::TryNew(areaName, basePolyMesh);
}
if (!faMeshPtr)
{
if (!isDefaultRegion)
{
// Report any area region specified but missing
// - silently ignore region0
Info<< "No area-mesh [" << polyMesh::regionName(areaName)
<< "] on volume-region ["
<< basePolyMesh.regionName() << "]" << endl;
}
continue;
}
const auto& areaMesh = faMeshPtr();
Info<< "Using area-mesh [" << polyMesh::regionName(areaName)
<< "] on volume-region ["
<< basePolyMesh.regionName() << "]" << endl;
reportFields::area(Info, faObjs);
reportFields::area(Info, faObjects);
const auto& pp = faMeshPtr->patch();
@ -91,12 +49,7 @@ if (doFiniteArea && !areaRegionNames.empty())
writeOpts,
(
outputDir/regionDir/"finite-area"
/ (
isDefaultRegion
? fileName("finiteArea")
: fileName(areaName/areaName)
)
+ timeDesc
/ "finiteArea" + timeDesc
),
UPstream::parRun()
);
@ -143,7 +96,7 @@ if (doFiniteArea && !areaRegionNames.empty())
(
writer,
areaMesh,
faObjs,
faObjects,
true // syncPar
);

View File

@ -135,13 +135,13 @@ Note
#include "emptyPolyPatch.H"
#include "volPointInterpolation.H"
#include "faceZoneMesh.H"
#include "faMesh.H"
#include "areaFields.H"
#include "fvMeshSubsetProxy.H"
#include "faceSet.H"
#include "pointSet.H"
#include "HashOps.H"
#include "regionProperties.H"
#include "stringListOps.H" // For stringListOps::findMatching()
#include "Cloud.H"
#include "readFields.H"
@ -434,7 +434,6 @@ int main(int argc, char *argv[])
argList::addOptionCompat("one-boundary", {"allPatches", 1806});
#include "addAllRegionOptions.H"
#include "addAllFaRegionOptions.H"
argList::addOption
(
@ -631,16 +630,6 @@ int main(int argc, char *argv[])
// Handle -allRegions, -regions, -region
#include "getAllRegionOptions.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
if (!doFiniteArea)
{
areaRegionNames.clear(); // For consistency
}
// ------------------------------------------------------------------------
// Names for sets and zones
word cellSelectionName;
word faceSetName;
@ -788,11 +777,8 @@ int main(int argc, char *argv[])
}
}
// fvMesh fields
IOobjectList objects;
// faMesh fields (multiple finite-area regions per volume region)
HashTable<IOobjectList> faObjects;
IOobjectList faObjects;
if (doConvertFields)
{
@ -800,51 +786,33 @@ int main(int argc, char *argv[])
objects =
IOobjectList(meshProxy.baseMesh(), runTime.timeName());
// List of area mesh objects (assuming single region)
faObjects =
IOobjectList
(
runTime,
runTime.timeName(),
faMesh::dbDir(meshProxy.baseMesh(), word::null),
IOobjectOption::NO_REGISTER
);
if (fieldSelector)
{
objects.filterObjects(fieldSelector);
faObjects.filterObjects(fieldSelector);
}
// Remove "*_0" restart fields
objects.prune_0();
faObjects.prune_0();
if (!doPointValues)
{
// Prune point fields if disabled
objects.filterClasses(Foam::fieldTypes::is_point, true);
}
// Lists of finite-area fields
faObjects.reserve(areaRegionNames.size());
for (const word& areaName : areaRegionNames)
{
// The finite-area objects for given area region.
// Add as method to faMeshesRegistry?
IOobjectList objs
(
faMesh::Registry(meshProxy.baseMesh()),
runTime.timeName(),
polyMesh::regionName(areaName),
IOobjectOption::NO_REGISTER
);
if (fieldSelector)
{
objs.filterObjects(fieldSelector);
}
// Remove "*_0" restart fields
objs.prune_0();
if (!objs.empty())
{
faObjects.emplace_set(areaName, std::move(objs));
}
}
}
if (processorFieldsOnly)
{
// Processor-patches only and continue

View File

@ -1,10 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/finiteArea/lnInclude
-I$(LIB_SRC)/finiteArea/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
EXE_LIBS = \
-lmeshTools \
-lfiniteVolume \
-lfiniteArea \
-lmeshTools \
-lgenericPatchFields

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2022-2025 OpenCFD Ltd.
Copyright (C) 2022-2024 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -47,7 +47,6 @@ Description
#include "areaFields.H"
#include "coupledFvPatch.H"
#include "coupledFaPatch.H"
#include "regionProperties.H"
using namespace Foam;
@ -87,18 +86,17 @@ bool consumeUnusedType(const fieldDescription& fieldDesc, Istream& is)
//? typedef GeometricField<Type, faePatchField, areaMesh> fieldType3;
//? typedef GeometricField<Type, fvsPatchField, volMesh> fieldType4;
const bool isExpectedType
if
(
fieldDesc.type() == fieldType1::typeName
|| fieldDesc.type() == fieldType2::typeName
);
if (isExpectedType)
)
{
(void) pTraits<Type>(is);
return true;
}
return isExpectedType;
return false;
}
@ -124,18 +122,13 @@ bool setCellFieldType
(
const fieldDescription& fieldDesc,
const fvMesh& mesh,
const labelUList& selectedCells,
const labelList& selectedCells,
Istream& is
)
{
typedef GeometricField<Type, fvPatchField, volMesh> fieldType;
const bool isExpectedType
(
fieldDesc.type() == fieldType::typeName
);
if (!isExpectedType)
if (fieldDesc.type() != fieldType::typeName)
{
return false;
}
@ -150,9 +143,7 @@ bool setCellFieldType
fieldDesc.name(),
mesh.thisDb().time().timeName(),
mesh.thisDb(),
IOobjectOption::MUST_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::NO_REGISTER
IOobject::MUST_READ
);
bool found = fieldHeader.typeHeaderOk<fieldType>(true);
@ -160,8 +151,13 @@ bool setCellFieldType
if (!found)
{
// Fallback to "constant" directory
fieldHeader.resetHeader();
fieldHeader.instance() = mesh.thisDb().time().constant();
fieldHeader = IOobject
(
fieldDesc.name(),
mesh.thisDb().time().constant(),
mesh.thisDb(),
IOobject::MUST_READ
);
found = fieldHeader.typeHeaderOk<fieldType>(true);
}
@ -175,7 +171,7 @@ bool setCellFieldType
fieldType field(fieldHeader, mesh, false);
if (isNull(selectedCells) || (selectedCells.size() == field.size()))
if (isNull(selectedCells) || selectedCells.size() == field.size())
{
field.primitiveFieldRef() = fieldValue;
}
@ -208,7 +204,7 @@ bool setCellFieldType
<< "Field " << fieldDesc.name() << " not found" << endl;
}
return isExpectedType;
return true;
}
@ -220,18 +216,13 @@ bool setAreaFieldType
(
const fieldDescription& fieldDesc,
const faMesh& mesh,
const labelUList& selectedFaces,
const labelList& selectedFaces,
Istream& is
)
{
typedef GeometricField<Type, faPatchField, areaMesh> fieldType;
const bool isExpectedType
(
fieldDesc.type() == fieldType::typeName
);
if (!isExpectedType)
if (fieldDesc.type() != fieldType::typeName)
{
return false;
}
@ -246,9 +237,7 @@ bool setAreaFieldType
fieldDesc.name(),
mesh.thisDb().time().timeName(),
mesh.thisDb(),
IOobjectOption::MUST_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::NO_REGISTER
IOobject::MUST_READ
);
bool found = fieldHeader.typeHeaderOk<fieldType>(true);
@ -256,8 +245,13 @@ bool setAreaFieldType
if (!found)
{
// Fallback to "constant" directory
fieldHeader.resetHeader();
fieldHeader.instance() = mesh.thisDb().time().constant(),
fieldHeader = IOobject
(
fieldDesc.name(),
mesh.thisDb().time().constant(),
mesh.thisDb(),
IOobject::MUST_READ
);
found = fieldHeader.typeHeaderOk<fieldType>(true);
}
@ -298,7 +292,7 @@ bool setAreaFieldType
<< "Field " << fieldDesc.name() << " not found" << endl;
}
return isExpectedType;
return true;
}
@ -310,18 +304,13 @@ bool setFaceFieldType
(
const fieldDescription& fieldDesc,
const fvMesh& mesh,
const labelUList& selectedFaces,
const labelList& selectedFaces,
Istream& is
)
{
typedef GeometricField<Type, fvPatchField, volMesh> fieldType;
const bool isExpectedType
(
fieldDesc.type() == fieldType::typeName
);
if (!isExpectedType)
if (fieldDesc.type() != fieldType::typeName)
{
return false;
}
@ -336,9 +325,7 @@ bool setFaceFieldType
fieldDesc.name(),
mesh.thisDb().time().timeName(),
mesh.thisDb(),
IOobjectOption::MUST_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::NO_REGISTER
IOobject::MUST_READ
);
bool found = fieldHeader.typeHeaderOk<fieldType>(true);
@ -346,8 +333,13 @@ bool setFaceFieldType
if (!found)
{
// Fallback to "constant" directory
fieldHeader.resetHeader();
fieldHeader.instance() = mesh.thisDb().time().constant();
fieldHeader = IOobject
(
fieldDesc.name(),
mesh.thisDb().time().constant(),
mesh.thisDb(),
IOobject::MUST_READ
);
found = fieldHeader.typeHeaderOk<fieldType>(true);
}
@ -363,14 +355,15 @@ bool setFaceFieldType
// Create flat list of selected faces and their value.
Field<Type> allBoundaryValues(mesh.nBoundaryFaces());
for (const auto& pfld : field.boundaryField())
forAll(field.boundaryField(), patchi)
{
SubField<Type>
(
allBoundaryValues,
pfld.size(),
pfld.patch().offset()
) = pfld;
field.boundaryField()[patchi].size(),
field.boundaryField()[patchi].patch().start()
- mesh.nInternalFaces()
) = field.boundaryField()[patchi];
}
// Override
@ -431,7 +424,8 @@ bool setFaceFieldType
(
allBoundaryValues,
fieldBf[patchi].size(),
fieldBf[patchi].patch().offset()
fieldBf[patchi].patch().start()
- mesh.nInternalFaces()
);
}
}
@ -451,7 +445,7 @@ bool setFaceFieldType
<< "Field " << fieldDesc.name() << " not found" << endl;
}
return isExpectedType;
return true;
}
@ -466,7 +460,7 @@ struct setCellField
(
const fieldDescription& fieldDesc,
const fvMesh& m,
const labelUList& selectedCells,
const labelList& selectedCells,
Istream& is
)
{
@ -483,11 +477,11 @@ struct setCellField
class iNew
{
const fvMesh& mesh_;
const labelUList& selected_;
const labelList& selected_;
public:
iNew(const fvMesh& mesh, const labelUList& selectedCells) noexcept
iNew(const fvMesh& mesh, const labelList& selectedCells)
:
mesh_(mesh),
selected_(selectedCells)
@ -534,7 +528,7 @@ struct setFaceField
(
const fieldDescription& fieldDesc,
const fvMesh& m,
const labelUList& selectedFaces,
const labelList& selectedFaces,
Istream& is
)
{
@ -551,11 +545,11 @@ struct setFaceField
class iNew
{
const fvMesh& mesh_;
const labelUList& selected_;
const labelList& selected_;
public:
iNew(const fvMesh& mesh, const labelUList& selectedFaces) noexcept
iNew(const fvMesh& mesh, const labelList& selectedFaces)
:
mesh_(mesh),
selected_(selectedFaces)
@ -602,7 +596,7 @@ struct setAreaField
(
const fieldDescription& fieldDesc,
const faMesh& m,
const labelUList& selectedFaces,
const labelList& selectedFaces,
Istream& is
)
{
@ -619,11 +613,11 @@ struct setAreaField
class iNew
{
const faMesh& mesh_;
const labelUList& selected_;
const labelList& selected_;
public:
iNew(const faMesh& mesh, const labelUList& selectedFaces) noexcept
iNew(const faMesh& mesh, const labelList& selectedFaces)
:
mesh_(mesh),
selected_(selectedFaces)
@ -681,7 +675,6 @@ int main(int argc, char *argv[])
);
#include "addRegionOption.H"
#include "addAllFaRegionOptions.H"
// -------------------------
@ -693,59 +686,15 @@ int main(int argc, char *argv[])
#include "createNamedMesh.H"
// Handle -all-area-regions, -area-regions, -area-region
#include "getAllFaRegionOptions.H"
autoPtr<faMesh> faMeshPtr;
if (args.found("no-finite-area"))
if (!args.found("no-finite-area"))
{
areaRegionNames.clear(); // For consistency
faMeshPtr = faMesh::TryNew(mesh);
}
// ------------------------------------------------------------------------
PtrList<faMesh> faMeshes;
// Setup all area meshes on this region
if (!areaRegionNames.empty()) // ie, !args.found("no-finite-area")
if (faMeshPtr)
{
faMeshes.resize(areaRegionNames.size());
label nGoodRegions(0);
for (const word& areaName : areaRegionNames)
{
autoPtr<faMesh> faMeshPtr = faMesh::TryNew(areaName, mesh);
if (faMeshPtr)
{
faMeshes.set(nGoodRegions++, std::move(faMeshPtr));
}
}
faMeshes.resize(nGoodRegions);
}
if (faMeshes.size() == 1)
{
Info<< "Using finite-area mesh";
if
(
const word& name = polyMesh::regionName(faMeshes[0].name());
!name.empty()
)
{
Info<< " [" << name << "]";
}
Info<< nl;
}
else if (faMeshes.size() > 1)
{
Info<< "Detected finite-area meshes:";
for (const faMesh& areaMesh : faMeshes)
{
Info<< " [" << areaMesh.name() << "]";
}
Info<< nl;
Info<< "Detected finite-area mesh" << nl;
}
const word dictName("setFieldsDict");
@ -763,45 +712,37 @@ int main(int argc, char *argv[])
// Default field values
if
(
const auto* eptr
= setFieldsDict.findEntry("defaultFieldValues", keyType::LITERAL)
)
{
ITstream& is = eptr->stream();
const entry* eptr =
setFieldsDict.findEntry("defaultFieldValues", keyType::LITERAL);
Info<< "Setting volume field default values" << endl;
PtrList<setCellField> defaultFieldValues
(
is,
setCellField::iNew(mesh, labelList::null())
);
for (const faMesh& areaMesh : faMeshes)
if (eptr)
{
Info<< "Setting area field default values";
ITstream& is = eptr->stream();
if
(
const word& name = polyMesh::regionName(areaMesh.name());
!name.empty()
)
{
Info<< " [" << name << "]";
}
Info<< endl;
Info<< "Setting volume field default values" << endl;
is.rewind();
PtrList<setAreaField> defaultAreaFieldValues
PtrList<setCellField> defaultFieldValues
(
is,
setAreaField::iNew(areaMesh, labelList::null())
setCellField::iNew(mesh, labelList::null())
);
if (faMeshPtr)
{
const faMesh& areaMesh = faMeshPtr();
is.rewind();
Info<< "Setting area field default values" << endl;
PtrList<setAreaField> defaultFieldValues
(
is,
setAreaField::iNew(areaMesh, labelList::null())
);
}
Info<< endl;
}
Info<< endl;
}
@ -827,15 +768,11 @@ int main(int argc, char *argv[])
labelList selectedCells(subset.sortedToc());
{
FixedList<label, 2> stats;
stats[0] = selectedCells.size();
stats[1] = mesh.nCells();
reduce(stats, sumOp<label>());
Info<< " Selected " << stats[0] << '/' << stats[1]
<< " cells" << nl;
}
Info<< " Selected "
<< returnReduce(selectedCells.size(), sumOp<label>())
<< '/'
<< returnReduce(mesh.nCells(), sumOp<label>())
<< " cells" << nl;
ITstream& is = region.dict().lookup("fieldValues");
@ -869,8 +806,10 @@ int main(int argc, char *argv[])
setFaceField::iNew(mesh, selectedFaces)
);
for (const faMesh& areaMesh : faMeshes)
if (faMeshPtr)
{
const faMesh& areaMesh = faMeshPtr();
const labelUList& faceLabels = areaMesh.faceLabels();
// Transcribe from mesh faces to finite-area addressing
@ -889,15 +828,11 @@ int main(int argc, char *argv[])
}
areaFaces.resize(nUsed);
{
FixedList<label, 2> stats;
stats[0] = areaFaces.size();
stats[1] = faceLabels.size();
reduce(stats, sumOp<label>());
Info<< " Selected " << stats[0] << '/' << stats[1]
<< " area faces for " << areaMesh.name() << endl;
}
Info<< " Selected "
<< returnReduce(areaFaces.size(), sumOp<label>())
<< '/'
<< returnReduce(faceLabels.size(), sumOp<label>())
<< " area faces" << nl;
is.rewind();

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2016-2020 OpenCFD Ltd.
Copyright (C) 2016-2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -189,6 +189,23 @@ public:
template<class TrackingData>
inline bool equal(const deltaData&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const deltaData& f0,
const label i1,
const deltaData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators
@ -296,6 +313,17 @@ public:
template<>
struct is_contiguous<LESModels::smoothDelta::deltaData> : std::true_type {};
//- Interpolation - used in e.g. cyclicAMI interpolation
inline LESModels::smoothDelta::deltaData lerp
(
const LESModels::smoothDelta::deltaData& a,
const LESModels::smoothDelta::deltaData& b,
const scalar t
)
{
return LESModels::smoothDelta::deltaData(lerp(a.delta(), b.delta(), t));
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2015 OpenFOAM Foundation
Copyright (C) 2019 OpenCFD Ltd.
Copyright (C) 2019,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -191,6 +191,40 @@ inline bool Foam::LESModels::smoothDelta::deltaData::equal
}
template<class TrackingData>
inline bool Foam::LESModels::smoothDelta::deltaData::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const deltaData& f0,
const label i1,
const deltaData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td) && f1.valid(td))
{
const deltaData w2(lerp(f0.delta(), f1.delta(), weight));
return update(w2, 1.0, tol, td);
}
else if (f0.valid(td))
{
return update(f0, 1.0, tol, td);
}
else if (f1.valid(td))
{
return update(f1, 1.0, tol, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::LESModels::smoothDelta::deltaData::operator==

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019-2020 OpenCFD Ltd.
Copyright (C) 2019-2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -242,6 +242,23 @@ public:
template<class TrackingData>
inline bool equal(const directionInfo&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const directionInfo& f0,
const label i1,
const directionInfo& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
Copyright (C) 2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -287,6 +287,35 @@ inline bool Foam::directionInfo::equal
}
template<class TrackingData>
inline bool Foam::directionInfo::interpolate
(
const polyMesh& mesh,
const point& pt,
const label i0,
const directionInfo& f0,
const label i1,
const directionInfo& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
return updateFace(mesh, -1, f0, tol, td);
}
if (f1.valid(td))
{
return updateFace(mesh, -1, f1, tol, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::directionInfo::operator==

View File

@ -186,6 +186,23 @@ public:
TrackingData& td
);
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const wallNormalInfo& f0,
const label i1,
const wallNormalInfo& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
//- Test for equality, with TrackingData
template<class TrackingData>
inline bool equal(const wallNormalInfo&, TrackingData& td) const;

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
Copyright (C) 2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -191,6 +191,35 @@ inline bool Foam::wallNormalInfo::equal
}
template<class TrackingData>
inline bool Foam::wallNormalInfo::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const wallNormalInfo& f0,
const label i1,
const wallNormalInfo& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
return update(f0, td);
}
if (f1.valid(td))
{
return update(f1, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::wallNormalInfo::operator==

View File

@ -197,6 +197,23 @@ public:
template<class TrackingData>
inline bool equal(const refinementData&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const refinementData& f0,
const label i1,
const refinementData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
Copyright (C) 2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -248,6 +248,35 @@ inline bool Foam::refinementData::equal
}
template<class TrackingData>
inline bool Foam::refinementData::interpolate
(
const polyMesh& mesh,
const point& pt,
const label i0,
const refinementData& f0,
const label i1,
const refinementData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
return updateFace(mesh, -1, f0, tol, td);
}
if (f1.valid(td))
{
return updateFace(mesh, -1, f1, tol, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::refinementData::operator==

View File

@ -229,6 +229,23 @@ public:
inline bool equal(const refinementDistanceData&, TrackingData&)
const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const refinementDistanceData& f0,
const label i1,
const refinementDistanceData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019-2020 OpenCFD Ltd.
Copyright (C) 2019-2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -278,6 +278,35 @@ inline bool Foam::refinementDistanceData::equal
}
template<class TrackingData>
inline bool Foam::refinementDistanceData::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const refinementDistanceData& f0,
const label i1,
const refinementDistanceData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
return update(pt, f0, tol, td);
}
if (f1.valid(td))
{
return update(pt, f1, tol, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::refinementDistanceData::operator==

View File

@ -9,7 +9,7 @@ faceSetOption/faceSetOption.C
derivedSources=sources/derived
$(derivedSources)/externalHeatFluxSource/externalHeatFluxSource.C
$(derivedSources)/jouleHeatingSource/jouleHeatingSource.cxx
$(derivedSources)/jouleHeatingSource/jouleHeatingSource.C
$(derivedSources)/contactHeatFluxSource/contactHeatFluxSource.C
$(derivedSources)/externalFileSource/externalFileSource.C

View File

@ -81,33 +81,6 @@ Foam::fa::option::option
active_(dict.getOrDefault("active", true)),
log(true)
{
// Suffix hint for variable names
if
(
coeffs_.readIfPresent("suffixing", suffixHint_)
|| dict.readIfPresent("suffixing", suffixHint_)
)
{
Switch sw = Switch::find(suffixHint_);
if (sw.good())
{
if (!sw) // No suffix
{
suffixHint_.clear();
}
}
else if (suffixHint_ == "default")
{
sw = true;
}
if (sw) // Default suffix
{
suffixHint_ = '_' + regionName_;
}
}
if (dict.readIfPresent("area", areaName_))
{
if (!sameRegionNames(areaName_, defaultAreaName))
@ -122,14 +95,9 @@ Foam::fa::option::option
}
}
Log << incrIndent << indent << "Source: " << name_;
if (!polyMesh::regionName(areaName_).empty())
{
Log<< " [" << areaName_ << ']';
}
Info<< decrIndent << endl;
Log << incrIndent << indent << "Source: " << name_
<< " [" << polyMesh::regionName(areaName_) << ']' << endl
<< decrIndent;
}
@ -149,20 +117,13 @@ Foam::autoPtr<Foam::fa::option> Foam::fa::option::New
coeffs.readIfPresent("area", areaName);
Info<< indent
<< "Selecting finite-area option, type " << modelType;
<< "Selecting finite-area option, type " << modelType
<< " [" << polyMesh::regionName(areaName) << ']';
if (sameRegionNames(areaName, defaultAreaName))
if (!sameRegionNames(areaName, defaultAreaName))
{
if (!polyMesh::regionName(areaName).empty())
{
Info<< " [" << areaName << ']';
}
Info<< " != " << defaultAreaName << nl;
}
else
{
Info<< " [" << areaName << "] != [" << defaultAreaName << ']';
}
Info<< endl;
mesh.time().libs().open

View File

@ -65,13 +65,10 @@ Usage
--> the selected faOption settings | dictionary | no | -
active | Flag to (de)activate faOption | bool | no | true
log | Flag to log faOption-related info | bool | no | true
suffixing | Suffix hint for option variables | word | no | -
\endtable
Note
- Source/sink options are to be added to the right-hand side of equations.
- Suffixing (true|false|none|default|...) currently selects between
no suffix and ('_'+region).
SourceFiles
faOption.C
@ -139,9 +136,6 @@ protected:
//- The model region name (finite-area)
word regionName_;
//- Suffix hint for variable names (default: "")
word suffixHint_;
// Protected Member Functions
@ -279,7 +273,7 @@ public:
//- The source name
const word& name() const noexcept { return name_; }
//- Return const access to the volume mesh
//- Return const access to the mesh database
const fvMesh& mesh() const noexcept { return mesh_; }
//- Return dictionary
@ -310,18 +304,6 @@ public:
inline bool active(bool on) noexcept;
// Helper Functions
//- The suffix hint when generating variable names
const word& suffixHint() const noexcept { return suffixHint_; }
//- Return the concatenation of \p base and the suffix hint
inline word suffixed(const char* base) const;
//- Return the concatenation of \p base and the suffix hint
inline word suffixed(const std::string& base) const;
// Checks
//- Is the source active?

View File

@ -61,18 +61,4 @@ inline const Foam::volSurfaceMapping& Foam::fa::option::vsm() const
}
inline Foam::word
Foam::fa::option::suffixed(const char* base) const
{
return word(base + suffixHint_);
}
inline Foam::word
Foam::fa::option::suffixed(const std::string& base) const
{
return word(base + suffixHint_);
}
// ************************************************************************* //

View File

@ -27,6 +27,7 @@ License
#include "faOptions.H"
#include "faMesh.H"
#include "faMeshesRegistry.H"
#include "Time.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@ -67,7 +68,7 @@ Foam::IOobject createIOobject
lookupName,
mesh.time().constant(),
// located under finite-area
faMesh::Registry(mesh),
faMeshesRegistry::New(mesh).thisDb(),
IOobjectOption::MUST_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::REGISTER
@ -187,7 +188,11 @@ Foam::fa::options& Foam::fa::options::New
);
// Registered under finite-area?
auto* ptr = faMesh::Registry(mesh).getObjectPtr<fa::options>(lookupName);
auto* ptr =
faMeshesRegistry::New(mesh).thisDb().getObjectPtr<fa::options>
(
lookupName
);
if (!ptr && polyMesh::regionName(defaultAreaName).empty())
{

View File

@ -30,6 +30,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "volFields.H"
#include "famSup.H"
#include "zeroGradientFaPatchFields.H"
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
@ -56,10 +57,13 @@ Foam::fa::contactHeatFluxSource::contactHeatFluxSource
:
fa::faceSetOption(sourceName, modelType, dict, mesh, defaultAreaName),
TName_(dict.getOrDefault<word>("T", "T")),
TprimaryName_(dict.getOrDefault<word>("Tprimary", "T")),
TprimaryName_(dict.get<word>("Tprimary")),
Tprimary_(mesh_.lookupObject<volScalarField>(TprimaryName_)),
thicknessLayers_(),
kappaLayers_(),
contactRes_(0),
curTimeIndex_(-1)
curTimeIndex_(-1),
coupling_()
{
fieldNames_.resize(1, TName_);
@ -127,30 +131,30 @@ void Foam::fa::contactHeatFluxSource::addSup
const label fieldi
)
{
if (!isActive())
if (isActive())
{
return;
}
DebugInfo
<< name() << ": applying source to "
<< eqn.psi().name() << endl;
DebugInfo<< name() << ": applying source to " << eqn.psi().name() << endl;
if (curTimeIndex_ != mesh().time().timeIndex())
{
tmp<DimensionedField<scalar, areaMesh>> htcw(htc());
if (curTimeIndex_ != mesh().time().timeIndex())
{
tmp<DimensionedField<scalar, areaMesh>> htcw(htc());
// Wall temperature - mapped from primary field to finite-area
auto Twall = DimensionedField<scalar, areaMesh>::New
(
"Tw_" + option::name(),
regionMesh(),
dimensionedScalar(dimTemperature, Zero)
);
// Wall temperature - mapped from primary field to finite-area
auto Twall = DimensionedField<scalar, areaMesh>::New
(
"Tw_" + option::name(),
regionMesh(),
dimensionedScalar(dimTemperature, Zero)
);
vsm().mapInternalToSurface<scalar>(Tprimary_, Twall.ref().field());
vsm().mapInternalToSurface<scalar>(Tprimary_, Twall.ref().field());
eqn += -fam::Sp(htcw(), eqn.psi()) + htcw()*Twall;
eqn += -fam::Sp(htcw(), eqn.psi()) + htcw()*Twall;
curTimeIndex_ = mesh().time().timeIndex();
curTimeIndex_ = mesh().time().timeIndex();
}
}
}
@ -186,14 +190,8 @@ bool Foam::fa::contactHeatFluxSource::read(const dictionary& dict)
const labelList& patches = regionMesh().whichPolyPatches();
if (patches.empty())
{
coupling_.clear();
}
else
{
coupling_.resize_null(patches.back()+1);
}
coupling_.clear();
coupling_.resize(patches.empty() ? 0 : (patches.last()+1));
for (const label patchi : patches)
{

View File

@ -40,8 +40,6 @@ Usage
{
// Mandatory entries (unmodifiable)
type contactHeatFluxSource;
// Optional entries (unmodifiable)
Tprimary <TprimaryFieldName>;
// Optional entries (runtime modifiable)
@ -62,19 +60,12 @@ Usage
\table
Property | Description | Type | Reqd | Dflt
type | Type name: contactHeatFluxSource | word | yes | -
Tprimary | Name of primary temperature field | word | yes | -
T | Name of operand temperature field | word | no | T
Tprimary | Name of primary temperature field | word | no | T
thicknessLayers | List of thicknesses of layers | scalarList | no | -
kappaLayers | List of conductivities of layers | scalarList | cndtnl | -
\endtable
Fields/variables used:
\table
Property | Description | Type | Deflt
T | Temperature | shell | T
Tprimary | Temperature | volume | T
\endtable
The inherited entries are elaborated in:
- \link faOption.H \endlink
- \link faceSetOption.H \endlink
@ -133,13 +124,13 @@ class contactHeatFluxSource
// Private Data
//- Name of shell temperature field (default: "T")
//- Name of temperature field
word TName_;
//- Name of volume temperature field (default: "T")
const word TprimaryName_;
//- Name of primary temperature field
word TprimaryName_;
//- Primary (volume) region temperature
//- Primary region temperature
const volScalarField& Tprimary_;
//- Thickness of layers

View File

@ -72,6 +72,10 @@ Foam::fa::externalHeatFluxSource::externalHeatFluxSource
fa::faceSetOption(sourceName, modelType, dict, m, defaultAreaName),
mode_(operationModeNames.get("mode", dict)),
TName_(dict.getOrDefault<word>("T", "T")),
Q_(nullptr),
q_(nullptr),
h_(nullptr),
Ta_(nullptr),
emissivity_(dict.getOrDefault<scalar>("emissivity", 0))
{
fieldNames_.resize(1, TName_);
@ -201,11 +205,6 @@ bool Foam::fa::externalHeatFluxSource::read(const dictionary& dict)
mode_ = operationModeNames.get("mode", dict);
Q_.reset(nullptr);
q_.reset(nullptr);
h_.reset(nullptr);
Ta_.reset(nullptr);
switch (mode_)
{
case fixedPower:

View File

@ -103,7 +103,7 @@ Usage
\verbatim
power | Use fixed power (supply Q)
flux | Use fixed heat flux (supply q)
coefficient | Use fixes heat transfer coefficient (supply h and Ta)
coefficient | Use fixes heat transfer coefficient (supply h and T)
\endverbatim
See also
@ -159,7 +159,7 @@ private:
//- Operation mode
operationMode mode_;
//- Name of shell temperature field (default: "T")
//- Name of temperature field
word TName_;
//- Heat power [W]
@ -174,7 +174,7 @@ private:
//- Ambient temperature [K]
autoPtr<Function1<scalar>> Ta_;
//- Surface emissivity for radiative transfer to ambient (default: 0)
//- Optional surface emissivity for radiative transfer to ambient
scalar emissivity_;

View File

@ -42,12 +42,6 @@ namespace fa
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Implementation
#include "jouleHeatingSourceImpl.cxx"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::fa::jouleHeatingSource::jouleHeatingSource
@ -65,7 +59,7 @@ Foam::fa::jouleHeatingSource::jouleHeatingSource
(
IOobject
(
suffixed(IOobject::scopedName(typeName, "V")),
IOobject::scopedName(typeName, "V_" + regionName_),
regionMesh().thisDb().time().timeName(),
regionMesh().thisDb(),
IOobject::MUST_READ,
@ -74,6 +68,8 @@ Foam::fa::jouleHeatingSource::jouleHeatingSource
),
regionMesh()
),
scalarSigmaVsTPtr_(nullptr),
tensorSigmaVsTPtr_(nullptr),
curTimeIndex_(-1),
nIter_(1),
anisotropicElectricalConductivity_(false)
@ -82,8 +78,6 @@ Foam::fa::jouleHeatingSource::jouleHeatingSource
fa::option::resetApplied();
read(dict);
if (anisotropicElectricalConductivity_)
{
Info<< " Using tensor electrical conductivity" << endl;
@ -96,6 +90,8 @@ Foam::fa::jouleHeatingSource::jouleHeatingSource
initialiseSigma(coeffs_, scalarSigmaVsTPtr_);
}
read(dict);
}
@ -109,14 +105,12 @@ void Foam::fa::jouleHeatingSource::addSup
const label fieldi
)
{
if (!isActive())
if (isActive())
{
return;
}
DebugInfo
<< name() << ": applying source to "
<< eqn.psi().name() << endl;
DebugInfo<< name() << ": applying source to " << eqn.psi().name() << endl;
{
if (curTimeIndex_ != mesh().time().timeIndex())
{
for (label i = 0; i < nIter_; ++i)
@ -124,7 +118,8 @@ void Foam::fa::jouleHeatingSource::addSup
if (anisotropicElectricalConductivity_)
{
// Update sigma as a function of T if required
const auto& sigma = updateSigma(tensorSigmaVsTPtr_);
const areaTensorField& sigma =
updateSigma(tensorSigmaVsTPtr_);
// Solve the electrical potential equation
faScalarMatrix VEqn(fam::laplacian(h*sigma, V_));
@ -134,7 +129,8 @@ void Foam::fa::jouleHeatingSource::addSup
else
{
// Update sigma as a function of T if required
const auto& sigma = updateSigma(scalarSigmaVsTPtr_);
const areaScalarField& sigma =
updateSigma(scalarSigmaVsTPtr_);
// Solve the electrical potential equation
faScalarMatrix VEqn(fam::laplacian(h*sigma, V_));
@ -149,7 +145,7 @@ void Foam::fa::jouleHeatingSource::addSup
// Add the Joule heating contribution
const word sigmaName
(
IOobject::scopedName(typeName, "sigma") + suffixHint()
IOobject::scopedName(typeName, "sigma_" + regionName_)
);
areaVectorField gradV("gradV", fac::grad(V_));
@ -166,13 +162,15 @@ void Foam::fa::jouleHeatingSource::addSup
if (anisotropicElectricalConductivity_)
{
const auto& sigma = obr.lookupObject<areaTensorField>(sigmaName);
const auto& sigma =
obr.lookupObject<areaTensorField>(sigmaName);
tsource = (h*sigma & gradV) & gradV;
}
else
{
const auto& sigma = obr.lookupObject<areaScalarField>(sigmaName);
const auto& sigma =
obr.lookupObject<areaScalarField>(sigmaName);
tsource = (h*sigma*gradV) & gradV;
}

View File

@ -115,30 +115,23 @@ Usage
- \link faceSetOption.H \endlink
Note
If the \c sigma entry is present, the electrical conductivity is specified
as a function of temperature using a \c Function1 type, otherwise
the \c sigma field will be read from file.
When the \c anisotropicElectricalConductivity flag is set to \c true,
\c sigma should be specified as a \em tensor quantity instead of as
an isotropic \em scalar quantity.
BREAKING Naming changes from 2056 and earlier for the fields:
\table
Field | Scoped names | Scoped names (old)
V | \<scope\>:V (suffix) | \<scope\>:V_ + regionName
sigma | \<scope\>:sigma (suffix) | \<scope\>:sigma_ + regionName
\endtable
It is possible to replicate the older naming by specifying
the \c suffixing to ('_' + regionName).
- \c anisotropicElectricalConductivity=true enables
anisotropic (tensorial) electrical conductivity.
- \c anisotropicElectricalConductivity=false enables
isotropic (scalar) electrical conductivity.
- The electrical conductivity can be specified using either:
- If the \c sigma entry is present the electrical conductivity is specified
as a function of temperature using a \c Function1 type.
- If not present the \c sigma field will be read from file.
- If the \c anisotropicElectricalConductivity flag is set to \c true,
\c sigma should be specified as a tensor quantity.
See also
- Foam::Function1
- Foam::fv::jouleHeatingSource
SourceFiles
jouleHeatingSource.cxx
jouleHeatingSourceImpl.cxx
jouleHeatingSource.C
jouleHeatingSourceTemplates.C
\*---------------------------------------------------------------------------*/
@ -195,13 +188,13 @@ class jouleHeatingSource
void initialiseSigma
(
const dictionary& dict,
autoPtr<Function1<Type>>& sigmaFunctionPtr
autoPtr<Function1<Type>>& sigmaVsTPtr
);
//- Update the electrical conductivity field
template<class Type>
const GeometricField<Type, faPatchField, areaMesh>&
updateSigma(const autoPtr<Function1<Type>>& sigmaFunctionPtr) const;
updateSigma(const autoPtr<Function1<Type>>& sigmaVsTPtr) const;
public:
@ -236,22 +229,22 @@ public:
// Member Functions
// Evaluation
// Evaluation
//- Add explicit contribution to energy equation
virtual void addSup
(
const areaScalarField& h,
const areaScalarField& rho,
faMatrix<scalar>& eqn,
const label fieldi
);
//- Add explicit contribution to compressible momentum equation
virtual void addSup
(
const areaScalarField& h,
const areaScalarField& rho,
faMatrix<scalar>& eqn,
const label fieldi
);
// IO
// IO
//- Read source dictionary
virtual bool read(const dictionary& dict);
//- Read source dictionary
virtual bool read(const dictionary& dict);
};
@ -262,6 +255,12 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "jouleHeatingSourceTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2019-2025 OpenCFD Ltd.
Copyright (C) 2019-2024 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -33,21 +33,16 @@ template<class Type>
void Foam::fa::jouleHeatingSource::initialiseSigma
(
const dictionary& dict,
autoPtr<Function1<Type>>& sigmaFunctionPtr
autoPtr<Function1<Type>>& sigmaVsTPtr
)
{
typedef GeometricField<Type, faPatchField, areaMesh> FieldType;
const word sigmaName
(
IOobject::scopedName(typeName, "sigma") + suffixHint()
);
auto& obr = regionMesh().thisDb();
IOobject io
(
sigmaName,
IOobject::scopedName(typeName, "sigma_" + regionName_),
obr.time().timeName(),
obr,
IOobject::NO_READ,
@ -57,11 +52,11 @@ void Foam::fa::jouleHeatingSource::initialiseSigma
autoPtr<FieldType> tsigma;
// Is sigma defined using a Function1 type?
sigmaFunctionPtr = Function1<Type>::NewIfPresent("sigma", dict, &mesh_);
if (sigmaFunctionPtr)
if (dict.found("sigma"))
{
// Sigma to be defined using a Function1 type
sigmaVsTPtr = Function1<Type>::New("sigma", dict, &mesh_);
tsigma.reset
(
new FieldType
@ -94,34 +89,31 @@ template<class Type>
const Foam::GeometricField<Type, Foam::faPatchField, Foam::areaMesh>&
Foam::fa::jouleHeatingSource::updateSigma
(
const autoPtr<Function1<Type>>& sigmaFunctionPtr
const autoPtr<Function1<Type>>& sigmaVsTPtr
) const
{
typedef GeometricField<Type, faPatchField, areaMesh> FieldType;
const word sigmaName
(
IOobject::scopedName(typeName, "sigma") + suffixHint()
);
const auto& obr = regionMesh().thisDb();
auto& sigma = obr.lookupObjectRef<FieldType>(sigmaName);
auto& sigma =
obr.lookupObjectRef<FieldType>
(
IOobject::scopedName(typeName, "sigma_" + regionName_)
);
if (!sigmaFunctionPtr)
if (!sigmaVsTPtr)
{
// Electrical conductivity field, sigma, was specified by the user
return sigma;
}
const auto& sigmaFunction = sigmaFunctionPtr();
const auto& T = obr.lookupObject<areaScalarField>(TName_);
// Internal field
forAll(sigma, i)
{
sigma[i] = sigmaFunction.value(T[i]);
sigma[i] = sigmaVsTPtr->value(T[i]);
}
@ -135,7 +127,7 @@ Foam::fa::jouleHeatingSource::updateSigma
const scalarField& Tbf = T.boundaryField()[patchi];
forAll(pf, facei)
{
pf[facei] = sigmaFunction.value(Tbf[facei]);
pf[facei] = sigmaVsTPtr->value(Tbf[facei]);
}
}
}

View File

@ -28,7 +28,6 @@ License
#include "faMesh.H"
#include "faMeshBoundaryHalo.H"
#include "faMeshesRegistry.H"
#include "faGlobalMeshData.H"
#include "Time.H"
#include "polyMesh.H"
@ -160,12 +159,6 @@ const Foam::objectRegistry* Foam::faMesh::registry(const polyMesh& pMesh)
// return obr.cfindObject<objectRegistry>(faMesh::prefix());
// }
// Forwarding
const Foam::objectRegistry& Foam::faMesh::Registry(const polyMesh& pMesh)
{
return faMeshesRegistry::Registry(pMesh);
}
const Foam::faMesh& Foam::faMesh::mesh
(

View File

@ -751,14 +751,10 @@ public:
// Database
//- Find the singleton parent registry (on the polyMesh)
//- that contains all objects related to finite-area.
//- The parent registry containing all finite-area meshes
//- on the polyMesh.
static const objectRegistry* registry(const polyMesh& pMesh);
//- Return the singleton parent registry (on the polyMesh)
//- that contains all objects related to finite-area.
static const objectRegistry& Registry(const polyMesh& pMesh);
//- The single-region finite-area region on the polyMesh.
//- Uses lookupObject semantics - Fatal if non-existent
static const faMesh& mesh(const polyMesh& pMesh);

View File

@ -27,8 +27,8 @@ License
#include "faMesh.H"
#include "faMeshesRegistry.H"
#include "polyMesh.H"
#include "Time.H"
#include "polyMesh.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
@ -43,7 +43,6 @@ Foam::faMeshRegistry::faMeshRegistry
IOobject
(
(areaName.empty() ? polyMesh::defaultRegion : areaName),
mesh.thisDb().time().timeName(),
faMeshesRegistry::New(mesh).thisDb(),
IOobjectOption::NO_READ,
IOobjectOption::AUTO_WRITE,

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2023-2025 OpenCFD Ltd.
Copyright (C) 2023-2024 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -86,8 +86,6 @@ bool Foam::faMeshesRegistry::writeObject
const bool writeOnProc
) const
{
// Could also restrict to faMesh only...
//
// for (const faMesh& m : objects_.csorted<faMesh>())
// {
// m.write(writeOnProc);

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2023-2025 OpenCFD Ltd.
Copyright (C) 2023-2024 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -144,34 +144,14 @@ public:
explicit faMeshesRegistry(const polyMesh& mesh);
// Factory Methods
//- Return the registry of objects on the singleton.
// Same as New(mesh).thisDb()
FOAM_NO_DANGLING_REFERENCE //< Reference stored in registry
static const objectRegistry& Registry(const polyMesh& mesh)
{
return MeshObject_type::New(mesh).thisDb();
}
// Database
// It redirects to the private objects but uses some
// objectRegistry method naming
//- The registry of the objects
//- Return the object registry
const objectRegistry& thisDb() const noexcept
{
return objects_;
}
//- Local relative to time
const fileName& dbDir() const
{
return objects_.dbDir();
}
//- The polyMesh reference
const polyMesh& mesh() const noexcept
{

View File

@ -227,9 +227,15 @@ bool Foam::cyclicACMIFvPatchField<Type>::all_ready() const
recvRequests_.start(),
recvRequests_.size()
)
&& UPstream::finishedRequests
(
recvRequests1_.start(),
recvRequests1_.size()
)
)
{
recvRequests_.clear();
recvRequests1_.clear();
++done;
}
@ -240,9 +246,15 @@ bool Foam::cyclicACMIFvPatchField<Type>::all_ready() const
sendRequests_.start(),
sendRequests_.size()
)
&& UPstream::finishedRequests
(
sendRequests1_.start(),
sendRequests1_.size()
)
)
{
sendRequests_.clear();
sendRequests1_.clear();
++done;
}
@ -260,9 +272,15 @@ bool Foam::cyclicACMIFvPatchField<Type>::ready() const
recvRequests_.start(),
recvRequests_.size()
)
&& UPstream::finishedRequests
(
recvRequests1_.start(),
recvRequests1_.size()
)
)
{
recvRequests_.clear();
recvRequests1_.clear();
if
(
@ -271,9 +289,15 @@ bool Foam::cyclicACMIFvPatchField<Type>::ready() const
sendRequests_.start(),
sendRequests_.size()
)
&& UPstream::finishedRequests
(
sendRequests1_.start(),
sendRequests1_.size()
)
)
{
sendRequests_.clear();
sendRequests1_.clear();
}
return true;
@ -483,7 +507,7 @@ void Foam::cyclicACMIFvPatchField<Type>::initEvaluate
const Field<Type> pnf(this->primitiveField(), nbrFaceCells);
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -494,14 +518,20 @@ void Foam::cyclicACMIFvPatchField<Type>::initEvaluate
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
cyclicACMIPatch_.initInterpolate
(
pnf,
sendRequests_,
sendBufs_,
recvRequests_,
recvBufs_
sendBufs_,
recvBufs_,
sendRequests1_,
recvRequests1_,
sendBufs1_,
recvBufs1_
);
}
}
@ -547,12 +577,15 @@ void Foam::cyclicACMIFvPatchField<Type>::evaluate
(
Field<Type>::null(), // Not used for distributed
recvRequests_,
recvBufs_
recvBufs_,
recvRequests1_,
recvBufs1_
).ptr()
);
// Receive requests all handled by last function call
recvRequests_.clear();
recvRequests1_.clear();
if (doTransform())
{
@ -608,7 +641,7 @@ void Foam::cyclicACMIFvPatchField<Type>::initInterfaceMatrixUpdate
transformCoupleField(pnf, cmpt);
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -619,14 +652,20 @@ void Foam::cyclicACMIFvPatchField<Type>::initInterfaceMatrixUpdate
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
cyclicACMIPatch_.initInterpolate
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_
scalarSendBufs_,
scalarRecvBufs_,
sendRequests1_,
recvRequests1_,
scalarSendBufs1_,
scalarRecvBufs1_
);
}
@ -681,11 +720,14 @@ void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix
(
solveScalarField::null(), // Not used for distributed
recvRequests_,
scalarRecvBufs_
scalarRecvBufs_,
recvRequests1_,
scalarRecvBufs1_
);
// Receive requests all handled by last function call
recvRequests_.clear();
recvRequests1_.clear();
}
else
{
@ -738,7 +780,7 @@ void Foam::cyclicACMIFvPatchField<Type>::initInterfaceMatrixUpdate
transformCoupleField(pnf);
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -749,14 +791,20 @@ void Foam::cyclicACMIFvPatchField<Type>::initInterfaceMatrixUpdate
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
cyclicACMIPatch_.initInterpolate
(
pnf,
sendRequests_,
sendBufs_,
recvRequests_,
recvBufs_
sendBufs_,
recvBufs_,
sendRequests1_,
recvRequests1_,
sendBufs1_,
recvBufs1_
);
}
@ -798,11 +846,14 @@ void Foam::cyclicACMIFvPatchField<Type>::updateInterfaceMatrix
(
Field<Type>::null(), // Not used for distributed
recvRequests_,
recvBufs_
recvBufs_,
recvRequests1_,
recvBufs1_
);
// Receive requests all handled by last function call
recvRequests_.clear();
recvRequests1_.clear();
}
else
{

View File

@ -102,6 +102,28 @@ class cyclicACMIFvPatchField
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs_;
// Only used for AMI caching
//- Current range of send requests (non-blocking)
mutable labelRange sendRequests1_;
//- Current range of recv requests (non-blocking)
mutable labelRange recvRequests1_;
//- Send buffers
mutable PtrList<List<Type>> sendBufs1_;
//- Receive buffers_
mutable PtrList<List<Type>> recvBufs1_;
//- Scalar send buffers
mutable PtrList<List<solveScalar>> scalarSendBufs1_;
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs1_;
//- Neighbour coupled internal cell data
mutable autoPtr<Field<Type>> patchNeighbourFieldPtr_;

View File

@ -207,9 +207,15 @@ bool Foam::cyclicAMIFvPatchField<Type>::all_ready() const
recvRequests_.start(),
recvRequests_.size()
)
&& UPstream::finishedRequests
(
recvRequests1_.start(),
recvRequests1_.size()
)
)
{
recvRequests_.clear();
recvRequests1_.clear();
++done;
}
@ -220,9 +226,15 @@ bool Foam::cyclicAMIFvPatchField<Type>::all_ready() const
sendRequests_.start(),
sendRequests_.size()
)
&& UPstream::finishedRequests
(
sendRequests1_.start(),
sendRequests1_.size()
)
)
{
sendRequests_.clear();
sendRequests1_.clear();
++done;
}
@ -240,9 +252,15 @@ bool Foam::cyclicAMIFvPatchField<Type>::ready() const
recvRequests_.start(),
recvRequests_.size()
)
&& UPstream::finishedRequests
(
recvRequests1_.start(),
recvRequests1_.size()
)
)
{
recvRequests_.clear();
recvRequests1_.clear();
if
(
@ -251,9 +269,15 @@ bool Foam::cyclicAMIFvPatchField<Type>::ready() const
sendRequests_.start(),
sendRequests_.size()
)
&& UPstream::finishedRequests
(
sendRequests1_.start(),
sendRequests1_.size()
)
)
{
sendRequests_.clear();
sendRequests1_.clear();
}
return true;
@ -319,9 +343,18 @@ Foam::cyclicAMIFvPatchField<Type>::getNeighbourField
template<class Type>
bool Foam::cyclicAMIFvPatchField<Type>::cacheNeighbourField()
bool Foam::cyclicAMIFvPatchField<Type>::cacheNeighbourField() const
{
return (FieldBase::localBoundaryConsistency() != 0);
// const auto& AMI = this->ownerAMI();
// if (AMI.cacheActive())
// {
// return false;
// }
// else
{
return (FieldBase::localBoundaryConsistency() != 0);
}
}
@ -350,11 +383,12 @@ Foam::cyclicAMIFvPatchField<Type>::getPatchNeighbourField
}
const auto& fvp = this->patch();
const auto& mesh = fvp.boundaryMesh().mesh();
if
(
patchNeighbourFieldPtr_
&& !fvp.boundaryMesh().mesh().upToDatePoints(this->internalField())
&& !mesh.upToDatePoints(this->internalField())
)
{
//DebugPout
@ -418,7 +452,8 @@ template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::cyclicAMIFvPatchField<Type>::patchNeighbourField() const
{
return this->getPatchNeighbourField(true); // checkCommunicator = true
// checkCommunicator = true
return this->getPatchNeighbourField(true);
}
@ -491,7 +526,7 @@ void Foam::cyclicAMIFvPatchField<Type>::initEvaluate
const cyclicAMIPolyPatch& cpp = cyclicAMIPatch_.cyclicAMIPatch();
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -502,14 +537,20 @@ void Foam::cyclicAMIFvPatchField<Type>::initEvaluate
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
cpp.initInterpolate
(
pnf,
sendRequests_,
sendBufs_,
recvRequests_,
recvBufs_
sendBufs_,
recvBufs_,
sendRequests1_,
recvRequests1_,
sendBufs1_,
recvBufs1_
);
}
}
@ -562,12 +603,15 @@ void Foam::cyclicAMIFvPatchField<Type>::evaluate
Field<Type>::null(), // Not used for distributed
recvRequests_,
recvBufs_,
recvRequests1_,
recvBufs1_,
defaultValues
).ptr()
);
// Receive requests all handled by last function call
recvRequests_.clear();
recvRequests1_.clear();
if (doTransform())
{
@ -618,7 +662,7 @@ void Foam::cyclicAMIFvPatchField<Type>::initInterfaceMatrixUpdate
const cyclicAMIPolyPatch& cpp = cyclicAMIPatch_.cyclicAMIPatch();
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -629,14 +673,20 @@ void Foam::cyclicAMIFvPatchField<Type>::initInterfaceMatrixUpdate
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
cpp.initInterpolate
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_
scalarSendBufs_,
scalarRecvBufs_,
sendRequests1_,
recvRequests1_,
scalarSendBufs1_,
scalarRecvBufs1_
);
}
@ -691,11 +741,14 @@ void Foam::cyclicAMIFvPatchField<Type>::updateInterfaceMatrix
solveScalarField::null(), // Not used for distributed
recvRequests_,
scalarRecvBufs_,
recvRequests1_,
scalarRecvBufs1_,
defaultValues
);
// Receive requests all handled by last function call
recvRequests_.clear();
recvRequests1_.clear();
}
else
{
@ -757,7 +810,7 @@ void Foam::cyclicAMIFvPatchField<Type>::initInterfaceMatrixUpdate
const cyclicAMIPolyPatch& cpp = cyclicAMIPatch_.cyclicAMIPatch();
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -768,14 +821,20 @@ void Foam::cyclicAMIFvPatchField<Type>::initInterfaceMatrixUpdate
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
cpp.initInterpolate
(
pnf,
sendRequests_,
sendBufs_,
recvRequests_,
recvBufs_
sendBufs_,
recvBufs_,
sendRequests1_,
recvRequests1_,
sendBufs1_,
recvBufs1_
);
}
@ -829,11 +888,14 @@ void Foam::cyclicAMIFvPatchField<Type>::updateInterfaceMatrix
Field<Type>::null(), // Not used for distributed
recvRequests_,
recvBufs_,
recvRequests1_,
recvBufs1_,
defaultValues
);
// Receive requests all handled by last function call
recvRequests_.clear();
recvRequests1_.clear();
}
else
{
@ -918,7 +980,7 @@ void Foam::cyclicAMIFvPatchField<Type>::manipulateMatrix
}
// Set internalCoeffs and boundaryCoeffs in the assembly matrix
// on clyclicAMI patches to be used in the individual matrix by
// on cyclicAMI patches to be used in the individual matrix by
// matrix.flux()
if (matrix.psi(mat).mesh().fluxRequired(this->internalField().name()))
{

View File

@ -113,6 +113,28 @@ class cyclicAMIFvPatchField
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs_;
// Only used for AMI caching
//- Current range of send requests (non-blocking)
mutable labelRange sendRequests1_;
//- Current range of recv requests (non-blocking)
mutable labelRange recvRequests1_;
//- Send buffers
mutable PtrList<List<Type>> sendBufs1_;
//- Receive buffers_
mutable PtrList<List<Type>> recvBufs1_;
//- Scalar send buffers
mutable PtrList<List<solveScalar>> scalarSendBufs1_;
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs1_;
//- Neighbour coupled internal cell data
mutable autoPtr<Field<Type>> patchNeighbourFieldPtr_;
@ -134,7 +156,7 @@ class cyclicAMIFvPatchField
virtual bool all_ready() const;
//- Use neighbour field caching
static bool cacheNeighbourField();
bool cacheNeighbourField() const;
//- Return neighbour coupled internal cell data
tmp<Field<Type>> getNeighbourField(const UList<Type>&) const;

View File

@ -209,6 +209,23 @@ public:
template<class TrackingData>
inline bool equal(const smoothData&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const smoothData& f0,
const label i1,
const smoothData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2013 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
Copyright (C) 2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -191,6 +191,45 @@ inline bool Foam::smoothData::equal
}
template<class TrackingData>
inline bool Foam::smoothData::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const smoothData& f0,
const label i1,
const smoothData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
// TBD. What to interpolate? Do we have a position? cell/face centre?
if (!valid(td))
{
if (f0.valid(td))
{
operator=(f0);
}
else
{
operator=(f1);
}
}
else if (f0.valid(td))
{
operator=(f0);
}
else
{
operator=(f1);
}
return true;
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::smoothData::operator==

View File

@ -207,6 +207,23 @@ public:
template<class TrackingData>
inline bool equal(const sweepData&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const sweepData& f0,
const label i1,
const sweepData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -210,6 +210,45 @@ inline bool Foam::sweepData::equal
}
template<class TrackingData>
inline bool Foam::sweepData::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const sweepData& f0,
const label i1,
const sweepData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
// TBD. What to interpolate? Do we have a position? cell/face centre?
if (!valid(td))
{
if (f0.valid(td))
{
operator=(f0);
}
else
{
operator=(f1);
}
}
else if (f0.valid(td))
{
operator=(f0);
}
else
{
operator=(f1);
}
return true;
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::sweepData::operator==

View File

@ -220,9 +220,14 @@ public:
(
const Field<Type>& fld,
labelRange& sendRequests,
PtrList<List<Type>>& sendBuffers,
labelRange& recvRequests,
PtrList<List<Type>>& recvBuffers
PtrList<List<Type>>& sendBuffers,
PtrList<List<Type>>& recvBuffers,
labelRange& sendRequests1,
labelRange& recvRequests1,
PtrList<List<Type>>& sendBuffers1,
PtrList<List<Type>>& recvBuffers1
) const
{
// Make sure areas are up-to-date
@ -232,9 +237,14 @@ public:
(
fld,
sendRequests,
sendBuffers,
recvRequests,
recvBuffers
sendBuffers,
recvBuffers,
sendRequests1,
recvRequests1,
sendBuffers1,
recvBuffers1
);
}
@ -243,7 +253,9 @@ public:
(
const Field<Type>& localFld,
const labelRange& requests, // The receive requests
const PtrList<List<Type>>& recvBuffers
const PtrList<List<Type>>& recvBuffers,
const labelRange& requests1, // The receive requests
const PtrList<List<Type>>& recvBuffers1
) const
{
return cyclicACMIPolyPatch_.interpolate
@ -251,6 +263,8 @@ public:
localFld,
requests,
recvBuffers,
requests1,
recvBuffers1,
UList<Type>()
);
}

View File

@ -14,7 +14,7 @@ $(derivedSources)/buoyancyEnergy/buoyancyEnergy.C
$(derivedSources)/buoyancyForce/buoyancyForce.C
$(derivedSources)/directionalPressureGradientExplicitSource/directionalPressureGradientExplicitSource.C
$(derivedSources)/explicitPorositySource/explicitPorositySource.C
$(derivedSources)/jouleHeatingSource/jouleHeatingSource.cxx
$(derivedSources)/jouleHeatingSource/jouleHeatingSource.C
$(derivedSources)/meanVelocityForce/meanVelocityForce.C
$(derivedSources)/meanVelocityForce/patchMeanVelocityForce/patchMeanVelocityForce.C
$(derivedSources)/multiphaseStabilizedTurbulence/multiphaseStabilizedTurbulence.C

View File

@ -29,6 +29,7 @@ License
#include "fvMatrices.H"
#include "fvmLaplacian.H"
#include "fvcGrad.H"
#include "zeroGradientFvPatchField.H"
#include "basicThermo.H"
#include "addToRunTimeSelectionTable.H"
@ -44,12 +45,6 @@ namespace fv
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Implementation
#include "jouleHeatingSourceImpl.cxx"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::tmp<Foam::volSymmTensorField>
@ -120,8 +115,11 @@ Foam::fv::jouleHeatingSource::jouleHeatingSource
),
mesh
),
curTimeIndex_(-1),
anisotropicElectricalConductivity_(false)
anisotropicElectricalConductivity_(false),
scalarSigmaVsTPtr_(nullptr),
vectorSigmaVsTPtr_(nullptr),
csysPtr_(nullptr),
curTimeIndex_(-1)
{
// Set the field name to that of the energy
// field from which the temperature is obtained
@ -164,7 +162,7 @@ void Foam::fv::jouleHeatingSource::addSup
else
{
// Update sigma as a function of T if required
const auto& sigma = updateSigma(scalarSigmaVsTPtr_);
const volScalarField& sigma = updateSigma(scalarSigmaVsTPtr_);
// Solve the electrical potential equation
fvScalarMatrix VEqn(fvm::laplacian(sigma, V_));
@ -228,7 +226,7 @@ bool Foam::fv::jouleHeatingSource::read(const dictionary& dict)
initialiseSigma(coeffs_, scalarSigmaVsTPtr_);
csysPtr_.reset(nullptr); // Do not need coordinate system
csysPtr_.clear(); // Do not need coordinate system
}
return true;

View File

@ -137,28 +137,25 @@ Note
anisotropic (vectorial) electrical conductivity.
- \c anisotropicElectricalConductivity=false enables
isotropic (scalar) electrical conductivity.
The electrical conductivity can be specified using either:
- The electrical conductivity can be specified using either:
- If the \c sigma entry is present the electrical conductivity is specified
as a function of temperature using a \c Function1 type.
- If not present the \c sigma field will be read from file.
- If the \c anisotropicElectricalConductivity flag is set to \c true,
\c sigma should be specified as a vector quantity.
.
See also
- Foam::Function1
- Foam::coordSystem
- Foam::fa::jouleHeatingSource
SourceFiles
jouleHeatingSource.cxx
jouleHeatingSourceImpl.cxx
jouleHeatingSource.C
jouleHeatingSourceTemplates.C
\*---------------------------------------------------------------------------*/
#ifndef Foam_fv_jouleHeatingSource_H
#define Foam_fv_jouleHeatingSource_H
#ifndef fv_jouleHeatingSource_H
#define fv_jouleHeatingSource_H
#include "fvOption.H"
#include "Function1.H"
@ -188,6 +185,9 @@ class jouleHeatingSource
//- Electrical potential field / [V]
volScalarField V_;
//- Flag to indicate that the electrical conductivity is anisotropic
bool anisotropicElectricalConductivity_;
//- Electrical conductivity as a scalar function of temperature
autoPtr<Function1<scalar>> scalarSigmaVsTPtr_;
@ -200,9 +200,6 @@ class jouleHeatingSource
//- Current time index (used for updating)
label curTimeIndex_;
//- Flag to indicate that the electrical conductivity is anisotropic
bool anisotropicElectricalConductivity_;
// Private Member Functions
@ -217,13 +214,13 @@ class jouleHeatingSource
void initialiseSigma
(
const dictionary& dict,
autoPtr<Function1<Type>>& sigmaFunctionPtr
autoPtr<Function1<Type>>& sigmaVsTPtr
);
//- Update the electrical conductivity field
template<class Type>
const GeometricField<Type, fvPatchField, volMesh>&
updateSigma(const autoPtr<Function1<Type>>& sigmaFunctionPtr) const;
updateSigma(const autoPtr<Function1<Type>>& sigmaVsTPtr) const;
public:
@ -256,21 +253,21 @@ public:
// Member Functions
// Evaluation
// Evaluation
//- Add explicit contribution to energy equation
virtual void addSup
(
const volScalarField& rho,
fvMatrix<scalar>& eqn,
const label fieldi
);
//- Add explicit contribution to compressible momentum equation
virtual void addSup
(
const volScalarField& rho,
fvMatrix<scalar>& eqn,
const label fieldi
);
// IO
// IO
//- Read source dictionary
virtual bool read(const dictionary& dict);
//- Read source dictionary
virtual bool read(const dictionary& dict);
};
@ -281,6 +278,12 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "jouleHeatingSourceTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2016-2025 OpenCFD Ltd.
Copyright (C) 2016-2024 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -27,25 +27,18 @@ License
#include "emptyFvPatchField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type>
void Foam::fv::jouleHeatingSource::initialiseSigma
(
const dictionary& dict,
autoPtr<Function1<Type>>& sigmaFunctionPtr
autoPtr<Function1<Type>>& sigmaVsTPtr
)
{
typedef GeometricField<Type, fvPatchField, volMesh> FieldType;
const word sigmaName
(
IOobject::scopedName(typeName, "sigma")
);
IOobject io
(
sigmaName,
IOobject::scopedName(typeName, "sigma"),
mesh_.time().timeName(),
mesh_.thisDb(),
IOobject::NO_READ,
@ -55,11 +48,11 @@ void Foam::fv::jouleHeatingSource::initialiseSigma
autoPtr<FieldType> tsigma;
// Is sigma defined using a Function1 ?
sigmaFunctionPtr = Function1<Type>::NewIfPresent("sigma", dict, &mesh_);
if (sigmaFunctionPtr)
if (dict.found("sigma"))
{
// Sigma to be defined using a Function1 type
sigmaVsTPtr = Function1<Type>::New("sigma", dict, &mesh_);
tsigma.reset
(
new FieldType
@ -92,32 +85,28 @@ template<class Type>
const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&
Foam::fv::jouleHeatingSource::updateSigma
(
const autoPtr<Function1<Type>>& sigmaFunctionPtr
const autoPtr<Function1<Type>>& sigmaVsTPtr
) const
{
typedef GeometricField<Type, fvPatchField, volMesh> FieldType;
const word sigmaName
auto& sigma = mesh_.lookupObjectRef<FieldType>
(
IOobject::scopedName(typeName, "sigma")
);
auto& sigma = mesh_.lookupObjectRef<FieldType>(sigmaName);
if (!sigmaFunctionPtr)
if (!sigmaVsTPtr)
{
// Electrical conductivity field, sigma, was specified by the user
return sigma;
}
const auto& sigmaFunction = sigmaFunctionPtr();
const auto& T = mesh_.lookupObject<volScalarField>(TName_);
// Internal field
forAll(sigma, i)
{
sigma[i] = sigmaFunction.value(T[i]);
sigma[i] = sigmaVsTPtr->value(T[i]);
}
@ -131,7 +120,7 @@ Foam::fv::jouleHeatingSource::updateSigma
const scalarField& Tbf = T.boundaryField()[patchi];
forAll(pf, facei)
{
pf[facei] = sigmaFunction.value(Tbf[facei]);
pf[facei] = sigmaVsTPtr->value(Tbf[facei]);
}
}
}

View File

@ -251,6 +251,23 @@ public:
template<class TrackingData>
inline bool equal(const wallPoints&, TrackingData&) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const wallPoints& f0,
const label i1,
const wallPoints& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2018-2023 OpenCFD Ltd.
Copyright (C) 2018-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -443,6 +443,41 @@ inline bool Foam::wallPoints::equal
}
template<class TrackingData>
inline bool Foam::wallPoints::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const wallPoints& f0,
const label i1,
const wallPoints& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (valid(td))
{
return false;
}
else if (f0.valid(td))
{
operator=(f0);
return true;
}
else if (f1.valid(td))
{
operator=(f1);
return true;
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::wallPoints::operator==

View File

@ -0,0 +1,651 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
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 "AMICache.H"
#include "AMIInterpolation.H"
#include "mathematicalConstants.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
Foam::scalar Foam::AMICache::cacheThetaTolerance_ = 1e-8;
int Foam::AMICache::debug = 0;
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::scalar Foam::AMICache::getRotationAngle(const point& globalPoint) const
{
if (!coordSysPtr_)
{
FatalErrorInFunction
<< "No co-ordinate system available for theta evaluation"
<< abort(FatalError);
}
scalar theta = coordSysPtr_->localPosition(globalPoint).y();
// Ensure 0 < theta < 2pi
if (mag(theta) < cacheThetaTolerance_)
{
theta = 0;
}
else if (theta < 0)
{
theta += constant::mathematical::twoPi;
}
return theta;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::AMICache::AMICache(const dictionary& dict, const bool toSource)
:
size_(dict.getOrDefault<label>("cacheSize", 0)),
rotationAxis_(dict.getOrDefault<vector>("rotationAxis", Zero)),
rotationCentre_(dict.getOrDefault<point>("rotationCentre", Zero)),
maxThetaDeg_(dict.getOrDefault<scalar>("maxThetaDeg", 5)),
complete_(false),
index0_(-1),
index1_(-1),
interpWeight_(0),
coordSysPtr_(nullptr),
theta_(),
cachedSrcAddress_(),
cachedSrcWeights_(),
cachedSrcWeightsSum_(),
cachedSrcMapPtr_(),
cachedTgtAddress_(),
cachedTgtWeights_(),
cachedTgtWeightsSum_(),
cachedTgtMapPtr_(),
toSource_(toSource)
{
if (size_ != 0)
{
theta_.resize(size_, GREAT);
cachedSrcAddress_.resize(size_);
cachedSrcWeights_.resize(size_);
cachedSrcWeightsSum_.resize(size_);
cachedSrcMapPtr_.resize(size_);
cachedTgtAddress_.resize(size_);
cachedTgtWeights_.resize(size_);
cachedTgtWeightsSum_.resize(size_);
cachedTgtMapPtr_.resize(size_);
}
}
Foam::AMICache::AMICache(const bool toSource)
:
size_(0),
rotationAxis_(Zero),
rotationCentre_(Zero),
complete_(false),
index0_(-1),
index1_(-1),
interpWeight_(0),
coordSysPtr_(nullptr),
theta_(),
cachedSrcAddress_(),
cachedSrcWeights_(),
cachedSrcWeightsSum_(),
cachedSrcMapPtr_(),
cachedTgtAddress_(),
cachedTgtWeights_(),
cachedTgtWeightsSum_(),
cachedTgtMapPtr_(),
toSource_(toSource)
{}
Foam::AMICache::AMICache(const AMICache& cache)
:
size_(cache.size_),
rotationAxis_(cache.rotationAxis_),
rotationCentre_(cache.rotationCentre_),
complete_(cache.complete_),
index0_(cache.index0_),
index1_(cache.index1_),
interpWeight_(cache.interpWeight_),
coordSysPtr_(nullptr), // Need to clone as cylindricalCS
theta_(cache.theta_),
cachedSrcAddress_(cache.cachedSrcAddress_),
cachedSrcWeights_(cache.cachedSrcWeights_),
cachedSrcWeightsSum_(cache.cachedSrcWeightsSum_),
cachedSrcMapPtr_(cache.cachedSrcMapPtr_.size()), // Need to clone
cachedTgtAddress_(cache.cachedTgtAddress_),
cachedTgtWeights_(cache.cachedTgtWeights_),
cachedTgtWeightsSum_(cache.cachedTgtWeightsSum_),
cachedTgtMapPtr_(cache.cachedTgtMapPtr_.size()), // Need to clone
toSource_(cache.toSource_)
{
if (cache.coordSysPtr_)
{
coordSysPtr_.reset(new coordSystem::cylindrical(cache.coordSysPtr_()));
}
forAll(cachedSrcMapPtr_, cachei)
{
cachedSrcMapPtr_[cachei].reset(cache.cachedSrcMapPtr_[cachei].clone());
}
forAll(cachedTgtMapPtr_, cachei)
{
cachedTgtMapPtr_[cachei].reset(cache.cachedTgtMapPtr_[cachei].clone());
}
}
Foam::AMICache::AMICache
(
const AMICache& cache,
const AMIInterpolation& fineAMI,
const labelList& sourceRestrictAddressing,
const labelList& targetRestrictAddressing
)
:
size_(cache.size_),
rotationAxis_(cache.rotationAxis_),
rotationCentre_(cache.rotationCentre_),
complete_(cache.complete_),
index0_(cache.index0_),
index1_(cache.index1_),
interpWeight_(cache.interpWeight_),
coordSysPtr_(nullptr),
theta_(cache.theta_),
cachedSrcAddress_(cache.size_),
cachedSrcWeights_(cache.size_),
cachedSrcWeightsSum_(cache.size_),
cachedSrcMapPtr_(cache.size_),
cachedTgtAddress_(cache.size_),
cachedTgtWeights_(cache.size_),
cachedTgtWeightsSum_(cache.size_),
cachedTgtMapPtr_(cache.size_),
toSource_(cache.toSource_)
{
if (size_ > 0 && fineAMI.comm() != -1)
{
for (label cachei : {index0_, index1_})
{
if (cachei == -1) continue;
scalarField dummySrcMagSf;
labelListList srcAddress;
scalarListList srcWeights;
scalarField srcWeightsSum;
autoPtr<mapDistribute> tgtMapPtr;
AMIInterpolation::agglomerate
(
cache.cachedTgtMapPtr()[cachei],
fineAMI.srcMagSf(),
cache.cachedSrcAddress()[cachei],
cache.cachedSrcWeights()[cachei],
sourceRestrictAddressing,
targetRestrictAddressing,
dummySrcMagSf,
srcAddress,
srcWeights,
srcWeightsSum,
tgtMapPtr,
fineAMI.comm()
);
scalarField dummyTgtMagSf;
labelListList tgtAddress;
scalarListList tgtWeights;
scalarField tgtWeightsSum;
autoPtr<mapDistribute> srcMapPtr;
AMIInterpolation::agglomerate
(
cache.cachedSrcMapPtr()[cachei],
fineAMI.tgtMagSf(),
cache.cachedTgtAddress()[cachei],
cache.cachedTgtWeights()[cachei],
targetRestrictAddressing,
sourceRestrictAddressing,
dummyTgtMagSf,
tgtAddress,
tgtWeights,
tgtWeightsSum,
srcMapPtr,
fineAMI.comm()
);
cachedSrcAddress_[cachei] = srcAddress;
cachedSrcWeights_[cachei] = srcWeights;
cachedSrcWeightsSum_[cachei] = srcWeightsSum;
cachedSrcMapPtr_[cachei] = srcMapPtr.clone();
cachedTgtAddress_[cachei] = tgtAddress;
cachedTgtWeights_[cachei] = tgtWeights;
cachedTgtWeightsSum_[cachei] = tgtWeightsSum;
cachedTgtMapPtr_[cachei] = tgtMapPtr.clone();
}
}
}
Foam::AMICache::AMICache(Istream& is)
:
size_(readLabel(is)),
rotationAxis_(is),
rotationCentre_(is),
maxThetaDeg_(readScalar(is)),
complete_(readBool(is)),
index0_(-1),
index1_(-1),
interpWeight_(0),
coordSysPtr_(nullptr),
theta_(),
cachedSrcAddress_(),
cachedSrcWeights_(),
cachedSrcWeightsSum_(),
cachedSrcMapPtr_(),
cachedTgtAddress_(),
cachedTgtWeights_(),
cachedTgtWeightsSum_(),
cachedTgtMapPtr_()
{
const bitSet goodMap(is);
if (goodMap.size())
{
is >> index0_
>> index1_
>> interpWeight_
>> theta_;
const bool goodCoord(readBool(is));
if (goodCoord)
{
coordSysPtr_.reset(new coordSystem::cylindrical(is));
}
is >> cachedSrcAddress_
>> cachedSrcWeights_
>> cachedSrcWeightsSum_;
cachedSrcMapPtr_.setSize(goodMap.size());
forAll(goodMap, cachei)
{
if (goodMap[cachei])
{
cachedSrcMapPtr_[cachei].reset(new mapDistribute(is));
}
}
is >> cachedTgtAddress_
>> cachedTgtWeights_
>> cachedTgtWeightsSum_;
cachedTgtMapPtr_.setSize(goodMap.size());
forAll(goodMap, cachei)
{
if (goodMap[cachei])
{
cachedTgtMapPtr_[cachei].reset(new mapDistribute(is));
}
}
}
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::AMICache::addToCache
(
const AMIInterpolation& ami,
const point& globalPoint
)
{
DebugPout<< "-- addToCache" << endl;
if (!active())
{
DebugInfo<< "-- addToCache - deactivated" << endl;
return;
}
if (!coordSysPtr_)
{
DebugInfo
<< "Creating rotation co-ordinate system:"
<< " rotationCentre:" << rotationCentre_
<< " rotationAxis:" << rotationAxis_
<< " p:" << globalPoint
<< endl;
coordSysPtr_.reset
(
new coordSystem::cylindrical(rotationCentre_, rotationAxis_)
);
DebugPout<< "Coord sys:" << coordSysPtr_() << endl;
}
// Check if cache is complete
if (!complete_)
{
for (const scalar angle : theta_)
{
// Check against null value; if any are present, cache is incomplete
if (angle > constant::mathematical::twoPi)
{
complete_ = false;
break;
}
}
}
if (!complete_)
{
const scalar theta = getRotationAngle(globalPoint);
const label bini = theta/constant::mathematical::twoPi*size_;
DebugPout<< " -- bini:" << bini << " for theta:" << theta << endl;
// Check if already have entry for this bin
if (theta_[bini] > constant::mathematical::twoPi)
{
DebugPout<< " -- setting cache at index " << bini << endl;
// New entry
theta_[bini] = theta;
cachedSrcAddress_[bini] = ami.srcAddress();
cachedSrcWeights_[bini] = ami.srcWeights();
cachedSrcWeightsSum_[bini] = ami.srcWeightsSum();
if (ami.hasSrcMap())
{
cachedSrcMapPtr_[bini] = ami.srcMap().clone();
}
cachedTgtAddress_[bini] = ami.tgtAddress();
cachedTgtWeights_[bini] = ami.tgtWeights();
cachedTgtWeightsSum_[bini] = ami.tgtWeightsSum();
if (ami.hasTgtMap())
{
cachedTgtMapPtr_[bini] = ami.tgtMap().clone();
}
}
}
}
bool Foam::AMICache::restoreCache(const point& globalPoint)
{
DebugPout<< "-- restoreCache" << endl;
index0_ = -1;
index1_ = -1;
interpWeight_ = -1;
if (!coordSysPtr_ || size_ == -1)
{
return false;
}
const scalar theta = getRotationAngle(globalPoint);
const scalar twoPi = constant::mathematical::twoPi;
const label bini = theta/twoPi*size_;
DebugPout<< " -- bini:" << bini << " for theta:" << theta << endl;
const auto validIndex = [&](const scalar bini)
{
return theta_[bini] < constant::mathematical::twoPi;
};
// Maximum angle in degrees for which to search for cached bins
const scalar maxThetaStencil = maxThetaDeg_*constant::mathematical::pi/180.0;
if
(
validIndex(bini)
&& (
mag(theta - theta_[bini]) < cacheThetaTolerance_
|| mag(theta - twoPi - theta_[bini]) < cacheThetaTolerance_
)
)
{
// Hit cached value - no interpolation needed
// index1_ = -1 indicates no interpolation
index0_ = bini;
index1_ = -1;
interpWeight_ = 0;
DebugInfo
<< " -- t0:" << theta_[index0_] << " theta:" << theta
<< " i0:" << index0_ << " i1:" << index1_
<< " w:" << interpWeight_ << endl;
return true;
}
else
{
// Find indices and values bracketing theta
const label nBin = theta_.size();
// Participating theta values and bin addresses
// - Note we add wrap-around values at start and end
DynamicList<scalar> thetap(nBin+2);
DynamicList<label> binAddresses(nBin+2);
// Initialise wrap-around values
thetap.push_back(0);
binAddresses.push_back(-1);
forAll(theta_, thetai)
{
if (validIndex(thetai))
{
thetap.push_back(theta_[thetai]);
binAddresses.push_back(thetai);
}
}
// Check that we have enough data points for interpolation
// - We added storage for lower wrap-around value, and we then need
// at least 2 additional values for the interpolation
if (thetap.size() < 3)
{
DebugPout<< " -- no cache available" << endl;
return false;
}
// Set wrap-around values if we have sufficient data
thetap[0] = thetap.last() - twoPi;
binAddresses[0] = binAddresses.last();
thetap.push_back(thetap[1] + twoPi);
binAddresses.push_back(binAddresses[1]);
// Find interpolation indices
label loweri = labelMax;
label i = 0;
while (i < thetap.size())
{
if (thetap[i] < theta)
{
loweri = i;
}
else
{
break;
}
++i;
}
if (loweri == labelMax)
{
DebugPout<< " -- no cache available" << endl;
return false;
}
label upperi = labelMax;
i = thetap.size() - 1;
while (i >= 0)
{
if (thetap[i] > theta)
{
upperi = i;
}
else
{
break;
}
--i;
}
if (upperi == labelMax)
{
DebugPout<< " -- no cache available" << endl;
return false;
}
// Ensure distances are valid
if (upperi == loweri)
{
DebugPout
<< " -- no cache available: theta:" << theta
<< " lower:" << loweri << " upper:" << upperi << endl;
return false;
}
if (mag(theta - thetap[loweri]) > maxThetaStencil)
{
DebugPout
<< " -- no cache available: theta:" << theta
<< " lower:" << thetap[loweri] << endl;
return false;
}
if (mag(theta - thetap[upperi]) > maxThetaStencil)
{
DebugPout
<< " -- no cache available: theta:" << theta
<< " upper:" << thetap[upperi] << endl;
return false;
}
index0_ = binAddresses[loweri];
index1_ = binAddresses[upperi];
interpWeight_ =
(theta - theta_[index0_])/(theta_[index1_] - theta_[index0_]);
DebugInfo
<< theta_.size()
<< " -- t0:" << theta_[index0_] << " theta:" << theta
<< " t1:" << theta_[index1_]
<< " i0:" << index0_ << " i1:" << index1_
<< " w:" << interpWeight_ << endl;
return true;
}
// If we get here then no valid cache found within stencil
DebugPout<< " -- no cache available" << endl;
return false;
}
void Foam::AMICache::write(Ostream& os) const
{
if (size_ > 0)
{
os.writeEntry("cacheSize", size_);
os.writeEntry("rotationAxis", rotationAxis_);
os.writeEntry("rotationCentre", rotationCentre_);
os.writeEntry("maxThetaDeg", maxThetaDeg_);
}
}
bool Foam::AMICache::writeData(Ostream& os) const
{
os << token::SPACE<< size_
<< token::SPACE<< rotationAxis_
<< token::SPACE<< rotationCentre_
<< token::SPACE<< maxThetaDeg_
<< token::SPACE<< complete_;
bitSet goodMap(cachedSrcMapPtr_.size());
forAll(goodMap, cachei)
{
goodMap.set(cachei, cachedSrcMapPtr_[cachei].good());
}
os << token::SPACE << goodMap;
if (goodMap.size())
{
os << token::SPACE << index0_
<< token::SPACE << index1_
<< token::SPACE << interpWeight_
<< token::SPACE << theta_;
os << token::SPACE << coordSysPtr_.good();
if (coordSysPtr_.good())
{
os << token::SPACE << coordSysPtr_();
}
os << token::SPACE << cachedSrcAddress_
<< token::SPACE << cachedSrcWeights_
<< token::SPACE << cachedSrcWeightsSum_;
for (const auto& index : goodMap)
{
os << token::SPACE << cachedSrcMapPtr_[index]();
}
os << token::SPACE << cachedTgtAddress_
<< token::SPACE << cachedTgtWeights_
<< token::SPACE << cachedTgtWeightsSum_;
for (const auto& index : goodMap)
{
os << token::SPACE << cachedTgtMapPtr_[index]();
}
}
return true;
}
// ************************************************************************* //

View File

@ -0,0 +1,388 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
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::AMICache
Description
Provides caching of weights and addressing to AMIInterpolation
SourceFiles
AMICache.C
SeeAlso
Foam::AMIInterpolation.H
\*---------------------------------------------------------------------------*/
#ifndef Foam_AMICache_H
#define Foam_AMICache_H
#include "cylindricalCS.H"
#include "mapDistribute.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class AMIInterpolation;
/*---------------------------------------------------------------------------*\
Class AMICache Declaration
\*---------------------------------------------------------------------------*/
class AMICache
{
public:
// Public Data Types
//- Tolerance used when caching the AMI to identify e.g. if the
//- current rotation angle has already been captured
static scalar cacheThetaTolerance_;
static int debug;
private:
// Private Data
//- Cache size
label size_;
//- Axis of rotation for rotational cyclics
vector rotationAxis_;
//- Point on axis of rotation for rotational cyclics
point rotationCentre_;
//- Maximum angle in degrees for which to search for cached bins
scalar maxThetaDeg_;
//- Flag to indicate that the cache is complete
bool complete_;
//- Cache index 0 in lists
label index0_;
//- Cache index 1 in lists
label index1_;
//- Interpolation weight between cache indices 0 and 1
scalar interpWeight_;
//- Local co-ordinate system
autoPtr<coordSystem::cylindrical> coordSysPtr_;
//- List of cached angle snapshots
List<scalar> theta_;
//- List of source addresses
List<labelListList> cachedSrcAddress_;
//- List of source weights
List<scalarListList> cachedSrcWeights_;
//- List of source weights sums
List<scalarField> cachedSrcWeightsSum_;
//- List of source parallel maps
List<autoPtr<mapDistribute>> cachedSrcMapPtr_;
//- List of target addresses
List<labelListList> cachedTgtAddress_;
//- List of target weights
List<scalarListList> cachedTgtWeights_;
//- List of target weights sums
List<scalarField> cachedTgtWeightsSum_;
//- List of target parallel maps
List<autoPtr<mapDistribute>> cachedTgtMapPtr_;
//- Flag to indicate interpolation direction
mutable bool toSource_;
// Private Member Functions
//- Get rotation angle for point using local co-ordinate system
scalar getRotationAngle(const point& globalPoint) const;
public:
// Constructors
//- Null constructor
AMICache(const bool toSource = true);
//- Construct from dictionary
AMICache(const dictionary& dict, const bool toSource = true);
//- Construct as copy
AMICache(const AMICache& cache);
//- Construct from agglomeration of AMIInterpolation. Agglomeration
//- passed in as new coarse size and addressing from fine from coarse
AMICache
(
const AMICache& cache,
const AMIInterpolation& fineAMI,
const labelList& sourceRestrictAddressing,
const labelList& targetRestrictAddressing
);
//- Construct from stream
AMICache(Istream& is);
// Member Functions
//- Return true if cache is active
constexpr bool active() const noexcept { return size_ > 0; }
//- Return cache size
constexpr label size() const noexcept { return size_; }
//- Return true if cache is complete
constexpr bool complete() const noexcept { return complete_; }
//- Return cache lower bound index
constexpr label index0() const noexcept { return index0_; }
//- Return cache upper bound index
constexpr label index1() const noexcept { return index1_; }
//- Return cache interpolation weight
constexpr label weight() const noexcept { return interpWeight_; }
//- Return list of cached rotation angles
const List<scalar>& theta() const noexcept { return theta_; }
//- Return List of source addresses
const List<labelListList>& cachedSrcAddress() const noexcept
{
return cachedSrcAddress_;
}
//- Return List of source weights
const List<scalarListList>& cachedSrcWeights() const noexcept
{
return cachedSrcWeights_;
}
//- Return List of source weights sums
const List<scalarField>& cachedSrcWeightsSum() const noexcept
{
return cachedSrcWeightsSum_;
}
//- Return List of source parallel maps
const List<autoPtr<mapDistribute>>& cachedSrcMapPtr() const noexcept
{
return cachedSrcMapPtr_;
}
//- Return List of target addresses
const List<labelListList>& cachedTgtAddress() const noexcept
{
return cachedTgtAddress_;
}
//- Return List of target weights
const List<scalarListList>& cachedTgtWeights() const noexcept
{
return cachedTgtWeights_;
}
//- Return List of target weights sums
const List<scalarField>& cachedTgtWeightsSum() const noexcept
{
return cachedTgtWeightsSum_;
}
//- Return List of target parallel maps
const List<autoPtr<mapDistribute>>& cachedTgtMapPtr() const noexcept
{
return cachedTgtMapPtr_;
}
//- Apply cached evaluation based on user supplied evaluation function
template<class Type, class EvalFunction>
bool apply(List<Type>& result, const EvalFunction& eval) const;
//- Flag that lower bound is applicable
constexpr bool applyLower() const noexcept
{
return index0_ != -1 && index1_ == -1;
}
//- Flag that upper bound is applicable
constexpr bool applyUpper() const noexcept
{
return index0_ == -1 && index1_ != -1;
}
//- Flag that interpolation is applicable
constexpr bool applyInterpolate() const noexcept
{
return index0_ != -1 && index1_ != -1;
}
//- Set the interpolation direction
constexpr bool setDirection(bool toSource) const noexcept
{
toSource_ = toSource;
return toSource_;
}
//- Restore AMI weights and addressing from the cache
bool restoreCache(const point& globalPoint);
//- Add AMI weights and addressing to the cache
void addToCache
(
const AMIInterpolation& ami,
const point& globalPoint
);
//- Check cache index is valid
void checkIndex(const label index) const
{
if ((index < 0) || (index > size_-1))
{
FatalErrorInFunction
<< "Supplied out of bounds: " << index << "/"
<< size_ << abort(FatalError);
}
}
// Helper functions to retrieve cached values at index0 and index1
// Note: uses interpolation direction toSource_
#undef defineMethods01
#define defineMethods01(Src, Tgt, idx) \
const labelListList& c##Src##Address##idx() const \
{ \
checkIndex(index##idx##_); \
return toSource_ ? \
cached##Src##Address_[index##idx##_] \
: cached##Tgt##Address_[index##idx##_]; \
} \
const scalarListList& c##Src##Weights##idx() const \
{ \
checkIndex(index##idx##_); \
return toSource_ ? \
cached##Src##Weights_[index##idx##_] \
: cached##Tgt##Weights_[index##idx##_]; \
} \
const scalarField& c##Src##WeightsSum##idx() const \
{ \
checkIndex(index##idx##_); \
return toSource_ ? \
cached##Src##WeightsSum_[index##idx##_] \
: cached##Tgt##WeightsSum_[index##idx##_]; \
} \
const autoPtr<mapDistribute>& c##Src##MapPtr##idx() const \
{ \
checkIndex(index##idx##_); \
return toSource_ ? \
cached##Src##MapPtr_[index##idx##_] \
: cached##Tgt##MapPtr_[index##idx##_]; \
}
defineMethods01(Src, Tgt, 0)
defineMethods01(Src, Tgt, 1)
defineMethods01(Tgt, Src, 0)
defineMethods01(Tgt, Src, 1)
// Helper functions to retrieve cached values at supplied index
// Note: uses interpolation direction toSource_
#undef defineMethodsIndex
#define defineMethodsIndex(Src, Tgt) \
const labelListList& c##Src##Address(const label index) const \
{ \
checkIndex(index); \
return toSource_ ? \
cached##Src##Address_[index] \
: cached##Tgt##Address_[index]; \
} \
const scalarListList& c##Src##Weights(const label index) const \
{ \
checkIndex(index); \
return toSource_ ? \
cached##Src##Weights_[index] \
: cached##Tgt##Weights_[index]; \
} \
const scalarField& c##Src##WeightsSum(const label index) const \
{ \
checkIndex(index); \
return toSource_ ? \
cached##Src##WeightsSum_[index] \
: cached##Tgt##WeightsSum_[index]; \
} \
const autoPtr<mapDistribute>& c##Src##MapPtr(const label index) \
const \
{ \
checkIndex(index); \
return toSource_ ? \
cached##Src##MapPtr_[index] \
: cached##Tgt##MapPtr_[index]; \
}
defineMethodsIndex(Src, Tgt)
defineMethodsIndex(Tgt, Src)
// I-O
//- Write AMI as a dictionary
void write(Ostream& os) const;
//- Write AMI raw
bool writeData(Ostream& os) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "AMICacheTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,63 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
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/>.
\*---------------------------------------------------------------------------*/
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type, class EvalFunction>
bool Foam::AMICache::apply(List<Type>& result, const EvalFunction& eval) const
{
if (applyLower())
{
eval(result, index0_);
return true;
}
else if (applyUpper())
{
eval(result, index1_);
return true;
}
else if (applyInterpolate())
{
List<Type> r0(result);
eval(r0, index0_);
List<Type> r1(result);
eval(r1, index1_);
//result = (r1 - r0)*interpWeight_ + r0;
forAll(result, i)
{
result[i] = lerp(r0[i], r1[i], interpWeight_);
}
return true;
}
return false;
}
// ************************************************************************* //

View File

@ -61,11 +61,26 @@ registerOptSwitch
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
void Foam::AMIInterpolation::addToCache(const point& refPt)
{
DebugInfo<< "-- addToCache" << endl;
cache_.addToCache(*this, refPt);
}
bool Foam::AMIInterpolation::restoreCache(const point& refPt)
{
DebugInfo<< "-- restoreCache" << endl;
upToDate_ = cache_.restoreCache(refPt);
return upToDate_;
}
Foam::autoPtr<Foam::indexedOctree<Foam::AMIInterpolation::treeType>>
Foam::AMIInterpolation::createTree
(
const primitivePatch& patch
) const
Foam::AMIInterpolation::createTree(const primitivePatch &patch) const
{
treeBoundBox bb(patch.points(), patch.meshPoints());
bb.inflate(0.01);
@ -700,7 +715,8 @@ Foam::AMIInterpolation::AMIInterpolation
tgtWeightsSum_(),
tgtCentroids_(),
tgtMapPtr_(nullptr),
upToDate_(false)
upToDate_(false),
cache_(dict)
{}
@ -730,7 +746,8 @@ Foam::AMIInterpolation::AMIInterpolation
tgtCentroids_(),
tgtPatchPts_(),
tgtMapPtr_(nullptr),
upToDate_(false)
upToDate_(false),
cache_()
{}
@ -743,7 +760,7 @@ Foam::AMIInterpolation::AMIInterpolation
:
requireMatch_(fineAMI.requireMatch_),
reverseTarget_(fineAMI.reverseTarget_),
lowWeightCorrection_(-1.0),
lowWeightCorrection_(-1.0), // Deactivated?
singlePatchProc_(fineAMI.singlePatchProc_),
comm_(fineAMI.comm()), // use fineAMI geomComm if present, comm otherwise
geomComm_(),
@ -759,16 +776,17 @@ Foam::AMIInterpolation::AMIInterpolation
tgtWeightsSum_(),
tgtPatchPts_(),
tgtMapPtr_(nullptr),
upToDate_(false)
upToDate_(false),
cache_(fineAMI.cache(), fineAMI, sourceRestrictAddressing, targetRestrictAddressing)
{
label sourceCoarseSize =
const label sourceCoarseSize =
(
sourceRestrictAddressing.size()
? max(sourceRestrictAddressing)+1
: 0
);
label neighbourCoarseSize =
const label neighbourCoarseSize =
(
targetRestrictAddressing.size()
? max(targetRestrictAddressing)+1
@ -895,14 +913,15 @@ Foam::AMIInterpolation::AMIInterpolation(const AMIInterpolation& ami)
srcWeights_(ami.srcWeights_),
srcWeightsSum_(ami.srcWeightsSum_),
srcCentroids_(ami.srcCentroids_),
srcMapPtr_(nullptr),
srcMapPtr_(ami.srcMapPtr_.clone()),
tgtMagSf_(ami.tgtMagSf_),
tgtAddress_(ami.tgtAddress_),
tgtWeights_(ami.tgtWeights_),
tgtWeightsSum_(ami.tgtWeightsSum_),
tgtCentroids_(ami.tgtCentroids_),
tgtMapPtr_(nullptr),
upToDate_(false)
tgtMapPtr_(ami.tgtMapPtr_.clone()),
upToDate_(ami.upToDate_),
cache_(ami.cache_)
{}
@ -930,7 +949,9 @@ Foam::AMIInterpolation::AMIInterpolation(Istream& is)
//tgtPatchPts_(is),
tgtMapPtr_(nullptr),
upToDate_(readBool(is))
upToDate_(readBool(is)),
cache_(is)
{
// Hopefully no need to stream geomComm_ since only used in processor
// agglomeration?
@ -1361,7 +1382,6 @@ const
}
else if (ray.distance() < nearest.distance())
{
nearest = ray;
nearestFacei = srcFacei;
}
@ -1539,6 +1559,8 @@ void Foam::AMIInterpolation::write(Ostream& os) const
{
os.writeEntry("lowWeightCorrection", lowWeightCorrection_);
}
cache_.write(os);
}
@ -1564,6 +1586,8 @@ bool Foam::AMIInterpolation::writeData(Ostream& os) const
<< token::SPACE<< upToDate();
cache_.writeData(os);
if (distributed() && comm() != -1)
{
os << token::SPACE<< srcMap()

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2016-2024 OpenCFD Ltd.
Copyright (C) 2016-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -61,7 +61,7 @@ SourceFiles
#include "indexedOctree.H"
#include "treeDataPrimitivePatch.H"
#include "runTimeSelectionTables.H"
#include "AMICache.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -78,6 +78,7 @@ public:
// Public Data Types
//- Flag to store face-to-face intersection triangles; default = false
static bool cacheIntersections_;
//- Control use of local communicator for AMI communication
@ -86,6 +87,10 @@ public:
// localComm : 2 : like 1 but always include master (for messages)
static int useLocalComm_;
//- Tolerance used when caching the AMI to identify e.g. if the
//- current rotation angle has already been captured
static scalar cacheThetaTolerance_;
protected:
@ -145,7 +150,6 @@ protected:
autoPtr<mapDistribute> srcMapPtr_;
// Target patch
//- Target face areas
@ -173,9 +177,13 @@ protected:
//- Target map pointer - parallel running only
autoPtr<mapDistribute> tgtMapPtr_;
//- Up-to-date flag
bool upToDate_;
//- Cache
AMICache cache_;
// Protected Member Functions
@ -212,10 +220,10 @@ protected:
// Access
//- Return the orginal src patch with optionally updated points
//- Return the original src patch with optionally updated points
inline const primitivePatch& srcPatch0() const;
//- Return the orginal tgt patch with optionally updated points
//- Return the original tgt patch with optionally updated points
inline const primitivePatch& tgtPatch0() const;
@ -240,27 +248,6 @@ protected:
);
// Constructor helpers
static void agglomerate
(
const autoPtr<mapDistribute>& targetMap,
const scalarList& fineSrcMagSf,
const labelListList& fineSrcAddress,
const scalarListList& fineSrcWeights,
const labelList& sourceRestrictAddressing,
const labelList& targetRestrictAddressing,
scalarList& srcMagSf,
labelListList& srcAddress,
scalarListList& srcWeights,
scalarField& srcWeightsSum,
autoPtr<mapDistribute>& tgtMap,
const label comm
);
public:
//- Runtime type information
@ -366,6 +353,26 @@ public:
// Member Functions
// Constructor helpers
static void agglomerate
(
const autoPtr<mapDistribute>& targetMap,
const scalarList& fineSrcMagSf,
const labelListList& fineSrcAddress,
const scalarListList& fineSrcWeights,
const labelList& sourceRestrictAddressing,
const labelList& targetRestrictAddressing,
scalarList& srcMagSf,
labelListList& srcAddress,
scalarListList& srcWeights,
scalarField& srcWeightsSum,
autoPtr<mapDistribute>& tgtMap,
const label comm
);
// Access
//- Is up-to-date?
@ -379,6 +386,10 @@ public:
return old;
}
//- Return a reference to the AMI cache
const AMICache& cache() const noexcept { return cache_; }
AMICache& cache() noexcept { return cache_; }
//- Distributed across processors (singlePatchProc == -1)
inline bool distributed() const noexcept;
@ -411,6 +422,8 @@ public:
// \returns old value
inline label comm(label communicator) noexcept;
//- Return true if caching is active
inline bool cacheActive() const;
// Source patch
@ -531,7 +544,8 @@ public:
// Low-level
//- Weighted sum of contributions
//- Weighted sum of contributions. Note: cop operates on
//- single Type only.
template<class Type, class CombineOp>
static void weightedSum
(
@ -545,20 +559,35 @@ public:
const UList<Type>& defaultValues
);
//- Weighted sum of contributions
//- Weighted sum of contributions (includes caching+interp)
//- (for primitive types only)
template<class Type>
void weightedSum
(
const bool interpolateToSource,
const bool toSource,
const UList<Type>& fld,
List<Type>& result,
const UList<Type>& defaultValues
) const;
//- Interpolate from target to source with supplied op
//- Weighted sum of (potentially distributed) contributions
//- and apply caching+interpolation. Note: iop
//- operates on single Type only
template<class Type, class CombineOp, class InterpolateOp>
void interpolate
(
const bool toSource,
const UList<Type>& fld,
const CombineOp& cop,
const InterpolateOp& iop,
List<Type>& result,
const UList<Type>& defaultValues
) const;
//- Interpolate from source to target with supplied op
//- to combine existing value with remote value and weight
template<class Type, class CombineOp>
void interpolateToSource
void interpolateToTarget
(
const UList<Type>& fld,
const CombineOp& cop,
@ -566,10 +595,10 @@ public:
const UList<Type>& defaultValues = UList<Type>::null()
) const;
//- Interpolate from source to target with supplied op
//- Interpolate from target to source with supplied op
//- to combine existing value with remote value and weight
template<class Type, class CombineOp>
void interpolateToTarget
void interpolateToSource
(
const UList<Type>& fld,
const CombineOp& cop,
@ -671,6 +700,12 @@ public:
)
const;
//- Restore AMI weights and addressing from the cache
bool restoreCache(const point& refPt);
//- Add AMI weights and addressing to the cache
void addToCache(const point& refPt);
// Checks

View File

@ -123,6 +123,12 @@ inline Foam::label Foam::AMIInterpolation::comm(label communicator) noexcept
}
inline bool Foam::AMIInterpolation::cacheActive() const
{
return cache_.active();
}
inline const Foam::List<Foam::scalar>& Foam::AMIInterpolation::srcMagSf() const
{
return srcMagSf_;

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2015-2023 OpenCFD Ltd.
Copyright (C) 2015-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -83,18 +83,19 @@ void Foam::AMIInterpolation::weightedSum
template<class Type>
void Foam::AMIInterpolation::weightedSum
(
const bool interpolateToSource,
const bool toSource,
const UList<Type>& fld,
List<Type>& result,
const UList<Type>& defaultValues
) const
{
// Note: using non-caching AMI
weightedSum
(
lowWeightCorrection_,
(interpolateToSource ? srcAddress_ : tgtAddress_),
(interpolateToSource ? srcWeights_ : tgtWeights_),
(interpolateToSource ? srcWeightsSum_ : tgtWeightsSum_),
(toSource ? srcAddress_ : tgtAddress_),
(toSource ? srcWeights_ : tgtWeights_),
(toSource ? srcWeightsSum_ : tgtWeightsSum_),
fld,
multiplyWeightedOp<Type, plusEqOp<Type>>(plusEqOp<Type>()),
result,
@ -103,6 +104,223 @@ void Foam::AMIInterpolation::weightedSum
}
template<class Type, class CombineOp, class InterpolateOp>
void Foam::AMIInterpolation::interpolate
(
const bool toSource,
const UList<Type>& fld,
const CombineOp& cop,
const InterpolateOp& iop,
List<Type>& result,
const UList<Type>& defaultValues
) const
{
// Note
// - behaves as old AMIInterpolation::interpolateToSource if toSource=true
// - `result` should be preallocated to correct size and initialized to
// an appropriate value (e.g. Zero)
// Get data locally and do a weighted sum
addProfiling(ami, "AMIInterpolation::interpolate");
cache_.setDirection(toSource);
auto checkSizes = [&](
const UList<Type>& fld,
const labelListList& srcAddr,
const labelListList& tgtAddr,
const UList<Type>& defVals
)
{
if (fld.size() != tgtAddr.size())
{
FatalErrorInFunction
<< "Supplied field size is not equal to "
<< (toSource ? "target" : "source") << " patch size" << nl
<< " source patch = " << srcAddr.size() << nl
<< " target patch = " << tgtAddr.size() << nl
<< " supplied field = " << fld.size()
<< abort(FatalError);
}
else if
(
(lowWeightCorrection_ > 0) && (defVals.size() != srcAddr.size())
)
{
FatalErrorInFunction
<< "Employing default values when sum of weights falls below "
<< lowWeightCorrection_
<< " but number of default values is not equal to "
<< (toSource ? "source" : "target") << " patch size" << nl
<< " default values = " << defVals.size() << nl
<< " source patch = " << srcAddr.size() << nl
<< abort(FatalError);
}
};
// Work space for if distributed
List<Type> work;
List<Type> result0;
if (cache_.index0() != -1)
{
result0 = result;
const auto& srcAddress = cache_.cSrcAddress0();
const auto& srcWeights = cache_.cSrcWeights0();
const auto& srcWeightsSum = cache_.cSrcWeightsSum0();
const auto& tgtAddress = cache_.cTgtAddress0();
checkSizes(fld, srcAddress, tgtAddress, defaultValues);
if (distributed() && cache_.cTgtMapPtr0())
{
const mapDistribute& map = cache_.cTgtMapPtr0()();
if (map.comm() == -1)
{
return;
}
work.resize_nocopy(map.constructSize());
SubList<Type>(work, fld.size()) = fld; // deep copy
map.distribute(work);
}
// if constexpr (is_contiguous_scalar<Type>::value)
// {
// result0 = Zero;
// }
weightedSum
(
lowWeightCorrection_,
srcAddress,
srcWeights,
srcWeightsSum,
(distributed() ? work : fld),
cop,
result0,
defaultValues
);
}
List<Type> result1;
if (cache_.index1() != -1)
{
result1 = result;
const auto& srcAddress = cache_.cSrcAddress1();
const auto& srcWeights = cache_.cSrcWeights1();
const auto& srcWeightsSum = cache_.cSrcWeightsSum1();
const auto& tgtAddress = cache_.cTgtAddress1();
checkSizes(fld, srcAddress, tgtAddress, defaultValues);
if (distributed() && cache_.cTgtMapPtr1())
{
const mapDistribute& map = cache_.cTgtMapPtr1()();
if (map.comm() == -1)
{
return;
}
work.resize_nocopy(map.constructSize());
SubList<Type>(work, fld.size()) = fld; // deep copy
map.distribute(work);
}
// if constexpr (is_contiguous_scalar<Type>::value)
// {
// result1 = Zero;
// }
weightedSum
(
lowWeightCorrection_,
srcAddress,
srcWeights,
srcWeightsSum,
(distributed() ? work : fld),
cop,
result1,
defaultValues
);
}
if (cache_.applyLower())
{
result = result0;
}
else if (cache_.applyUpper())
{
result = result1;
}
else if (cache_.applyInterpolate())
{
forAll(result, i)
{
iop(result[i], i, i, result0[i], i, result1[i], cache_.weight());
}
}
else
{
// No cache - evaluate the AMI
const auto& srcAddress = (toSource ? srcAddress_ : tgtAddress_);
const auto& srcWeights = (toSource ? srcWeights_ : tgtWeights_);
const auto& srcWeightsSum =
(toSource ? srcWeightsSum_ : tgtWeightsSum_);
const auto& tgtAddress = (toSource ? tgtAddress_ : srcAddress_);
checkSizes(fld, srcAddress, tgtAddress, defaultValues);
if (distributed() && tgtMapPtr_)
{
const mapDistribute& map =
(
toSource
? tgtMapPtr_()
: srcMapPtr_()
);
if (map.comm() == -1)
{
return;
}
work.resize_nocopy(map.constructSize());
SubList<Type>(work, fld.size()) = fld; // deep copy
map.distribute(work);
}
// result.resize_nocopy(srcAddress.size());
// if constexpr (is_contiguous_scalar<Type>::value)
// {
// result = Zero;
// }
weightedSum
(
lowWeightCorrection_,
srcAddress,
srcWeights,
srcWeightsSum,
(distributed() ? work : fld),
cop,
result,
defaultValues
);
}
}
// Leave API intact below!
template<class Type, class CombineOp>
void Foam::AMIInterpolation::interpolateToTarget
(
@ -112,58 +330,31 @@ void Foam::AMIInterpolation::interpolateToTarget
const UList<Type>& defaultValues
) const
{
// In-place interpolation
addProfiling(ami, "AMIInterpolation::interpolateToTarget");
if (fld.size() != srcAddress_.size())
{
FatalErrorInFunction
<< "Supplied field size is not equal to source patch size" << nl
<< " source patch = " << srcAddress_.size() << nl
<< " target patch = " << tgtAddress_.size() << nl
<< " supplied field = " << fld.size()
<< abort(FatalError);
}
else if
// Wrap lerp operator to operate inplace
auto iop = [&]
(
(lowWeightCorrection_ > 0)
&& (defaultValues.size() != tgtAddress_.size())
Type& res,
const label i,
const label ia,
const Type& a,
const label ib,
const Type& b,
const scalar w
)
{
FatalErrorInFunction
<< "Employing default values when sum of weights falls below "
<< lowWeightCorrection_
<< " but supplied default field size is not equal to target "
<< "patch size" << nl
<< " default values = " << defaultValues.size() << nl
<< " target patch = " << tgtAddress_.size() << nl
<< abort(FatalError);
}
res = lerp(a, b, w);
};
result.setSize(tgtAddress_.size());
List<Type> work;
if (distributed() && srcMapPtr_)
{
const mapDistribute& map = srcMapPtr_();
if (map.comm() == -1)
{
return;
}
work.resize_nocopy(map.constructSize());
SubList<Type>(work, fld.size()) = fld; // deep copy
map.distribute(work);
}
weightedSum
interpolate
(
lowWeightCorrection_,
tgtAddress_,
tgtWeights_,
tgtWeightsSum_,
(distributed() ? work : fld),
false, // interpolate to target
fld,
cop,
iop,
result,
defaultValues
);
@ -179,58 +370,32 @@ void Foam::AMIInterpolation::interpolateToSource
const UList<Type>& defaultValues
) const
{
// In-place interpolation
addProfiling(ami, "AMIInterpolation::interpolateToSource");
if (fld.size() != tgtAddress_.size())
{
FatalErrorInFunction
<< "Supplied field size is not equal to target patch size" << nl
<< " source patch = " << srcAddress_.size() << nl
<< " target patch = " << tgtAddress_.size() << nl
<< " supplied field = " << fld.size()
<< abort(FatalError);
}
else if
// Wrap lerp operator to operate inplace
auto iop = [&]
(
(lowWeightCorrection_ > 0)
&& (defaultValues.size() != srcAddress_.size())
Type& res,
const label i,
const label ia,
const Type& a,
const label ib,
const Type& b,
const scalar w
)
{
FatalErrorInFunction
<< "Employing default values when sum of weights falls below "
<< lowWeightCorrection_
<< " but number of default values is not equal to source "
<< "patch size" << nl
<< " default values = " << defaultValues.size() << nl
<< " source patch = " << srcAddress_.size() << nl
<< abort(FatalError);
}
res = lerp(a, b, w);
};
result.setSize(srcAddress_.size());
List<Type> work;
if (distributed() && tgtMapPtr_)
{
const mapDistribute& map = tgtMapPtr_();
if (map.comm() == -1)
{
return;
}
work.resize_nocopy(map.constructSize());
SubList<Type>(work, fld.size()) = fld; // deep copy
map.distribute(work);
}
weightedSum
interpolate
(
lowWeightCorrection_,
srcAddress_,
srcWeights_,
srcWeightsSum_,
(distributed() ? work : fld),
true, // toSource,
fld,
cop,
iop,
result,
defaultValues
);

View File

@ -67,9 +67,7 @@ Foam::cyclicACMIGAMGInterfaceField::cyclicACMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, fineInterface),
cyclicACMIInterface_(refCast<const cyclicACMIGAMGInterface>(GAMGCp)),
doTransform_(false),
rank_(0),
sendRequests_(),
recvRequests_()
rank_(0)
{
const cyclicAMILduInterfaceField& p =
refCast<const cyclicAMILduInterfaceField>(fineInterface);
@ -89,9 +87,7 @@ Foam::cyclicACMIGAMGInterfaceField::cyclicACMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, doTransform, rank),
cyclicACMIInterface_(refCast<const cyclicACMIGAMGInterface>(GAMGCp)),
doTransform_(doTransform),
rank_(rank),
sendRequests_(),
recvRequests_()
rank_(rank)
{}
@ -104,9 +100,7 @@ Foam::cyclicACMIGAMGInterfaceField::cyclicACMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, is),
cyclicACMIInterface_(refCast<const cyclicACMIGAMGInterface>(GAMGCp)),
doTransform_(readBool(is)),
rank_(readLabel(is)),
sendRequests_(),
recvRequests_()
rank_(readLabel(is))
{}
@ -120,9 +114,7 @@ Foam::cyclicACMIGAMGInterfaceField::cyclicACMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, local),
cyclicACMIInterface_(refCast<const cyclicACMIGAMGInterface>(GAMGCp)),
doTransform_(false),
rank_(0),
sendRequests_(),
recvRequests_()
rank_(0)
{
const auto& p = refCast<const cyclicACMILduInterfaceField>(local);
@ -142,9 +134,15 @@ bool Foam::cyclicACMIGAMGInterfaceField::ready() const
recvRequests_.start(),
recvRequests_.size()
)
&& UPstream::finishedRequests
(
recvRequests1_.start(),
recvRequests1_.size()
)
)
{
recvRequests_.clear();
recvRequests1_.clear();
if
(
@ -153,9 +151,15 @@ bool Foam::cyclicACMIGAMGInterfaceField::ready() const
sendRequests_.start(),
sendRequests_.size()
)
&& UPstream::finishedRequests
(
sendRequests1_.start(),
sendRequests1_.size()
)
)
{
sendRequests_.clear();
sendRequests1_.clear();
}
return true;
@ -210,15 +214,9 @@ void Foam::cyclicACMIGAMGInterfaceField::initInterfaceMatrixUpdate
// Transform according to the transformation tensors
transformCoupleField(pnf, cmpt);
const auto& map =
(
cyclicACMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -226,22 +224,63 @@ void Foam::cyclicACMIGAMGInterfaceField::initInterfaceMatrixUpdate
<< abort(FatalError);
}
// Assume that sends are also OK
sendRequests_.clear();
const auto& cache = AMI.cache();
// Insert send/receive requests (non-blocking). See e.g.
// cyclicAMIPolyPatchTemplates.C
const label oldWarnComm = UPstream::commWarn(AMI.comm());
map.send
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_,
19462+cyclicACMIInterface_.index() // unique offset + patch index
);
UPstream::commWarn(oldWarnComm);
if (cache.index0() == -1 && cache.index1() == -1)
{
const auto& map =
(
cyclicACMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
// Insert send/receive requests (non-blocking). See e.g.
// cyclicAMIPolyPatchTemplates.C
const label oldWarnComm = UPstream::commWarn(AMI.comm());
map.send
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_,
19462+cyclicACMIInterface_.index() // unique offset + patch index
);
UPstream::commWarn(oldWarnComm);
}
else
{
cache.setDirection(cyclicACMIInterface_.owner());
if (cache.index0() != -1)
{
const auto& map0 = cache.cTgtMapPtr0()();
map0.send
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_,
19462+cyclicACMIInterface_.index() // unique offset + patch index
);
}
if (cache.index1() != -1)
{
const auto& map1 = cache.cTgtMapPtr1()();
map1.send
(
pnf,
sendRequests1_,
scalarSendBufs1_,
recvRequests1_,
scalarRecvBufs1_,
19463+cyclicACMIInterface_.index() // unique offset + patch index
);
}
}
}
this->updatedMatrix(false);
@ -260,6 +299,8 @@ void Foam::cyclicACMIGAMGInterfaceField::updateInterfaceMatrix
const Pstream::commsTypes
) const
{
typedef multiplyWeightedOp<scalar, plusEqOp<scalar>> opType;
const labelUList& faceCells = lduAddr.patchAddr(patchId);
const auto& AMI =
@ -269,48 +310,122 @@ void Foam::cyclicACMIGAMGInterfaceField::updateInterfaceMatrix
: cyclicACMIInterface_.neighbPatch().AMI()
);
DebugPout<< "cyclicACMIGAMGInterfaceField::updateInterfaceMatrix() :"
<< " interface:" << cyclicACMIInterface_.index()
<< " size:" << cyclicACMIInterface_.size()
<< " owner:" << cyclicACMIInterface_.owner()
<< " AMI distributed:" << AMI.distributed()
<< endl;
const auto& cache = AMI.cache();
if (AMI.distributed() && AMI.comm() != -1)
{
const auto& map =
(
cyclicACMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
if (cache.index0() == -1 && cache.index1() == -1)
{
const auto& map =
(
cyclicACMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
solveScalarField work;
map.receive
(
recvRequests_,
scalarRecvBufs_,
work,
19462+cyclicACMIInterface_.index() // unique offset + patch index
);
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
solveScalarField work;
map.receive
(
recvRequests_,
scalarRecvBufs_,
work,
19462+cyclicACMIInterface_.index() // unique offset + patch index
);
// Receive requests all handled by last function call
recvRequests_.clear();
// Receive requests all handled by last function call
recvRequests_.clear();
solveScalarField pnf(faceCells.size(), Zero);
AMI.weightedSum
(
cyclicACMIInterface_.owner(),
work,
pnf, // result
solveScalarField::null()
);
solveScalarField pnf(faceCells.size(), Zero);
AMI.weightedSum
(
cyclicACMIInterface_.owner(),
work,
pnf, // result
solveScalarField::null()
);
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
else
{
cache.setDirection(cyclicACMIInterface_.owner());
solveScalarField pnf(faceCells.size());
solveScalarField work(faceCells.size());
if (cache.index0() != -1)
{
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
work = Zero;
cache.cTgtMapPtr0()().receive
(
recvRequests_,
scalarRecvBufs_,
work,
19462+cyclicACMIInterface_.index() // unique offset + patch index
);
// Receive requests all handled by last function call
recvRequests_.clear();
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress0(),
cache.cSrcWeights0(),
cache.cSrcWeightsSum0(),
work,
opType(plusEqOp<scalar>()),
pnf,
solveScalarField::null()
);
pnf *= (1-cache.weight());
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
if (cache.index1() != -1)
{
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
work = Zero;
cache.cTgtMapPtr1()().receive
(
recvRequests1_,
scalarRecvBufs1_,
work,
19463+cyclicACMIInterface_.index() // unique offset + patch index
);
// Receive requests all handled by last function call
recvRequests1_.clear();
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress1(),
cache.cSrcWeights1(),
cache.cSrcWeightsSum1(),
work,
opType(plusEqOp<scalar>()),
pnf,
solveScalarField::null()
);
pnf *= cache.weight();
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
}
}
else
{
@ -318,23 +433,73 @@ void Foam::cyclicACMIGAMGInterfaceField::updateInterfaceMatrix
const labelUList& nbrFaceCells =
lduAddr.patchAddr(cyclicACMIInterface_.neighbPatchID());
solveScalarField pnf(psiInternal, nbrFaceCells);
solveScalarField work(psiInternal, nbrFaceCells);
// Transform according to the transformation tensors
transformCoupleField(pnf, cmpt);
transformCoupleField(work, cmpt);
if (cyclicACMIInterface_.owner())
if (cache.index0() == -1 && cache.index1() == -1)
{
pnf = AMI.interpolateToSource(pnf);
solveScalarField pnf(faceCells.size(), Zero);
AMI.weightedSum
(
cyclicACMIInterface_.owner(),
work,
pnf, // result
solveScalarField::null()
);
const labelUList& faceCells = lduAddr.patchAddr(patchId);
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
else
{
pnf = AMI.interpolateToTarget(pnf);
cache.setDirection(cyclicACMIInterface_.owner());
solveScalarField pnf(faceCells.size());
if (cache.index0() != -1)
{
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress0(),
cache.cSrcWeights0(),
cache.cSrcWeightsSum0(),
work,
opType(plusEqOp<scalar>()),
pnf,
solveScalarField::null()
);
pnf *= (1 - cache.weight());
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
if (cache.index1() != -1)
{
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress1(),
cache.cSrcWeights1(),
cache.cSrcWeightsSum1(),
work,
opType(plusEqOp<scalar>()),
pnf,
solveScalarField::null()
);
pnf *= cache.weight();
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
}
const labelUList& faceCells = lduAddr.patchAddr(patchId);
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
this->updatedMatrix(true);

View File

@ -83,6 +83,20 @@ class cyclicACMIGAMGInterfaceField
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs_;
// Only used for AMI caching
//- Current range of send requests (non-blocking)
mutable labelRange sendRequests1_;
//- Current range of recv requests (non-blocking)
mutable labelRange recvRequests1_;
//- Scalar send buffers
mutable PtrList<List<solveScalar>> scalarSendBufs1_;
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs1_;
// Private Member Functions

View File

@ -67,9 +67,7 @@ Foam::cyclicAMIGAMGInterfaceField::cyclicAMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, fineInterface),
cyclicAMIInterface_(refCast<const cyclicAMIGAMGInterface>(GAMGCp)),
doTransform_(false),
rank_(0),
sendRequests_(),
recvRequests_()
rank_(0)
{
const cyclicAMILduInterfaceField& p =
refCast<const cyclicAMILduInterfaceField>(fineInterface);
@ -89,9 +87,7 @@ Foam::cyclicAMIGAMGInterfaceField::cyclicAMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, doTransform, rank),
cyclicAMIInterface_(refCast<const cyclicAMIGAMGInterface>(GAMGCp)),
doTransform_(doTransform),
rank_(rank),
sendRequests_(),
recvRequests_()
rank_(rank)
{}
@ -104,9 +100,7 @@ Foam::cyclicAMIGAMGInterfaceField::cyclicAMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, is),
cyclicAMIInterface_(refCast<const cyclicAMIGAMGInterface>(GAMGCp)),
doTransform_(readBool(is)),
rank_(readLabel(is)),
sendRequests_(),
recvRequests_()
rank_(readLabel(is))
{}
@ -120,9 +114,7 @@ Foam::cyclicAMIGAMGInterfaceField::cyclicAMIGAMGInterfaceField
GAMGInterfaceField(GAMGCp, local),
cyclicAMIInterface_(refCast<const cyclicAMIGAMGInterface>(GAMGCp)),
doTransform_(false),
rank_(0),
sendRequests_(), // assume no requests in flight for input field
recvRequests_()
rank_(0)
{
const auto& p = refCast<const cyclicAMILduInterfaceField>(local);
@ -142,9 +134,15 @@ bool Foam::cyclicAMIGAMGInterfaceField::ready() const
recvRequests_.start(),
recvRequests_.size()
)
&& UPstream::finishedRequests
(
recvRequests1_.start(),
recvRequests1_.size()
)
)
{
recvRequests_.clear();
recvRequests1_.clear();
if
(
@ -153,9 +151,15 @@ bool Foam::cyclicAMIGAMGInterfaceField::ready() const
sendRequests_.start(),
sendRequests_.size()
)
&& UPstream::finishedRequests
(
sendRequests1_.start(),
sendRequests1_.size()
)
)
{
sendRequests_.clear();
sendRequests1_.clear();
}
return true;
@ -183,7 +187,6 @@ void Foam::cyclicAMIGAMGInterfaceField::initInterfaceMatrixUpdate
? cyclicAMIInterface_.AMI()
: cyclicAMIInterface_.neighbPatch().AMI()
);
if (AMI.distributed() && AMI.comm() != -1)
{
//DebugPout<< "cyclicAMIFvPatchField::initInterfaceMatrixUpdate() :"
@ -211,15 +214,8 @@ void Foam::cyclicAMIGAMGInterfaceField::initInterfaceMatrixUpdate
// Transform according to the transformation tensors
transformCoupleField(pnf, cmpt);
const auto& map =
(
cyclicAMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
// Assert that all receives are known to have finished
if (!recvRequests_.empty())
if (!recvRequests_.empty() || !recvRequests1_.empty())
{
FatalErrorInFunction
<< "Outstanding recv request(s) on patch "
@ -227,22 +223,63 @@ void Foam::cyclicAMIGAMGInterfaceField::initInterfaceMatrixUpdate
<< abort(FatalError);
}
// Assume that sends are also OK
sendRequests_.clear();
const auto& cache = AMI.cache();
// Insert send/receive requests (non-blocking). See e.g.
// cyclicAMIPolyPatchTemplates.C
const label oldWarnComm = UPstream::commWarn(AMI.comm());
map.send
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_,
19462+cyclicAMIInterface_.index() // unique offset + patch index
);
UPstream::commWarn(oldWarnComm);
if (cache.index0() == -1 && cache.index1() == -1)
{
const auto& map =
(
cyclicAMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
// Insert send/receive requests (non-blocking). See e.g.
// cyclicAMIPolyPatchTemplates.C
const label oldWarnComm = UPstream::commWarn(AMI.comm());
map.send
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_,
19462+cyclicAMIInterface_.index() // unique offset + patch index
);
UPstream::commWarn(oldWarnComm);
}
else
{
cache.setDirection(cyclicAMIInterface_.owner());
if (cache.index0() != -1)
{
const auto& map0 = cache.cTgtMapPtr0()();
map0.send
(
pnf,
sendRequests_,
scalarSendBufs_,
recvRequests_,
scalarRecvBufs_,
19462+cyclicAMIInterface_.index() // unique offset + patch index
);
}
if (cache.index1() != -1)
{
const auto& map1 = cache.cTgtMapPtr1()();
map1.send
(
pnf,
sendRequests1_,
scalarSendBufs1_,
recvRequests1_,
scalarRecvBufs1_,
19463+cyclicAMIInterface_.index() // unique offset + patch index
);
}
}
}
this->updatedMatrix(false);
@ -261,6 +298,8 @@ void Foam::cyclicAMIGAMGInterfaceField::updateInterfaceMatrix
const Pstream::commsTypes commsType
) const
{
typedef multiplyWeightedOp<scalar, plusEqOp<scalar>> opType;
const labelUList& faceCells = lduAddr.patchAddr(patchId);
const auto& AMI =
@ -284,6 +323,12 @@ void Foam::cyclicAMIGAMGInterfaceField::updateInterfaceMatrix
// << " AMI low-weight:" << AMI.applyLowWeightCorrection()
// << endl;
const auto& cache = AMI.cache();
// Assume that sends are also OK
sendRequests_.clear();
sendRequests1_.clear();
if (AMI.distributed() && AMI.comm() != -1)
{
if (commsType != UPstream::commsTypes::nonBlocking)
@ -293,38 +338,118 @@ void Foam::cyclicAMIGAMGInterfaceField::updateInterfaceMatrix
<< exit(FatalError);
}
const auto& map =
(
cyclicAMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
if (cache.index0() == -1 && cache.index1() == -1)
{
const auto& map =
(
cyclicAMIInterface_.owner()
? AMI.tgtMap()
: AMI.srcMap()
);
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
solveScalarField work;
map.receive
(
recvRequests_,
scalarRecvBufs_,
work,
19462+cyclicAMIInterface_.index() // unique offset + patch index
);
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
solveScalarField work;
map.receive
(
recvRequests_,
scalarRecvBufs_,
work,
19462+cyclicAMIInterface_.index() // unique offset + patch index
);
// Receive requests all handled by last function call
recvRequests_.clear();
// Receive requests all handled by last function call
recvRequests_.clear();
solveScalarField pnf(faceCells.size(), Zero);
AMI.weightedSum
(
cyclicAMIInterface_.owner(),
work,
pnf, // result
defaultValues
);
solveScalarField pnf(faceCells.size(), Zero);
AMI.weightedSum
(
cyclicAMIInterface_.owner(),
work,
pnf, // result
defaultValues
);
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
else
{
cache.setDirection(cyclicAMIInterface_.owner());
solveScalarField pnf(faceCells.size());
solveScalarField work(faceCells.size());
if (cache.index0() != -1)
{
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
work = Zero;
cache.cTgtMapPtr0()().receive
(
recvRequests_,
scalarRecvBufs_,
work,
19462+cyclicAMIInterface_.index() // unique offset + patch index
);
// Receive requests all handled by last function call
recvRequests_.clear();
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress0(),
cache.cSrcWeights0(),
cache.cSrcWeightsSum0(),
work,
opType(plusEqOp<scalar>()),
pnf,
defaultValues
);
pnf *= (1-cache.weight());
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
if (cache.index1() != -1)
{
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
work = Zero;
cache.cTgtMapPtr1()().receive
(
recvRequests1_,
scalarRecvBufs1_,
work,
19463+cyclicAMIInterface_.index() // unique offset + patch index
);
// Receive requests all handled by last function call
recvRequests1_.clear();
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress1(),
cache.cSrcWeights1(),
cache.cSrcWeightsSum1(),
work,
opType(plusEqOp<scalar>()),
pnf,
defaultValues
);
pnf *= cache.weight();
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
}
}
else
{
@ -337,17 +462,68 @@ void Foam::cyclicAMIGAMGInterfaceField::updateInterfaceMatrix
// Transform according to the transformation tensors
transformCoupleField(work, cmpt);
solveScalarField pnf(faceCells.size(), Zero);
AMI.weightedSum
(
cyclicAMIInterface_.owner(),
work,
pnf, // result
defaultValues
);
solveScalarField pnf(faceCells.size());
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
if (cache.index0() == -1 && cache.index1() == -1)
{
pnf = Zero;
AMI.weightedSum
(
cyclicAMIInterface_.owner(),
work,
pnf, // result
defaultValues
);
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
else
{
cache.setDirection(cyclicAMIInterface_.owner());
if (cache.index0() != -1)
{
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress0(),
cache.cSrcWeights0(),
cache.cSrcWeightsSum0(),
work,
opType(plusEqOp<scalar>()),
pnf,
defaultValues
);
pnf *= (1 - cache.weight());
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
if (cache.index1() != -1)
{
pnf = Zero;
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress1(),
cache.cSrcWeights1(),
cache.cSrcWeightsSum1(),
work,
opType(plusEqOp<scalar>()),
pnf,
defaultValues
);
pnf *= cache.weight();
// Add result using coefficients
this->addToInternalField(result, !add, faceCells, coeffs, pnf);
}
}
}
this->updatedMatrix(true);

View File

@ -83,6 +83,21 @@ class cyclicAMIGAMGInterfaceField
mutable PtrList<List<solveScalar>> scalarRecvBufs_;
// Only used for AMI caching
//- Current range of send requests (non-blocking)
mutable labelRange sendRequests1_;
//- Current range of recv requests (non-blocking)
mutable labelRange recvRequests1_;
//- Scalar send buffers
mutable PtrList<List<solveScalar>> scalarSendBufs1_;
//- Scalar receive buffers
mutable PtrList<List<solveScalar>> scalarRecvBufs1_;
// Private Member Functions
//- No copy construct

View File

@ -381,6 +381,46 @@ void Foam::cyclicAMIPolyPatch::resetAMI(const UList<point>& points) const
return;
}
point refPt = point::max;
bool restoredFromCache = false;
if (AMIPtr_->cacheActive())
{
const label comm = boundaryMesh().mesh().comm();
if (UPstream::parRun())
{
label refProci = -1;
if (size() > 0)
{
refProci = UPstream::myProcNo(comm);
}
reduce(refProci, maxOp<label>(), UPstream::msgType(), comm);
if (refProci == UPstream::myProcNo(comm))
{
refPt = points[meshPoints()[0]];
}
reduce(refPt, minOp<point>(), UPstream::msgType(), comm);
}
else
{
refPt = points[meshPoints()[0]];
}
// Sets cache indices to use and time interpolation weight
restoredFromCache = AMIPtr_->restoreCache(refPt);
if (returnReduceOr(restoredFromCache, comm))
{
// Restored AMI weight and addressing from cache - all done
Info<< "AMI: weights and addresses restored from cache" << endl;
return;
}
}
const cyclicAMIPolyPatch& nbr = neighbPatch();
const pointField srcPoints(points, meshPoints());
pointField nbrPoints(points, nbr.meshPoints());
@ -432,7 +472,7 @@ void Foam::cyclicAMIPolyPatch::resetAMI(const UList<point>& points) const
meshTools::writeOBJ(osN, nbrPatch0.localFaces(), nbrPoints);
OFstream osO(t.path()/name() + "_ownerPatch.obj");
meshTools::writeOBJ(osO, this->localFaces(), localPoints());
meshTools::writeOBJ(osO, this->localFaces(), srcPoints);
}
// Construct/apply AMI interpolation to determine addressing and weights
@ -446,6 +486,8 @@ void Foam::cyclicAMIPolyPatch::resetAMI(const UList<point>& points) const
{
AMIPtr_->checkSymmetricWeights(true);
}
AMIPtr_->addToCache(refPt);
}

View File

@ -99,7 +99,7 @@ protected:
//- Index of other half
mutable label nbrPatchID_;
//- Particle displacement fraction accross AMI
//- Particle displacement fraction across AMI
const scalar fraction_;
@ -166,6 +166,9 @@ protected:
//- Temporary storage for AMI face centres
mutable vectorField faceCentres0_;
// Cache
bool cacheAMI_;
// Protected Member Functions
@ -520,9 +523,14 @@ public:
(
const Field<Type>& fld,
labelRange& sendRequests,
PtrList<List<Type>>& sendBuffers,
labelRange& recvRequests,
PtrList<List<Type>>& recvBuffers
PtrList<List<Type>>& sendBuffers,
PtrList<List<Type>>& recvBuffers,
labelRange& sendRequests1,
labelRange& recvRequests1,
PtrList<List<Type>>& sendBuffers1,
PtrList<List<Type>>& recvBuffers1
) const;
template<class Type>
@ -530,9 +538,14 @@ public:
(
const Field<Type>& fld,
labelRange& sendRequests,
PtrList<List<Type>>& sendBuffers,
labelRange& recvRequests,
PtrList<List<Type>>& recvBuffers
PtrList<List<Type>>& sendBuffers,
PtrList<List<Type>>& recvBuffers,
labelRange& sendRequests1,
labelRange& recvRequests1,
PtrList<List<Type>>& sendBuffers1,
PtrList<List<Type>>& recvBuffers1
) const;
template<class Type>
@ -541,6 +554,8 @@ public:
const Field<Type>& localFld,
const labelRange& requests, // The receive requests
const PtrList<List<Type>>& recvBuffers,
const labelRange& requests1, // The receive requests
const PtrList<List<Type>>& recvBuffers1,
const UList<Type>& defaultValues
) const;

View File

@ -63,6 +63,7 @@ Foam::tmp<Foam::Field<Type>> Foam::cyclicAMIPolyPatch::interpolate
if constexpr (transform_supported)
{
// Only creates the co-ord system if using periodic AMI
cs.reset(cylindricalCS());
}
@ -130,7 +131,7 @@ Foam::tmp<Foam::Field<Type>> Foam::cyclicAMIPolyPatch::interpolate
const tensorField ownT(cs().R(this->faceCentres()));
Field<Type> localDeflt(defaultValues.size());
if (defaultValues.size() == size())
if (defaultValues.size() != 0 && defaultValues.size() == size())
{
// Transform default values into cylindrical coords (using
// *this faceCentres)
@ -176,27 +177,77 @@ void Foam::cyclicAMIPolyPatch::initInterpolateUntransformed
(
const Field<Type>& fld,
labelRange& sendRequests,
PtrList<List<Type>>& sendBuffers,
labelRange& recvRequests,
PtrList<List<Type>>& recvBuffers
PtrList<List<Type>>& sendBuffers,
PtrList<List<Type>>& recvBuffers,
labelRange& sendRequests1,
labelRange& recvRequests1,
PtrList<List<Type>>& sendBuffers1,
PtrList<List<Type>>& recvBuffers1
) const
{
const auto& AMI = (owner() ? this->AMI() : neighbPatch().AMI());
if (AMI.distributed() && AMI.comm() != -1)
{
const auto& map = (owner() ? AMI.tgtMap() : AMI.srcMap());
const auto& cache = AMI.cache();
// Insert send/receive requests (non-blocking)
map.send
(
fld,
sendRequests,
sendBuffers,
recvRequests,
recvBuffers,
3894+this->index() // unique offset + patch index
);
if (cache.index0() == -1 && cache.index1() == -1)
{
// No caching
const auto& map = (owner() ? AMI.tgtMap() : AMI.srcMap());
// Insert send/receive requests (non-blocking)
map.send
(
fld,
sendRequests,
sendBuffers,
recvRequests,
recvBuffers,
3894+this->index() // unique offset + patch index
);
}
else
{
// Caching is active
cache.setDirection(owner());
if (cache.index0() != -1)
{
const auto& map0 = cache.cTgtMapPtr0()();
// Insert send/receive requests (non-blocking)
map0.send
(
fld,
sendRequests,
sendBuffers,
recvRequests,
recvBuffers,
3894+this->index() // unique offset + patch index
);
}
if (cache.index1() != -1)
{
const auto& map1 = cache.cTgtMapPtr1()();
// Insert send/receive requests (non-blocking)
map1.send
(
fld,
sendRequests1,
sendBuffers1,
recvRequests1,
recvBuffers1,
3895+this->index() // unique offset + patch index
);
}
}
}
}
@ -206,9 +257,14 @@ void Foam::cyclicAMIPolyPatch::initInterpolate
(
const Field<Type>& fld,
labelRange& sendRequests,
PtrList<List<Type>>& sendBuffers,
labelRange& recvRequests,
PtrList<List<Type>>& recvBuffers
PtrList<List<Type>>& sendBuffers,
PtrList<List<Type>>& recvBuffers,
labelRange& sendRequests1,
labelRange& recvRequests1,
PtrList<List<Type>>& sendBuffers1,
PtrList<List<Type>>& recvBuffers1
) const
{
const auto& AMI = (owner() ? this->AMI() : neighbPatch().AMI());
@ -221,46 +277,50 @@ void Foam::cyclicAMIPolyPatch::initInterpolate
// Can rotate fields (vector and non-spherical tensors)
constexpr bool transform_supported = is_rotational_vectorspace_v<Type>;
[[maybe_unused]]
autoPtr<coordSystem::cylindrical> cs;
if constexpr (transform_supported)
{
cs.reset(cylindricalCS());
}
// Only creates the co-ord system if using periodic AMI
// - convert to cylindrical coordinate system
auto cs = cylindricalCS();
if (!cs)
{
initInterpolateUntransformed
(
fld,
sendRequests,
sendBuffers,
recvRequests,
recvBuffers
);
}
else if constexpr (transform_supported)
{
const cyclicAMIPolyPatch& nbrPp = this->neighbPatch();
Field<Type> localFld(fld.size());
// Transform to cylindrical coords
if (cs)
{
Field<Type> localFld(fld.size());
const cyclicAMIPolyPatch& nbrPp = this->neighbPatch();
const tensorField nbrT(cs().R(nbrPp.faceCentres()));
Foam::invTransform(localFld, nbrT, fld);
}
initInterpolateUntransformed
(
localFld,
sendRequests,
sendBuffers,
recvRequests,
recvBuffers
);
initInterpolateUntransformed
(
localFld,
sendRequests,
recvRequests,
sendBuffers,
recvBuffers,
sendRequests1,
recvRequests1,
sendBuffers1,
recvBuffers1
);
return;
}
}
initInterpolateUntransformed
(
fld,
sendRequests,
recvRequests,
sendBuffers,
recvBuffers,
sendRequests1,
recvRequests1,
sendBuffers1,
recvBuffers1
);
}
@ -270,14 +330,21 @@ Foam::tmp<Foam::Field<Type>> Foam::cyclicAMIPolyPatch::interpolate
const Field<Type>& localFld,
const labelRange& requests,
const PtrList<List<Type>>& recvBuffers,
const labelRange& requests1,
const PtrList<List<Type>>& recvBuffers1,
const UList<Type>& defaultValues
) const
{
// Note: cannot be localFld.size() -> is set to null for distributed AMI
auto tresult = tmp<Field<Type>>::New(this->size(), Zero);
const auto& AMI = (owner() ? this->AMI() : neighbPatch().AMI());
const auto& cache = AMI.cache();
cache.setDirection(owner());
Field<Type> work;
Field<Type> work1;
if (AMI.distributed())
{
if (AMI.comm() == -1)
@ -285,71 +352,157 @@ Foam::tmp<Foam::Field<Type>> Foam::cyclicAMIPolyPatch::interpolate
return tresult;
}
const auto& map = (owner() ? AMI.tgtMap() : AMI.srcMap());
if (cache.index0() == -1 && cache.index1() == -1)
{
// No caching
const auto& map = (owner() ? AMI.tgtMap() : AMI.srcMap());
// Receive (= copy) data from buffers into work. TBD: receive directly
// into slices of work.
map.receive
(
requests,
recvBuffers,
work,
3894+this->index() // unique offset + patch index
);
// Receive (= copy) data from buffers into work. TBD: receive
// directly into slices of work.
map.receive
(
requests,
recvBuffers,
work,
3894+this->index() // unique offset + patch index
);
}
else
{
// Using AMI cache
if (cache.index0() != -1)
{
cache.cTgtMapPtr0()().receive
(
requests,
recvBuffers,
work,
3894+this->index() // unique offset + patch index
);
}
if (cache.index1() != -1)
{
cache.cTgtMapPtr1()().receive
(
requests1,
recvBuffers1,
work1,
3895+this->index() // unique offset + patch index
);
}
}
}
const Field<Type>& fld = (AMI.distributed() ? work : localFld);
const Field<Type>& fld1 = (AMI.distributed() ? work1 : localFld);
// Rotate fields (vector and non-spherical tensors)
constexpr bool transform_supported = is_rotational_vectorspace_v<Type>;
[[maybe_unused]]
autoPtr<coordSystem::cylindrical> cs;
// Rotate fields (vector and non-spherical tensors) for periodic AMI
tensorField ownTransform;
Field<Type> localDeflt;
// Rotate fields (vector and non-spherical tensors)
if constexpr (transform_supported)
{
cs.reset(cylindricalCS());
}
// Only creates the co-ord system if using periodic AMI
// - convert to cylindrical coordinate system
auto cs = cylindricalCS();
if (!cs)
{
AMI.weightedSum
(
owner(),
fld,
tresult.ref(),
defaultValues
);
}
else if constexpr (transform_supported)
{
const tensorField ownT(cs().R(this->faceCentres()));
Field<Type> localDeflt(defaultValues.size());
if (defaultValues.size() == size())
if (cs)
{
// Transform default values into cylindrical coords (using
// *this faceCentres)
// We get in UList (why? Copied from cyclicAMI). Convert to
// Field so we can use transformField routines.
const SubField<Type> defaultSubFld(defaultValues);
const Field<Type>& defaultFld(defaultSubFld);
Foam::invTransform(localDeflt, ownT, defaultFld);
}
ownTransform = cs().R(this->faceCentres());
localDeflt = defaultValues;
if (defaultValues.size() == size())
{
// Transform default values into cylindrical coords (using
// *this faceCentres)
// We get in UList (why? Copied from cyclicAMI). Convert to
// Field so we can use transformField routines.
const SubField<Type> defaultSubFld(defaultValues);
const Field<Type>& defaultFld(defaultSubFld);
Foam::invTransform(localDeflt, ownTransform, defaultFld);
}
}
}
const auto& localDefaultValues =
localDeflt.size() ? localDeflt : defaultValues;
if (cache.index0() == -1 && cache.index1() == -1)
{
// No caching
AMI.weightedSum
(
owner(),
fld,
tresult.ref(),
localDeflt
localDefaultValues
);
// Transform back
Foam::transform(tresult.ref(), ownT, tresult());
}
if (ownTransform.size())
{
Foam::transform(tresult.ref(), ownTransform, tresult());
}
return tresult;
return tresult;
}
else
{
if (cache.index0() != -1)
{
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress0(),
cache.cSrcWeights0(),
cache.cSrcWeightsSum0(),
fld,
multiplyWeightedOp<Type, plusEqOp<Type>>(plusEqOp<Type>()),
tresult.ref(),
localDefaultValues
);
if (ownTransform.size())
{
Foam::transform(tresult.ref(), ownTransform, tresult());
}
// Assuming cache weight is zero when index1 is inactive (==-1)
tresult.ref() *= (1 - cache.weight());
}
if (cache.index1() != -1)
{
auto tresult1 = tmp<Field<Type>>::New(this->size(), Zero);
AMIInterpolation::weightedSum
(
AMI.lowWeightCorrection(),
cache.cSrcAddress1(),
cache.cSrcWeights1(),
cache.cSrcWeightsSum1(),
fld1,
multiplyWeightedOp<Type, plusEqOp<Type>>(plusEqOp<Type>()),
tresult1.ref(),
localDefaultValues
);
if (ownTransform.size())
{
Foam::transform(tresult1.ref(), ownTransform, tresult1());
}
tresult1.ref() *= cache.weight();
tresult.ref() += tresult1();
}
return tresult;
}
}

View File

@ -281,6 +281,7 @@ processorLOD/faceBox/faceBox.C
AMI=AMIInterpolation
$(AMI)/AMIInterpolation/AMIInterpolation.C
$(AMI)/AMIInterpolation/AMIInterpolationNew.C
$(AMI)/AMIInterpolation/AMICache.C
$(AMI)/AMIInterpolation/advancingFrontAMI/advancingFrontAMI.C
$(AMI)/AMIInterpolation/advancingFrontAMI/advancingFrontAMIParallelOps.C
$(AMI)/AMIInterpolation/faceAreaWeightAMI/faceAreaWeightAMI.C

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2018-2024 OpenCFD Ltd.
Copyright (C) 2018-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -95,6 +95,102 @@ namespace Foam
}
}
};
template<class Type, class TrackingData>
class interpolate
{
//- Combine operator for AMIInterpolation
FaceCellWave<Type, TrackingData>& solver_;
const cyclicAMIPolyPatch& patch_;
public:
interpolate
(
FaceCellWave<Type, TrackingData>& solver,
const cyclicAMIPolyPatch& patch
)
:
solver_(solver),
patch_(patch)
{}
void operator()
(
Type& x,
const label localFacei,
const label f0i,
const Type& y0,
const label f1i,
const Type& y1,
const scalar weight
) const
{
if (y0.valid(solver_.data()))
{
if (y1.valid(solver_.data()))
{
x.interpolate
(
solver_.mesh(),
patch_.faceCentres()[localFacei],
f0i,
y0,
f1i,
y1,
weight,
solver_.propagationTol(),
solver_.data()
);
}
else
{
// Convert patch face into mesh face
label meshFacei = -1;
if (patch_.owner())
{
meshFacei = patch_.start() + f0i;
}
else
{
meshFacei = patch_.neighbPatch().start() + f0i;
}
x.updateFace
(
solver_.mesh(),
meshFacei,
y0,
solver_.propagationTol(),
solver_.data()
);
}
}
else if (y1.valid(solver_.data()))
{
// Convert patch face into mesh face
label meshFacei = -1;
if (patch_.owner())
{
meshFacei = patch_.start() + f1i;
}
else
{
meshFacei = patch_.neighbPatch().start() + f1i;
}
x.updateFace
(
solver_.mesh(),
meshFacei,
y1,
solver_.propagationTol(),
solver_.data()
);
}
}
};
}
@ -784,18 +880,43 @@ void Foam::FaceCellWave<Type, TrackingData>::handleAMICyclicPatches()
// Transfer sendInfo to cycPatch
combine<Type, TrackingData> cmb(*this, cycPatch);
// Linear interpolation
interpolate<Type, TrackingData> interp(*this, cycPatch);
const auto& AMI =
(
cycPatch.owner()
? cycPatch.AMI()
: cycPatch.neighbPatch().AMI()
);
if (cycPatch.applyLowWeightCorrection())
{
List<Type> defVals
const List<Type> defVals
(
cycPatch.patchInternalList(allCellInfo_)
);
cycPatch.interpolate(sendInfo, cmb, receiveInfo, defVals);
AMI.interpolate
(
cycPatch.owner(),
sendInfo,
cmb,
interp,
receiveInfo,
defVals
);
}
else
{
cycPatch.interpolate(sendInfo, cmb, receiveInfo);
AMI.interpolate
(
cycPatch.owner(),
sendInfo,
cmb,
interp,
receiveInfo,
UList<Type>::null() // no default values needed
);
}
}

View File

@ -197,6 +197,23 @@ public:
template<class TrackingData>
inline bool equal(const cellInfo&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const cellInfo& f0,
const label i1,
const cellInfo& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
Copyright (C) 2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -249,6 +249,31 @@ inline bool Foam::cellInfo::equal
}
template<class TrackingData>
inline bool Foam::cellInfo::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const cellInfo& f0,
const label i1,
const cellInfo& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (!f0.valid(td))
{
return update(f1, -1, -1, -1, -1, td);
}
else
{
return update(f0, -1, -1, -1, -1, td);
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::cellInfo::operator==

View File

@ -76,6 +76,16 @@ class wallPoint
// Private Member Functions
// Update this with w2 if distance less
template<class TrackingData>
inline bool update
(
const wallPoint& w2,
const scalar dist2,
const scalar tol,
TrackingData& td
);
//- Evaluate distance to point.
// Update distSqr, origin from whomever is nearer pt.
// \return true if w2 is closer to point, false otherwise.
@ -206,10 +216,41 @@ public:
TrackingData& td
);
//- Interpolate different values
template<class TrackingData>
wallPoint
(
const polyMesh&,
const scalar weight,
const label face0,
const wallPoint& info0,
const label face1,
const wallPoint& info1,
const scalar tol,
TrackingData& td
);
//- Test for equality, with TrackingData
template<class TrackingData>
inline bool equal(const wallPoint&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const wallPoint& f0,
const label i1,
const wallPoint& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2020 OpenCFD Ltd.
Copyright (C) 2020,2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -35,21 +35,12 @@ License
template<class TrackingData>
inline bool Foam::wallPoint::update
(
const point& pt,
const wallPoint& w2,
const scalar dist2,
const scalar tol,
TrackingData& td
)
{
//Already done in calling algorithm
//if (w2.origin() == origin_)
//{
// // Shortcut. Same input so same distance.
// return false;
//}
const scalar dist2 = magSqr(pt - w2.origin());
if (!valid(td))
{
// current not yet set so use any value
@ -83,6 +74,26 @@ inline bool Foam::wallPoint::update
}
// Update this with w2 if w2 nearer to pt.
template<class TrackingData>
inline bool Foam::wallPoint::update
(
const point& pt,
const wallPoint& w2,
const scalar tol,
TrackingData& td
)
{
return update
(
w2,
magSqr(pt - w2.origin()),
tol,
td
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
inline Foam::wallPoint::wallPoint()
@ -260,6 +271,58 @@ inline bool Foam::wallPoint::equal
}
template<class TrackingData>
inline bool Foam::wallPoint::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const wallPoint& f0,
const label i1,
const wallPoint& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
if (f1.valid(td))
{
// What do we do about the origin? Should be consistent with the
// distance? But then interpolating between a point 'on the left'
// and a point 'on the right' the interpolate might just be on
// top of us so suddenly setting the distance to be zero...
const scalar dist2 =
sqr(lerp(sqrt(f0.distSqr()), sqrt(f1.distSqr()), weight));
if (f0.distSqr() < f1.distSqr())
{
// Update with interpolated distance and f0 origin
return update(f0, dist2, tol, td);
}
else
{
return update(f1, dist2, tol, td);
}
}
else
{
return update(f0, f0.distSqr(), tol, td);
}
}
else if (f1.valid(td))
{
return update(f1, f1.distSqr(), tol, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::wallPoint::operator==

View File

@ -196,6 +196,23 @@ public:
TrackingData& td
) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const topoDistanceData<Type>& f0,
const label i1,
const topoDistanceData<Type>& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2016 OpenFOAM Foundation
Copyright (C) 2019-2020 OpenCFD Ltd.
Copyright (C) 2019-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -199,6 +199,38 @@ inline bool Foam::topoDistanceData<Type>::equal
}
template<class Type>
template<class TrackingData>
inline bool Foam::topoDistanceData<Type>::interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const topoDistanceData<Type>& f0,
const label i1,
const topoDistanceData<Type>& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
// Topological distance - how to interpolate?
if (!valid(td))
{
if (f0.valid(td))
{
this->operator=(f0);
}
else
{
this->operator=(f1);
}
return true;
}
return false;
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
template<class Type>

View File

@ -174,6 +174,23 @@ public:
template<class TrackingData>
inline bool equal(const minData&, TrackingData& td) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const minData& f0,
const label i1,
const minData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2014-2016 OpenFOAM Foundation
Copyright (C) 2015-2020 OpenCFD Ltd.
Copyright (C) 2015-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -173,6 +173,35 @@ inline bool Foam::minData::equal
}
template<class TrackingData>
inline bool Foam::minData::interpolate
(
const polyMesh& mesh,
const point& pt,
const label i0,
const minData& f0,
const label i1,
const minData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
return updateFace(mesh, -1, f0, tol, td);
}
if (f1.valid(td))
{
return updateFace(mesh, -1, f1, tol, td);
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::minData::operator==

View File

@ -200,6 +200,23 @@ public:
TrackingData&
) const;
//- Interpolate between two values (lerp). Returns true if
//- causes changes. Not sure if needs to be specialised between
//- face and cell and what index is needed...
template<class TrackingData>
inline bool interpolate
(
const polyMesh&,
const point& pt,
const label i0,
const meshToMeshData& f0,
const label i1,
const meshToMeshData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
);
// Member Operators

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2017-2020 OpenCFD Ltd.
Copyright (C) 2017-2025 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -205,6 +205,35 @@ inline bool Foam::meshToMeshData::equal
}
template<class TrackingData>
inline bool Foam::meshToMeshData::interpolate
(
const polyMesh& mesh,
const point& pt,
const label i0,
const meshToMeshData& f0,
const label i1,
const meshToMeshData& f1,
const scalar weight,
const scalar tol,
TrackingData& td
)
{
if (f0.valid(td))
{
return updateFace(mesh, -1, f0, tol, td);
}
if (f1.valid(td))
{
return updateFace(mesh, -1, f1, tol, td);;
}
else
{
return false;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline bool Foam::meshToMeshData::operator==

View File

@ -480,7 +480,6 @@ void Foam::faMeshReconstructor::createMesh()
(
new faMesh
(
procMesh_.name(),
*serialVolMesh_,
identity(singlePatchFaces_.size()) // faceLabels
)

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2019-2025 OpenCFD Ltd.
Copyright (C) 2019-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -144,23 +144,11 @@ KirchhoffShell::KirchhoffShell
)
:
vibrationShellModel(modelType, mesh, dict),
h_
(
IOobject
(
dict.getOrDefault<word>("h", suffixed("hs")),
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
regionMesh()
),
ps_
(
IOobject
(
dict.getOrDefault<word>("ps", suffixed("ps")),
"ps_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::READ_IF_PRESENT,
@ -169,11 +157,23 @@ KirchhoffShell::KirchhoffShell
regionMesh(),
dimensionedScalar(dimPressure, Zero)
),
h_
(
IOobject
(
"h_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
regionMesh()
),
laplaceW_
(
IOobject
(
suffixed("laplaceW"),
"laplaceW_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::NO_READ,
@ -186,7 +186,7 @@ KirchhoffShell::KirchhoffShell
(
IOobject
(
suffixed("laplace2W"),
"laplace2W_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::NO_READ,
@ -199,7 +199,7 @@ KirchhoffShell::KirchhoffShell
(
IOobject
(
suffixed("w0"),
"w0_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::NO_READ,
@ -212,7 +212,7 @@ KirchhoffShell::KirchhoffShell
(
IOobject
(
suffixed("w00"),
"w00_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::NO_READ,
@ -225,7 +225,7 @@ KirchhoffShell::KirchhoffShell
(
IOobject
(
suffixed("laplaceW0"),
"laplaceW0_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::NO_READ,
@ -238,7 +238,7 @@ KirchhoffShell::KirchhoffShell
(
IOobject
(
suffixed("laplace2W0"),
"laplace2W0_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::NO_READ,
@ -256,7 +256,6 @@ KirchhoffShell::KirchhoffShell
init(dict);
}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void KirchhoffShell::preEvolveRegion()

View File

@ -24,7 +24,7 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::regionModels::KirchhoffShell (Foam::regionFaModels)
Foam::regionFaModels::KirchhoffShell
Description
Vibration-shell finite-area model.
@ -51,25 +51,9 @@ Usage
\table
Property | Description | Type | Reqd | Deflt
vibrationShellModel | Type name: KirchhoffShell | word | yes | -
f0 | Damping coefficient [1/s] | scalar | yes | -
f1 | Damping coefficient [1/s] | scalar | yes | -
f2 | Damping coefficient [1/s] | scalar | yes | -
h | Name of thickness field | word | no | hs (suffix)
ps | Name of pressure field | word | no | ps (suffix)
\endtable
Fields/variables used:
\table
Property | Description | Type | Deflt
h | Thickness | shell | hs (suffix)
ps | Pressure on shell | shell | ps (suffix)
\endtable
Naming changes from 2056 and earlier:
\table
Keyword | Description | Keyword (old) | Deflt (old)
h | Thickness (shell) | - | "h_" + regionName
ps | Pressure (shell) | - | "ps_" + regionName
f0 | Damping coefficient [1/s] | scalar | yes | -
f1 | Damping coefficient [1/s] | scalar | yes | -
f2 | Damping coefficient [1/s] | scalar | yes | -
\endtable
The inherited entries are elaborated in:
@ -106,12 +90,12 @@ class KirchhoffShell
// Source term fields
//- Shell thickness [m]
areaScalarField h_;
//- External surface source [Pa]
const areaScalarField ps_;
//- Thickness [m]
areaScalarField h_;
//- Laplace of the displacement
areaScalarField laplaceW_;

View File

@ -27,6 +27,7 @@ License
#include "regionFaModel.H"
#include "faMesh.H"
#include "faMeshesRegistry.H"
#include "Time.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@ -73,7 +74,7 @@ Foam::IOobject createModelIOobject
(
objName,
mesh.time().constant(),
faMesh::Registry(mesh),
faMeshesRegistry::New(mesh).thisDb(),
IOobjectOption::NO_READ,
IOobjectOption::NO_WRITE,
IOobjectOption::REGISTER
@ -140,71 +141,12 @@ Foam::IOobject createPropertiesIOobject
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
void Foam::regionModels::regionFaModel::constructMeshObjects
(
// Just for error reference
const dictionary& dict
)
void Foam::regionModels::regionFaModel::constructMeshObjects()
{
regionMeshPtr_.reset(nullptr);
#if 1
regionMeshPtr_.reset
(
new faMesh(areaName_, primaryMesh_)
);
#else
// With try/catch and error messages
// DIY
// regionMeshPtr_ = faMesh::TryNew(areaName_, primaryMesh_);
// More heavy handed, but gives a better chance of locating
// the source of the error.
{
const bool oldThrowingError = FatalError.throwing(true);
const bool oldThrowingIOerr = FatalIOError.throwing(true);
try
{
regionMeshPtr_.reset
(
new faMesh(areaName_, primaryMesh_)
);
}
catch (const Foam::error& err)
{
Warning << err << nl << endl;
// Trickery to get original message
err.write(Warning, false);
}
catch (const Foam::IOerror& err)
{
Warning << err << nl << endl;
// Trickery to get original message
err.write(Warning, false);
}
FatalError.throwing(oldThrowingError);
FatalIOError.throwing(oldThrowingIOerr);
}
if (!regionMeshPtr_)
{
FatalError
<< "Failed to create finite-area mesh [" << areaName_
<< "] for model: "<< modelName_ << nl
<< "A common cause is an incorrect or"
" missing 'area' entry in the setup" << nl
<< ">>>>" << nl
<< dict.relativeName() << dict << "<<<<" << endl
<< exit(FatalError);
}
#endif
}
@ -292,35 +234,7 @@ Foam::regionModels::regionFaModel::regionFaModel
regionName_(dict.get<word>("region")),
coeffs_(dict.subOrEmptyDict(modelName + "Coeffs"))
{
// Suffix hint for variable names
if
(
coeffs_.readIfPresent("suffixing", suffixHint_)
|| dict.readIfPresent("suffixing", suffixHint_)
)
{
Switch sw = Switch::find(suffixHint_);
if (sw.good())
{
if (!sw) // No suffix
{
suffixHint_.clear();
}
}
else if (suffixHint_ == "default")
{
// Use default suffix
sw = true;
}
if (sw) // Default suffix
{
suffixHint_ = '_' + regionName_;
}
}
constructMeshObjects(dict);
constructMeshObjects();
initialise();
if (readFields)
@ -336,14 +250,9 @@ void Foam::regionModels::regionFaModel::evolve()
{
if (active_)
{
Info<< "\nEvolving " << modelName_
<< " for region " << regionMesh().name();
if (!polyMesh::regionName(areaName_).empty())
{
Info<< " [" << areaName_ << "]";
}
Info<< endl;
Info<< "\nEvolving " << modelName_ << " for region "
<< regionMesh().name() << " : "
<< polyMesh::regionName(areaName_) << endl;
preEvolveRegion();
@ -356,7 +265,7 @@ void Foam::regionModels::regionFaModel::evolve()
{
Info<< incrIndent;
info();
Info<< decrIndent << endl;
Info<< endl << decrIndent;
}
}
}

View File

@ -57,20 +57,9 @@ Usage
region | Name of operand region | word | yes | -
area | Name of the finite-area mesh | word | no | region0
active | Flag to activate the model | bool | yes | -
suffixing | Suffix hint for model variables | word | no | -
infoOutput | Flag to activate information output | bool | no | false
\endtable
Note
The \c suffixing parameter is a hint that individual models may use
when automatically generating variable names internally. For example,
a model \em may provide ("T" + suffixing) and ("q" + suffixing)
as its default names temperature and flux fields, respectively.
If the user specifies \c suffixing = "_foo", those \em default names
would then become "T_foo" and "q_foo", respectively.
Suffixing (true|false|none|default|...) currently selects between
no suffix and ('_'+region).
SourceFiles
regionFaModelI.H
regionFaModel.C
@ -101,7 +90,7 @@ class regionFaModel
// Private Member Functions
//- Construct region mesh and fields
void constructMeshObjects(const dictionary&);
void constructMeshObjects();
//- Initialise the region
void initialise();
@ -129,21 +118,18 @@ protected:
//- Model name
const word modelName_;
//- Suffix hint for automatic model variable names (default: "")
word suffixHint_;
//- The finite-area mesh name (default: region0)
//- The finite-area mesh name
word areaName_;
//- Region name
word regionName_;
//- Model coefficients dictionary
dictionary coeffs_;
//- Pointer to the region mesh database
autoPtr<faMesh> regionMeshPtr_;
//- Model coefficients dictionary
dictionary coeffs_;
//- Dictionary of output properties
autoPtr<IOdictionary> outputPropertiesPtr_;
@ -238,15 +224,6 @@ public:
//- Return mapping between surface and volume fields
const volSurfaceMapping& vsm() const;
//- The suffix hint for automatic model variable names
const word& suffixHint() const noexcept { return suffixHint_; }
//- Return the concatenation of \p base and the suffix hint
inline word suffixed(const char* base) const;
//- Return the concatenation of \p base and the suffix hint
inline word suffixed(const std::string& base) const;
// Evolution

View File

@ -27,7 +27,7 @@ License
#include "regionFaModel.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
inline const Foam::faMesh& Foam::regionModels::regionFaModel::regionMesh() const
{
@ -106,18 +106,4 @@ inline bool Foam::regionModels::regionFaModel::isRegionPatch
}
inline Foam::word
Foam::regionModels::regionFaModel::suffixed(const char* base) const
{
return word(base + suffixHint_);
}
inline Foam::word
Foam::regionModels::regionFaModel::suffixed(const std::string& base) const
{
return word(base + suffixHint_);
}
// ************************************************************************* //

View File

@ -123,29 +123,14 @@ thermalShell::thermalShell
)
:
thermalShellModel(modelType, mesh, dict),
hName_(dict.getOrDefault<word>("h", suffixed("hs"))),
qsName_(dict.getOrDefault<word>("qs", suffixed("qs"))),
qrName_(dict.getOrDefault<word>("qr", "none")),
nNonOrthCorr_(1),
// Only need/want thermal solid properties
thermo_(dict.subDict("thermo"), solidProperties::THERMAL),
h_
(
IOobject
(
hName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
regionMesh()
),
qs_
(
IOobject
(
qsName_,
"qs_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::READ_IF_PRESENT,
@ -154,6 +139,19 @@ thermalShell::thermalShell
regionMesh(),
dimensionedScalar(dimPower/dimArea, Zero)
),
h_
(
IOobject
(
"h_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
regionMesh()
),
qrName_(dict.getOrDefault<word>("qr", "none")),
thickness_(dict.getOrDefault<scalar>("thickness", 0))
{
init(dict);

View File

@ -5,7 +5,7 @@
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2019-2025 OpenCFD Ltd.
Copyright (C) 2019-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@ -56,28 +56,11 @@ Usage
where the entries mean:
\table
Property | Description | Type | Reqd | Deflt
thermalShellModel | Type name: thermalShell | word | yes | -
h | Name of thickness field | word | no | hs (suffix)
qs | Name of source field | word | no | qs (suffix)
qr | Name of radiative heat flux field | word | no | none
thermo | Solid thermal properties | dict | yes | -
thickness | Uniform shell thickness [m] | scalar | choice | -
\endtable
Fields/variables used:
\table
Property | Description | Type | Deflt
h | Thickness | shell | hs (suffix)
qs | Source field | shell | qs (suffix)
qr | Radiative heat flux field | volume | none
\endtable
Note the following naming changes from 2056 and earlier:
\table
Keyword | Description | Keyword (old) | Deflt (old)
h | Temperature (volume) | h | "h_" + regionName
qs | Temperature (shell) | qs | "qs_" + regionName
Property | Description | Type | Reqd | Deflt
thermalShellModel | Type name: thermalShell | word | yes | -
thermo | Solid thermal properties | dictionary | yes | -
qr | Name of radiative heat flux field | word | no | none
thickness | Uniform film thickness [m] | scalar | choice | -
\endtable
The inherited entries are elaborated in:
@ -85,6 +68,7 @@ Usage
SourceFiles
thermalShell.C
thermalShellI.H
\*---------------------------------------------------------------------------*/
@ -124,38 +108,31 @@ protected:
// Protected Data
//- Name of shell thickness [height] field (default: "hs" + suffix)
const word hName_;
// Solution parameters
//- Name of surface energy source (default: "qs" + suffix)
const word qsName_;
//- Name of (volume) radiative flux field (default: none)
const word qrName_;
//- Number of non orthogonal correctors
label nNonOrthCorr_;
// Solution Parameters
// Thermo properties
//- Number of non orthogonal correctors
label nNonOrthCorr_;
//- Solid properties
solidProperties thermo_;
// Thermo properties
// Source term fields
//- Solid properties
solidProperties thermo_;
//- External surface energy source [J/m2/s]
areaScalarField qs_;
//- Film thickness [m]
areaScalarField h_;
// Source term fields
//- Name of the primary region radiative flux
const word qrName_;
//- Shell thickness field [m]
areaScalarField h_;
//- External surface energy source [J/m2/s]
areaScalarField qs_;
//- Uniform shell thickness [m]
scalar thickness_;
//- Uniform film thickness [m]
scalar thickness_;
// Protected Member Functions
@ -197,7 +174,7 @@ public:
// Fields
//- Return the shell specific heat capacity [J/kg/K]
//- Return the film specific heat capacity [J/kg/K]
const tmp<areaScalarField> Cp() const;
//- Return density [kg/m3]

View File

@ -50,14 +50,13 @@ thermalShellModel::thermalShellModel
)
:
regionFaModel(mesh, "thermalShell", modelType, dict, true),
TName_(dict.getOrDefault<word>("T", suffixed("Ts"))),
TprimaryName_(dict.getOrDefault<word>("Tprimary", "T")),
Tp_(mesh.lookupObject<volScalarField>(TprimaryName_)),
TName_(dict.get<word>("T")),
Tp_(mesh.lookupObject<volScalarField>(TName_)),
T_
(
IOobject
(
TName_,
"Ts_" + regionName_,
regionMesh().time().timeName(),
regionMesh().thisDb(),
IOobject::MUST_READ,
@ -72,8 +71,8 @@ thermalShellModel::thermalShellModel
{
if (faOptions_.optionList::empty())
{
Info<< "No finite-area options present for area:"
<< regionFaModel::areaName() << endl;
Info<< "No finite area options present for area : "
<< polyMesh::regionName(regionFaModel::areaName()) << endl;
}
}

View File

@ -24,7 +24,7 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::regionModels::thermalShellModel
Foam::regionModels::thermalShellModels::thermalShellModel
Description
Intermediate class for thermal-shell finite-area models.
@ -34,12 +34,12 @@ Usage
\verbatim
<patchName>
{
// Mandatory entries
T <word>;
// Optional entries
thermalShellModel <word>;
T <word>;
Tprimary <word>;
// Inherited entries
...
}
@ -47,27 +47,12 @@ Usage
where the entries mean:
\table
Property | Description | Type | Reqd | Deflt
T | Name of operand temperature field | word | no | Ts (suffix)
Tprimary | Name of primary temperature field | word | no | T
thermalShellModel | Name of thermalShellModel thermal-shell model <!--
Property | Description | Type | Reqd | Deflt
T | Name of operand temperature field | word | yes | -
thermalShellModel | Name of thermal-shell model <!--
--> | word | choice | -
\endtable
Fields/variables used:
\table
Property | Description | Type | Deflt
T | Temperature | shell | Ts (suffix)
Tprimary | Temperature | volume | T
\endtable
\b BREAKING Naming changes from 2056 and earlier:
\table
Keyword | Description | Keyword (old) | Deflt (old)
T | Temperature (shell) | - | "Ts_" + regionName
Tprimary | Temperature (volume) | T | -
\endtable
The inherited entries are elaborated in:
- \link regionFaModel.H \endlink
@ -104,11 +89,8 @@ protected:
// Protected Data
//- Name of shell temperature field (default: "Ts" + suffix)
const word TName_;
//- Name of volume temperature field (default: "T")
const word TprimaryName_;
//- Name of the temperature field
word TName_;
//- Primary (volume) region temperature
const volScalarField& Tp_;

Some files were not shown because too many files have changed in this diff Show More