openfoam/src/lagrangian/solidParticle/solidParticleIO.C
Mark Olesen b6ca9887b0 ENH: add parcel property types (fixes #109)
- assists in decoding what the binary IO content means

ENH: use label instead of bool for the KinematicParcel active state (fixes #111)

- avoids internal padding of the data structure and simplifies
  downstream use.

ENH: make particle sizeofFields public (fixes #110)

Also fixes #108 (missing faceI and stepFraction entries) which had
also been fixed in the upstream as well.
2016-04-21 08:47:50 +02:00

145 lines
3.6 KiB
C

/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "solidParticle.H"
#include "IOstreams.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const std::size_t Foam::solidParticle::sizeofFields
(
sizeof(solidParticle) - sizeof(particle)
);
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::solidParticle::solidParticle
(
const polyMesh& mesh,
Istream& is,
bool readFields
)
:
particle(mesh, is, readFields)
{
if (readFields)
{
if (is.format() == IOstream::ASCII)
{
d_ = readScalar(is);
is >> U_;
}
else
{
is.read(reinterpret_cast<char*>(&d_), sizeofFields);
}
}
// Check state of Istream
is.check("solidParticle::solidParticle(Istream&)");
}
void Foam::solidParticle::readFields(Cloud<solidParticle>& c)
{
if (!c.size())
{
return;
}
particle::readFields(c);
IOField<scalar> d(c.fieldIOobject("d", IOobject::MUST_READ));
c.checkFieldIOobject(c, d);
IOField<vector> U(c.fieldIOobject("U", IOobject::MUST_READ));
c.checkFieldIOobject(c, U);
label i = 0;
forAllIter(Cloud<solidParticle>, c, iter)
{
solidParticle& p = iter();
p.d_ = d[i];
p.U_ = U[i];
i++;
}
}
void Foam::solidParticle::writeFields(const Cloud<solidParticle>& c)
{
particle::writeFields(c);
label np = c.size();
IOField<scalar> d(c.fieldIOobject("d", IOobject::NO_READ), np);
IOField<vector> U(c.fieldIOobject("U", IOobject::NO_READ), np);
label i = 0;
forAllConstIter(Cloud<solidParticle>, c, iter)
{
const solidParticle& p = iter();
d[i] = p.d_;
U[i] = p.U_;
i++;
}
d.write();
U.write();
}
// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
Foam::Ostream& Foam::operator<<(Ostream& os, const solidParticle& p)
{
if (os.format() == IOstream::ASCII)
{
os << static_cast<const particle&>(p)
<< token::SPACE << p.d_
<< token::SPACE << p.U_;
}
else
{
os << static_cast<const particle&>(p);
os.write
(
reinterpret_cast<const char*>(&p.d_),
solidParticle::sizeofFields
);
}
// Check state of Ostream
os.check("Ostream& operator<<(Ostream&, const solidParticle&)");
return os;
}
// ************************************************************************* //