Commit 0e1ed9ee authored by Yuanxin Liu's avatar Yuanxin Liu
Browse files

Fix parallel streakline filter bugs

vtkParticleTracerBase::PreviousOutputPointData needs to shallow copy
OutputPointData instead of pointer-copy. This is fixed. It is also
renamed to something better: ParticlePointData, with comments added.

The "<" operator for StreakParticle needs to return false for
equality. Otherwise, stl sort will crap out with memory corruption.

Fixed a spurious warning

Change-Id: I2e2c6aa1defab135fc6f9b5767ed7f03ba1bcc92
parent b48917a6
......@@ -259,6 +259,17 @@ int vtkParticleTracerBase::RequestInformation(
{
vtkWarningMacro(<<"Not enough input time steps for particle integration");
}
//clamp the default start time to be within the data time range
if(this->StartTime < this->InputTimeValues[0])
{
this->SetStartTime(this->InputTimeValues[0]);
}
else if (this->StartTime > this->InputTimeValues.back())
{
this->SetStartTime(this->InputTimeValues.back());
}
}
else
{
......@@ -730,8 +741,6 @@ vtkPolyData* vtkParticleTracerBase::Execute(vtkInformationVector** inputVector)
this->CurrentTimeStep==this->StartTimeStep? StartTime:
(this->CurrentTimeStep==this->TerminationTimeStep? this->TerminationTime : this->GetCacheDataTime(1));
this->PreviousOutputPointData = this->OutputPointData;
//set up the output
vtkPolyData *output = vtkPolyData::New();
//
......@@ -769,7 +778,6 @@ vtkPolyData* vtkParticleTracerBase::Execute(vtkInformationVector** inputVector)
this->OutputPointData->Initialize();
vtkDebugMacro(<< "About to Interpolate allocate space");
this->OutputPointData->InterpolateAllocate(this->DataReferenceT[0]->GetPointData());
//
this->ParticleAge->SetName("ParticleAge");
this->ParticleIds->SetName("ParticleId");
this->ParticleSourceIds->SetName("ParticleSourceId");
......@@ -973,6 +981,9 @@ vtkPolyData* vtkParticleTracerBase::Execute(vtkInformationVector** inputVector)
this->OutputPointData->AddArray(this->ParticleAngularVel);
}
this->ParticlePointData = vtkSmartPointer<vtkPointData>::New();
this->ParticlePointData->ShallowCopy(this->OutputPointData);
// save some locator building, by re-using them as time progresses
this->Interpolator->AdvanceOneTimeStep();
......@@ -993,6 +1004,10 @@ int vtkParticleTracerBase::RequestData(
vtkInformationVector *outputVector)
{
PRINT("RD start: "<<this->CurrentTimeStep<<" "<<this->CurrentTime<<" "<<this->StartTimeStep<<" "<<this->TerminationTimeStep);
if(this->StartTimeStep<0)
{
return 0;
}
vtkInformation *outInfo = outputVector->GetInformationObject(0);
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
......@@ -1120,8 +1135,14 @@ void vtkParticleTracerBase::IntegrateParticle(
if (!this->RetryWithPush(info, point1, delT, substeps))
{
Assert(previous.PointId>=0); //the particle must have already been added;
this->SendParticleToAnotherProcess(info,previous, this->PreviousOutputPointData);
if(previous.PointId <0)
{
vtkWarningMacro("the particle should have been added");
}
else
{
this->SendParticleToAnotherProcess(info,previous, this->ParticlePointData);
}
this->ParticleHistories.erase(it);
particle_good = false;
break;
......@@ -1367,7 +1388,7 @@ void vtkParticleTracerBase::SetTerminationTime(double t)
this->ResetCache();
}
if( t <= this->StartTime)
if( t < this->StartTime)
{
vtkWarningMacro("Can't go backward");
t = this->StartTime;
......@@ -1592,6 +1613,16 @@ vtkFloatArray* vtkParticleTracerBase::GetParticleAngularVel(vtkPointData* pd)
return vtkFloatArray::SafeDownCast(pd->GetArray("AngularVelocity"));
}
void vtkParticleTracerBase::PrintParticleHistories()
{
cout<<"Particle id, ages: "<<endl;
for(ParticleListIterator itr = this->ParticleHistories.begin(); itr!=this->ParticleHistories.end();itr++)
{
ParticleInformation& info(*itr);
cout<<info.InjectedPointId<<" "<<info.age<<" "<<endl;
}
cout<<endl;
}
......
......@@ -102,6 +102,7 @@ public:
vtkTypeMacro(vtkParticleTracerBase,vtkPolyDataAlgorithm)
void PrintSelf(ostream& os, vtkIndent indent);
void PrintParticleHistories();
// Description
// Turn on/off vorticity computation at streamline points
......@@ -222,6 +223,8 @@ public:
vtkSmartPointer<vtkPointData> ProtoPD;
vtkIdType UniqueIdCounter;// global Id counter used to give particles a stamp
vtkParticleTracerBaseNamespace::ParticleDataList ParticleHistories;
vtkSmartPointer<vtkPointData> ParticlePointData; //the current particle point data consistent
//with particle history
int ReinjectionCounter;
//Everything related to time
......@@ -466,8 +469,7 @@ private:
vtkSmartPointer<vtkFloatArray> ParticleAngularVel;
vtkSmartPointer<vtkDoubleArray> CellVectors;
vtkSmartPointer<vtkPointData> OutputPointData;
vtkSmartPointer<vtkPointData> PreviousOutputPointData;
vtkSmartPointer<vtkDataSet> DataReferenceT[2];
vtkSmartPointer<vtkDataSet> DataReferenceT[2];
vtkSmartPointer<vtkCellArray> ParticleCells;
vtkParticleTracerBase(const vtkParticleTracerBase&); // Not implemented.
......
......@@ -29,6 +29,7 @@ PURPOSE. See the above copyright notice for more information.
#include <assert.h>
#include <algorithm>
#define DEBUGSTREAKLINEFILTER 1
#ifdef DEBUGSTREAKLINEFILTER
#define Assert(a) assert(a)
#else
......@@ -48,7 +49,7 @@ public:
bool operator<(const StreakParticle& other) const
{
return this->Age >= other.Age;
return this->Age > other.Age;
}
};
......@@ -99,7 +100,6 @@ void StreaklineFilterInternal::Finalize()
Streak& streak = streaks[streakId];
float age = particleAge->GetValue(i);
streak.push_back(StreakParticle(i,age));
}
//sort streaks based on age
......@@ -109,14 +109,13 @@ void StreaklineFilterInternal::Finalize()
std::sort(streak.begin(),streak.end());
}
this->Filter->Output->SetLines(vtkSmartPointer<vtkCellArray>::New());
vtkCellArray* outLines = this->Filter->Output->GetLines();
Assert(outLines->GetNumberOfCells()==0);
Assert(outLines);
for(unsigned int i=0; i<streaks.size();i++)
{
Streak& streak(streaks[i]);
const Streak& streak(streaks[i]);
vtkNew<vtkIdList> ids;
for(unsigned int j=0; j<streak.size();j++)
......
......@@ -69,9 +69,8 @@ void vtkPStreaklineFilter::Finalize()
}
append->Update();
vtkPolyData* appoutput = append->GetOutput();
this->Output->CopyStructure(appoutput);
this->Output->GetPointData()->PassData(appoutput->GetPointData());
this->Output->GetCellData()->PassData(appoutput->GetCellData());
this->Output->Initialize();
this->Output->ShallowCopy(appoutput);
assert(this->Output->GetNumberOfPoints()==totalNumPts);
this->It.Finalize();
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment