Changeset f57182 for src/base/interface.cpp

Timestamp:

Mar 30, 2013, 2:44:52 PM (13 years ago)

Author:

Julian Iseringhausen <isering@…>

Children:

8180d8

Parents:

d13e27

git-author:

Julian Iseringhausen <isering@…> (06/11/12 14:02:16)

git-committer:

Julian Iseringhausen <isering@…> (03/30/13 14:44:52)

Message:

Open boundary conditions.

Conflicts:

lib/vmg/src/Makefile.am
lib/vmg/src/base/factory.cpp
lib/vmg/test/unit_test/library/dirichlet_fas_lr_mpi.cpp

File:

• src/base/interface.cpp (modified) (2 diffs)

src/base/interface.cpp

@@ -38 +38 @@
 using namespace VMG;
 
+static Index GetGlobalIndex(const Vector& pos, const SpatialExtent& extent, const BT& bt)
+{
+  const Index index = (pos - extent.Begin()) / extent.MeshWidth() + 0.5;
+  return index + (bt == LocallyRefined ? 1 : 0);
+}
+
 void Interface::InitInterface(const Vector& box_offset, const vmg_float& box_size,
                               const int& coarseningSteps, const vmg_float& alpha)
 {
   int i;
-  Index num_cells = Helper::intpow(2, levelMax);
   Index size_factor;
 
   const Vector box_center = box_offset + 0.5 * box_size;
 
+  /*
+   * Get Extents
+   */
   for (i=0; i<coarseningSteps; ++i) {
-
-    global.push_back(GlobalIndices());
     extent.push_back(SpatialExtent());
-
     for (int j=0; j<3; ++j)
       size_factor[j] = (bc[j] == Periodic ? 1 : Helper::intpow(2, static_cast<int>(log(pow(alpha, i+1)) / log(2.0) + 1.0)));
@@ -58 +63 @@
     extent.back().Begin() = box_center - 0.5 * extent.back().Size();
     extent.back().End() = extent.back().Begin() + extent.back().Size();
-    extent.back().MeshWidth() = pow(2.0, i-levelMax);
-
-    num_cells = Helper::intpow(2,levelMax-i) * size_factor;
-
-    global.back().GlobalSize() = num_cells + 1;
-    global.back().LocalSize() = num_cells + 1;
-    global.back().LocalBegin() = 0;
-    global.back().LocalEnd() = num_cells + 1;
-
-    global.back().FinestAbsSize() = Helper::intpow(2,i) * num_cells + 1;
-    global.back().FinestAbsBegin() = ((global.back().FinestAbsSize()-1) * (1-size_factor)) / (2*size_factor);
-    global.back().FinestAbsEnd() = global.back().FinestAbsBegin() + global.back().FinestAbsSize();
-
-    if (i==0)
-      global.back().GlobalFinerBegin() = (num_cells - Helper::intpow(2, levelMax-i))/2;
-    else
-      global.back().GlobalFinerBegin() = (global.back().FinestAbsSize() - (++global.rbegin())->FinestAbsSize()) / Helper::intpow(2,i+1);
-
-    global.back().GlobalFinerEnd() = global.back().GlobalSize() - global.back().GlobalFinerBegin();
-    global.back().GlobalFinerSize() = global.back().GlobalFinerEnd() - global.back().GlobalFinerBegin();
-
-    global.back().LocalFinerBegin() = global.back().GlobalFinerBegin();
-    global.back().LocalFinerEnd() = global.back().GlobalFinerEnd();
-    global.back().LocalFinerSize() = global.back().GlobalFinerSize();
+    extent.back().MeshWidth() = pow(2.0, i-levelMax) * extent.back().Size() / size_factor;
   }
 
-  while (global.size() == 0 || global.back().GlobalSize().Min() > Helper::intpow(2, levelMin)+1) {
+  if (extent.size() == 0) {
+    extent.push_back(SpatialExtent());
+    extent.back().Size() = box_size;
+    extent.back().Begin() = box_offset;
+    extent.back().End() = box_offset + box_size;
+    extent.back().MeshWidth() = box_size / Helper::intpow(2,levelMax);
+  }
 
-    if (global.size() > 0)
-      num_cells /= 2;
+  while ((extent.back().Size() / extent.back().MeshWidth()).Min() > Helper::intpow(2,levelMin) + 1) {
+    extent.push_back(SpatialExtent(extent.back()));
+    extent.back().MeshWidth() *= 2;
+  }
 
+  /*
+   * Find GlobalMax
+   */
+  for (std::vector<SpatialExtent>::const_iterator iter=extent.begin(); iter!=extent.end(); ++iter) {
     global.push_back(GlobalIndices());
-    extent.push_back(SpatialExtent());
+    global.back().BoundaryType() = GlobalCoarsened;
+  }
 
-    if (global.size() == 1) {
-      extent.back().Size() = box_size;
-      extent.back().Begin() = box_offset;
-      extent.back().End() = box_offset + box_size;
-      extent.back().MeshWidth() = box_size / static_cast<Vector>(num_cells);
-    }else {
-      extent.back().Size() = (++extent.rbegin())->Size();
-      extent.back().Begin() = (++extent.rbegin())->Begin();
-      extent.back().End() = (++extent.rbegin())->End();
-      extent.back().MeshWidth() = 2.0 * (++extent.rbegin())->MeshWidth();
+  for (i=global.size()-1; i>=0; --i) {
+    if (i == 0 || extent[i-1].Size() != extent[i].Size()) {
+      global[i].BoundaryType() = GlobalMax;
+      break;
+    }
+  }
+  for (--i; i>=0; --i)
+    global[i].BoundaryType() = LocallyRefined;
+
+  /*
+   * Compute global grid values
+   */
+  for (i=0; i<global.size(); ++i) {
+
+    const Index num_cells = extent[i].Size() / extent[i].MeshWidth() + 0.5;
+
+    for (int j=0; j<3; ++j) {
+
+      if (bc[j] == Dirichlet || bc[j] == Open) {
+        /*
+         * Dirichlet
+         */
+
+        if (global[i].BoundaryType() == LocallyRefined) {
+          /*
+           * Locally Refined
+           */
+
+          global[i].GlobalSize()[j] = num_cells[j] + 3;
+
+          if (i == 0) {
+            global[i].GlobalFinerBegin()[j] = GetGlobalIndex(box_offset, extent[i], global[i].BoundaryType())[j];
+            global[i].GlobalFinerEnd()[j] = GetGlobalIndex(box_offset + box_size, extent[i], global[i].BoundaryType())[j] + 1;
+            global[i].GlobalFinerSize()[j] = global[i].GlobalFinerEnd()[j] - global[i].GlobalFinerBegin()[j];
+          } else {
+            global[i].GlobalFinerBegin()[j] = GetGlobalIndex(extent[i-1].Begin(), extent[i], global[i].BoundaryType())[j];
+            global[i].GlobalFinerEnd()[j] = GetGlobalIndex(extent[i-1].End(), extent[i], global[i].BoundaryType())[j] + 1;
+            global[i].GlobalFinerSize()[j] = global[i].GlobalFinerEnd()[j] - global[i].GlobalFinerBegin()[j];
+          }
+
+        } else if (global[i].BoundaryType() == GlobalMax) {
+          /*
+           * Global Max
+           */
+
+          global[i].GlobalSize()[j] = num_cells[j] + 1;
+
+          if (i == 0) {
+            global[i].GlobalFinerBegin()[j] = 0;
+            global[i].GlobalFinerEnd()[j] = num_cells[j] + 1;
+            global[i].GlobalFinerSize()[j] = num_cells[j] + 1;
+          } else {
+            global[i].GlobalFinerBegin()[j] = GetGlobalIndex(extent[i-1].Begin(), extent[i], global[i].BoundaryType())[j];
+            global[i].GlobalFinerEnd()[j] = GetGlobalIndex(extent[i-1].End(), extent[i], global[i].BoundaryType())[j] + 1;
+            global[i].GlobalFinerSize()[j] = global[i].GlobalFinerEnd()[j] - global[i].GlobalFinerBegin()[j];
+          }
+
+        } else {
+          /*
+           * Global Coarsened
+           */
+
+          global[i].GlobalSize()[j] = num_cells[j] + 1;
+
+          global[i].GlobalFinerBegin()[j] = 0;
+          global[i].GlobalFinerEnd()[j] = num_cells[j] + 1;
+          global[i].GlobalFinerSize()[j] = num_cells[j] + 1;
+
+        }
+
+      } else if (bc[j] == Periodic) {
+        /*
+         * Periodic
+         */
+
+        global[i].GlobalSize()[j] = num_cells[j];
+
+        global[i].GlobalFinerBegin()[j] = 0;
+        global[i].GlobalFinerEnd()[j] = num_cells[j];
+        global[i].GlobalFinerSize()[j] = num_cells[j];
+
+      }
     }
 
-    global.back().GlobalSize() = num_cells + 1;
-    global.back().LocalSize() = num_cells + 1;
-    global.back().LocalBegin() = 0;
-    global.back().LocalEnd() = num_cells + 1;
+    global[i].LocalBegin() = 0;
+    global[i].LocalEnd() = global[i].GlobalSize();
+    global[i].LocalSize() = global[i].GlobalSize();
 
-    if (global.size() == 1) {
-      global.back().FinestAbsBegin() = 0;
-      global.back().FinestAbsEnd() = global.back().GlobalSize();
-      global.back().FinestAbsSize() = global.back().GlobalSize();
-    }else {
-      global.back().FinestAbsBegin() = (++global.rbegin())->FinestAbsBegin();
-      global.back().FinestAbsEnd() = (++global.rbegin())->FinestAbsEnd();
-      global.back().FinestAbsSize() = (++global.rbegin())->FinestAbsSize();
-    }
-
-    global.back().GlobalFinerBegin() = 0;
-    global.back().GlobalFinerEnd() = global.back().GlobalSize();
-    global.back().GlobalFinerSize() = global.back().GlobalSize();
-
-    global.back().LocalFinerBegin() = global.back().GlobalFinerBegin();
-    global.back().LocalFinerEnd() = global.back().GlobalFinerEnd();
-    global.back().LocalFinerSize() = global.back().GlobalFinerSize();
+    global[i].LocalFinerBegin() = global[i].GlobalFinerBegin();
+    global[i].LocalFinerEnd() = global[i].GlobalFinerEnd();
+    global[i].LocalFinerSize() = global[i].GlobalFinerSize();
 
   }
 
-  for (i=0; i<3; ++i)
-    if (bc[i] == Periodic)
-      for (unsigned int j=0; j<global.size(); ++j) {
-        global[j].GlobalSize()[i] -= 1;
-        global[j].FinestAbsSize()[i] -= Helper::intpow(2, j);
-        global[j].FinestAbsEnd()[i] -= Helper::intpow(2, j);
-        global[j].LocalSize()[i] -= 1;
-        global[j].LocalEnd()[i] -= 1;
-        global[j].GlobalFinerSize()[i] -= 1;
-        global[j].GlobalFinerEnd()[i] -= 1;
-        global[j].LocalFinerSize()[i] -= 1;
-        global[j].LocalFinerEnd()[i] -= 1;
-      }
-
   levelMin = levelMax - global.size() + 1;
 
-  global.back().BoundaryType() = GlobalCoarsened;
-
-  for (i=global.size()-2; i>=0; --i) {
-    if (global[i].FinestAbsSize().Product() >= global[i+1].FinestAbsSize().Product()) {
-      global[i].BoundaryType() = GlobalCoarsened;
-    }else {
-      global[i].BoundaryType() = LocallyRefined;
-      global[i+1].BoundaryType() = GlobalMax;
-      break;
-    }
-  }
-
-  for (; i>=0; --i)
-    global[i].BoundaryType() = LocallyRefined;
-
-  if (global.front().BoundaryType() != LocallyRefined &&
-      global.front().BoundaryType() != GlobalMax)
-    global.front().BoundaryType() = GlobalMax;
-
 }