Commit 7ab57cc5 authored by Mark OLESEN's avatar Mark OLESEN
Browse files

DEFEATURE: remove CloudToVTK function object (issue #985)

- superseded by the vtkCloud function object, which provides
  significantly more functionality and more versatile output
  (using the vtp format).
parent 0c0bc572
......@@ -95,41 +95,14 @@ addGeometryToScene
return;
}
inputFileName_.clear();
// The CloudToVTK functionObject from
//
// lagrangian/intermediate/submodels/CloudFunctionObjects/
//
// stores file state on cloud OutputProperties itself,
//
// whereas the vtkCloud functionObject treats it like other
// output and stores via the stateFunctionObject.
// Since it uses VTP format, there is only a single file with all fields
// - lookup by cloudName.
const dictionary& cloudDict =
geometryBase::parent_.mesh().lookupObject<IOdictionary>
(
cloudName_ + "OutputProperties"
);
if (cloudDict.found("cloudFunctionObject"))
{
const dictionary& foDict = cloudDict.subDict("cloudFunctionObject");
if (foDict.found(functionObjectName_))
{
foDict.subDict(functionObjectName_)
.readIfPresent("file", inputFileName_);
}
}
else
{
inputFileName_ = getFileName("file", cloudName_);
}
// The vtkCloud stores 'file' via the stateFunctionObject
// (lookup by cloudName).
// It only generates VTP format, which means there is a single file
// containing all fields.
inputFileName_ = getFileName("file", cloudName_);
inputFileName_.expand();
if (inputFileName_.empty())
{
WarningInFunction
......@@ -150,17 +123,11 @@ addGeometryToScene
dataset = reader->GetOutput();
}
else if (inputFileName_.hasExt("vtk"))
{
auto reader = vtkSmartPointer<vtkPolyDataReader>::New();
reader->SetFileName(inputFileName_.c_str());
reader->Update();
dataset = reader->GetOutput();
}
else
{
// Invalid name - ignore
// Invalid name - ignore.
// Don't support VTK legacy format at all - it is too wasteful
// and cumbersome.
inputFileName_.clear();
}
......
submodels/CloudFunctionObjects/CloudToVTK/vtkTools.C
PARCELS=parcels
BASEPARCELS=$(PARCELS)/baseClasses
DERIVEDPARCELS=$(PARCELS)/derived
......
......@@ -28,7 +28,6 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "CloudToVTK.H"
#include "FacePostProcessing.H"
#include "ParticleCollector.H"
#include "ParticleErosion.H"
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
\\/ M anipulation | Copyright (C) 2015-2016 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 "CloudToVTK.H"
#include "vtkTools.H"
#include "floatScalar.H"
// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
template<class CloudType>
void Foam::CloudToVTK<CloudType>::writeData
(
std::ostream& vtkOs,
const bool binary,
const List<floatScalar>& data
) const
{
const label procI = Pstream::myProcNo();
List<List<floatScalar>> allProcData(Pstream::nProcs());
allProcData[procI] = data;
Pstream::gatherList(allProcData);
List<floatScalar> allData =
ListListOps::combine<List<floatScalar>>
(
allProcData,
accessOp<List<floatScalar>>()
);
vtkTools::write(vtkOs, binary, allData);
}
template<class CloudType>
template<class Type>
void Foam::CloudToVTK<CloudType>::writeFieldData
(
std::ostream& vtkOs,
const bool binary,
const List<floatScalar>& data,
const word& title,
const label nParcels
) const
{
vtkOs
<< title << ' '
<< int(pTraits<Type>::nComponents) << ' '
<< nParcels << " float" << std::endl;
writeData(vtkOs, binary, data);
}
template<class CloudType>
void Foam::CloudToVTK<CloudType>::write()
{
label nParcels = this->owner().size();
DynamicList<floatScalar> position(3*nParcels);
DynamicList<floatScalar> U(3*nParcels);
DynamicList<floatScalar> d(nParcels);
DynamicList<floatScalar> age(nParcels);
DynamicList<floatScalar> rho(nParcels);
forAllConstIter(typename CloudType, this->owner(), iter)
{
vtkTools::insert(iter().position(), position);
vtkTools::insert(iter().U(), U);
vtkTools::insert(iter().d(), d);
vtkTools::insert(iter().age(), age);
vtkTools::insert(iter().rho(), rho);
}
reduce(nParcels, sumOp<label>());
binary_ = false;
if (Pstream::master())
{
// Create directory if does not exist
mkDir(this->outputTimeDir());
// Open new file at start up
const fileName fName = this->outputTimeDir()/(type() + ".vtk");
this->setModelProperty("file", fName);
OFstream os(fName, binary_ ? IOstream::BINARY : IOstream::ASCII);
std::ostream& vtkOs = os.stdStream();
vtkTools::writeHeader(vtkOs, binary_, this->modelName().c_str());
vtkOs
<< "DATASET POLYDATA" << std::endl
<< "POINTS " << nParcels << " float" << std::endl;
writeData(vtkOs, binary_, position);
vtkOs
<< "POINT_DATA " << nParcels << std::endl
<< "FIELD attributes " << 4
<< std::endl;
writeFieldData<vector>(vtkOs, binary_, U, "U", nParcels);
writeFieldData<scalar>(vtkOs, binary_, d, "d", nParcels);
writeFieldData<scalar>(vtkOs, binary_, age, "age", nParcels);
writeFieldData<scalar>(vtkOs, binary_, rho, "rho", nParcels);
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class CloudType>
Foam::CloudToVTK<CloudType>::CloudToVTK
(
const dictionary& dict,
CloudType& owner,
const word& modelName
)
:
CloudFunctionObject<CloudType>(dict, owner, modelName, typeName),
binary_(dict.lookupOrDefault("binary", true))
{}
template<class CloudType>
Foam::CloudToVTK<CloudType>::CloudToVTK
(
const CloudToVTK<CloudType>& c
)
:
CloudFunctionObject<CloudType>(c),
binary_(c.binary_)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
template<class CloudType>
Foam::CloudToVTK<CloudType>::~CloudToVTK()
{}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
\\/ M anipulation | Copyright (C) 2015-2016 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::CloudToVTK
Group
grpLagrangianIntermediateFunctionObjects
Description
Write cloud data in VTK format
SourceFiles
CloudToVTK.C
\*---------------------------------------------------------------------------*/
#ifndef CloudToVTK_H
#define CloudToVTK_H
#include "CloudFunctionObject.H"
#include "volFields.H"
#include "OFstream.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class CloudToVTK Declaration
\*---------------------------------------------------------------------------*/
template<class CloudType>
class CloudToVTK
:
public CloudFunctionObject<CloudType>
{
protected:
// Protected data
//- Ascii/binary output flag
bool binary_;
// Protected Member Functions
//- Write post-processing info
virtual void write();
//- Helper function to write VTK data
void writeData
(
std::ostream& vtkOs,
const bool binary,
const List<floatScalar>& data
) const;
//- Helper function to write VTK field data
template<class Type>
void writeFieldData
(
std::ostream& vtkOs,
const bool binary,
const List<floatScalar>& data,
const word& title,
const label nParcels
) const;
public:
//- Runtime type information
TypeName("cloudToVTK");
// Constructors
//- Construct from dictionary
CloudToVTK
(
const dictionary& dict,
CloudType& owner,
const word& modelName
);
//- Construct copy
CloudToVTK(const CloudToVTK<CloudType>& c);
//- Construct and return a clone
virtual autoPtr<CloudFunctionObject<CloudType>> clone() const
{
return autoPtr<CloudFunctionObject<CloudType>>
(
new CloudToVTK<CloudType>(*this)
);
}
//- Destructor
virtual ~CloudToVTK();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "CloudToVTK.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2014 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 "vtkTools.H"
#include "endian.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
void Foam::vtkTools::swapWord(label& word32)
{
char* mem = reinterpret_cast<char*>(&word32);
char a = mem[0];
mem[0] = mem[3];
mem[3] = a;
a = mem[1];
mem[1] = mem[2];
mem[2] = a;
}
void Foam::vtkTools::swapWords(const label nWords, label* words32)
{
for (label i = 0; i < nWords; i++)
{
swapWord(words32[i]);
}
}
void Foam::vtkTools::write
(
std::ostream& os,
const bool binary,
List<floatScalar>& fField
)
{
if (binary)
{
#ifdef WM_LITTLE_ENDIAN
swapWords(fField.size(), reinterpret_cast<label*>(fField.begin()));
#endif
os.write
(
reinterpret_cast<char*>(fField.begin()),
fField.size()*sizeof(float)
);
os << std::endl;
}
else
{
forAll(fField, i)
{
os << fField[i];
if (i > 0 && (i % 10) == 0)
{
os << std::endl;
}
else
{
os << ' ';
}
}
os << std::endl;
}
}
void Foam::vtkTools::write
(
std::ostream& os,
const bool binary,
DynamicList<floatScalar>& fField
)
{
List<floatScalar>& fld = fField.shrink();
write(os, binary, fld);
}
void Foam::vtkTools::write
(
std::ostream& os,
const bool binary,
labelList& elems
)
{
if (binary)
{
#ifdef WM_LITTLE_ENDIAN
swapWords(elems.size(), reinterpret_cast<label*>(elems.begin()));
#endif
os.write
(
reinterpret_cast<char*>(elems.begin()),
elems.size()*sizeof(label)
);
os << std::endl;
}
else
{
forAll(elems, i)
{
os << elems[i];
if (i > 0 && (i % 10) == 0)
{
os << std::endl;
}
else
{
os << ' ';
}
}
os << std::endl;
}
}
void Foam::vtkTools::write
(
std::ostream& os,
const bool binary,
DynamicList<label>& elems
)
{
labelList& fld = elems.shrink();
write(os, binary, fld);
}
void Foam::vtkTools::writeHeader
(
std::ostream& os,
const bool binary,
const std::string& title
)
{
os << "# vtk DataFile Version 2.0" << std::endl
<< title << std::endl;
if (binary)
{
os << "BINARY" << std::endl;
}
else
{
os << "ASCII" << std::endl;
}
}
void Foam::vtkTools::writeCellDataHeader
(
std::ostream& os,
const label nCells,
const label nFields
)
{
os << "CELL_DATA " << nCells << std::endl
<< "FIELD attributes " << nFields << std::endl;
}
void Foam::vtkTools::writePointDataHeader
(
std::ostream& os,
const label nPoints,
const label nFields
)