/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011-2017 OpenFOAM Foundation
Copyright (C) 2020-2023 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 .
\*---------------------------------------------------------------------------*/
#include "processorCyclicPointPatchField.H"
#include "transformField.H"
#include "processorPolyPatch.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template
Foam::processorCyclicPointPatchField::processorCyclicPointPatchField
(
const pointPatch& p,
const DimensionedField& iF
)
:
coupledPointPatchField(p, iF),
procPatch_(refCast(p))
{}
template
Foam::processorCyclicPointPatchField::processorCyclicPointPatchField
(
const pointPatch& p,
const DimensionedField& iF,
const dictionary& dict
)
:
coupledPointPatchField(p, iF, dict),
procPatch_(refCast(p, dict))
{}
template
Foam::processorCyclicPointPatchField::processorCyclicPointPatchField
(
const processorCyclicPointPatchField& ptf,
const pointPatch& p,
const DimensionedField& iF,
const pointPatchFieldMapper& mapper
)
:
coupledPointPatchField(ptf, p, iF, mapper),
procPatch_(refCast(ptf.patch()))
{}
template
Foam::processorCyclicPointPatchField::processorCyclicPointPatchField
(
const processorCyclicPointPatchField& ptf,
const DimensionedField& iF
)
:
coupledPointPatchField(ptf, iF),
procPatch_(refCast(ptf.patch()))
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template
void Foam::processorCyclicPointPatchField::initSwapAddSeparated
(
const Pstream::commsTypes commsType,
Field& pField
) const
{
if (Pstream::parRun())
{
// Get internal field into correct order for opposite side. Note use
// of member data buffer since using non-blocking. Could be optimised
// out if not using non-blocking...
sendBuf_ = this->patchInternalField
(
pField,
procPatch_.reverseMeshPoints()
);
if (commsType == Pstream::commsTypes::nonBlocking)
{
recvBuf_.resize_nocopy(sendBuf_.size());
UIPstream::read
(
commsType,
procPatch_.neighbProcNo(),
recvBuf_.data_bytes(),
recvBuf_.size_bytes(),
procPatch_.tag(),
procPatch_.comm()
);
}
UOPstream::write
(
commsType,
procPatch_.neighbProcNo(),
sendBuf_.cdata_bytes(),
sendBuf_.size_bytes(),
procPatch_.tag(),
procPatch_.comm()
);
}
}
template
void Foam::processorCyclicPointPatchField::swapAddSeparated
(
const Pstream::commsTypes commsType,
Field& pField
) const
{
if (Pstream::parRun())
{
// If nonblocking, data is already in receive buffer
if (commsType != Pstream::commsTypes::nonBlocking)
{
recvBuf_.resize_nocopy(this->size());
UIPstream::read
(
commsType,
procPatch_.neighbProcNo(),
recvBuf_.data_bytes(),
recvBuf_.size_bytes(),
procPatch_.tag(),
procPatch_.comm()
);
}
if (doTransform())
{
const processorCyclicPolyPatch& ppp =
procPatch_.procCyclicPolyPatch();
const tensor& forwardT = ppp.forwardT()[0];
transform(recvBuf_, forwardT, recvBuf_);
}
// All points are separated
this->addToInternalField(pField, recvBuf_);
}
}
// ************************************************************************* //