]> https://gitweb.dealii.org/ - dealii.git/commitdiff
SparsityPattern::matrix_position
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 13 Feb 2002 09:27:49 +0000 (09:27 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 13 Feb 2002 09:27:49 +0000 (09:27 +0000)
git-svn-id: https://svn.dealii.org/trunk@5513 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/2002/c-3-3.html
deal.II/lac/include/lac/sparsity_pattern.h
deal.II/lac/source/sparsity_pattern.cc
tests/lac/sparsity_pattern.cc

index 0ca59afb35aa8a6f9208d650287d59bc35a8699a..ae2c8335c9362d116d5b9c7a267e364b96e1990d 100644 (file)
@@ -26,51 +26,57 @@ documentation, etc</a>.
 <h3>General</h3>
 
 <ol>
-  <li> <p> Changed: The build system was entirely revised. Object
-  files in debug more now have the suffix <code>.g.o</code> instead of
-  <code>.go</code>. All object files from the subdirectories are now
-  placed into the <code>/lib</code> top-level directory, rather than
-  in library directories in the individual subdirs.
-  <br>
-  (WB 2002/02/11)
-  </p>
+  <li> <p> 
+       Changed: The build system was entirely revised. Object
+       files in debug more now have the suffix <code>.g.o</code>
+       instead of <code>.go</code>. All object files from the
+       subdirectories are now placed into the <code>/lib</code>
+       top-level directory, rather than in library directories in the
+       individual subdirs.
+       <br>
+       (WB 2002/02/11)
+       </p>
 </ol>
 
 
 
 <a name="base"></a>
 <h3>base</h3>
-
-  <li> <p> Changed: The quite logorrhoeic function name <code
-  class="class">TensorProductPolynomials</code>::<code
-  class="member">n_tensor_product_polynomials</code> was changed to
-  <code class="member">n</code> to be compliant wth the new class <code
-  class="class">PolynomialSpace</code>.
-  <br>
-  (GK 2002/02/11)
-  </p>
-
-  <li> <p> New: The class <code class="class">PolynomialSpace</code>
-  implements the space of polynomials at most a certain degree in
-  arbitrary space dimensions.
-  <br>
-  (GK 2002/02/11)
-  </p>
-
-  <li> <p> New: The function <code class="class">DataOutBase</code>::
-  <code class="member">write_tecplot_binary</code> has been added.  This
-  function will write Tecplot binary files if the Tecplot API is detected
-  by ./configure.  To use this feature be sure that the environment
-  variable TECHOME points to a valid Tecplot installation and that the files
-  $TECHOME/include/TECIO.h and $TECHOME/lib/tecio.a exist.  The name of the
-  file to be written is specified through the <code class="class">DataOutBase</code>
-  ::<code class="member">TecplotFlags</code>.
-  <code class="member">tecplot_binary_file_name</code> variable. If the API is not
-  available this code simply calls the existing ASCII output function.
-  <br>
-  (BK 2002/02/11)
-  </p>
-
+  <li> <p> 
+       Changed: The quite logorrhoeic function name <code
+       class="class">TensorProductPolynomials</code>::<code
+       class="member">n_tensor_product_polynomials</code> was changed to
+       <code class="member">n</code> to be compliant wth the new class <code
+       class="class">PolynomialSpace</code>.
+       <br>
+       (GK 2002/02/11)
+       </p>
+
+  <li> <p> 
+       New: The class <code class="class">PolynomialSpace</code>
+       implements the space of polynomials at most a certain degree in
+       arbitrary space dimensions.
+       <br>
+       (GK 2002/02/11)
+       </p>
+
+  <li> <p> 
+       New: The function <code class="class">DataOutBase</code>::
+       <code class="member">write_tecplot_binary</code> has been
+       added.  This function will write Tecplot binary files if the
+       Tecplot API is detected by ./configure.  To use this feature be
+       sure that the environment variable TECHOME points to a valid
+       Tecplot installation and that the files
+       $TECHOME/include/TECIO.h and $TECHOME/lib/tecio.a exist.  The
+       name of the file to be written is specified through the <code
+       class="class">DataOutBase</code> ::<code
+       class="member">TecplotFlags</code>.  <code
+       class="member">tecplot_binary_file_name</code> variable. If the
+       API is not available this code simply calls the existing ASCII
+       output function.
+       <br>
+       (BK 2002/02/11)
+       </p>
 <ol>
 </ol>
 
@@ -80,26 +86,35 @@ documentation, etc</a>.
 <h3>lac</h3>
 
 <ol>
-  <li> <p> New: Functions <code
-  class="member">SparsityPattern::copy_from</code> and <code
-  class="member">SparseMatrix::copy_from</code> allow to copy a full
-  matrix into a sparse matrix.
-  <br>
-  (WB 2002/02/06)
-  </p>
+  <li> <p> 
+       New: Function <code
+       class="member">SparsityPattern::matrix_position</code> is the
+       inverse function for <code
+       class="member">SparsityPattern::operator()</code>.
+       <br>
+       (WB 2002/02/13)
+       </p>
+
+  <li> <p> 
+       New: Functions <code
+       class="member">SparsityPattern::copy_from</code> and <code
+       class="member">SparseMatrix::copy_from</code> allow to copy a full
+       matrix into a sparse matrix.
+       <br>
+       (WB 2002/02/06)
+       </p>
 </ol>
 
 
 
 <a name="deal.II"></a>
 <h3>deal.II</h3>
-
-  <li> <p> New: Finite element family with complete polynomial spaces
-  for discontinuous Galerkin: <code class="class">FE_DGP</code>
-  <br>
-  (GK 2002/02/11)
-  </p>
-
+  <li> <p> 
+       New: Finite element family with complete polynomial spaces
+       for discontinuous Galerkin: <code class="class">FE_DGP</code>
+       <br> 
+       (GK 2002/02/11)
+       </p>
 <ol>
 </ol>
 
index 0155c938e52bca09341ee0b208c1c9e25a86403d..f917b9ca7eeb46ca042aa8c78079788eb56dcbbd 100644 (file)
@@ -256,10 +256,12 @@ class SparsityPattern : public Subscriptor
     ~SparsityPattern ();
 
                                     /**
-                                     * 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.
+                                     * 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.
                                      */
     SparsityPattern & operator = (const SparsityPattern &);
     
@@ -510,10 +512,11 @@ class SparsityPattern : public Subscriptor
 
                                     /**
                                      * Return the index of the matrix
-                                     * element with row number @p{i} and
-                                     * column number @p{j}. If the matrix
-                                     * element is not a nonzero one,
-                                     * return @p{SparsityPattern::invalid_entry}.
+                                     * element with row number @p{i}
+                                     * and column number @p{j}. If
+                                     * the matrix element is not a
+                                     * nonzero one, return
+                                     * @p{SparsityPattern::invalid_entry}.
                                      *
                                      * This function is usually called
                                      * by the @p{operator()} of the
@@ -525,10 +528,39 @@ class SparsityPattern : public Subscriptor
                                      * with a binary sort algorithm
                                      * because the column numbers are
                                      * sorted.
+                                     *
+                                     * If @p{m} is the number of
+                                     * entries in @p{row}, then the
+                                     * complexity of this function is
+                                     * @p{log(m)} if the sparsity
+                                     * pattern is compressed.
                                      */
     unsigned int operator() (const unsigned int i, 
                             const unsigned int j) const;
 
+                                    /**
+                                     * This is the inverse operation
+                                     * to @p{operator()}: given a
+                                     * global index, find out row and
+                                     * column of the matrix entry to
+                                     * which it belongs. The returned
+                                     * value is the pair composed of
+                                     * row and column index.
+                                     *
+                                     * This function may only be
+                                     * called if the sparsity pattern
+                                     * is closed. The global index
+                                     * must then be between zero and
+                                     * @p{n_nonzero_elements}.
+                                     *
+                                     * If @p{N} is the number of
+                                     * rows of this matrix, then the
+                                     * complexity of this function is
+                                     * @p{log(N)}.
+                                     */
+    std::pair<unsigned int, unsigned int>
+    matrix_position (const unsigned int global_index) const;
+    
                                     /**
                                      * Add a nonzero entry to the matrix.
                                      * This function may only be called
index 0d26610f400fb8ade89adb428530f849a6ec6fcf..5edda629201d0f2b0b5b657bdee80a6b025ebbe2 100644 (file)
@@ -661,6 +661,38 @@ SparsityPattern::exists (const unsigned int i,
 }
 
 
+
+std::pair<unsigned int, unsigned int>
+SparsityPattern::matrix_position (const unsigned int global_index) const
+{
+  Assert (compressed == true, ExcNotCompressed());
+  Assert (global_index < n_nonzero_elements(),
+         ExcIndexRange (global_index, 0, n_nonzero_elements()));
+
+                                  // first find the row in which the
+                                  // entry is located. for this note
+                                  // that the rowstart array indexes
+                                  // the global indices at which each
+                                  // row starts. since it is sorted,
+                                  // and since there is an element
+                                  // for the one-past-last row, we
+                                  // can simply use a bisection
+                                  // search on it
+  const unsigned int row
+    = (std::upper_bound (&rowstart[0], &rowstart[rows], global_index)
+       - &rowstart[0] - 1);
+
+                                  // now, the column index is simple
+                                  // since that is what the colnums
+                                  // array stores:
+  const unsigned int col = colnums[global_index];
+
+                                  // so return the respective pair
+  return std::make_pair (row,col);
+};
+
+
+
 void
 SparsityPattern::symmetrize () 
 {
index 023540907b1b295bb870c0f30f8d454f98edf4b3..0de153d6794c250a7ded7be2f9bf358ec5bf9ad7 100644 (file)
@@ -96,6 +96,31 @@ main ()
       for (; sp3_p != sp3.get_column_numbers()+sp3.get_rowstart_indices()[row+1]; ++sp3_p, ++sp4_p)
        Assert (*sp3_p == *sp4_p, ExcInternalError());
     };
+
+
+                                  // check the matrix_position
+                                  // function with sparsity patterns
+                                  // sp1 through sp4. the checked
+                                  // function should be the inverse
+                                  // of operator()
+                                  //
+                                  // check inverseness property first
+                                  // forward, then backward
+  for (unsigned int loop=1; loop<=4; ++loop)
+    {
+      const SparsityPattern &
+       sp = (loop==1 ? sp1 : (loop==2 ? sp2 : (loop==3 ? sp3 : sp4)));
+      for (unsigned int i=0; i<sp.n_nonzero_elements(); ++i)
+       Assert (sp(sp.matrix_position(i).first,
+                   sp.matrix_position(i).second) == i,
+               ExcInternalError());
+      for (unsigned int row=0; row<sp.n_rows(); ++row)
+       for (unsigned int col=0; col<sp.n_cols(); ++col)
+         if (sp(row,col) != SparsityPattern::invalid_entry)
+           Assert (sp.matrix_position(sp(row,col)) ==
+                   std::make_pair(row,col),
+                   ExcInternalError());
+    };
 };
 
   

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.