From 0558ae80a37f9a0048c4322c0bb4c3f2961fc6fb Mon Sep 17 00:00:00 2001
From: wolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Fri, 24 Nov 2000 15:42:32 +0000
Subject: [PATCH] Remove block numbers as templates arguments of
 BlockSparsityPattern.

git-svn-id: https://svn.dealii.org/trunk@3501 0785d39b-7218-0410-832d-ea1e28bc413d
---
 .../deal.II/include/dofs/dof_constraints.h    |   5 +-
 .../include/dofs/dof_constraints.templates.h  | 143 +-----------
 .../deal.II/source/dofs/dof_constraints.cc    | 151 ++++++++++++-
 deal.II/deal.II/source/dofs/dof_tools.cc      |  12 +-
 deal.II/deal.II/source/numerics/matrices.cc   |   2 +-
 deal.II/lac/include/lac/block_indices.h       |  19 ++
 deal.II/lac/include/lac/block_sparse_matrix.h |  30 +--
 .../lac/block_sparse_matrix.templates.h       |   6 +-
 .../lac/include/lac/block_sparsity_pattern.h  | 161 ++++++++++----
 .../lac/block_sparsity_pattern.templates.h    | 184 ----------------
 deal.II/lac/include/lac/block_vector.h        |  22 +-
 deal.II/lac/source/block_sparsity_pattern.cc  | 207 +++++++++++++++++-
 tests/deal.II/block_matrices.cc               |  32 ++-
 13 files changed, 540 insertions(+), 434 deletions(-)
 delete mode 100644 deal.II/lac/include/lac/block_sparsity_pattern.templates.h

diff --git a/deal.II/deal.II/include/dofs/dof_constraints.h b/deal.II/deal.II/include/dofs/dof_constraints.h
index 316495496a..526304a923 100644
--- a/deal.II/deal.II/include/dofs/dof_constraints.h
+++ b/deal.II/deal.II/include/dofs/dof_constraints.h
@@ -20,7 +20,7 @@
 #include <base/subscriptor.h>
 
 class SparsityPattern;
-template <int rows, int columns> class BlockSparsityPattern;
+class BlockSparsityPattern;
 template <typename number> class SparseMatrix;
 template <typename number, int rows, int columns> class BlockSparseMatrix;
 class BlockIndices;
@@ -255,8 +255,7 @@ class ConstraintMatrix : public Subscriptor
 				      * condenses square block sparsity
 				      * patterns.
 				      */
-    template <int blocks>
-    void condense (BlockSparsityPattern<blocks,blocks> &sparsity) const;
+    void condense (BlockSparsityPattern &sparsity) const;
     
 				     /**
 				      * Condense a given matrix. The associated
diff --git a/deal.II/deal.II/include/dofs/dof_constraints.templates.h b/deal.II/deal.II/include/dofs/dof_constraints.templates.h
index 6deb241407..20dcee477b 100644
--- a/deal.II/deal.II/include/dofs/dof_constraints.templates.h
+++ b/deal.II/deal.II/include/dofs/dof_constraints.templates.h
@@ -18,145 +18,6 @@
 
 
 
-template <int blocks>
-void ConstraintMatrix::condense (BlockSparsityPattern<blocks,blocks> &sparsity) const
-{
-  Assert (sorted == true, ExcMatrixNotClosed());
-  Assert (sparsity.is_compressed() == false, ExcMatrixIsClosed());
-  Assert (sparsity.n_rows() == sparsity.n_cols(),
-	  ExcMatrixNotSquare());
-  Assert (sparsity.get_column_indices() == sparsity.get_row_indices(),
-	  ExcMatrixNotSquare());
-  
-  const BlockIndices &
-    index_mapping = sparsity.get_column_indices();
-  
-				   // store for each index whether it
-				   // must be distributed or not. If entry
-				   // is -1, no distribution is necessary.
-				   // otherwise, the number states which
-				   // line in the constraint matrix handles
-				   // this index
-  vector<int> distribute (sparsity.n_rows(), -1);
-  
-  for (unsigned int c=0; c<lines.size(); ++c)
-    distribute[lines[c].line] = static_cast<signed int>(c);
-
-  const unsigned int n_rows = sparsity.n_rows();
-  for (unsigned int row=0; row<n_rows; ++row)
-    {
-				       // get index of this row
-				       // within the blocks
-      const pair<unsigned int,unsigned int>
-	block_index = index_mapping.global_to_local(row);
-      const unsigned int block_row = block_index.first;
-      
-      if (distribute[row] == -1)
-					 // regular line. loop over
-					 // all columns and see
-					 // whether this column must
-					 // be distributed
-	{
-
-					   // to loop over all entries
-					   // in this row, we have to
-					   // loop over all blocks in
-					   // this blockrow and the
-					   // corresponding row
-					   // therein
-	  for (unsigned int block_col=0; block_col<blocks; ++block_col)
-	    {
-	      const SparsityPattern &
-		block_sparsity = sparsity.block(block_row, block_col);
-	      
-	      const unsigned int
-		first = block_sparsity.get_rowstart_indices()[block_index.second],
-		last  = block_sparsity.get_rowstart_indices()[block_index.second+1];
-	      for (unsigned int j=first; j<last; ++j)
-						 // end of row reached?
-		if (block_sparsity.get_column_numbers()[j] == SparsityPattern::invalid_entry)
-		  {
-						     // nothing more
-						     // to do
-		    break;
-		  }
-		else
-		  {
-		    const unsigned int global_col
-		      = index_mapping.local_to_global(block_col,
-						      block_sparsity.get_column_numbers()[j]);
-		    
-		    if (distribute[global_col] != -1)
-						       // distribute entry at regular
-						       // row @p{row} and irregular column
-						       // global_col
-		      {
-			for (unsigned int q=0;
-			     q!=lines[distribute[global_col]]
-					    .entries.size(); ++q)
-			  sparsity.add (row,
-					lines[distribute[global_col]].entries[q].first);
-		      };
-		  };
-	    };
-	}
-      else
-	{
-					   // row must be
-					   // distributed. split the
-					   // whole row into the
-					   // chunks defined by the
-					   // blocks
-	  for (unsigned int block_col=0; block_col<blocks; ++block_col)
-	    {
-	      const SparsityPattern &
-		block_sparsity = sparsity.block(block_row,block_col);
-	      
-	      const unsigned int
-		first = block_sparsity.get_rowstart_indices()[block_index.second],
-		last  = block_sparsity.get_rowstart_indices()[block_index.second+1];
-      
-	      for (unsigned int j=first; j<last; ++j)
-						 // end of row reached?
-		if (block_sparsity.get_column_numbers()[j] == SparsityPattern::invalid_entry)
-		  {
-						     // nothing more to do
-		    break;
-		  }
-		else
-		  {
-		    const unsigned int global_col
-		      = index_mapping.local_to_global (block_col,
-						       block_sparsity.get_column_numbers()[j]);
-		    
-		    if (distribute[global_col] == -1)
-						       // distribute entry at irregular
-						       // row @p{row} and regular column
-						       // global_col.
-		      {
-			for (unsigned int q=0; q!=lines[distribute[row]].entries.size(); ++q) 
-			  sparsity.add (lines[distribute[row]].entries[q].first,
-					global_col);
-		      }
-		    else
-						       // distribute entry at irregular
-						       // row @p{row} and irregular column
-						       // @p{global_col}
-		      {
-			for (unsigned int p=0; p!=lines[distribute[row]].entries.size(); ++p)
-			  for (unsigned int q=0; q!=lines[distribute[global_col]].entries.size(); ++q)
-			    sparsity.add (lines[distribute[row]].entries[p].first,
-					  lines[distribute[global_col]].entries[q].first);
-		      };
-		  };
-	    };
-	};
-    };
-  
-  sparsity.compress();
-};
-
-
 
 template<typename number>
 void
@@ -405,13 +266,15 @@ template <typename number, int blocks>
 void
 ConstraintMatrix::condense (BlockSparseMatrix<number,blocks,blocks> &uncondensed) const
 {
-  const BlockSparsityPattern<blocks,blocks> &
+  const BlockSparsityPattern &
     sparsity = uncondensed.get_sparsity_pattern ();
 
   Assert (sorted == true, ExcMatrixNotClosed());
   Assert (sparsity.is_compressed() == true, ExcMatrixNotClosed());
   Assert (sparsity.n_rows() == sparsity.n_cols(),
 	  ExcMatrixNotSquare());
+  Assert (sparsity.n_block_rows() == sparsity.n_block_cols(),
+	  ExcMatrixNotSquare());
   Assert (sparsity.get_column_indices() == sparsity.get_row_indices(),
 	  ExcMatrixNotSquare());
 
diff --git a/deal.II/deal.II/source/dofs/dof_constraints.cc b/deal.II/deal.II/source/dofs/dof_constraints.cc
index d0bbdbbd15..d6975d3bdb 100644
--- a/deal.II/deal.II/source/dofs/dof_constraints.cc
+++ b/deal.II/deal.II/source/dofs/dof_constraints.cc
@@ -410,6 +410,149 @@ void ConstraintMatrix::condense (SparsityPattern &sparsity) const
 
 
 
+void ConstraintMatrix::condense (BlockSparsityPattern &sparsity) const
+{
+  Assert (sorted == true, ExcMatrixNotClosed());
+  Assert (sparsity.is_compressed() == false, ExcMatrixIsClosed());
+  Assert (sparsity.n_rows() == sparsity.n_cols(),
+	  ExcMatrixNotSquare());
+  Assert (sparsity.n_block_rows() == sparsity.n_block_cols(),
+	  ExcMatrixNotSquare());
+  Assert (sparsity.get_column_indices() == sparsity.get_row_indices(),
+	  ExcMatrixNotSquare());
+  
+  const BlockIndices &
+    index_mapping = sparsity.get_column_indices();
+
+  const unsigned int n_blocks = sparsity.n_block_rows();
+  
+				   // store for each index whether it
+				   // must be distributed or not. If entry
+				   // is -1, no distribution is necessary.
+				   // otherwise, the number states which
+				   // line in the constraint matrix handles
+				   // this index
+  vector<int> distribute (sparsity.n_rows(), -1);
+  
+  for (unsigned int c=0; c<lines.size(); ++c)
+    distribute[lines[c].line] = static_cast<signed int>(c);
+
+  const unsigned int n_rows = sparsity.n_rows();
+  for (unsigned int row=0; row<n_rows; ++row)
+    {
+				       // get index of this row
+				       // within the blocks
+      const pair<unsigned int,unsigned int>
+	block_index = index_mapping.global_to_local(row);
+      const unsigned int block_row = block_index.first;
+      
+      if (distribute[row] == -1)
+					 // regular line. loop over
+					 // all columns and see
+					 // whether this column must
+					 // be distributed
+	{
+
+					   // to loop over all entries
+					   // in this row, we have to
+					   // loop over all blocks in
+					   // this blockrow and the
+					   // corresponding row
+					   // therein
+	  for (unsigned int block_col=0; block_col<n_blocks; ++block_col)
+	    {
+	      const SparsityPattern &
+		block_sparsity = sparsity.block(block_row, block_col);
+	      
+	      const unsigned int
+		first = block_sparsity.get_rowstart_indices()[block_index.second],
+		last  = block_sparsity.get_rowstart_indices()[block_index.second+1];
+	      for (unsigned int j=first; j<last; ++j)
+						 // end of row reached?
+		if (block_sparsity.get_column_numbers()[j] == SparsityPattern::invalid_entry)
+		  {
+						     // nothing more
+						     // to do
+		    break;
+		  }
+		else
+		  {
+		    const unsigned int global_col
+		      = index_mapping.local_to_global(block_col,
+						      block_sparsity.get_column_numbers()[j]);
+		    
+		    if (distribute[global_col] != -1)
+						       // distribute entry at regular
+						       // row @p{row} and irregular column
+						       // global_col
+		      {
+			for (unsigned int q=0;
+			     q!=lines[distribute[global_col]]
+					    .entries.size(); ++q)
+			  sparsity.add (row,
+					lines[distribute[global_col]].entries[q].first);
+		      };
+		  };
+	    };
+	}
+      else
+	{
+					   // row must be
+					   // distributed. split the
+					   // whole row into the
+					   // chunks defined by the
+					   // blocks
+	  for (unsigned int block_col=0; block_col<n_blocks; ++block_col)
+	    {
+	      const SparsityPattern &
+		block_sparsity = sparsity.block(block_row,block_col);
+	      
+	      const unsigned int
+		first = block_sparsity.get_rowstart_indices()[block_index.second],
+		last  = block_sparsity.get_rowstart_indices()[block_index.second+1];
+      
+	      for (unsigned int j=first; j<last; ++j)
+						 // end of row reached?
+		if (block_sparsity.get_column_numbers()[j] == SparsityPattern::invalid_entry)
+		  {
+						     // nothing more to do
+		    break;
+		  }
+		else
+		  {
+		    const unsigned int global_col
+		      = index_mapping.local_to_global (block_col,
+						       block_sparsity.get_column_numbers()[j]);
+		    
+		    if (distribute[global_col] == -1)
+						       // distribute entry at irregular
+						       // row @p{row} and regular column
+						       // global_col.
+		      {
+			for (unsigned int q=0; q!=lines[distribute[row]].entries.size(); ++q) 
+			  sparsity.add (lines[distribute[row]].entries[q].first,
+					global_col);
+		      }
+		    else
+						       // distribute entry at irregular
+						       // row @p{row} and irregular column
+						       // @p{global_col}
+		      {
+			for (unsigned int p=0; p!=lines[distribute[row]].entries.size(); ++p)
+			  for (unsigned int q=0; q!=lines[distribute[global_col]].entries.size(); ++q)
+			    sparsity.add (lines[distribute[row]].entries[p].first,
+					  lines[distribute[global_col]].entries[q].first);
+		      };
+		  };
+	    };
+	};
+    };
+  
+  sparsity.compress();
+};
+
+
+
 unsigned int ConstraintMatrix::n_constraints () const
 {
   return lines.size();
@@ -535,11 +678,3 @@ matrix_functions_2;
 #undef MatrixType
 
 
-#define MatrixType BlockSparsityPattern<2,2>
-matrix_functions_2;
-#undef MatrixType
-
-#define MatrixType BlockSparsityPattern<3,3>
-matrix_functions_2;
-#undef MatrixType
-
diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc
index 9c6aa0b71d..46ab22e16c 100644
--- a/deal.II/deal.II/source/dofs/dof_tools.cc
+++ b/deal.II/deal.II/source/dofs/dof_tools.cc
@@ -1692,10 +1692,7 @@ DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
 
 template void
 DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension> &dof,
-				 BlockSparsityPattern<2,2>    &sparsity);
-template void
-DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension> &dof,
-				 BlockSparsityPattern<3,3>    &sparsity);
+				 BlockSparsityPattern                &sparsity);
 
 template void 
 DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
@@ -1705,12 +1702,7 @@ DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
 template void 
 DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
 				 const vector<vector<bool> > &mask,
-				 BlockSparsityPattern<2,2>   &sparsity);
-
-template void 
-DoFTools::make_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
-				 const vector<vector<bool> > &mask,
-				 BlockSparsityPattern<3,3>   &sparsity);
+				 BlockSparsityPattern        &sparsity);
 
 template
 void
diff --git a/deal.II/deal.II/source/numerics/matrices.cc b/deal.II/deal.II/source/numerics/matrices.cc
index d0dd1a726e..9b9172c178 100644
--- a/deal.II/deal.II/source/numerics/matrices.cc
+++ b/deal.II/deal.II/source/numerics/matrices.cc
@@ -719,7 +719,7 @@ MatrixTools<dim>::apply_boundary_values (const map<unsigned int,double> &boundar
   map<unsigned int,double>::const_iterator  dof  = boundary_values.begin(),
 					    endd = boundary_values.end();
   const unsigned int n_dofs = matrix.m();
-  const BlockSparsityPattern<blocks,blocks> &
+  const BlockSparsityPattern &
     sparsity_pattern = matrix.get_sparsity_pattern();
 
 				   // if a diagonal entry is zero
diff --git a/deal.II/lac/include/lac/block_indices.h b/deal.II/lac/include/lac/block_indices.h
index 2baeaff2ff..a48c81800f 100644
--- a/deal.II/lac/include/lac/block_indices.h
+++ b/deal.II/lac/include/lac/block_indices.h
@@ -47,6 +47,14 @@ class BlockIndices
 				      * size of the vector
 				      */
     BlockIndices (const vector<unsigned int> &n);
+
+				     /**
+				      * Reinitialize the number of
+				      * blocks and assign each block
+				      * the same number of elements.
+				      */
+    void reinit (const unsigned int n_blocks,
+		 const unsigned int n_elements_per_block);
     
 				     /**
 				      * Reinitialize the number of
@@ -152,6 +160,17 @@ BlockIndices::BlockIndices (const vector<unsigned int> &n)
 
 
 
+inline
+void
+BlockIndices::reinit (const unsigned int n_blocks,
+		      const unsigned int n_elements_per_block)
+{
+  const vector<unsigned int> v(n_blocks, n_elements_per_block);
+  reinit (v);
+};
+
+
+
 inline
 void
 BlockIndices::reinit (const vector<unsigned int> &n)
diff --git a/deal.II/lac/include/lac/block_sparse_matrix.h b/deal.II/lac/include/lac/block_sparse_matrix.h
index 94962632c8..a1b70e92a2 100644
--- a/deal.II/lac/include/lac/block_sparse_matrix.h
+++ b/deal.II/lac/include/lac/block_sparse_matrix.h
@@ -97,7 +97,7 @@ class BlockSparseMatrix : public Subscriptor
 				      * long as @p{reinit} is not called
 				      * with a new sparsity structure.
 				      */
-    BlockSparseMatrix (const BlockSparsityPattern<rows,columns> &sparsity);
+    BlockSparseMatrix (const BlockSparsityPattern &sparsity);
 
     
 
@@ -144,7 +144,7 @@ class BlockSparseMatrix : public Subscriptor
 				      * The elements of the matrix are
 				      * set to zero by this function.
 				      */
-    virtual void reinit (const BlockSparsityPattern<rows,columns> &sparsity);
+    virtual void reinit (const BlockSparsityPattern &sparsity);
 
     
 				     /**
@@ -178,17 +178,17 @@ class BlockSparseMatrix : public Subscriptor
 				      */
     virtual void clear ();
     
-				   /**
-				    * Return the number of blocks in a
-				    * column.
-				    */
-  unsigned int n_block_rows () const;
-  
-				   /**
-				    * Return the number of blocks in a
-				    * row.
-				    */
-  unsigned int n_block_cols () const;
+				     /**
+				      * Return the number of blocks in a
+				      * column.
+				      */
+    unsigned int n_block_rows () const;
+    
+				     /**
+				      * Return the number of blocks in a
+				      * row.
+				      */
+    unsigned int n_block_cols () const;
   
 				     /**
 				      * Return whether the object is
@@ -434,7 +434,7 @@ class BlockSparseMatrix : public Subscriptor
 				      * function of objects which
 				      * operate on it.
 				      */
-    const BlockSparsityPattern<rows,columns> &
+    const BlockSparsityPattern &
     get_sparsity_pattern () const;
 
 				     /**
@@ -452,7 +452,7 @@ class BlockSparseMatrix : public Subscriptor
 				      * it using the @p{SmartPointer}
 				      * class.
 				      */
-    SmartPointer<const BlockSparsityPattern<rows,columns> > sparsity_pattern;
+    SmartPointer<const BlockSparsityPattern> sparsity_pattern;
 
 				     /**
 				      * Array of sub-matrices.
diff --git a/deal.II/lac/include/lac/block_sparse_matrix.templates.h b/deal.II/lac/include/lac/block_sparse_matrix.templates.h
index b5dfc66abe..e7cb56603c 100644
--- a/deal.II/lac/include/lac/block_sparse_matrix.templates.h
+++ b/deal.II/lac/include/lac/block_sparse_matrix.templates.h
@@ -27,7 +27,7 @@ BlockSparseMatrix<number,rows,columns>::BlockSparseMatrix () :
 
 template <typename number, int  rows, int columns>
 BlockSparseMatrix<number,rows,columns>::
-BlockSparseMatrix (const BlockSparsityPattern<rows,columns> &sparsity)
+BlockSparseMatrix (const BlockSparsityPattern &sparsity)
 {
   reinit (sparsity);
 };
@@ -65,7 +65,7 @@ BlockSparseMatrix<number,rows,columns>::reinit ()
 template <typename number, int  rows, int columns>
 void
 BlockSparseMatrix<number,rows,columns>::
-reinit (const BlockSparsityPattern<rows,columns> &sparsity)
+reinit (const BlockSparsityPattern &sparsity)
 {
   sparsity_pattern = &sparsity;
   
@@ -112,7 +112,7 @@ BlockSparseMatrix<number,rows,columns>::n_nonzero_elements () const
 
 
 template <typename number, int  rows, int columns>
-const BlockSparsityPattern<rows,columns> &
+const BlockSparsityPattern &
 BlockSparseMatrix<number,rows,columns>::get_sparsity_pattern () const
 {
   return *sparsity_pattern;
diff --git a/deal.II/lac/include/lac/block_sparsity_pattern.h b/deal.II/lac/include/lac/block_sparsity_pattern.h
index e072ff269e..029cbbb2a8 100644
--- a/deal.II/lac/include/lac/block_sparsity_pattern.h
+++ b/deal.II/lac/include/lac/block_sparsity_pattern.h
@@ -16,6 +16,7 @@
 
 #include <base/exceptions.h>
 #include <base/subscriptor.h>
+#include <base/smartpointer.h>
 #include <lac/sparsity_pattern.h>
 #include <lac/block_indices.h>
 
@@ -28,12 +29,11 @@
  * calls to its member functions to calls to the respective member
  * functions of the member sparsity patterns.
  *
- * The largest difference between the @p{SparsityPattern} class and this
- * class is probably the absence of several of the constructors as
- * well as the absenceof the @p{reinit} functions. The reason is that
- * mostly, the matrices have different properties and you will want to
- * initialize the blocks making up the matrix separately. You can
- * access the different blocks using the @p{block(row,col)} function.
+ * The largest difference between the @p{SparsityPattern} class and
+ * this class is that mostly, the matrices have different properties
+ * and you will want to work on the blocks making up the matrix rather
+ * than the whole matrix. You can access the different blocks using
+ * the @p{block(row,col)} function.
  *
  * Attention: this object is not automatically notified if the size of
  * one of its subobjects' size is changed. After you initialize the
@@ -57,7 +57,6 @@
  *
  * @author Wolfgang Bangerth, 2000
  */
-template <int rows, int columns=rows>
 class BlockSparsityPattern : public Subscriptor
 {
   public:
@@ -88,10 +87,70 @@ class BlockSparsityPattern : public Subscriptor
     BlockSparsityPattern ();
 
 				     /**
-				      * Copy operator. For this the same holds
-				      * as for the copy constructor: it is
-				      * declared, defined and fine to be called,
-				      * but the latter only for empty objects.
+				      * Initialize the matrix with the
+				      * given number of block rows and
+				      * columns. The blocks themselves
+				      * are still empty, and you have
+				      * to call @p{collect_args} after
+				      * you assign them sizes.
+				      */
+    BlockSparsityPattern (const unsigned int n_rows,
+			  const unsigned int n_columns);
+
+				     /**
+				      * Copy constructor. This
+				      * constructor is only allowed to
+				      * be called if the sparsity pattern to be
+				      * copied is empty, i.e. there
+				      * are no block allocated at
+				      * present. This is for the same
+				      * reason as for the
+				      * @p{SparsityPattern}, see there
+				      * for the details.
+				      */
+    BlockSparsityPattern (const BlockSparsityPattern &bsp);
+
+				     /**
+				      * Destructor.
+				      */
+    ~BlockSparsityPattern ();
+    
+				     /**
+				      * Resize the matrix. This
+				      * deletes all previously
+				      * existant blocks and replaces
+				      * them by unitialized ones,
+				      * i.e. ones for which also the
+				      * sizes are not yet set. You
+				      * have to do that by calling the
+				      * @p{reinit} functions of the
+				      * blocks themselves. Do not
+				      * forget to call
+				      * @p{collect_size} after that on
+				      * this object.
+				      *
+				      * The reason that you have to
+				      * set sizes of the blocks
+				      * yourself is that the sizes may
+				      * be varying, the maximum number
+				      * of elements per row may be
+				      * varying, etc. It is simpler
+				      * not to reproduce the interface
+				      * of the @p{SparsityPattern}
+				      * class here but rather let the
+				      * user call whatever function
+				      * she desires.
+				      */
+    void reinit (const unsigned int n_rows,
+		 const unsigned int n_columns);
+    
+				     /**
+				      * Copy operator. For this the
+				      * same holds as for the copy
+				      * constructor: it is declared,
+				      * defined and fine to be called,
+				      * but the latter only for empty
+				      * objects.
 				      */
     BlockSparsityPattern & operator = (const BlockSparsityPattern &);
 
@@ -154,17 +213,17 @@ class BlockSparsityPattern : public Subscriptor
 				      */
     void compress ();
 
-				   /**
-				    * Return the number of blocks in a
-				    * column.
-				    */
-  unsigned int n_block_rows () const;
-  
-				   /**
-				    * Return the number of blocks in a
-				    * row.
-				    */
-  unsigned int n_block_cols () const;
+				     /**
+				      * Return the number of blocks in a
+				      * column.
+				      */
+    unsigned int n_block_rows () const;
+    
+				     /**
+				      * Return the number of blocks in a
+				      * row.
+				      */
+    unsigned int n_block_cols () const;
   
 				     /**
 				      * Return whether the object is
@@ -202,7 +261,7 @@ class BlockSparsityPattern : public Subscriptor
 				      * and then relays to that block.
 				      */
     void add (const unsigned int i, const unsigned int j);
-	
+    
 				     /**
 				      * Return number of rows of this
 				      * matrix, which equals the
@@ -272,13 +331,30 @@ class BlockSparsityPattern : public Subscriptor
 		    int, int, int, int,
 		    << "The blocks [" << arg1 << ',' << arg2 << "] and ["
 		    << arg3 << ',' << arg4 << "] have differing column numbers.");
-    
+				     /**
+				      * Exception
+				      */
+    DeclException2 (ExcIncompatibleSizes,
+		    int, int,
+		    << "The number of blocks " << arg1 << " and " << arg2
+		    << " are different.");
     
   private:
+
+				     /**
+				      * Number of block rows.
+				      */
+    unsigned int rows;
+
+				     /**
+				      * Number of block columns.
+				      */
+    unsigned int columns;
+    
 				     /**
 				      * Array of sparsity patterns.
 				      */
-    SparsityPattern sub_objects[rows][columns];
+    vector<vector<SmartPointer<SparsityPattern> > > sub_objects;
 
 				     /**
 				      * Object storing and managing
@@ -312,53 +388,48 @@ class BlockSparsityPattern : public Subscriptor
 
 
 
-template <int rows, int columns>
 inline
 SparsityPattern &
-BlockSparsityPattern<rows,columns>::block (const unsigned int row,
-					   const unsigned int column)
+BlockSparsityPattern::block (const unsigned int row,
+			     const unsigned int column)
 {
-  return sub_objects[row][column];
+  return *sub_objects[row][column];
 };
 
 
 
-template <int rows, int columns>
 inline
 const SparsityPattern &
-BlockSparsityPattern<rows,columns>::block (const unsigned int row,
-					   const unsigned int column) const
+BlockSparsityPattern::block (const unsigned int row,
+			     const unsigned int column) const
 {
-  return sub_objects[row][column];
+  return *sub_objects[row][column];
 };
 
 
 
-template <int rows, int columns>
 inline
 const BlockIndices &
-BlockSparsityPattern<rows,columns>::get_row_indices () const
+BlockSparsityPattern::get_row_indices () const
 {
   return row_indices;
 };
 
 
 
-template <int rows, int columns>
 inline
 const BlockIndices &
-BlockSparsityPattern<rows,columns>::get_column_indices () const
+BlockSparsityPattern::get_column_indices () const
 {
   return column_indices;
 };
 
 
 
-template <int rows, int columns>
 inline
 void
-BlockSparsityPattern<rows,columns>::add (const unsigned int i,
-					 const unsigned int j)
+BlockSparsityPattern::add (const unsigned int i,
+			   const unsigned int j)
 {
 				   // if you get an error here, are
 				   // you sure you called
@@ -366,26 +437,24 @@ BlockSparsityPattern<rows,columns>::add (const unsigned int i,
   const pair<unsigned int,unsigned int>
     row_index = row_indices.global_to_local (i),
     col_index = column_indices.global_to_local (j);
-  sub_objects[row_index.first][col_index.first].add (row_index.second,
-						     col_index.second);
+  sub_objects[row_index.first][col_index.first]->add (row_index.second,
+						      col_index.second);
 };
 
 
 
-template <int rows, int columns>
 inline
 unsigned int
-BlockSparsityPattern<rows,columns>::n_block_cols () const
+BlockSparsityPattern::n_block_cols () const
 {
   return columns;
 }
 
 
 
-template <int rows, int columns>
 inline
 unsigned int
-BlockSparsityPattern<rows,columns>::n_block_rows () const
+BlockSparsityPattern::n_block_rows () const
 {
   return rows;
 }
diff --git a/deal.II/lac/include/lac/block_sparsity_pattern.templates.h b/deal.II/lac/include/lac/block_sparsity_pattern.templates.h
deleted file mode 100644
index c878cd3ccf..0000000000
--- a/deal.II/lac/include/lac/block_sparsity_pattern.templates.h
+++ /dev/null
@@ -1,184 +0,0 @@
-//----------------------------  block_sparsity_pattern.templates.h  ---------------------------
-//    $Id$
-//    Version: $Name$
-//
-//    Copyright (C) 1998, 1999, 2000 by the deal.II authors
-//
-//    This file is subject to QPL and may not be  distributed
-//    without copyright and license information. Please refer
-//    to the file deal.II/doc/license.html for the  text  and
-//    further information on this license.
-//
-//----------------------------  block_sparsity_pattern.templates.h  ---------------------------
-#ifndef __deal2__block_sparsity_pattern_templates_h
-#define __deal2__block_sparsity_pattern_templates_h
-
-
-#include <lac/block_sparsity_pattern.h>
-
-
-template <int rows, int columns>
-BlockSparsityPattern<rows,columns>::BlockSparsityPattern () 
-		: row_indices (rows), column_indices(columns)
-{
-  Assert (rows>0,    ExcInvalidSize (rows));
-  Assert (columns>0, ExcInvalidSize (columns));
-};
-
-  
-
-template <int rows, int columns>
-BlockSparsityPattern<rows,columns> &
-BlockSparsityPattern<rows,columns>::operator = (const BlockSparsityPattern<rows,columns> &bsp)
-{
-				   // copy objects
-  for (unsigned int i=0; i<rows; ++i)
-    for (unsigned int j=0; j<columns; ++j)
-      sub_objects[i][j] = bsp.sub_objects[i][j];
-				   // update index objects
-  collect_sizes ();
-
-  return *this;
-};
-
-  
-
-
-template <int rows, int columns>
-void
-BlockSparsityPattern<rows,columns>::collect_sizes ()
-{
-  vector<unsigned int> row_sizes (rows);
-  vector<unsigned int> col_sizes (columns);
-
-				   // first find out the row sizes
-				   // from the first block column
-  for (unsigned int r=0; r<rows; ++r)
-    row_sizes[r] = sub_objects[r][0].n_rows();
-				   // then check that the following
-				   // block columns have the same
-				   // sizes
-  for (unsigned int c=1; c<columns; ++c)
-    for (unsigned int r=0; r<rows; ++r)
-      Assert (row_sizes[r] == sub_objects[r][c].n_rows(),
-	      ExcIncompatibleRowNumbers (r,0,r,c));
-
-				   // finally initialize the row
-				   // indices with this array
-  row_indices.reinit (row_sizes);
-  
-  
-				   // then do the same with the columns
-  for (unsigned int c=0; c<columns; ++c)
-    col_sizes[c] = sub_objects[0][c].n_cols();
-  for (unsigned int r=1; r<rows; ++r)
-    for (unsigned int c=0; c<columns; ++c)
-      Assert (col_sizes[c] == sub_objects[r][c].n_cols(),
-	      ExcIncompatibleRowNumbers (0,c,r,c));
-
-				   // finally initialize the row
-				   // indices with this array
-  column_indices.reinit (col_sizes);
-};
-
-
-
-template <int rows, int columns>
-void
-BlockSparsityPattern<rows,columns>::compress ()
-{
-  for (unsigned int i=0; i<rows; ++i)
-    for (unsigned int j=0; j<columns; ++j)
-      sub_objects[i][j].compress ();
-};
-
-
-
-template <int rows, int columns>
-bool
-BlockSparsityPattern<rows,columns>::empty () const
-{
-  for (unsigned int i=0; i<rows; ++i)
-    for (unsigned int j=0; j<columns; ++j)
-      if (sub_objects[i][j].empty () == false)
-	return false;
-  return true;
-};
-
-
-
-template <int rows, int columns>
-unsigned int
-BlockSparsityPattern<rows,columns>::max_entries_per_row () const
-{
-  unsigned int max_entries = 0;
-  for (unsigned int block_row=0; block_row<rows; ++block_row)
-    {
-      unsigned int this_row = 0;
-      for (unsigned int c=0; c<columns; ++c)
-	this_row += sub_objects[block_row][c].max_entries_per_row ();
-
-      if (this_row > max_entries)
-	max_entries = this_row;
-    };
-  return max_entries;
-};
-
-
-template <int rows, int columns>
-unsigned int
-BlockSparsityPattern<rows,columns>::n_rows () const
-{
-				   // only count in first column, since
-				   // all rows should be equivalent
-  unsigned int count = 0;
-  for (unsigned int r=0; r<rows; ++r)
-    count += sub_objects[r][0].n_rows();
-  return count;
-};
-
-
-
-template <int rows, int columns>
-unsigned int
-BlockSparsityPattern<rows,columns>::n_cols () const
-{
-				   // only count in first row, since
-				   // all rows should be equivalent
-  unsigned int count = 0;
-  for (unsigned int c=0; c<columns; ++c)
-    count += sub_objects[0][c].n_cols();
-  return count;
-};
-
-
-
-
-
-template <int rows, int columns>
-unsigned int
-BlockSparsityPattern<rows,columns>::n_nonzero_elements () const
-{
-  unsigned int count = 0;
-  for (unsigned int i=0; i<rows; ++i)
-    for (unsigned int j=0; j<columns; ++j)
-      count += sub_objects[i][j].n_nonzero_elements ();
-  return count;
-};
-
-
-
-template <int rows, int columns>
-bool
-BlockSparsityPattern<rows,columns>::is_compressed () const
-{
-  for (unsigned int i=0; i<rows; ++i)
-    for (unsigned int j=0; j<columns; ++j)
-      if (sub_objects[i][j].is_compressed () == false)
-	return false;
-  return true;
-};
-
-
-
-#endif // ifdef block_sparsity_pattern_templates_h
diff --git a/deal.II/lac/include/lac/block_vector.h b/deal.II/lac/include/lac/block_vector.h
index 4f24700c39..c7fce8a4c4 100644
--- a/deal.II/lac/include/lac/block_vector.h
+++ b/deal.II/lac/include/lac/block_vector.h
@@ -130,17 +130,17 @@ class BlockVector
 				      */
     ~BlockVector ();
 
-				   /**
-				    * Reinitialize the BlockVector to
-				    * contain @p{num_blocks} blocks of
-				    * size @p{block_size} each.
-				    *
-				    * If @p{fast==false}, the vector
-				    * is filled with zeros.
-				    */
-  void reinit (const unsigned int num_blocks,
-	       const unsigned int block_size,
-	       const bool fast = false);
+				     /**
+				      * Reinitialize the BlockVector to
+				      * contain @p{num_blocks} blocks of
+				      * size @p{block_size} each.
+				      *
+				      * If @p{fast==false}, the vector
+				      * is filled with zeros.
+				      */
+    void reinit (const unsigned int num_blocks,
+		 const unsigned int block_size,
+		 const bool fast = false);
   
 				     /**
 				      * Reinitialize the BlockVector
diff --git a/deal.II/lac/source/block_sparsity_pattern.cc b/deal.II/lac/source/block_sparsity_pattern.cc
index e886af0c12..452e12e1e1 100644
--- a/deal.II/lac/source/block_sparsity_pattern.cc
+++ b/deal.II/lac/source/block_sparsity_pattern.cc
@@ -11,11 +11,206 @@
 //
 //----------------------------  block_block_sparsity_pattern.cc  ---------------------------
 
-#include <lac/block_sparsity_pattern.templates.h>
+#include <lac/block_sparsity_pattern.h>
 
-// explicit instantiations
-template class BlockSparsityPattern<2,2>;
-template class BlockSparsityPattern<2,3>;
 
-template class BlockSparsityPattern<3,2>;
-template class BlockSparsityPattern<3,3>;
+
+BlockSparsityPattern::BlockSparsityPattern ()
+		:
+		rows (0),
+		columns (0)
+{};
+
+
+
+BlockSparsityPattern::BlockSparsityPattern (const unsigned int r,
+					    const unsigned int c)
+		:
+		rows (0),
+		columns (0)
+{
+  reinit (r,c);
+};
+
+
+
+BlockSparsityPattern::~BlockSparsityPattern ()
+{
+				   // clear all memory
+  reinit (0,0);
+};
+
+
+
+void
+BlockSparsityPattern::reinit (const unsigned int r,
+			      const unsigned int c)
+{
+				   // delete previous content
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      {
+	SparsityPattern * sp = sub_objects[i][j];
+	sub_objects[i][j] = 0;
+	delete sp;
+      };
+  sub_objects.clear ();
+
+				   // set new sizes
+  rows = r;
+  columns = c;
+  sub_objects.resize (rows, vector<SmartPointer<SparsityPattern> > (columns));
+
+				   // allocate new objects
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      sub_objects[i][j] = new SparsityPattern;
+};
+
+
+
+BlockSparsityPattern &
+BlockSparsityPattern::operator = (const BlockSparsityPattern &bsp)
+{
+  Assert (rows == bsp.rows, ExcIncompatibleSizes(rows, bsp.rows));
+  Assert (columns == bsp.columns, ExcIncompatibleSizes(columns, bsp.columns));
+				   // copy objects
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      *sub_objects[i][j] = *bsp.sub_objects[i][j];
+				   // update index objects
+  collect_sizes ();
+
+  return *this;
+};
+
+  
+
+
+void
+BlockSparsityPattern::collect_sizes ()
+{
+  vector<unsigned int> row_sizes (rows);
+  vector<unsigned int> col_sizes (columns);
+
+				   // first find out the row sizes
+				   // from the first block column
+  for (unsigned int r=0; r<rows; ++r)
+    row_sizes[r] = sub_objects[r][0]->n_rows();
+				   // then check that the following
+				   // block columns have the same
+				   // sizes
+  for (unsigned int c=1; c<columns; ++c)
+    for (unsigned int r=0; r<rows; ++r)
+      Assert (row_sizes[r] == sub_objects[r][c]->n_rows(),
+	      ExcIncompatibleRowNumbers (r,0,r,c));
+
+				   // finally initialize the row
+				   // indices with this array
+  row_indices.reinit (row_sizes);
+  
+  
+				   // then do the same with the columns
+  for (unsigned int c=0; c<columns; ++c)
+    col_sizes[c] = sub_objects[0][c]->n_cols();
+  for (unsigned int r=1; r<rows; ++r)
+    for (unsigned int c=0; c<columns; ++c)
+      Assert (col_sizes[c] == sub_objects[r][c]->n_cols(),
+	      ExcIncompatibleRowNumbers (0,c,r,c));
+
+				   // finally initialize the row
+				   // indices with this array
+  column_indices.reinit (col_sizes);
+};
+
+
+
+void
+BlockSparsityPattern::compress ()
+{
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      sub_objects[i][j]->compress ();
+};
+
+
+
+bool
+BlockSparsityPattern::empty () const
+{
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      if (sub_objects[i][j]->empty () == false)
+	return false;
+  return true;
+};
+
+
+
+unsigned int
+BlockSparsityPattern::max_entries_per_row () const
+{
+  unsigned int max_entries = 0;
+  for (unsigned int block_row=0; block_row<rows; ++block_row)
+    {
+      unsigned int this_row = 0;
+      for (unsigned int c=0; c<columns; ++c)
+	this_row += sub_objects[block_row][c]->max_entries_per_row ();
+
+      if (this_row > max_entries)
+	max_entries = this_row;
+    };
+  return max_entries;
+};
+
+
+unsigned int
+BlockSparsityPattern::n_rows () const
+{
+				   // only count in first column, since
+				   // all rows should be equivalent
+  unsigned int count = 0;
+  for (unsigned int r=0; r<rows; ++r)
+    count += sub_objects[r][0]->n_rows();
+  return count;
+};
+
+
+
+unsigned int
+BlockSparsityPattern::n_cols () const
+{
+				   // only count in first row, since
+				   // all rows should be equivalent
+  unsigned int count = 0;
+  for (unsigned int c=0; c<columns; ++c)
+    count += sub_objects[0][c]->n_cols();
+  return count;
+};
+
+
+
+
+
+unsigned int
+BlockSparsityPattern::n_nonzero_elements () const
+{
+  unsigned int count = 0;
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      count += sub_objects[i][j]->n_nonzero_elements ();
+  return count;
+};
+
+
+
+bool
+BlockSparsityPattern::is_compressed () const
+{
+  for (unsigned int i=0; i<rows; ++i)
+    for (unsigned int j=0; j<columns; ++j)
+      if (sub_objects[i][j]->is_compressed () == false)
+	return false;
+  return true;
+};
+
diff --git a/tests/deal.II/block_matrices.cc b/tests/deal.II/block_matrices.cc
index 3f1852197f..94aabdf21a 100644
--- a/tests/deal.II/block_matrices.cc
+++ b/tests/deal.II/block_matrices.cc
@@ -60,7 +60,7 @@ class LaplaceProblem
     DoFHandler<2>        dof_handler;
 
     Sparsity      sparsity_pattern;
-    Matrix system_matrix;
+    Matrix        system_matrix;
 
     Vector       system_rhs;
 };
@@ -69,9 +69,27 @@ class LaplaceProblem
 template <class Vector, class Matrix, class Sparsity>
 LaplaceProblem<Vector,Matrix,Sparsity>::LaplaceProblem () :
 		dof_handler (triangulation)
+{
+  sparsity_pattern.reinit (system_matrix.n_block_rows(),
+			   system_matrix.n_block_cols());
+};
+
+
+
+template <>
+LaplaceProblem<Vector<float>,SparseMatrix<float>,SparsityPattern>::LaplaceProblem () :
+		dof_handler (triangulation)
 {};
 
 
+
+template <>
+LaplaceProblem<Vector<double>,SparseMatrix<double>,SparsityPattern>::LaplaceProblem () :
+		dof_handler (triangulation)
+{};
+
+
+
 template <class Vector, class Matrix, class Sparsity>
 void LaplaceProblem<Vector,Matrix,Sparsity>::make_grid_and_dofs ()
 {
@@ -138,7 +156,7 @@ void LaplaceProblem<Vector<float>,SparseMatrix<float>,SparsityPattern>::reinit_v
 
 
 template <>
-void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSparsityPattern<2,2> >::reinit_sparsity () 
+void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSparsityPattern>::reinit_sparsity () 
 {
   const unsigned int n_dofs = dof_handler.n_dofs();
   const unsigned int block_size[2] = { n_dofs/3, n_dofs - n_dofs/3 };
@@ -153,7 +171,7 @@ void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSpars
 
 
 template <>
-void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSparsityPattern<2,2> >::reinit_vectors () 
+void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSparsityPattern>::reinit_vectors () 
 {
   const unsigned int n_dofs = dof_handler.n_dofs();
   const unsigned int block_size_[2] = { n_dofs/3, n_dofs - n_dofs/3 };
@@ -167,7 +185,7 @@ void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSpars
 
 
 template <>
-void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSparsityPattern<3,3> >::reinit_sparsity () 
+void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSparsityPattern>::reinit_sparsity () 
 {
   const unsigned int n_dofs = dof_handler.n_dofs();
   const unsigned int block_size[3] = { n_dofs/5, n_dofs/7, n_dofs - n_dofs/5 - n_dofs/7 };
@@ -182,7 +200,7 @@ void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSpars
 
 
 template <>
-void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSparsityPattern<3,3> >::reinit_vectors () 
+void LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSparsityPattern>::reinit_vectors () 
 {
   const unsigned int n_dofs = dof_handler.n_dofs();
   const unsigned int block_size_[3] = { n_dofs/5, n_dofs/7, n_dofs - n_dofs/5 - n_dofs/7 };
@@ -330,7 +348,7 @@ int main ()
   
   if (true)
     {
-      LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSparsityPattern<2,2> >
+      LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,2,2>,BlockSparsityPattern>
 	laplace_problem;  
       laplace_problem.run ();
       
@@ -342,7 +360,7 @@ int main ()
 
   if (true)
     {
-      LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSparsityPattern<3,3> >
+      LaplaceProblem<BlockVector<double>,BlockSparseMatrix<double,3,3>,BlockSparsityPattern>
 	laplace_problem;  
       laplace_problem.run ();
       
-- 
2.39.5