Commits

Matthew Turk committed 2d119c0

Incomplete change to correct centering types and implicitly extended fields.

Comments (0)

Files changed (6)

field_objects/FieldDefinitions.h

 
 enum CenteringType {
   CellCentered,
-  VertexCentered,
-  FaceCentered
+  CornerCentered,
+  FaceCenteredX,
+  FaceCenteredY,
+  FaceCenteredZ,
+  EdgeCenteredX,
+  EdgeCenteredY,
+  EdgeCenteredZ
 };
 
 enum InterpolationType {

field_objects/FieldDescriptor.C

 
 FieldDescriptor::FieldDescriptor(
     CenteringType ValueCentering, int Rank,
-    int Dimensions[MAX_DIMENSIONS],
+    int CellDimensions[MAX_DIMENSIONS],
     long long LeftEdge[MAX_DIMENSIONS],
     InterpolationType InterpolationMethod,
     const char* Name,
   this->ValueCentering = ValueCentering;
   this->Rank = Rank;
   for (dim = 0; dim < this->Rank; dim++) {
-    this->Dimensions[dim] = Dimensions[dim];
+    this->CellDimensions[dim] = CellDimensions[dim];
     this->LeftEdge[dim] = LeftEdge[dim];
   }
   this->InterpolationMethod = InterpolationMethod;
   if (NewUnitsName == NULL) NewUnitsName = this->UnitsName;
   // Note that this will allocate and set the bits to deallocate later
   FieldDescriptor *fd = new FieldDescriptor(this->ValueCentering, this->Rank,
-        this->Dimensions, this->LeftEdge, this->InterpolationMethod,
+        this->CellDimensions, this->LeftEdge, this->InterpolationMethod,
         NewName, NewUnitsName, NULL);
   if(this->FieldPointer != NULL) {
     // We copy if there are any FieldPointer
 }
 
 int FieldDescriptor::Index(int i, int j, int k) {
-  return ((k*this->Dimensions[1]) + j)*this->Dimensions[0] + i;
+  /* This should be inlined or sped up or something. */
+  int FieldDimensions[MAX_DIMENSIONS];
+  this->GetFieldDimensions(FieldDimensions);
+  return ((k*FieldDimensions[1]) + j)*FieldDimensions[0] + i;
+}
+
+CenteringType FieldDescriptor::GetValueCentering() {
+  return this->ValueCentering;
 }
 
 char* FieldDescriptor::GetName() {
 }
 
 int FieldDescriptor::GetSize() {
-  int size = 1;
-  int dim;
+  int dim, Dimensions[MAX_DIMENSIONS], size = 1;
+  this->GetFieldDimensions(Dimensions);
   for (dim = 0; dim < this->Rank; dim++) {
-    size *= this->Dimensions[dim];
+    size *= Dimensions[dim];
   }
   return size;
 }
 
-int FieldDescriptor::GetDimensions(int Dimensions[MAX_DIMENSIONS]) {
+void FieldDescriptor::GetCellDimensions(int Dimensions[MAX_DIMENSIONS]) {
   int dim;
   for (dim = 0; dim < this->Rank; dim++) {
-    Dimensions[dim] = this->Dimensions[dim];
+    Dimensions[dim] = this->CellDimensions[dim];
   }
 }
 
-int FieldDescriptor::GetFieldDimensions(int Dimensions[MAX_DIMENSIONS]) {
+void FieldDescriptor::GetFieldDimensions(int Dimensions[MAX_DIMENSIONS]) {
   int dim;
+  int Extension[MAX_DIMENSIONS];
+  this->GetCellDimensions(Dimensions);
+  this->GetFieldExtension(Extension);
+  for (dim = 0; dim < this->Rank; dim++) {
+    Dimensions[dim] += Extension[dim];
+  }
+}
+
+void FieldDescriptor::GetFieldExtension(int Extension[MAX_DIMENSIONS]) {
   switch (this->ValueCentering) {
     case CellCentered:
-      for (dim = 0; dim < this->Rank; dim++) {
-        Dimensions[dim] = this->Dimensions[dim];
-      }
+      Extension[0] = Extension[1] = Extension[2] = 0;
       break;
-    case VertexCentered:
-      for (dim = 0; dim < this->Rank; dim++) {
-        Dimensions[dim] = this->Dimensions[dim] + 1;
-      }
+    case CornerCentered:
+      Extension[0] = Extension[1] = Extension[2] = 1;
       break;
-    case FaceCentered:
-      for (dim = 0; dim < this->Rank; dim++) {
-        Dimensions[dim] = this->Dimensions[dim] + 2; // INCORRECT
-      }
-    default:
+    case FaceCenteredX:
+      Extension[0] = 1;
+      Extension[1] = Extension[2] = 0;
+      break;
+    case FaceCenteredY:
+      Extension[1] = 1;
+      Extension[0] = Extension[2] = 0;
+      break;
+    case FaceCenteredZ:
+      Extension[2] = 1;
+      Extension[0] = Extension[1] = 0;
+      break;
+    case EdgeCenteredX:
+      Extension[0] = 0;
+      Extension[1] = Extension[2] = 1;
+      break;
+    case EdgeCenteredY:
+      Extension[1] = 0;
+      Extension[0] = Extension[2] = 1;
+      break;
+    case EdgeCenteredZ:
+      Extension[2] = 0;
+      Extension[0] = Extension[1] = 1;
       break;
   }
 }
       Other->GetFieldDimensions(OtherDims);
       for (dim = 0; dim < this->Rank; dim++) {
         OverlapStart = MAX(GlobalLeftEdgeThis[dim], GlobalLeftEdgeOther[dim]);
-        OverlapEnd = MIN(GlobalLeftEdgeThis[dim] + this->Dimensions[dim],
+        OverlapEnd = MIN(GlobalLeftEdgeThis[dim] + this->CellDimensions[dim],
                          GlobalLeftEdgeOther[dim] + OtherDims[dim]);
         LeftEdgeThis[dim] = OverlapStart - GlobalLeftEdgeThis[dim];
         LeftEdgeOther[dim] = OverlapStart - GlobalLeftEdgeOther[dim];
 
 void FieldDescriptor::CopyFrom(double val) {
   static int Zero[MAX_DIMENSIONS] = {0, 0, 0};
-  this->InPlaceBinaryOperation<CopyVal>(val, Zero, this->Dimensions);
+  this->InPlaceBinaryOperation<CopyVal>(val, Zero, this->CellDimensions);
 }
 
 void FieldDescriptor::Add(FieldDescriptor *Other) {
 
 void FieldDescriptor::Add(double val) {
   static int Zero[MAX_DIMENSIONS] = {0, 0, 0};
-  this->InPlaceBinaryOperation<AddVal>(val, Zero, this->Dimensions);
+  this->InPlaceBinaryOperation<AddVal>(val, Zero, this->CellDimensions);
 }
 
 void FieldDescriptor::Subtract(FieldDescriptor *Other) {
 
 void FieldDescriptor::Subtract(double val) {
   static int Zero[MAX_DIMENSIONS] = {0, 0, 0};
-  this->InPlaceBinaryOperation<SubVal>(val, Zero, this->Dimensions);
+  this->InPlaceBinaryOperation<SubVal>(val, Zero, this->CellDimensions);
 }
 
 void FieldDescriptor::Multiply(FieldDescriptor *Other) {
 
 void FieldDescriptor::Multiply(double val) {
   static int Zero[MAX_DIMENSIONS] = {0, 0, 0};
-  this->InPlaceBinaryOperation<MultVal>(val, Zero, this->Dimensions);
+  this->InPlaceBinaryOperation<MultVal>(val, Zero, this->CellDimensions);
 }
 
 void FieldDescriptor::Divide(FieldDescriptor *Other) {
 
 void FieldDescriptor::Divide(double val) {
   static int Zero[MAX_DIMENSIONS] = {0, 0, 0};
-  this->InPlaceBinaryOperation<DivVal>(val, Zero, this->Dimensions);
+  this->InPlaceBinaryOperation<DivVal>(val, Zero, this->CellDimensions);
 }
 
 // Operator Overloading
 template <MathFunction function>
 double FieldDescriptor::UnaryAccumulator(
     int *LeftEdge, int *RightEdge, double InitialValue) {
+  // NOTE: This takes a RightEdge in *cell* values.  This will *do the right
+  // thing* for face and vertex centered fields.
   int i, j, k;
   double val = InitialValue;
   int ind;
   static int Zero[MAX_DIMENSIONS] = {0, 0, 0};
+  int Extension[MAX_DIMENSIONS];
   if (LeftEdge == NULL) {
     LeftEdge = Zero;
   }
   if (RightEdge == NULL) {
-    RightEdge = this->Dimensions;
+    RightEdge = this->CellDimensions;
   }
+  this->GetFieldExtension(Extension);
   double *vals = this->GetValues();
 
-  for (i = LeftEdge[0]; i < RightEdge[0]; i++) {
-    for (j = LeftEdge[1]; j < RightEdge[1]; j++) {
-      for (k = LeftEdge[2]; k < RightEdge[2]; k++) {
+  for (i = LeftEdge[0]; i < RightEdge[0] + Extension[0]; i++) {
+    for (j = LeftEdge[1]; j < RightEdge[1] + Extension[1]; j++) {
+      for (k = LeftEdge[2]; k < RightEdge[2] + Extension[2]; k++) {
         ind = this->Index(i, j, k);
         val = function(val, vals[ind]);
       }
     throw FieldsIncompatible("FieldDescriptors are not compatible.");
   }
 
+  int Extension[MAX_DIMENSIONS];
+  this->GetFieldExtension(Extension);
+
   int dim, i, j, k, i1, i2, j1, j2, k1, k2, ind1, ind2;
 
   double *v1 = this->GetValues();
   double *v2 = Other->GetValues();
 
   for (i = 0, i1 = LeftEdgeThis[0], i2 = LeftEdgeOther[0];
-       i < CopyDims[0];
+       i < CopyDims[0] + Extension[0];
        i++, i1++, i2++) {
     for (j = 0, j1 = LeftEdgeThis[1], j2 = LeftEdgeOther[1];
-        j < CopyDims[1];
+        j < CopyDims[1] + Extension[1];
         j++, j1++, j2++) {
       for (k = 0, k1 = LeftEdgeThis[2], k2 = LeftEdgeOther[2];
-          k < CopyDims[2];
+          k < CopyDims[2] + Extension[2];
           k++, k1++, k2++) {
         ind1 = this->Index(i1, j1, k1);
         ind2 = Other->Index(i2, j2, k2);
 
   double *v1 = this->GetValues();
 
-  for (i = 0, i1 = LeftEdge[0]; i < CopyDims[0]; i++, i1++) {
-    for (j = 0, j1 = LeftEdge[1]; j < CopyDims[1]; j++, j1++) {
-      for (k = 0, k1 = LeftEdge[2]; k < CopyDims[2]; k++, k1++) {
+  int Extension[MAX_DIMENSIONS];
+  this->GetFieldExtension(Extension);
+
+  for (i = 0, i1 = LeftEdge[0]; i < CopyDims[0] + Extension[0]; i++, i1++) {
+    for (j = 0, j1 = LeftEdge[1]; j < CopyDims[1] + Extension[0]; j++, j1++) {
+      for (k = 0, k1 = LeftEdge[2]; k < CopyDims[2] + Extension[0]; k++, k1++) {
         ind1 = this->Index(i1, j1, k1);
         v1[ind1] = function(v1[ind1], OtherValue);
       }

field_objects/FieldDescriptor.h

       FieldDescriptor *Duplicate(char *NewName = NULL, char *NewUnitsName = NULL);
 
       int Index(int i, int j, int k);
+      CenteringType GetValueCentering();
       char* GetName();
       void SetName(std::string NewName);
       void SetName(const char* NewName);
       double *GetValues();
       int GetSize();
-      int GetDimensions(int Dimensions[MAX_DIMENSIONS]);
-      int GetFieldDimensions(int Dimensions[MAX_DIMENSIONS]);
+      void GetCellDimensions(int Dimensions[MAX_DIMENSIONS]);
+      void GetFieldDimensions(int Dimensions[MAX_DIMENSIONS]);
       void GetLeftEdge(long long LeftEdge[MAX_DIMENSIONS]);
 
       void GetOverlapRegion(FieldDescriptor *Other,
       CenteringType ValueCentering;
       int Rank;
       InterpolationType InterpolationMethod;
+      void GetFieldExtension(int Extension[MAX_DIMENSIONS]);
 
     private:
       char* Name;
       char* UnitsName;
-      int Dimensions[MAX_DIMENSIONS];
+      int CellDimensions[MAX_DIMENSIONS];
       long long LeftEdge[MAX_DIMENSIONS];
       double **FieldPointer;
       int DeallocateFieldPointer;

tests/TestFieldCenterings.C

       virtual void SetUp() {
         int dims[MAX_DIMENSIONS] = {8, 8, 8};
         long long le[MAX_DIMENSIONS] = {0, 0, 0};
-        this->fd_cell = new FieldDescriptor(
+        int fdi;
+        this->fds[fdi++] = new FieldDescriptor(
             CellCentered, 3,
             dims, le, FirstOrderA,
             "Density", "g/cc", NULL);
-        this->fd_vert = new FieldDescriptor(
-            VertexCentered, 3,
+        this->fds[fdi++] = new FieldDescriptor(
+            CornerCentered, 3,
             dims, le, FirstOrderA,
             "Density", "g/cc", NULL);
-        this->fd_face = new FieldDescriptor(
-            FaceCentered, 3,
+
+        this->fds[fdi++] = new FieldDescriptor(
+            FaceCenteredX, 3,
+            dims, le, FirstOrderA,
+            "Density", "g/cc", NULL);
+        this->fds[fdi++] = new FieldDescriptor(
+            FaceCenteredY, 3,
+            dims, le, FirstOrderA,
+            "Density", "g/cc", NULL);
+        this->fds[fdi++] = new FieldDescriptor(
+            FaceCenteredZ, 3,
+            dims, le, FirstOrderA,
+            "Density", "g/cc", NULL);
+
+        this->fds[fdi++] = new FieldDescriptor(
+            EdgeCenteredX, 3,
+            dims, le, FirstOrderA,
+            "Density", "g/cc", NULL);
+        this->fds[fdi++] = new FieldDescriptor(
+            EdgeCenteredY, 3,
+            dims, le, FirstOrderA,
+            "Density", "g/cc", NULL);
+        this->fds[fdi++] = new FieldDescriptor(
+            EdgeCenteredZ, 3,
             dims, le, FirstOrderA,
             "Density", "g/cc", NULL);
       }
 
       virtual void TearDown() {
-        delete this->fd_cell;
-        delete this->fd_vert;
-        delete this->fd_face;
+
       }
 
-      FieldDescriptor *fd_cell;
-      FieldDescriptor *fd_vert;
-      FieldDescriptor *fd_face;
+      FieldDescriptor *fds[8];
   };
 
 }
 
 TEST_F(FieldCenteringSimpleTest, TestCanCombine) {
-  ASSERT_EQ(this->fd_cell->CanCombine(this->fd_cell), 1);
-  ASSERT_EQ(this->fd_cell->CanCombine(this->fd_vert), 0);
-  ASSERT_EQ(this->fd_cell->CanCombine(this->fd_face), 0);
-
-  ASSERT_EQ(this->fd_vert->CanCombine(this->fd_cell), 0);
-  ASSERT_EQ(this->fd_vert->CanCombine(this->fd_vert), 1);
-  ASSERT_EQ(this->fd_vert->CanCombine(this->fd_face), 0);
-
-  ASSERT_EQ(this->fd_face->CanCombine(this->fd_cell), 0);
-  ASSERT_EQ(this->fd_face->CanCombine(this->fd_vert), 0);
-  ASSERT_EQ(this->fd_face->CanCombine(this->fd_face), 1);
+  int i, j, cc;
+  FieldDescriptor *fdi, *fdj;
+  for (i = 0; i < 8; i++) {
+    fdi = this->fds[i];
+    for (j = 0; j < 8; j++) {
+      fdj = this->fds[j];
+      cc = (i == j) ? 1 : 0;
+      ASSERT_EQ(fdi->CanCombine(fdj), cc);
+    }
+  }
 }
 
 TEST_F(FieldCenteringSimpleTest, TestAdd) {
-  ASSERT_NO_THROW(this->fd_cell->Add(this->fd_cell));
-  ASSERT_THROW(this->fd_cell->Add(this->fd_vert), FieldsIncompatible);
-  ASSERT_THROW(this->fd_cell->Add(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_vert->Add(this->fd_cell), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_vert->Add(this->fd_vert));
-  ASSERT_THROW(this->fd_vert->Add(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_face->Add(this->fd_cell), FieldsIncompatible);
-  ASSERT_THROW(this->fd_face->Add(this->fd_vert), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_face->Add(this->fd_face));
+  int i, j, cc;
+  FieldDescriptor *fdi, *fdj;
+  for (i = 0; i < 8; i++) {
+    fdi = this->fds[i];
+    for (j = 0; j < 8; j++) {
+      fdj = this->fds[j];
+      if (i == j) {
+        ASSERT_NO_THROW(fdi->Add(fdj));
+      } else {
+        ASSERT_THROW(fdi->Add(fdj), FieldsIncompatible);
+      }
+    }
+  }
 }
 
 TEST_F(FieldCenteringSimpleTest, TestSubtract) {
-  ASSERT_NO_THROW(this->fd_cell->Subtract(this->fd_cell));
-  ASSERT_THROW(this->fd_cell->Subtract(this->fd_vert), FieldsIncompatible);
-  ASSERT_THROW(this->fd_cell->Subtract(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_vert->Subtract(this->fd_cell), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_vert->Subtract(this->fd_vert));
-  ASSERT_THROW(this->fd_vert->Subtract(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_face->Subtract(this->fd_cell), FieldsIncompatible);
-  ASSERT_THROW(this->fd_face->Subtract(this->fd_vert), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_face->Subtract(this->fd_face));
+  int i, j, cc;
+  FieldDescriptor *fdi, *fdj;
+  for (i = 0; i < 8; i++) {
+    fdi = this->fds[i];
+    for (j = 0; j < 8; j++) {
+      fdj = this->fds[j];
+      if (i == j) {
+        ASSERT_NO_THROW(fdi->Subtract(fdj));
+      } else {
+        ASSERT_THROW(fdi->Subtract(fdj), FieldsIncompatible);
+      }
+    }
+  }
 }
 
 TEST_F(FieldCenteringSimpleTest, TestMultiply) {
-  ASSERT_NO_THROW(this->fd_cell->Multiply(this->fd_cell));
-  ASSERT_THROW(this->fd_cell->Multiply(this->fd_vert), FieldsIncompatible);
-  ASSERT_THROW(this->fd_cell->Multiply(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_vert->Multiply(this->fd_cell), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_vert->Multiply(this->fd_vert));
-  ASSERT_THROW(this->fd_vert->Multiply(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_face->Multiply(this->fd_cell), FieldsIncompatible);
-  ASSERT_THROW(this->fd_face->Multiply(this->fd_vert), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_face->Multiply(this->fd_face));
+  int i, j, cc;
+  FieldDescriptor *fdi, *fdj;
+  for (i = 0; i < 8; i++) {
+    fdi = this->fds[i];
+    for (j = 0; j < 8; j++) {
+      fdj = this->fds[j];
+      if (i == j) {
+        ASSERT_NO_THROW(fdi->Multiply(fdj));
+      } else {
+        ASSERT_THROW(fdi->Multiply(fdj), FieldsIncompatible);
+      }
+    }
+  }
 }
 
 TEST_F(FieldCenteringSimpleTest, TestDivide) {
-  ASSERT_NO_THROW(this->fd_cell->Divide(this->fd_cell));
-  ASSERT_THROW(this->fd_cell->Divide(this->fd_vert), FieldsIncompatible);
-  ASSERT_THROW(this->fd_cell->Divide(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_vert->Divide(this->fd_cell), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_vert->Divide(this->fd_vert));
-  ASSERT_THROW(this->fd_vert->Divide(this->fd_face), FieldsIncompatible);
-
-  ASSERT_THROW(this->fd_face->Divide(this->fd_cell), FieldsIncompatible);
-  ASSERT_THROW(this->fd_face->Divide(this->fd_vert), FieldsIncompatible);
-  ASSERT_NO_THROW(this->fd_face->Divide(this->fd_face));
+  int i, j, cc;
+  FieldDescriptor *fdi, *fdj;
+  for (i = 0; i < 8; i++) {
+    fdi = this->fds[i];
+    for (j = 0; j < 8; j++) {
+      fdj = this->fds[j];
+      if (i == j) {
+        ASSERT_NO_THROW(fdi->Divide(fdj));
+      } else {
+        ASSERT_THROW(fdi->Divide(fdj), FieldsIncompatible);
+      }
+    }
+  }
 }

tests/TestFieldDescriptors.C

     ASSERT_STREQ("DensityThree", fd->GetName());
 
     int Dimensions[MAX_DIMENSIONS];
-    fd->GetDimensions(Dimensions);
+    fd->GetCellDimensions(Dimensions);
     ASSERT_EQ(Dimensions[0], 16);
     ASSERT_EQ(Dimensions[1], 24);
     ASSERT_EQ(Dimensions[2], 32);
   ASSERT_EQ(Position[1], -2);
   ASSERT_EQ(Position[2], -2);
   int Dimensions[MAX_DIMENSIONS];
-  fd->GetDimensions(Dimensions);
+  fd->GetCellDimensions(Dimensions);
   ASSERT_EQ(Dimensions[0], 20);
   ASSERT_EQ(Dimensions[1], 20);
   ASSERT_EQ(Dimensions[2], 20);