Changeset ac6d04 for src/solver/solver_singular.cpp

Timestamp:

Apr 10, 2012, 1:55:49 PM (14 years ago)

Author:

Julian Iseringhausen <isering@…>

Children:

a40eea

Parents:

d24c2f

Message:

Merge recent changes of the vmg library into ScaFaCos.

Includes a fix for the communication problems on Jugene.

git-svn-id: ​https://svn.version.fz-juelich.de/scafacos/trunk@1666 5161e1c8-67bf-11de-9fd5-51895aff932f

File:

• src/solver/solver_singular.cpp (modified) (5 diffs)

src/solver/solver_singular.cpp

@@ -42 +42 @@
 
   // Make sure that arrays are big enough to hold expanded system of equations
-  this->Realloc(rhs.Global().SizeGlobal().Product() + 1);
+  this->Realloc(rhs.Global().GlobalSize().Product() + 1);
 
   for (grid_iter = rhs.Iterators().Local().Begin(); grid_iter != rhs.Iterators().Local().End(); ++grid_iter) {
 
     // Compute 1-dimensional index from 3-dimensional grid
-    index = rhs.GlobalLinearIndex(*grid_iter - rhs.Local().Begin() + rhs.Global().BeginLocal());
+    index = rhs.GlobalLinearIndex(*grid_iter - rhs.Local().Begin() + rhs.Global().LocalBegin());
 
     // Check if we computed the index correctly
@@ -64 +64 @@
     for (stencil_iter = A.begin(); stencil_iter != A.end(); ++stencil_iter) {
 
-      i = *grid_iter - rhs.Local().Begin() + rhs.Global().BeginLocal() + stencil_iter->Disp();
+      i = *grid_iter - rhs.Local().Begin() + rhs.Global().LocalBegin() + stencil_iter->Disp();
 
       for (int j=0; j<3; ++j)
         if (comm->BoundaryConditions()[j] == Periodic) {
           if (i[j] < 0)
-            i[j] += rhs.Global().SizeGlobal()[j];
-          else if (i[j] >= rhs.Global().SizeGlobal()[j])
-            i[j] -= rhs.Global().SizeGlobal()[j];
+            i[j] += rhs.Global().GlobalSize()[j];
+          else if (i[j] >= rhs.Global().GlobalSize()[j])
+            i[j] -= rhs.Global().GlobalSize()[j];
         }
 
@@ -87 +87 @@
     row_sum += iter->Val();
 
-  if (fabs(row_sum) <= std::numeric_limits<vmg_float>::epsilon()) {
+  if (std::abs(row_sum) <= (A.size()+1) * std::numeric_limits<vmg_float>::epsilon()) {
 
     // Expand equation system in order to make the system regular.
@@ -108 +108 @@
 {
   int index;
+  const vmg_float correction = this->Sol(this->Size()-1);
 
   for (int i=0; i<sol.Local().Size().X(); i++)
@@ -114 +115 @@
 
         // Compute global 1-dimensional index
-        index = sol.GlobalLinearIndex(sol.Global().BeginLocal().X()+i,
-                                      sol.Global().BeginLocal().Y()+j,
-                                      sol.Global().BeginLocal().Z()+k);
+        index = sol.GlobalLinearIndex(sol.Global().LocalBegin().X()+i,
+                                      sol.Global().LocalBegin().Y()+j,
+                                      sol.Global().LocalBegin().Z()+k);
 
         // Set solution
-        sol(sol.Local().Begin().X()+i, sol.Local().Begin().Y()+j, sol.Local().Begin().Z()+k) = this->Sol(index);
+        sol(sol.Local().Begin().X()+i, sol.Local().Begin().Y()+j, sol.Local().Begin().Z()+k) = this->Sol(index) - correction;
 
       }