]> https://gitweb.dealii.org/ - dealii.git/commitdiff
avoid compiler warning -unused-variable
authorBenjamin Brands <benjamin.brands@fau.de>
Tue, 6 Feb 2018 15:53:53 +0000 (16:53 +0100)
committerBenjamin Brands <benjamin.brands@fau.de>
Wed, 14 Feb 2018 20:10:16 +0000 (21:10 +0100)
doc/news/changes/minor/20180207BenjaminBrands [new file with mode: 0644]
include/deal.II/lac/scalapack.h
include/deal.II/lac/scalapack.templates.h
source/lac/scalapack.cc
tests/scalapack/scalapack_12.cc
tests/scalapack/scalapack_12.mpirun=1.output
tests/scalapack/scalapack_12.mpirun=11.output
tests/scalapack/scalapack_12.mpirun=9.output

diff --git a/doc/news/changes/minor/20180207BenjaminBrands b/doc/news/changes/minor/20180207BenjaminBrands
new file mode 100644 (file)
index 0000000..bda0349
--- /dev/null
@@ -0,0 +1,3 @@
+New: Added routines to perform addition, multiplication and scaling for ScaLAPACKMatrix
+<br>
+(Benjamin Brands, 2018/02/07)
index 788355e4a11352d7e3913f0616a03c4d89d04f4a..0cc03e9b87f570de4d5f5152ad0f69d2a288a919 100644 (file)
@@ -191,31 +191,23 @@ public:
                const std::pair<unsigned int,unsigned int> &submatrix_size) const;
 
   /**
-   * Transposing assignment: <i>A = B<sup>T</sup></i>
+   * Transposing assignment: $\mathbf{A} = \mathbf{B}^T$
    *
-   * The matrices <tt>A</tt> and <tt>B</tt> must have the same process grid
+   * The matrices $\mathbf{A}$ and $\mathbf{B}$ must have the same process grid.
    *
-   * Following alignment conditions have to be fulfilled: MB_A = NB_B and NB_A = MB_B
+   * The following alignment conditions have to be fulfilled: $MB_A=NB_B$ and $NB_A=MB_B$.
    */
   void copy_transposed(const ScaLAPACKMatrix<NumberType> &B);
 
   /**
-   * Matrix-addition:
-   *
-   * <i>A = a A + b op(B)</i>
-   *
-   * if(transpose_B)  <i>op(B) = B<sup>T</sup></i>
+   * The operations based on the input parameter @p transpose_B and the alignment conditions are summarized in the following table:
    *
-   * else  <i>op(B) = B</i>
+   * | transpose_B |          Block Sizes         |                    Operation                  |
+   * | :---------: | :--------------------------: | :-------------------------------------------: |
+   * |   false     | $MB_A=MB_B$ <br> $NB_A=NB_B$ |  $\mathbf{A} = a \mathbf{A} + b \mathbf{B}$   |
+   * |   true      | $MB_A=NB_B$ <br> $NB_A=MB_B$ | $\mathbf{A} = a \mathbf{A} + b \mathbf{B}^T$  |
    *
-   * The matrices <tt>A</tt> and <tt>B</tt> must have the same process grid
-   *
-   * Following alignment conditions have to be fulfilled:
-   *
-   * | transpose_B |          Block Sizes         |
-   * | :---------: | :--------------------------: |
-   * |   false     | MB_A = MB_B <br> NB_A = NB_B |
-   * |   true      | MB_A = NB_B <br> NB_A = MB_B |
+   * The matrices $\mathbf{A}$ and $\mathbf{B}$ must have the same process grid.
    */
   void add(const ScaLAPACKMatrix<NumberType> &B,
            const NumberType a=0.,
@@ -224,24 +216,22 @@ public:
 
   /**
    * Matrix-addition:
+   * $\mathbf{A} = \mathbf{A} + b \mathbf{B}$
    *
-   * <i>A += b B</i>
-   *
-   * The matrices <tt>A</tt> and <tt>B</tt> must have the same process grid
+   * The matrices $\mathbf{A}$ and $\mathbf{B}$ must have the same process grid.
    *
-   * Following alignment conditions have to be fulfilled: MB_A = MB_B and NB_A = NB_B
+   * The following alignment conditions have to be fulfilled: $MB_A=MB_B$ and $NB_A=NB_B$.
    */
   void add(const NumberType b,
            const ScaLAPACKMatrix<NumberType> &B);
 
   /**
    * Matrix-addition:
+   * $\mathbf{A} = \mathbf{A} + b \mathbf{B}^T$
    *
-   * <i>A += b B<sup>T</sup></i>
+   * The matrices $\mathbf{A}$ and $\mathbf{B}$ must have the same process grid.
    *
-   * The matrices <tt>A</tt> and <tt>B</tt> must have the same process grid
-   *
-   * Following alignment conditions have to be fulfilled: MB_A = NB_B and NB_A = MB_B
+   * The following alignment conditions have to be fulfilled: $MB_A=NB_B$ and $NB_A=MB_B$.
    */
   void Tadd(const NumberType b,
             const ScaLAPACKMatrix<NumberType> &B);
@@ -251,97 +241,97 @@ public:
    *
    * The operations based on the input parameters and the alignment conditions are summarized in the following table:
    *
-   * | transpose_A | transpose_B |                   Block Sizes                 |                     Operation                    |
-   * | :---------: | :---------: | :-------------------------------------------: | :----------------------------------------------: |
-   * | false       |   false     | MB_A = MB_C <br> NB_A = MB_B <br> NB_B = NB_C |             <i>C = a A * B + b C</i>             |
-   * | false       |   true      | MB_A = MB_C <br> NB_A = NB_B <br> MB_B = NB_C |      <i>C = a A * B<sup>T</sup> + b C</i>        |
-   * | true        |   false     | MB_A = MB_B <br> NB_A = MB_C <br> NB_B = NB_C |        <i>C = a A<sup>T</sup> * B + b C</i>      |
-   * | true        |   true      | MB_A = NB_B <br> NB_A = MB_C <br> MB_B = NB_C | <i>C = a A<sup>T</sup> * B<sup>T</sup> + b C</i> |
+   * | transpose_A | transpose_B |                  Block Sizes                  |                             Operation                           |
+   * | :---------: | :---------: | :-------------------------------------------: | :-------------------------------------------------------------: |
+   * | false       |   false     | $MB_A=MB_C$ <br> $NB_A=MB_B$ <br> $NB_B=NB_C$ |   $\mathbf{C} = b \mathbf{A} \cdot \mathbf{B} + c \mathbf{C}$   |
+   * | false       |   true      | $MB_A=MB_C$ <br> $NB_A=NB_B$ <br> $MB_B=NB_C$ |  $\mathbf{C} = b \mathbf{A} \cdot \mathbf{B}^T + c \mathbf{C}$  |
+   * | true        |   false     | $MB_A=MB_B$ <br> $NB_A=MB_C$ <br> $NB_B=NB_C$ | $\mathbf{C} = b \mathbf{A}^T \cdot \mathbf{B} + c \mathbf{C}$   |
+   * | true        |   true      | $MB_A=NB_B$ <br> $NB_A=MB_C$ <br> $MB_B=NB_C$ | $\mathbf{C} = b \mathbf{A}^T \cdot \mathbf{B}^T + c \mathbf{C}$ |
    *
-   * It is assumed that <tt>A</tt> and <tt>B</tt> have compatible sizes and that
-   * <tt>C</tt> already has the right size.
+   * It is assumed that $\mathbf{A}$ and $\mathbf{B}$ have compatible sizes and that
+   * $\mathbf{C}$ already has the right size.
    *
-   * The matrices <tt>A</tt>, <tt>B</tt> and <tt>C</tt> must have the same process grid.
+   * The matrices $\mathbf{A}$, $\mathbf{B}$ and $\mathbf{C}$ must have the same process grid.
    */
-  void mult(ScaLAPACKMatrix<NumberType> &C,
+  void mult(const NumberType b,
             const ScaLAPACKMatrix<NumberType> &B,
-            const NumberType a=1,
-            const NumberType b=0,
+            const NumberType c,
+            ScaLAPACKMatrix<NumberType> &C,
             const bool transpose_A=false,
             const bool transpose_B=false) const;
 
   /**
    * Matrix-matrix-multiplication.
    *
-   * The optional parameter <tt>adding</tt> determines, whether the result is
-   * stored in <tt>C</tt> or added to <tt>C</tt>.
+   * The optional parameter @p adding determines whether the result is
+   * stored in $\mathbf{C}$ or added to $\mathbf{C}$.
    *
-   * if (adding) <i>C += A*B</i>
+   * if (@p adding) $\mathbf{C} = \mathbf{C} + \mathbf{A} \cdot \mathbf{B}$
    *
-   * else <i>C = A*B</i>
+   * else $\mathbf{C} = \mathbf{A} \cdot \mathbf{B}$
    *
-   * Assumes that <tt>A</tt> and <tt>B</tt> have compatible sizes and that
-   * <tt>C</tt> already has the right size.
+   * It is assumed that $\mathbf{A}$ and $\mathbf{B}$ have compatible sizes and that
+   * $\mathbf{C}$ already has the right size.
    *
-   * Following alignment conditions have to be fulfilled: MB_A = MB_C, NB_A = MB_B and NB_B = NB_C
+   * The following alignment conditions have to be fulfilled: $MB_A=MB_C$, $NB_A=MB_B$ and $NB_B=NB_C$.
    */
   void mmult(ScaLAPACKMatrix<NumberType> &C,
              const ScaLAPACKMatrix<NumberType> &B,
              const bool adding=false) const;
 
   /**
-   * Matrix-matrix-multiplication using transpose of <tt>this</tt>.
+   * Matrix-matrix-multiplication using transpose of $\mathbf{A}$.
    *
-   * The optional parameter <tt>adding</tt> determines, whether the result is
-   * stored in <tt>C</tt> or added to <tt>C</tt>.
+   * The optional parameter @p adding determines whether the result is
+   * stored in $\mathbf{C}$ or added to $\mathbf{C}$.
    *
-   * if (adding) <i>C += A<sup>T</sup>*B</i>
+   * if (@p adding) $\mathbf{C} = \mathbf{C} + \mathbf{A}^T \cdot \mathbf{B}$
    *
-   * else <i>C = A<sup>T</sup>*B</i>
+   * else $\mathbf{C} = \mathbf{A}^T \cdot \mathbf{B}$
    *
-   * Assumes that <tt>A</tt> and <tt>B</tt> have compatible sizes and that
-   * <tt>C</tt> already has the right size.
+   * It is assumed that $\mathbf{A}$ and $\mathbf{B}$ have compatible sizes and that
+   * $\mathbf{C}$ already has the right size.
    *
-   * Following alignment conditions have to be fulfilled: MB_A = MB_B, NB_A = MB_C and NB_B = NB_C
+   * The following alignment conditions have to be fulfilled: $MB_A=MB_B$, $NB_A=MB_C$ and $NB_B=NB_C$.
    */
   void Tmmult (ScaLAPACKMatrix<NumberType> &C,
                const ScaLAPACKMatrix<NumberType> &B,
                const bool adding=false) const;
 
   /**
-   * Matrix-matrix-multiplication using transpose of <tt>B</tt>.
+   * Matrix-matrix-multiplication using the transpose of $\mathbf{B}$.
    *
-   * The optional parameter <tt>adding</tt> determines, whether the result is
-   * stored in <tt>C</tt> or added to <tt>C</tt>.
+   * The optional parameter @p adding determines whether the result is
+   * stored in $\mathbf{C}$ or added to $\mathbf{C}$.
    *
-   * if (adding) <i>C += A*B<sup>T</sup></i>
+   * if (@p adding) $\mathbf{C} = \mathbf{C} + \mathbf{A} \cdot \mathbf{B}^T$
    *
-   * else <i>C = A*B<sup>T</sup></i>
+   * else $\mathbf{C} = \mathbf{A} \cdot \mathbf{B}^T$
    *
-   * Assumes that <tt>A</tt> and <tt>B</tt> have compatible sizes and that
-   * <tt>C</tt> already has the right size.
+   * It is assumed that $\mathbf{A}$ and $\mathbf{B}$ have compatible sizes and that
+   * $\mathbf{C}$ already has the right size.
    *
-   * Following alignment conditions have to be fulfilled: MB_A = MB_C, NB_A = NB_B and MB_B = NB_C
+   * The following alignment conditions have to be fulfilled: $MB_A=MB_C$, $NB_A=NB_B$ and $MB_B=NB_C$.
    */
   void mTmult (ScaLAPACKMatrix<NumberType> &C,
                const ScaLAPACKMatrix<NumberType> &B,
                const bool adding=false) const;
 
   /**
-   * Matrix-matrix-multiplication using transpose of <tt>this</tt> and
-   * <tt>B</tt>.
+   * Matrix-matrix-multiplication using transpose of $\mathbf{A}$ and
+   * $\mathbf{B}$.
    *
-   * The optional parameter <tt>adding</tt> determines, whether the result is
-   * stored in <tt>C</tt> or added to <tt>C</tt>.
+   * The optional parameter @p adding determines whether the result is
+   * stored in $\mathbf{C}$ or added to $\mathbf{C}$.
    *
-   * if (adding) <i>C += A<sup>T</sup>*B<sup>T</sup></i>
+   * if (@p adding) $\mathbf{C} = \mathbf{C} + \mathbf{A}^T \cdot \mathbf{B}^T$
    *
-   * else <i>C = A<sup>T</sup>*B<sup>T</sup></i>
+   * else $\mathbf{C} = \mathbf{A}^T \cdot \mathbf{B}^T$
    *
-   * Assumes that <tt>A</tt> and <tt>B</tt> have compatible sizes and that
-   * <tt>C</tt> already has the right size.
+   * It is assumed that $\mathbf{A}$ and $\mathbf{B}$ have compatible sizes and that
+   * $\mathbf{C}$ already has the right size.
    *
-   * Following alignment conditions have to be fulfilled: MB_A = NB_B, NB_A = MB_C and MB_B = NB_C
+   * The following alignment conditions have to be fulfilled: $MB_A=NB_B$, $NB_A=MB_C$ and $MB_B=NB_C$.
    */
   void TmTmult (ScaLAPACKMatrix<NumberType> &C,
                 const ScaLAPACKMatrix<NumberType> &B,
@@ -548,12 +538,20 @@ public:
   NumberType &local_el(const unsigned int loc_row, const unsigned int loc_column);
 
   /**
-   * scaling the columns of the distributed matrix by scalars in array @p factors
+   * Scale the columns of the distributed matrix by the scalars provided in the array @p factors.
+   *
+   * The array @p factors must have as many entries as the matrix columns.
+   *
+   * Copies of @p factors have to be available on all processes of the underlying MPI communicator.
    */
   void scale_columns(const ArrayView<const NumberType> &factors);
 
   /**
-   * scaling the rows of the distributed matrix by scalars in array @p factors
+   * Scale the rows of the distributed matrix by the scalars provided in the array @p factors.
+   *
+   * The array @p factors must have as many entries as the matrix rows.
+   *
+   * Copies of @p factors have to be available on all processes of the underlying MPI communicator.
    */
   void scale_rows(const ArrayView<const NumberType> &factors);
 
index bd02083a49f26bbbee43d4f2346953c9cfcf2407..0f06675cc04ac6f99a8ddd09659f9a36ad065268 100644 (file)
@@ -1551,19 +1551,19 @@ inline void pgels(const char *trans,
 
 
 template <typename number>
-inline void pgeadd(const char *transa,
-                   const int *m,
-                   const int *n,
-                   const number *alpha,
-                   const number *A,
-                   const int *IA,
-                   const int *JA,
-                   const int *DESCA,
-                   const number *beta,
-                   number *C,
-                   const int *IC,
-                   const int *JC,
-                   const int *DESCC)
+inline void pgeadd(const char * /*transa*/,
+                   const int * /*m*/,
+                   const int * /*n*/,
+                   const number * /*alpha*/,
+                   const number * /*A*/,
+                   const int * /*IA*/,
+                   const int * /*JA*/,
+                   const int * /*DESCA*/,
+                   const number * /*beta*/,
+                   number * /*C*/,
+                   const int * /*IC*/,
+                   const int * /*JC*/,
+                   const int * /*DESCC*/)
 {
   Assert (false, dealii::ExcNotImplemented());
 }
@@ -1604,18 +1604,18 @@ inline void pgeadd(const char *transa,
 
 
 template <typename number>
-inline void ptran(const int *m,
-                  const int *n,
-                  const number *alpha,
-                  const number *A,
-                  const int *IA,
-                  const int *JA,
-                  const int *DESCA,
-                  const number *beta,
-                  number *C,
-                  const int *IC,
-                  const int *JC,
-                  const int *DESCC)
+inline void ptran(const int * /*m*/,
+                  const int * /*n*/,
+                  const number * /*alpha*/,
+                  const number * /*A*/,
+                  const int * /*IA*/,
+                  const int * /*JA*/,
+                  const int * /*DESCA*/,
+                  const number * /*beta*/,
+                  number * /*C*/,
+                  const int * /*IC*/,
+                  const int * /*JC*/,
+                  const int * /*DESCC*/)
 {
   Assert (false, dealii::ExcNotImplemented());
 }
index 1898d38c92c79e69edb6c40812623599638d4508..864962d2df355b474123fe5a95775926d4618c33 100644 (file)
@@ -503,10 +503,10 @@ void ScaLAPACKMatrix<NumberType>::Tadd(const NumberType a,
 
 
 template <typename NumberType>
-void ScaLAPACKMatrix<NumberType>::mult(ScaLAPACKMatrix<NumberType> &C,
+void ScaLAPACKMatrix<NumberType>::mult(const NumberType b,
                                        const ScaLAPACKMatrix<NumberType> &B,
-                                       const NumberType alpha,
-                                       const NumberType beta,
+                                       const NumberType c,
+                                       ScaLAPACKMatrix<NumberType> &C,
                                        const bool transpose_A,
                                        const bool transpose_B) const
 {
@@ -547,12 +547,10 @@ void ScaLAPACKMatrix<NumberType>::mult(ScaLAPACKMatrix<NumberType> &C,
       Assert(this->n_rows==B.n_columns,ExcDimensionMismatch(this->n_rows,B.n_columns));
       Assert(this->n_columns==C.n_rows,ExcDimensionMismatch(this->n_columns,C.n_rows));
       Assert(B.n_rows==C.n_columns,ExcDimensionMismatch(B.n_rows,C.n_columns));
-
       Assert(this->column_block_size==C.row_block_size,ExcDimensionMismatch(this->row_block_size,C.row_block_size));
       Assert(this->row_block_size==B.column_block_size,ExcDimensionMismatch(this->column_block_size,B.row_block_size));
       Assert(B.row_block_size==C.column_block_size,ExcDimensionMismatch(B.column_block_size,C.column_block_size));
     }
-  Threads::Mutex::ScopedLock lock (mutex);
 
   if (this->grid->mpi_process_is_active)
     {
@@ -567,9 +565,9 @@ void ScaLAPACKMatrix<NumberType>::mult(ScaLAPACKMatrix<NumberType> &C,
       int k = transpose_A ? this->n_rows : this->n_columns;
 
       pgemm(&trans_a,&trans_b,&m,&n,&k,
-            &alpha,A_loc,&(this->submatrix_row),&(this->submatrix_column),this->descriptor,
+            &b,A_loc,&(this->submatrix_row),&(this->submatrix_column),this->descriptor,
             B_loc,&B.submatrix_row,&B.submatrix_column,B.descriptor,
-            &beta,C_loc,&C.submatrix_row,&C.submatrix_column,C.descriptor);
+            &c,C_loc,&C.submatrix_row,&C.submatrix_column,C.descriptor);
     }
 }
 
@@ -581,9 +579,9 @@ void ScaLAPACKMatrix<NumberType>::mmult(ScaLAPACKMatrix<NumberType> &C,
                                         const bool adding) const
 {
   if (adding)
-    mult(C,B,1.,1.,false,false);
+    mult(1.,B,1.,C,false,false);
   else
-    mult(C,B,1.,0.,false,false);
+    mult(1.,B,0,C,false,false);
 }
 
 
@@ -594,9 +592,9 @@ void ScaLAPACKMatrix<NumberType>::Tmmult(ScaLAPACKMatrix<NumberType> &C,
                                          const bool adding) const
 {
   if (adding)
-    mult(C,B,1.,1.,true,false);
+    mult(1.,B,1.,C,true,false);
   else
-    mult(C,B,1.,0.,true,false);
+    mult(1.,B,0,C,true,false);
 }
 
 
@@ -607,9 +605,9 @@ void ScaLAPACKMatrix<NumberType>::mTmult(ScaLAPACKMatrix<NumberType> &C,
                                          const bool adding) const
 {
   if (adding)
-    mult(C,B,1.,1.,false,true);
+    mult(1.,B,1.,C,false,true);
   else
-    mult(C,B,1.,0.,false,true);
+    mult(1.,B,0,C,false,true);
 }
 
 
@@ -620,9 +618,9 @@ void ScaLAPACKMatrix<NumberType>::TmTmult(ScaLAPACKMatrix<NumberType> &C,
                                           const bool adding) const
 {
   if (adding)
-    mult(C,B,1.,1.,true,true);
+    mult(1.,B,1.,C,true,true);
   else
-    mult(C,B,1.,0.,true,true);
+    mult(1.,B,0,C,true,true);
 }
 
 
@@ -1514,16 +1512,13 @@ void ScaLAPACKMatrix<NumberType>::scale_columns(const ArrayView<const NumberType
   Assert(n_columns==(int)factors.size(),ExcDimensionMismatch(n_columns,factors.size()));
 
   if (this->grid->mpi_process_is_active)
-    {
-      for (int i=0; i<n_local_rows; ++i)
-        {
-          for (int j=0; j<n_local_columns; ++j)
-            {
-              const int glob_j = global_column(j);
-              local_el(i,j) *= factors[glob_j];
-            }
-        }
-    }
+    for (int i=0; i<n_local_columns; ++i)
+      {
+        const NumberType s = factors[global_column(i)];
+
+        for (int j=0; j<n_local_rows; ++j)
+          local_el(j,i) *= s;
+      }
 }
 
 
@@ -1534,16 +1529,13 @@ void ScaLAPACKMatrix<NumberType>::scale_rows(const ArrayView<const NumberType> &
   Assert(n_rows==(int)factors.size(),ExcDimensionMismatch(n_rows,factors.size()));
 
   if (this->grid->mpi_process_is_active)
-    {
-      for (int i=0; i<n_local_rows; ++i)
-        {
-          const int glob_i = global_row(i);
-          for (int j=0; j<n_local_columns; ++j)
-            {
-              local_el(i,j) *= factors[glob_i];
-            }
-        }
-    }
+    for (int i=0; i<n_local_rows; ++i)
+      {
+        const NumberType s = factors[global_row(i)];
+
+        for (int j=0; j<n_local_columns; ++j)
+          local_el(i,j) *= s;
+      }
 }
 
 
index a832eea9a5c6d2ab1776331837fffbd5e669e370..2518e5efbbdbb279b0d643d6f61a51054b175de1 100644 (file)
@@ -49,7 +49,7 @@ void test()
 
   const std::vector<unsigned int> sizes = {{300,400,500}};
 
-  // test C = alpha A*B + beta C
+  // test C = b A*B + c C
   {
     FullMatrix<NumberType> full_A(sizes[0],sizes[2]);
     FullMatrix<NumberType> full_B(sizes[2],sizes[1]);
@@ -70,17 +70,17 @@ void test()
     scalapack_B = full_B;
     scalapack_C = full_C;
 
-    const NumberType alpha = 1.4, beta = 0.1;
+    const NumberType b=1.4, c=0.1;
 
-    full_A *= alpha;
-    full_C *= beta;
+    full_A *= b;
+    full_C *= c;
     full_A.mmult(full_C,full_B,true);
 
-    scalapack_A.mult(scalapack_C,scalapack_B,alpha,beta,false,false);
+    scalapack_A.mult(b,scalapack_B,c,scalapack_C,false,false);
     FullMatrix<NumberType> tmp_full_C(full_C.m(),full_C.n());
     scalapack_C.copy_to(tmp_full_C);
 
-    pcout << "   computing C = alpha A * B + beta C with"
+    pcout << "   computing C = b A * B + c C with"
           << " A in R^(" << scalapack_A.m() << "x" << scalapack_A.n() << "),"
           << " B in R^(" << scalapack_B.m() << "x" << scalapack_B.n() << ") and"
           << " C in R^(" << scalapack_C.m() << "x" << scalapack_C.n() << ")" << std::endl;
@@ -109,17 +109,17 @@ void test()
     scalapack_B = full_B;
     scalapack_C = full_C;
 
-    const NumberType alpha = 1.4, beta = 0.1;
+    const NumberType b=1.4, c=0.1;
 
-    full_A *= alpha;
-    full_C *= beta;
+    full_A *= b;
+    full_C *= c;
     full_A.Tmmult(full_C,full_B,true);
 
-    scalapack_A.mult(scalapack_C,scalapack_B,alpha,beta,true,false);
+    scalapack_A.mult(b,scalapack_B,c,scalapack_C,true,false);
     FullMatrix<NumberType> tmp_full_C(full_C.m(),full_C.n());
     scalapack_C.copy_to(tmp_full_C);
 
-    pcout << "   computing C = alpha A^T * B + beta C with"
+    pcout << "   computing C = b A^T * B + c C with"
           << " A in R^(" << scalapack_A.m() << "x" << scalapack_A.n() << "),"
           << " B in R^(" << scalapack_B.m() << "x" << scalapack_B.n() << ") and"
           << " C in R^(" << scalapack_C.m() << "x" << scalapack_C.n() << ")" << std::endl;
@@ -148,17 +148,17 @@ void test()
     scalapack_B = full_B;
     scalapack_C = full_C;
 
-    const NumberType alpha = 1.4, beta = 0.1;
+    const NumberType b=1.4, c=0.1;
 
-    full_A *= alpha;
-    full_C *= beta;
+    full_A *= b;
+    full_C *= c;
     full_A.mTmult(full_C,full_B,true);
 
-    scalapack_A.mult(scalapack_C,scalapack_B,alpha,beta,false,true);
+    scalapack_A.mult(b,scalapack_B,c,scalapack_C,false,true);
     FullMatrix<NumberType> tmp_full_C(full_C.m(),full_C.n());
     scalapack_C.copy_to(tmp_full_C);
 
-    pcout << "   computing C = alpha A * B^T + beta C with"
+    pcout << "   computing C = b A * B^T + c C with"
           << " A in R^(" << scalapack_A.m() << "x" << scalapack_A.n() << "),"
           << " B in R^(" << scalapack_B.m() << "x" << scalapack_B.n() << ") and"
           << " C in R^(" << scalapack_C.m() << "x" << scalapack_C.n() << ")" << std::endl;
@@ -187,17 +187,17 @@ void test()
     scalapack_B = full_B;
     scalapack_C = full_C;
 
-    const NumberType alpha = 1.4, beta = 0.1;
+    const NumberType b=1.4, c=0.1;
 
-    full_A *= alpha;
-    full_C *= beta;
+    full_A *= b;
+    full_C *= c;
     full_A.TmTmult(full_C,full_B,true);
 
-    scalapack_A.mult(scalapack_C,scalapack_B,alpha,beta,true,true);
+    scalapack_A.mult(b,scalapack_B,c,scalapack_C,true,true);
     FullMatrix<NumberType> tmp_full_C(full_C.m(),full_C.n());
     scalapack_C.copy_to(tmp_full_C);
 
-    pcout << "   computing C = alpha A^T * B^T + beta C with"
+    pcout << "   computing C = b A^T * B^T + c C with"
           << " A in R^(" << scalapack_A.m() << "x" << scalapack_A.n() << "),"
           << " B in R^(" << scalapack_B.m() << "x" << scalapack_B.n() << ") and"
           << " C in R^(" << scalapack_C.m() << "x" << scalapack_C.n() << ")" << std::endl;
index a65b49530501fb8622fb907f6eecb1d44841d501..a90c724898078f352ed7312180a781d2936f9ca6 100644 (file)
@@ -1,15 +1,15 @@
 2D process grid: 1x1
 
-   computing C = alpha A * B + beta C with A in R^(300x500), B in R^(500x400) and C in R^(300x400)
+   computing C = b A * B + c C with A in R^(300x500), B in R^(500x400) and C in R^(300x400)
    norms: 60787.8449 & 60787.8449  for d
 
-   computing C = alpha A^T * B + beta C with A in R^(500x300), B in R^(500x400) and C in R^(300x400)
+   computing C = b A^T * B + c C with A in R^(500x300), B in R^(500x400) and C in R^(300x400)
    norms: 60655.07764 & 60655.07764  for d
 
-   computing C = alpha A * B^T + beta C with A in R^(300x500), B in R^(400x500) and C in R^(300x400)
+   computing C = b A * B^T + c C with A in R^(300x500), B in R^(400x500) and C in R^(300x400)
    norms: 60707.53954 & 60707.53954  for d
 
-   computing C = alpha A^T * B^T + beta C with A in R^(500x300), B in R^(400x500) and C in R^(300x400)
+   computing C = b A^T * B^T + c C with A in R^(500x300), B in R^(400x500) and C in R^(300x400)
    norms: 60757.09659 & 60757.09659  for d
 
 
index 9646ef5461a660a23e0dd68a990d87874168abd5..d71715a1755d4655f49233934b9b216e8fa6731b 100644 (file)
@@ -1,15 +1,15 @@
 2D process grid: 3x3
 
-   computing C = alpha A * B + beta C with A in R^(300x500), B in R^(500x400) and C in R^(300x400)
+   computing C = b A * B + c C with A in R^(300x500), B in R^(500x400) and C in R^(300x400)
    norms: 60787.8449 & 60787.8449  for d
 
-   computing C = alpha A^T * B + beta C with A in R^(500x300), B in R^(500x400) and C in R^(300x400)
+   computing C = b A^T * B + c C with A in R^(500x300), B in R^(500x400) and C in R^(300x400)
    norms: 60655.07764 & 60655.07764  for d
 
-   computing C = alpha A * B^T + beta C with A in R^(300x500), B in R^(400x500) and C in R^(300x400)
+   computing C = b A * B^T + c C with A in R^(300x500), B in R^(400x500) and C in R^(300x400)
    norms: 60707.53954 & 60707.53954  for d
 
-   computing C = alpha A^T * B^T + beta C with A in R^(500x300), B in R^(400x500) and C in R^(300x400)
+   computing C = b A^T * B^T + c C with A in R^(500x300), B in R^(400x500) and C in R^(300x400)
    norms: 60757.09659 & 60757.09659  for d
 
 
index 9646ef5461a660a23e0dd68a990d87874168abd5..d71715a1755d4655f49233934b9b216e8fa6731b 100644 (file)
@@ -1,15 +1,15 @@
 2D process grid: 3x3
 
-   computing C = alpha A * B + beta C with A in R^(300x500), B in R^(500x400) and C in R^(300x400)
+   computing C = b A * B + c C with A in R^(300x500), B in R^(500x400) and C in R^(300x400)
    norms: 60787.8449 & 60787.8449  for d
 
-   computing C = alpha A^T * B + beta C with A in R^(500x300), B in R^(500x400) and C in R^(300x400)
+   computing C = b A^T * B + c C with A in R^(500x300), B in R^(500x400) and C in R^(300x400)
    norms: 60655.07764 & 60655.07764  for d
 
-   computing C = alpha A * B^T + beta C with A in R^(300x500), B in R^(400x500) and C in R^(300x400)
+   computing C = b A * B^T + c C with A in R^(300x500), B in R^(400x500) and C in R^(300x400)
    norms: 60707.53954 & 60707.53954  for d
 
-   computing C = alpha A^T * B^T + beta C with A in R^(500x300), B in R^(400x500) and C in R^(300x400)
+   computing C = b A^T * B^T + c C with A in R^(500x300), B in R^(400x500) and C in R^(300x400)
    norms: 60757.09659 & 60757.09659  for d
 
 

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.