From 3c9efbebfb5c6a87d53132323c75e982e8026a09 Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Thu, 15 Jul 2010 13:45:32 +0000
Subject: [PATCH] Fixed 0-diagonal problem in parallel computations.

---
 AMDiS/libtool                            |  6 +++---
 AMDiS/src/parallel/GlobalMatrixSolver.cc | 23 ++++++++++++-----------
 AMDiS/src/parallel/MeshDistributor.cc    |  3 +++
 3 files changed, 18 insertions(+), 14 deletions(-)

diff --git a/AMDiS/libtool b/AMDiS/libtool
index 02151259..1a37a1c4 100755
--- a/AMDiS/libtool
+++ b/AMDiS/libtool
@@ -44,7 +44,7 @@ available_tags=" CXX F77"
 
 # ### BEGIN LIBTOOL CONFIG
 
-# Libtool was configured on host deimos103:
+# Libtool was configured on host deimos102:
 
 # Shell to use when invoking shell scripts.
 SHELL="/bin/sh"
@@ -6760,7 +6760,7 @@ build_old_libs=`case $build_libtool_libs in yes) $echo no;; *) $echo yes;; esac`
 # End:
 # ### BEGIN LIBTOOL TAG CONFIG: CXX
 
-# Libtool was configured on host deimos103:
+# Libtool was configured on host deimos102:
 
 # Shell to use when invoking shell scripts.
 SHELL="/bin/sh"
@@ -7065,7 +7065,7 @@ include_expsyms=""
 
 # ### BEGIN LIBTOOL TAG CONFIG: F77
 
-# Libtool was configured on host deimos103:
+# Libtool was configured on host deimos102:
 
 # Shell to use when invoking shell scripts.
 SHELL="/bin/sh"
diff --git a/AMDiS/src/parallel/GlobalMatrixSolver.cc b/AMDiS/src/parallel/GlobalMatrixSolver.cc
index 05bdb6e6..22823ca9 100644
--- a/AMDiS/src/parallel/GlobalMatrixSolver.cc
+++ b/AMDiS/src/parallel/GlobalMatrixSolver.cc
@@ -4,7 +4,7 @@
 #include "DOFVector.h"
 #include "Debug.h"
 #include "SystemVector.h"
-
+#include "VtkWriter.h"
 #include "petscksp.h"
 
 namespace AMDiS {
@@ -94,6 +94,8 @@ namespace AMDiS {
       int globalRowDof = meshDistributor->mapLocalToGlobal(*cursor);
       // Test if the current row dof is a periodic dof.
       bool periodicRow = meshDistributor->isPeriodicDof(globalRowDof);
+      // Calculate petsc row index.
+      int rowIndex = globalRowDof * dispMult + dispAddRow;
 
 
       // === Traverse all non zero entries of the row and produce vector cols ===
@@ -102,13 +104,14 @@ namespace AMDiS {
 
       for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); 
 	   icursor != icend; ++icursor) {
+	
+	// Global index of the current column index.
+	int globalColDof = meshDistributor->mapLocalToGlobal(col(*icursor));
+	// Calculate the exact position of the column index in the petsc matrix.
+	int colIndex = globalColDof * dispMult + dispAddCol;
 
 	// Set only non null values.
-	if (value(*icursor) != 0.0) {
-	  // Global index of the current column index.
-	  int globalColDof = meshDistributor->mapLocalToGlobal(col(*icursor));
-	  // Calculate the exact position of the column index in the petsc matrix.
-	  int colIndex = globalColDof * dispMult + dispAddCol;
+	if (value(*icursor) != 0.0 || rowIndex == colIndex) {
 
 	  // If the current row is not periodic, but the current dof index is periodic,
 	  // we have to duplicate the value to the other corresponding periodic columns.
@@ -140,10 +143,7 @@ namespace AMDiS {
 
       // === Up to now we have assembled on row. Now, the row must be send to the ===
       // === corresponding rows to the petsc matrix.                              ===
-
-      // Calculate petsc row index.
-      int rowIndex = globalRowDof * dispMult + dispAddRow;
-      
+     
       if (periodicRow) {
 	// The row dof is periodic, so send dof to all the corresponding rows.
 	std::set<int>& perAsc = meshDistributor->getPerDofAssociations(globalRowDof);
@@ -151,7 +151,7 @@ namespace AMDiS {
 	double scalFactor = 1.0 / (perAsc.size() + 1.0);
 
 	for (unsigned int i = 0; i < values.size(); i++)
-	  values[i] *= scalFactor;
+	  values[i] *= scalFactor;	
 
 	// Send the main row to the petsc matrix.
 	MatSetValues(petscMatrix, 1, &rowIndex, cols.size(), 
@@ -174,6 +174,7 @@ namespace AMDiS {
 	  }
 
 	  int perRowIndex = (meshDistributor->getPeriodicMapping(*perIt, globalRowDof) * dispMult) + dispAddRow;
+
 	  MatSetValues(petscMatrix, 1, &perRowIndex, perCols.size(), 
 		       &(perCols[0]), &(perValues[0]), ADD_VALUES);
 	}
diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc
index a12a6d16..1afe17ff 100644
--- a/AMDiS/src/parallel/MeshDistributor.cc
+++ b/AMDiS/src/parallel/MeshDistributor.cc
@@ -1815,10 +1815,13 @@ namespace AMDiS {
     debug::getAllDofs(feSpace, testDofs);
     for (std::set<const DegreeOfFreedom*>::iterator it = testDofs.begin();
 	 it != testDofs.end(); ++it) {
+//       WorldVector<double> coords;
+//       mesh->getDofIndexCoords(**it, feSpace, coords);
       MSG("DOF %d:   mapLocalGlobalDofs = %d   vertexDof = %d    isRankDof = %d\n", **it, 
 	  mapLocalGlobalDofs[**it], 
 	  vertexDof[*it],
 	  isRankDof[**it]);
+      //      MSG("coord = %f %f\n", coords[0], coords[1]);
     }
 #endif
 #endif
-- 
GitLab