]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add the additional instantiations for std::complex<float> and make the necessary...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 18 Nov 2007 23:58:36 +0000 (23:58 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 18 Nov 2007 23:58:36 +0000 (23:58 +0000)
git-svn-id: https://svn.dealii.org/trunk@15521 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/full_matrix.templates.h
deal.II/lac/source/full_matrix.float.cc

index 54ab3a673beb5becb0963cf1317b6bae178e3fb7..ba2341b72354f0c22a259502129ca082b72e7220 100644 (file)
@@ -116,7 +116,7 @@ FullMatrix<number>::all_zero () const
   const number* p = this->data();
   const number* const e = this->data() + this->n_elements();
   while (p!=e)
-    if (*p++ != 0.0)
+    if (*p++ != number(0.0))
       return false;
 
   return true;
@@ -153,7 +153,7 @@ FullMatrix<number>::operator /= (const number factor)
   number       *p = &this->el(0,0);
   const number *e = &this->el(0,0) + n()*m();
 
-  const number factor_inv = 1./factor;
+  const number factor_inv = number(1.)/factor;
 
   Assert (numbers::is_finite(factor_inv), 
           ExcMessage("The given value is not finite but either infinite or Not A Number (NaN)"));
@@ -265,7 +265,7 @@ number FullMatrix<number>::residual (Vector<number2>& dst,
                     size_n = n();
   for (unsigned int i=0; i<size_n; ++i)
     {
-      number s = right(i);
+      number s = number(right(i));
       for (unsigned int j=0; j<size_m; ++j)
        s -= number(src(j)) * this->el(i,j);
       dst(i) = s;
@@ -290,7 +290,7 @@ void FullMatrix<number>::forward (Vector<number2>       &dst,
   unsigned int nu = ( (m()<n()) ? m() : n());
   for (i=0; i<nu; ++i)
     {
-      number s = src(i);
+      number s = number(src(i));
       for (j=0; j<i; ++j)
        s -= number(dst(j)) * this->el(i,j);
       dst(i) = s/this->el(i,i);
@@ -889,7 +889,7 @@ FullMatrix<number>::symmetrize ()
   for (unsigned int i=0; i<N; ++i)
     for (unsigned int j=i+1; j<N; ++j)
       {
-       const number t = (this->el(i,j) + this->el(j,i)) / 2.;
+       const number t = (this->el(i,j) + this->el(j,i)) / number(2.);
        this->el(i,j) = this->el(j,i) = t;
       };
 }
@@ -1889,18 +1889,18 @@ FullMatrix<number>::invert (const FullMatrix<number2> &M)
   switch (this->n_cols()) 
     {
       case 1:
-           this->el(0,0) = 1.0/M.el(0,0);
+           this->el(0,0) = number(1.0)/M.el(0,0);
            return;
       case 2:
                                             // this is Maple output,
                                             // thus a bit unstructured
       {
-           const number t4 = 1.0/(M.el(0,0)*M.el(1,1)-M.el(0,1)*M.el(1,0));
-           this->el(0,0) = M.el(1,1)*t4;
-           this->el(0,1) = -M.el(0,1)*t4;
-           this->el(1,0) = -M.el(1,0)*t4;
-           this->el(1,1) = M.el(0,0)*t4;
-           return;
+       const number t4 = number(1.0)/(M.el(0,0)*M.el(1,1)-M.el(0,1)*M.el(1,0));
+       this->el(0,0) = M.el(1,1)*t4;
+       this->el(0,1) = -M.el(0,1)*t4;
+       this->el(1,0) = -M.el(1,0)*t4;
+       this->el(1,1) = M.el(0,0)*t4;
+       return;
       };
       
       case 3:
@@ -1911,7 +1911,7 @@ FullMatrix<number>::invert (const FullMatrix<number2> &M)
                        t00 = M.el(0,2)*M.el(1,0),
                        t01 = M.el(0,1)*M.el(2,0),
                        t04 = M.el(0,2)*M.el(2,0),
-                       t07 = 1.0/(t4*M.el(2,2)-t6*M.el(2,1)-t8*M.el(2,2)+
+                       t07 = number(1.0)/(t4*M.el(2,2)-t6*M.el(2,1)-t8*M.el(2,2)+
                                   t00*M.el(2,1)+t01*M.el(1,2)-t04*M.el(1,1));
            this->el(0,0) = (M.el(1,1)*M.el(2,2)-M.el(1,2)*M.el(2,1))*t07;
            this->el(0,1) = -(M.el(0,1)*M.el(2,2)-M.el(0,2)*M.el(2,1))*t07;
@@ -1964,7 +1964,7 @@ FullMatrix<number>::invert (const FullMatrix<number2> &M)
        const number t63 = t43*t20-t43*t22-t46*t33+t46*t35+
                           t49*t50-t49*t52-t54*t25+t54*t27+
                           t57*t38-t57*t40-t60*t50+t60*t52;
-       const number t65 = 1./(t42+t63);
+       const number t65 = number(1.)/(t42+t63);
        const number t71 = M.el(0,2)*M.el(2,1);
        const number t73 = M.el(0,3)*M.el(2,1);
        const number t75 = M.el(0,2)*M.el(3,1);
@@ -2095,7 +2095,7 @@ FullMatrix<number>::print_formatted (
       for (unsigned int j=0; j<n(); ++j)
        if (std::abs(this->el(i,j)) > threshold)
          out << std::setw(width)
-             << this->el(i,j) * denominator << ' ';
+             << this->el(i,j) * number(denominator) << ' ';
        else
          out << std::setw(width) << zero_string << ' ';
       out << std::endl;
@@ -2170,7 +2170,7 @@ FullMatrix<number>::gauss_jordan ()
        }
 
                                       // transformation
-      const number hr = 1./this->el(j,j);
+      const number hr = number(1.)/this->el(j,j);
       this->el(j,j) = hr;
       for (unsigned int k=0; k<N; ++k)
        {
index 2ad16b567d2d348cea04a03aef40363bf4a6719a..75b77eef07c95b2a725ffe345ba9a8970bd7796e 100644 (file)
@@ -114,4 +114,118 @@ void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
 template TYPEMAT FullMatrix<TYPEMAT>::residual<TYPEVEC,TYPERES>(
   Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
 
+
+
+
+
+
+
+
+
+
+
+
+
+#undef TYPEMAT
+#undef TYPEMAT2
+#undef TYPEVEC
+#undef TYPERES
+
+#define TYPEMAT std::complex<float>
+
+template class FullMatrix<TYPEMAT>;
+
+template void FullMatrix<TYPEMAT>::print(
+  LogStream&, const unsigned int, const unsigned int) const;
+template void FullMatrix<TYPEMAT>::print(
+  std::ostream&, const unsigned int, const unsigned int) const;
+
+template FullMatrix<TYPEMAT>& FullMatrix<TYPEMAT>::operator =(
+  const FullMatrix<std::complex<double> >&);
+
+#define TYPEMAT2 std::complex<float>
+
+//template FullMatrix<TYPEMAT>& FullMatrix<TYPEMAT>::operator =<>(const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::fill<TYPEMAT2> (
+  const FullMatrix<TYPEMAT2>&, unsigned, unsigned, unsigned, unsigned);
+template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&,
+                                                 const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&,
+                                                 const TYPEMAT, const FullMatrix<TYPEMAT2>&,
+                                                 const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (
+  const FullMatrix<TYPEMAT2>&, TYPEMAT, unsigned, unsigned, unsigned, unsigned);
+template void FullMatrix<TYPEMAT>::Tadd<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::Tadd<TYPEMAT2> (
+  const FullMatrix<TYPEMAT2>&, TYPEMAT, unsigned, unsigned, unsigned, unsigned);
+template void FullMatrix<TYPEMAT>::equ<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::equ<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&,
+                                                 const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::equ<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&,
+                                                 const TYPEMAT, const FullMatrix<TYPEMAT2>&,
+                                                 const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::mmult<TYPEMAT2> (FullMatrix<TYPEMAT2>&, const FullMatrix<TYPEMAT2>&, const bool) const;
+template void FullMatrix<TYPEMAT>::Tmmult<TYPEMAT2> (FullMatrix<TYPEMAT2>&, const FullMatrix<TYPEMAT2>&, const bool) const;
+template void FullMatrix<TYPEMAT>::add_diag<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::invert<TYPEMAT2> (const FullMatrix<TYPEMAT2>&);
+
+#define TYPEVEC std::complex<double>
+#define TYPERES std::complex<double>
+
+template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (
+  const FullMatrix<TYPEVEC>&,
+  const std::vector<unsigned int>&,
+  const std::vector<unsigned int>&);
+template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template TYPEMAT FullMatrix<TYPEMAT>::residual<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (
+  const Vector<TYPEVEC> &) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(
+  const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+
+template
+void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
+  Vector<TYPEVEC> &, const Vector<TYPEVEC> &, const TYPEMAT) const;
+
+#undef TYPEVEC
+#define TYPEVEC std::complex<float>
+
+template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (
+  const FullMatrix<TYPEVEC>&,
+  const std::vector<unsigned int>&,
+  const std::vector<unsigned int>&);
+template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template TYPEMAT FullMatrix<TYPEMAT>::residual<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (
+  const Vector<TYPEVEC> &) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(
+  const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template
+void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
+  Vector<TYPEVEC> &, const Vector<TYPEVEC> &, const TYPEMAT) const;
+
+
+#undef TYPERES
+#define TYPERES std::complex<float>
+
+template TYPEMAT FullMatrix<TYPEMAT>::residual<TYPEVEC,TYPERES>(
+  Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+
 DEAL_II_NAMESPACE_CLOSE

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.