Commit 0d9e5b52 authored by whitlocb's avatar whitlocb

Enhanced support for double precision.



git-svn-id: http://visit.ilight.com/svn/visit/trunk/src@17967 18c085ea-50e0-402c-830e-de6fd14e8384
parent 1c5ad1aa
This diff is collapsed.
......@@ -96,6 +96,10 @@ class avtSourceFromDatabase;
// Mark C. Miller, Mon Nov 9 10:40:34 PST 2009
// Changed interface to main transform method to operate on a single
// dataset instead of a dataset collection.
//
// Brad Whitlock, Sun Apr 22 00:01:35 PDT 2012
// I added some methods that test for excess precision.
//
// ****************************************************************************
class DATABASE_API avtTransformManager
......@@ -118,6 +122,11 @@ class DATABASE_API avtTransformManager
void ClearTimestep(int ts) { cache.ClearTimestep(ts); };
private:
bool CoordinatesHaveExcessPrecision(vtkDataSet *ds,
bool needNativePrecision) const;
bool DataHasExcessPrecision(vtkDataArray *da,
bool needNativePrecision) const;
vtkDataSet *NativeToFloat(const avtDatabaseMetaData *const md,
const avtDataRequest_p &spec,
vtkDataSet *ds, int dom);
......
......@@ -91,6 +91,9 @@ class avtMaterial;
// Hank Childs, Thu Feb 14 17:12:38 PST 2008
// Add virtual method "CanOnlyCreateGhostNodes".
//
// Brad Whitlock, Sun Apr 22 10:32:28 PDT 2012
// Added ExchangeDoubleVector.
//
// ****************************************************************************
class DATABASE_API avtDomainBoundaries
......@@ -110,6 +113,10 @@ class DATABASE_API avtDomainBoundaries
bool isPointData,
std::vector<vtkDataArray*> vectors) =0;
virtual std::vector<vtkDataArray*> ExchangeDoubleVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors) =0;
virtual std::vector<vtkDataArray*> ExchangeIntVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors) =0;
......
......@@ -992,6 +992,14 @@ avtNekDomainBoundaries::ExchangeFloatVector(vector<int> domainNum,
EXCEPTION0(ImproperUseException);
}
vector<vtkDataArray*>
avtNekDomainBoundaries::ExchangeDoubleVector(vector<int> domainNum,
bool isPointData,
vector<vtkDataArray*> vectors)
{
EXCEPTION0(ImproperUseException);
}
vector<vtkDataArray*>
avtNekDomainBoundaries::ExchangeIntVector(vector<int> domainNum,
......
......@@ -83,6 +83,9 @@ class avtMaterial;
// Add data member for multiple blocks. This is for the case where multiple
// domains have been put in the same VTK data set.
//
// Brad Whitlock, Sun Apr 22 10:29:43 PDT 2012
// Double support.
//
// ****************************************************************************
class DATABASE_API avtNekDomainBoundaries : public avtDomainBoundaries
......@@ -105,6 +108,10 @@ class DATABASE_API avtNekDomainBoundaries : public avtDomainBoundaries
bool isPointData,
std::vector<vtkDataArray*> vectors);
virtual std::vector<vtkDataArray*> ExchangeDoubleVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors);
virtual std::vector<vtkDataArray*> ExchangeIntVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors);
......
......@@ -178,6 +178,8 @@ class BoundaryHelperFunctions
private:
avtStructuredDomainBoundaries *sdb;
public:
typedef T Storage;
BoundaryHelperFunctions(avtStructuredDomainBoundaries *sdb_) : sdb(sdb_) { }
T ***InitializeBoundaryData();
......@@ -271,6 +273,9 @@ class BoundaryHelperFunctions
// REFINED_ZONE_IN_AMR_GRID to be properly represented in zones that are
// also DUPLICATED_ZONE_INTERNAL_TO_PROBLEM.
//
// Brad Whitlock, Sun Apr 22 09:59:55 PDT 2012
// Support for double.
//
// ****************************************************************************
class DATABASE_API avtStructuredDomainBoundaries : public avtDomainBoundaries
......@@ -309,6 +314,10 @@ class DATABASE_API avtStructuredDomainBoundaries : public avtDomainBoundaries
bool isPointData,
std::vector<vtkDataArray*> vectors);
virtual std::vector<vtkDataArray*> ExchangeDoubleVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors);
virtual std::vector<vtkDataArray*> ExchangeIntVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors);
......@@ -334,6 +343,10 @@ class DATABASE_API avtStructuredDomainBoundaries : public avtDomainBoundaries
bool isPointData,
std::vector<vtkDataArray*> scalars);
virtual std::vector<vtkDataArray*> ExchangeDoubleScalar(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> scalars);
virtual std::vector<vtkDataArray*> ExchangeIntScalar(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> scalars);
......@@ -360,9 +373,11 @@ class DATABASE_API avtStructuredDomainBoundaries : public avtDomainBoundaries
friend class BoundaryHelperFunctions<int>;
friend class BoundaryHelperFunctions<float>;
friend class BoundaryHelperFunctions<double>;
friend class BoundaryHelperFunctions<unsigned char>;
BoundaryHelperFunctions<int> *bhf_int;
BoundaryHelperFunctions<float> *bhf_float;
BoundaryHelperFunctions<double> *bhf_double;
BoundaryHelperFunctions<unsigned char> *bhf_uchar;
// helper methods
......@@ -386,6 +401,12 @@ class DATABASE_API avtCurvilinearDomainBoundaries
virtual std::vector<vtkDataSet*> ExchangeMesh(std::vector<int> domainNum,
std::vector<vtkDataSet*> meshes);
protected:
template <typename Helper>
void ExchangeMesh(Helper *bhf, int vtktype,
std::vector<int> domainNum,
std::vector<vtkDataSet*> meshes,
std::vector<vtkDataSet*> &out);
};
......
......@@ -79,6 +79,7 @@ namespace
template <> MPI_Datatype GetMPIDataType<int>() { return MPI_INT; }
template <> MPI_Datatype GetMPIDataType<char>() { return MPI_CHAR; }
template <> MPI_Datatype GetMPIDataType<float>() { return MPI_FLOAT;}
template <> MPI_Datatype GetMPIDataType<double>() { return MPI_DOUBLE;}
template <> MPI_Datatype GetMPIDataType<unsigned int>()
{ return MPI_UNSIGNED; }
......@@ -92,25 +93,12 @@ namespace
// Use the preprocessor to help ensure that the right template ExchangeData
// function is instantiated.
//
#if defined(_MSC_VER) && (_MSC_VER <= 1200)
// MSVC 6
static float hack_float;
static char hack_char;
static unsigned char hack_unsigned_char;
static int hack_int;
static unsigned int hack_unsigned_int;
#define ExchangeData_float(A,B,C) ExchangeData(A,B,C,hack_float);
#define ExchangeData_char(A,B,C) ExchangeData(A,B,C,hack_char);
#define ExchangeData_unsigned_char(A,B,C) ExchangeData(A,B,C,hack_unsigned_char);
#define ExchangeData_int(A,B,C) ExchangeData(A,B,C,hack_int);
#define ExchangeData_unsigned_int(A,B,C) ExchangeData(A,B,C,hack_unsigned_int);
#else
#define ExchangeData_float ExchangeData<float>
#define ExchangeData_double ExchangeData<double>
#define ExchangeData_char ExchangeData<char>
#define ExchangeData_unsigned_char ExchangeData<unsigned char>
#define ExchangeData_int ExchangeData<int>
#define ExchangeData_unsigned_int ExchangeData<unsigned int>
#endif
// ****************************************************************************
......@@ -599,6 +587,9 @@ avtUnstructuredDomainBoundaries::ExchangeMesh(vector<int> domainNum,
// depending on the platform to work around a problem with templates
// using the MSVC6.0 compiler on Windows.
//
// Brad Whitlock, Sun Apr 22 10:36:38 PDT 2012
// Double support.
//
// ****************************************************************************
vector<vtkDataArray*>
......@@ -626,6 +617,8 @@ avtUnstructuredDomainBoundaries::ExchangeScalar(vector<int> domainNum,
return ExchangeData_char(domainNum, isPointData, scalars);
case VTK_FLOAT:
return ExchangeData_float(domainNum, isPointData, scalars);
case VTK_DOUBLE:
return ExchangeData_double(domainNum, isPointData, scalars);
case VTK_UNSIGNED_CHAR:
return ExchangeData_unsigned_char(domainNum, isPointData, scalars);
case VTK_UNSIGNED_INT:
......@@ -673,6 +666,32 @@ avtUnstructuredDomainBoundaries::ExchangeFloatVector(vector<int> domainNum,
return ExchangeData_float(domainNum, isPointData, vectors);
}
// ****************************************************************************
// Method: avtUnstructuredDomainBoundaries::ExchangeDoubleVector
//
// Purpose:
// Exchange the ghost zone information for some vectors,
// returning the new ones.
//
// Arguments:
// domainNum an array of domain numbers for each mesh
// isPointData true if this is node-centered, false if cell-centered
// vectors an array of vectors
//
// Programmer: Akira Haddox
// Creation: August 11, 2003
//
// Modifications:
//
// ****************************************************************************
vector<vtkDataArray*>
avtUnstructuredDomainBoundaries::ExchangeDoubleVector(vector<int> domainNum,
bool isPointData,
vector<vtkDataArray*> vectors)
{
return ExchangeData_double(domainNum, isPointData, vectors);
}
// ****************************************************************************
// Method: avtUnstructuredDomainBoundaries::ExchangeIntVector
......@@ -1282,101 +1301,26 @@ avtUnstructuredDomainBoundaries::ConfirmMesh(vector<int> domainNum,
// Hank Childs, Wed Feb 14 15:48:00 PST 2007
// Fix bug where last entry in the array was being overwritten.
//
// Brad Whitlock, Sun Apr 22 10:38:45 PDT 2012
// Remove MSVC 6 code.
//
// ****************************************************************************
template <class T>
template <typename T>
vector<vtkDataArray*>
avtUnstructuredDomainBoundaries::ExchangeData(vector<int> &domainNum,
bool isPointData,
vector<vtkDataArray*> &data
#if defined(_MSC_VER) && (_MSC_VER <= 1200)
, T signature
#endif
)
vector<vtkDataArray*> &data)
{
// Gather the needed information
vector<int> domain2proc = CreateDomainToProcessorMap(domainNum);
T ***gainedData;
int **nGainedTuples;
#if defined(_MSC_VER) && (_MSC_VER <= 1200)
//
// This code is an "inline" copy of the CommunicateDataInformation method
// without the various parallel ifdefs. The MSVC 6.0 compiler refused to
// instantiate CommunicateDataInformation so I inlined it.
//
// Get the processor rank
int rank = 0;
int nComponents = 0;
if (data.size())
nComponents = data[0]->GetNumberOfComponents();
gainedData = new T**[nTotalDomains];
nGainedTuples = new int*[nTotalDomains];
int sendDom;
for (sendDom = 0; sendDom < nTotalDomains; ++sendDom)
{
gainedData[sendDom] = new T*[nTotalDomains];
nGainedTuples[sendDom] = new int[nTotalDomains];
int recvDom;
for (recvDom = 0; recvDom < nTotalDomains; ++recvDom)
{
gainedData[sendDom][recvDom] = NULL;
nGainedTuples[sendDom][recvDom] = 0;
// Cases where no computation is required.
if (sendDom == recvDom)
continue;
if (domain2proc[sendDom] == -1 || domain2proc[recvDom] == -1)
continue;
// If this process owns both of the domains, it's an internal
// calculation: no communication needed
if (domain2proc[sendDom] == rank && domain2proc[recvDom] == rank)
{
int i;
for (i = 0; i < domainNum.size(); ++i)
if (domainNum[i] == sendDom)
break;
int index = GetGivenIndex(sendDom, recvDom);
// If no domain boundary, then there's no work to do.
if (index < 0)
continue;
vector<int> &mapRef = isPointData ? givenPoints[index]
: givenCells[index];
int nTuples = mapRef.size();
nGainedTuples[sendDom][recvDom] = nTuples;
gainedData[sendDom][recvDom] = new T[nTuples * nComponents];
T * origPtr = (T*)(data[i]->GetVoidPointer(0));
T * dataPtr = gainedData[sendDom][recvDom];
for (i = 0; i < nTuples; ++i)
{
T *ptr = origPtr + mapRef[i] * nComponents;
int j;
for (j = 0; j < nComponents; ++j)
{
*(dataPtr++) = *(ptr++);
}
}
}
}
}
nComponents = 0;
#else
CommunicateDataInformation<T> (domain2proc, domainNum, data, isPointData,
gainedData, nGainedTuples);
int nComponents = 0;
#endif
vector<vtkDataArray*> out(data.size(), NULL);
if (data.size())
nComponents = data[0]->GetNumberOfComponents();
......
......@@ -118,6 +118,10 @@ class DATABASE_API avtUnstructuredDomainBoundaries : public avtDomainBoundaries
bool isPointData,
std::vector<vtkDataArray*> vectors);
virtual std::vector<vtkDataArray*> ExchangeDoubleVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors);
virtual std::vector<vtkDataArray*> ExchangeIntVector(std::vector<int> domainNum,
bool isPointData,
std::vector<vtkDataArray*> vectors);
......@@ -163,12 +167,7 @@ class DATABASE_API avtUnstructuredDomainBoundaries : public avtDomainBoundaries
template <class T>
std::vector<vtkDataArray*> ExchangeData(std::vector<int> &domainNum,
bool isPointData,
std::vector<vtkDataArray*> &data
#if defined(_MSC_VER) && (_MSC_VER <= 1200)
// Extra argument to help the compiler instantiate the right function.
, T signature
#endif
);
std::vector<vtkDataArray*> &data);
// Communication methods
std::vector<int> CreateDomainToProcessorMap(const std::vector<int> &domainNum);
......
......@@ -37,7 +37,7 @@
*****************************************************************************/
// ************************************************************************* //
// avtBinaryMathExpression.C //
// avtBinaryMathExpression.C //
// ************************************************************************* //
#include <avtBinaryMathExpression.h>
......@@ -218,7 +218,7 @@ avtBinaryMathExpression::DeriveVariable(vtkDataSet *in_ds)
nvals = data2->GetNumberOfTuples();
vtkDataArray *dv = CreateArray(data1);
vtkDataArray *dv = CreateArray(data1, data2);
dv->SetNumberOfComponents(ncomps);
dv->SetNumberOfTuples(nvals);
......@@ -244,6 +244,54 @@ avtBinaryMathExpression::DeriveVariable(vtkDataSet *in_ds)
return dv;
}
// ****************************************************************************
// Method: PrecisionScore
//
// Purpose:
// Assign a "score" to the precision of VTK types so we can come up with
// a winner when comparing 2 types.
//
// Arguments:
// dt : The VTK data type.
//
// Returns: A score for the input data type.
//
// Note:
//
// Programmer: Brad Whitlock
// Creation: Wed Apr 18 14:00:14 PDT 2012
//
// Modifications:
//
// ****************************************************************************
static int
PrecisionScore(int dt)
{
// Assign a score to an input precision so we can compare the precisions
// of the input data arrays to see which wins. This helps us return the
// "higher" type so operations like (int,double) would return double.
int score = 0;
switch(dt)
{
default:
case VTK_VOID: score = 0; break;
case VTK_BIT: score = 1; break;
case VTK_CHAR: score = 2; break;
case VTK_UNSIGNED_CHAR: score = 2; break;
case VTK_SHORT: score = 3; break;
case VTK_UNSIGNED_SHORT: score = 3; break;
case VTK_INT: score = 4; break;
case VTK_UNSIGNED_INT: score = 4; break;
case VTK_LONG: score = 5; break;
case VTK_UNSIGNED_LONG: score = 5; break;
case VTK_FLOAT: score = 6; break;
case VTK_DOUBLE: score = 7; break;
case VTK_ID_TYPE: score = 5; break;
case VTK_SIGNED_CHAR: score = 2; break;
}
return score;
}
// ****************************************************************************
// Method: avtBinaryMathExpression::CreateArray
......@@ -256,12 +304,27 @@ avtBinaryMathExpression::DeriveVariable(vtkDataSet *in_ds)
// Programmer: Hank Childs
// Creation: August 20, 2003
//
// Modifications:
// Brad Whitlock, Wed Apr 18 14:01:10 PDT 2012
// Score the input types so we promote to the best one. This helps us when
// we get inputs like (int, double), in which case we'd want to promote to
// double.
//
// ****************************************************************************
vtkDataArray *
avtBinaryMathExpression::CreateArray(vtkDataArray *in1)
avtBinaryMathExpression::CreateArray(vtkDataArray *in1, vtkDataArray *in2)
{
return in1->NewInstance();
if(in1->GetDataType() == in2->GetDataType())
return in1->NewInstance();
int s1 = PrecisionScore(in1->GetDataType());
int s2 = PrecisionScore(in2->GetDataType());
if(s1 >= s2)
return in1->NewInstance();
return in2->NewInstance();
}
......
......@@ -78,6 +78,9 @@ class vtkDataArray;
// Added support for rectilinear grids with an inherent transform.
// Binary math filters can handle these with no modifications.
//
// Kathleen Biagas, Thu Apr 5 08:55:32 PDT 2012
// Added second array to CreateArray method.
//
// ****************************************************************************
class EXPRESSION_API avtBinaryMathExpression
......@@ -93,7 +96,7 @@ class EXPRESSION_API avtBinaryMathExpression
protected:
virtual vtkDataArray *DeriveVariable(vtkDataSet *);
virtual vtkDataArray *CreateArray(vtkDataArray *);
virtual vtkDataArray *CreateArray(vtkDataArray *, vtkDataArray *);
virtual void DoOperation(vtkDataArray *in1, vtkDataArray *in2,
vtkDataArray *out, int, int) = 0;
virtual int GetNumberOfComponentsInOutput(int ncompsIn1,
......
......@@ -37,7 +37,7 @@
*****************************************************************************/
// ************************************************************************* //
// avtComparisonExpression.C //
// avtComparisonExpression.C //
// ************************************************************************* //
#include <avtComparisonExpression.h>
......
......@@ -37,7 +37,7 @@
*****************************************************************************/
// ************************************************************************* //
// avtComparisonExpression.h //
// avtComparisonExpression.h //
// ************************************************************************* //
#ifndef AVT_COMPARISON_FILTER_H
......@@ -72,8 +72,8 @@ class avtComparisonExpression : public avtBinaryMathExpression
virtual ~avtComparisonExpression();
protected:
virtual int GetNumberOfComponentsInOutput() { return 1; };
virtual vtkDataArray *CreateArray(vtkDataArray *)
virtual int GetNumberOfComponentsInOutput(int, int) { return 1; };
virtual vtkDataArray *CreateArray(vtkDataArray *, vtkDataArray*)
{ return vtkUnsignedCharArray::New(); };
};
......
......@@ -53,6 +53,8 @@
#include <vtkDataSet.h>
#include <vtkPointData.h>
#include <vtkPointDataToCellData.h>
#include <vtkPointSet.h>
#include <vtkRectilinearGrid.h>
#include <vtkUnsignedCharArray.h>
#include <avtExprNode.h>
......@@ -146,7 +148,6 @@ avtExpressionFilter::ProcessArguments(ArgsExpr *args, ExprPipelineState *state)
std::vector<ArgExpr*>::iterator i;
for (i=arguments->begin(); i != arguments->end(); i++)
{
ExprParseTreeNode *n = (*i)->GetExpr();
avtExprNode *expr_node = dynamic_cast<avtExprNode*>((*i)->GetExpr());
if (expr_node == NULL)
{
......@@ -913,3 +914,43 @@ avtExpressionFilter::GetNumericVal(ExprNode *node, double &val)
return ok;
}
// ****************************************************************************
// Method: avtExpressionFilter::CreateArrayFromMesh
//
// Purpose:
//
// Arguments:
// ds The mesh the variable lays on.
//