]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
number templates improved for long double
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 5 May 1999 11:31:31 +0000 (11:31 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 5 May 1999 11:31:31 +0000 (11:31 +0000)
git-svn-id: https://svn.dealii.org/trunk@1272 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/full_matrix.h
deal.II/lac/include/lac/full_matrix.templates.h
deal.II/lac/include/lac/sparse_matrix.templates.h
deal.II/lac/source/sparse_matrix.cc
deal.II/lac/source/sparse_matrix.double.cc [new file with mode: 0644]
deal.II/lac/source/sparse_matrix.float.cc [new file with mode: 0644]
deal.II/lac/source/sparse_matrix.long_double.cc [new file with mode: 0644]
deal.II/lac/source/vector.long_double.cc [new file with mode: 0644]

index 88635a9d41b41e22c3bd4b8ad3c965b20932af8e..732efc0b345c3b29627097deeed507927a63f2d7 100644 (file)
@@ -36,6 +36,9 @@ class iVector;
  * data types, the implementation of non-inline functions is in
  * "fullmatrix.templates.h". Driver files are in the source tree.
  *
+ * Internal calculations are usually done with the accuracy of the vector argument to
+ * functions. If there is no argument with a number type, the matrix number type is used.
+ *
  * <TABLE BORDER=1>
  * <TR><TH ALIGN=CENTER><B>this</B><TH ALIGN=CENTER><B>other
  * matrix</B><TH ALIGN=CENTER><B>vector</B><TH ALIGN=CENTER><B>rhs in
index dff65e6cb2ac5c840fd4e324107889bb8e721f8c..f6193a393d328948d08833c4cc48844db911b15a 100644 (file)
@@ -115,10 +115,10 @@ FullMatrix<number>::vmult (Vector<number2>& dst,
   Assert(dst.size() == m(), ExcDimensionMismatch(dst.size(), m()));
   Assert(src.size() == n(), ExcDimensionMismatch(src.size(), n()));
 
-  double s;
+  number2 s;
   if ((n()==3) && (m()==3))
   {
-    double s0,s1,s2;
+    number2 s0,s1,s2;
     s   = src(0);
     s0  = s*val[0]; s1  = s*val[3]; s2  = s*val[6]; 
     s   = src(1);
@@ -141,7 +141,7 @@ FullMatrix<number>::vmult (Vector<number2>& dst,
   }
   else if ((n()==4) && (m()==4))
   {
-    double s0,s1,s2,s3;
+    number2 s0,s1,s2,s3;
     s = src(0);
     s0  = s*val[0]; s1  = s*val[4]; s2  = s*val[8];  s3  = s*val[12];
     s = src(1);
@@ -168,7 +168,7 @@ FullMatrix<number>::vmult (Vector<number2>& dst,
   }
   else if ((n()==8) && (m()==8))
   {
-    double s0,s1,s2,s3,s4,s5,s6,s7;
+    number2 s0,s1,s2,s3,s4,s5,s6,s7;
     s = src(0);
     s0 = s*val[0]; s1 = s*val[8]; s2 = s*val[16]; s3 = s*val[24];
     s4 = s*val[32]; s5 = s*val[40]; s6 = s*val[48]; s7 = s*val[56];
@@ -244,10 +244,10 @@ void FullMatrix<number>::gsmult (Vector<number2>& dst, const Vector<number2>& sr
   Assert(src.size() == n(), ExcDimensionMismatch(src.size(), n()));
   Assert(gl.n() == n(), ExcDimensionMismatch(gl.n(), n()));
 
-  double s;
+  number2 s;
   if ((n()==3) && (m()==3))
   {
-    double s0=0.,s1=0.,s2=0.;
+    number2 s0=0.,s1=0.,s2=0.;
     s   = src(0);
     if(gl(1)<gl(0)) s1  = s*val[3]; if(gl(2)<gl(0))  s2  = s*val[6]; 
     s   = src(1);
@@ -261,7 +261,7 @@ void FullMatrix<number>::gsmult (Vector<number2>& dst, const Vector<number2>& sr
   }
   else if ((n()==4) && (m()==4))
   {
-    double s0=0.,s1=0.,s2=0.,s3=0.;
+    number2 s0=0.,s1=0.,s2=0.,s3=0.;
     s = src(0);
     if(gl(1)<gl(0)) s1 = s*val[4];  if(gl(2)<gl(0)) s2 = s*val[8]; if(gl(3)<gl(0)) s3 = s*val[12];
     s = src(1);
@@ -278,7 +278,7 @@ void FullMatrix<number>::gsmult (Vector<number2>& dst, const Vector<number2>& sr
   }
   else if ((n()==8) && (m()==8))
   {
-    double s0=0.,s1=0.,s2=0.,s3=0.,s4=0.,s5=0.,s6=0.,s7=0.;
+    number2 s0=0.,s1=0.,s2=0.,s3=0.,s4=0.,s5=0.,s6=0.,s7=0.;
     s = src(0);
     if(gl(1)<gl(0)) s1 = s*val[8]; 
     if(gl(2)<gl(0)) s2 = s*val[16]; 
@@ -380,7 +380,7 @@ void FullMatrix<number>::Tvmult (Vector<number2>& dst,
   Assert(src.size() == m(), ExcDimensionMismatch(src.size(), m()));
 
   unsigned int i,j;
-  double s;
+  number2 s;
   const unsigned int size_m = m(),
                     size_n = n();
   for (i=0; i<size_m; ++i)
@@ -406,7 +406,7 @@ double FullMatrix<number>::residual (Vector<number2>& dst,
   Assert(right.size() == m(), ExcDimensionMismatch(right.size(), m()));
 
   unsigned int i,j;
-  double s, res = 0.;
+  number2 s, res = 0.;
   const unsigned int size_m = m(),
                     size_n = n();
   for (i=0; i<size_n; ++i)
@@ -432,7 +432,7 @@ void FullMatrix<number>::forward (Vector<number2>& dst,
 
   unsigned int i,j;
   unsigned int nu = ( (m()<n()) ? m() : n());
-  double s;
+  number2 s;
   for (i=0; i<nu; ++i)
     {
       s = src(i);
@@ -450,7 +450,7 @@ void FullMatrix<number>::backward (Vector<number2>& dst,
 {
   unsigned int j;
   unsigned int nu = (m()<n() ? m() : n());
-  double s;
+  number2 s;
   for (int i=nu-1; i>=0; --i)
     {
       s = src(i);
@@ -593,7 +593,7 @@ void FullMatrix<number>::mmult (FullMatrix<number2>& dst,
 {
   Assert (n() == src.m(), ExcDimensionMismatch(n(), src.m()));
   unsigned int i,j,k;
-  double s = 1.;
+  number2 s = 1.;
   dst.reinit(m(), src.n());
 
   for (i=0;i<m();i++)
@@ -612,7 +612,7 @@ void FullMatrix<number>::mmult (FullMatrix<number2>& dst,
   Assert (m() == src.n(), ExcDimensionMismatch(m(), src.n()));
 
   unsigned int i,j,k;
-  double s = 1.;
+  number2 s = 1.;
 
   dst.reinit(n(), src.m());
 
@@ -634,7 +634,7 @@ void FullMatrix<number>::Tmmult (FullMatrix<number2>& dst, const FullMatrix<numb
   Assert (n() == src.n(), ExcDimensionMismatch(n(), src.n()));
 
   unsigned int i,j,k;
-  double s = 1.;
+  number2 s = 1.;
   dst.reinit(m(), src.m());
 
   for (i=0;i<m();i++)
@@ -653,7 +653,7 @@ void FullMatrix<number>::Tmmult (FullMatrix<number2>& dst, const FullMatrix<numb
   Assert (m() == src.n(), ExcDimensionMismatch(m(), src.n()));
 
   unsigned int i,j,k;
-  double s = 1.;
+  number2 s = 1.;
   
   dst.reinit(n(), src.m());
 
@@ -675,14 +675,14 @@ double FullMatrix<number>::matrix_norm (const Vector<number2> &v) const
   Assert(m() == v.size(), ExcDimensionMismatch(m(),v.size()));
   Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size()));
 
-  double sum = 0.;
+  number2 sum = 0.;
   const unsigned int n_rows = m();
   const number *val_ptr = &val[0];
   const number2 *v_ptr;
   
   for (unsigned int row=0; row<n_rows; ++row)
     {
-      double s = 0.;
+      number2 s = 0.;
       const number * const val_end_of_row = val_ptr+n_rows;
       v_ptr = v.begin();
       while (val_ptr != val_end_of_row)
@@ -703,7 +703,7 @@ double FullMatrix<number>::matrix_scalar_product (const Vector<number2> &u, cons
   Assert(m() == u.size(), ExcDimensionMismatch(m(),v.size()));
   Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size()));
 
-  double sum = 0.;
+  number2 sum = 0.;
   const unsigned int n_rows = m();
   const unsigned int n_cols = n();
   const number *val_ptr = &val[0];
@@ -711,7 +711,7 @@ double FullMatrix<number>::matrix_scalar_product (const Vector<number2> &u, cons
   
   for (unsigned int row=0; row<n_rows; ++row)
     {
-      double s = 0.;
+      number2 s = 0.;
       const number * const val_end_of_row = val_ptr+n_cols;
       v_ptr = v.begin();
       while (val_ptr != val_end_of_row)
@@ -1188,7 +1188,7 @@ template <typename number>
 double
 FullMatrix<number>::norm2 () const
 {
-  double s = 0.;
+  number s = 0.;
   for (unsigned int i=0;i<dim_image*dim_range;++i)
     s += val[i]*val[i];
   return s;
@@ -1223,7 +1223,7 @@ FullMatrix<number>::invert (const FullMatrix<number> &M)
                                             // this is Maple output,
                                             // thus a bit unstructured
       {
-           const double t4 = 1.0/(M.el(0,0)*M.el(1,1)-M.el(0,1)*M.el(1,0));
+           const number t4 = 1.0/(M.el(0,0)*M.el(1,1)-M.el(0,1)*M.el(1,0));
            el(0,0) = M.el(1,1)*t4;
            el(0,1) = -M.el(0,1)*t4;
            el(1,0) = -M.el(1,0)*t4;
@@ -1233,7 +1233,7 @@ FullMatrix<number>::invert (const FullMatrix<number> &M)
       
       case 3:
       {
-           const double t4 = M.el(0,0)*M.el(1,1),
+           const number t4 = M.el(0,0)*M.el(1,1),
                         t6 = M.el(0,0)*M.el(1,2),
                         t8 = M.el(0,1)*M.el(1,0),
                        t00 = M.el(0,2)*M.el(1,0),
@@ -1262,61 +1262,61 @@ FullMatrix<number>::invert (const FullMatrix<number> &M)
                                         // readlib(C);
                                         // C(ai,optimized,filename=x4);
 
-       const double t14 = M.el(0,0)*M.el(1,1);
-       const double t15 = M.el(2,2)*M.el(3,3);
-       const double t17 = M.el(2,3)*M.el(3,2);
-       const double t19 = M.el(0,0)*M.el(2,1);
-       const double t20 = M.el(1,2)*M.el(3,3);
-       const double t22 = M.el(1,3)*M.el(3,2);
-       const double t24 = M.el(0,0)*M.el(3,1);
-       const double t25 = M.el(1,2)*M.el(2,3);
-       const double t27 = M.el(1,3)*M.el(2,2);
-       const double t29 = M.el(1,0)*M.el(0,1);
-       const double t32 = M.el(1,0)*M.el(2,1);
-       const double t33 = M.el(0,2)*M.el(3,3);
-       const double t35 = M.el(0,3)*M.el(3,2);
-       const double t37 = M.el(1,0)*M.el(3,1);
-       const double t38 = M.el(0,2)*M.el(2,3);
-       const double t40 = M.el(0,3)*M.el(2,2);
-       const double t42 = t14*t15-t14*t17-t19*t20+t19*t22+
+       const number t14 = M.el(0,0)*M.el(1,1);
+       const number t15 = M.el(2,2)*M.el(3,3);
+       const number t17 = M.el(2,3)*M.el(3,2);
+       const number t19 = M.el(0,0)*M.el(2,1);
+       const number t20 = M.el(1,2)*M.el(3,3);
+       const number t22 = M.el(1,3)*M.el(3,2);
+       const number t24 = M.el(0,0)*M.el(3,1);
+       const number t25 = M.el(1,2)*M.el(2,3);
+       const number t27 = M.el(1,3)*M.el(2,2);
+       const number t29 = M.el(1,0)*M.el(0,1);
+       const number t32 = M.el(1,0)*M.el(2,1);
+       const number t33 = M.el(0,2)*M.el(3,3);
+       const number t35 = M.el(0,3)*M.el(3,2);
+       const number t37 = M.el(1,0)*M.el(3,1);
+       const number t38 = M.el(0,2)*M.el(2,3);
+       const number t40 = M.el(0,3)*M.el(2,2);
+       const number t42 = t14*t15-t14*t17-t19*t20+t19*t22+
                           t24*t25-t24*t27-t29*t15+t29*t17+
                           t32*t33-t32*t35-t37*t38+t37*t40;
-       const double t43 = M.el(2,0)*M.el(0,1);
-       const double t46 = M.el(2,0)*M.el(1,1);
-       const double t49 = M.el(2,0)*M.el(3,1);
-       const double t50 = M.el(0,2)*M.el(1,3);
-       const double t52 = M.el(0,3)*M.el(1,2);
-       const double t54 = M.el(3,0)*M.el(0,1);
-       const double t57 = M.el(3,0)*M.el(1,1);
-       const double t60 = M.el(3,0)*M.el(2,1);
-       const double t63 = t43*t20-t43*t22-t46*t33+t46*t35+
+       const number t43 = M.el(2,0)*M.el(0,1);
+       const number t46 = M.el(2,0)*M.el(1,1);
+       const number t49 = M.el(2,0)*M.el(3,1);
+       const number t50 = M.el(0,2)*M.el(1,3);
+       const number t52 = M.el(0,3)*M.el(1,2);
+       const number t54 = M.el(3,0)*M.el(0,1);
+       const number t57 = M.el(3,0)*M.el(1,1);
+       const number t60 = M.el(3,0)*M.el(2,1);
+       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 double t65 = 1/(t42+t63);
-       const double t71 = M.el(0,2)*M.el(2,1);
-       const double t73 = M.el(0,3)*M.el(2,1);
-       const double t75 = M.el(0,2)*M.el(3,1);
-       const double t77 = M.el(0,3)*M.el(3,1);
-       const double t81 = M.el(0,1)*M.el(1,2);
-       const double t83 = M.el(0,1)*M.el(1,3);
-       const double t85 = M.el(0,2)*M.el(1,1);
-       const double t87 = M.el(0,3)*M.el(1,1);
-       const double t101 = M.el(1,0)*M.el(2,2);
-       const double t103 = M.el(1,0)*M.el(2,3);
-       const double t105 = M.el(2,0)*M.el(1,2);
-       const double t107 = M.el(2,0)*M.el(1,3);
-       const double t109 = M.el(3,0)*M.el(1,2);
-       const double t111 = M.el(3,0)*M.el(1,3);
-       const double t115 = M.el(0,0)*M.el(2,2);
-       const double t117 = M.el(0,0)*M.el(2,3);
-       const double t119 = M.el(2,0)*M.el(0,2);
-       const double t121 = M.el(2,0)*M.el(0,3);
-       const double t123 = M.el(3,0)*M.el(0,2);
-       const double t125 = M.el(3,0)*M.el(0,3);
-       const double t129 = M.el(0,0)*M.el(1,2);
-       const double t131 = M.el(0,0)*M.el(1,3);
-       const double t133 = M.el(1,0)*M.el(0,2);
-       const double t135 = M.el(1,0)*M.el(0,3);
+       const number t65 = 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);
+       const number t77 = M.el(0,3)*M.el(3,1);
+       const number t81 = M.el(0,1)*M.el(1,2);
+       const number t83 = M.el(0,1)*M.el(1,3);
+       const number t85 = M.el(0,2)*M.el(1,1);
+       const number t87 = M.el(0,3)*M.el(1,1);
+       const number t101 = M.el(1,0)*M.el(2,2);
+       const number t103 = M.el(1,0)*M.el(2,3);
+       const number t105 = M.el(2,0)*M.el(1,2);
+       const number t107 = M.el(2,0)*M.el(1,3);
+       const number t109 = M.el(3,0)*M.el(1,2);
+       const number t111 = M.el(3,0)*M.el(1,3);
+       const number t115 = M.el(0,0)*M.el(2,2);
+       const number t117 = M.el(0,0)*M.el(2,3);
+       const number t119 = M.el(2,0)*M.el(0,2);
+       const number t121 = M.el(2,0)*M.el(0,3);
+       const number t123 = M.el(3,0)*M.el(0,2);
+       const number t125 = M.el(3,0)*M.el(0,3);
+       const number t129 = M.el(0,0)*M.el(1,2);
+       const number t131 = M.el(0,0)*M.el(1,3);
+       const number t133 = M.el(1,0)*M.el(0,2);
+       const number t135 = M.el(1,0)*M.el(0,3);
        el(0,0) = (M.el(1,1)*M.el(2,2)*M.el(3,3)-M.el(1,1)*M.el(2,3)*M.el(3,2)-
                   M.el(2,1)*M.el(1,2)*M.el(3,3)+M.el(2,1)*M.el(1,3)*M.el(3,2)+
                   M.el(3,1)*M.el(1,2)*M.el(2,3)-M.el(3,1)*M.el(1,3)*M.el(2,2))*t65;
@@ -1396,7 +1396,7 @@ FullMatrix<number>::gauss_jordan()
   iVector p(n());
 
   unsigned int i,j,k,r;
-  double max, hr;
+  number max, hr;
 
   for (i=0; i<n(); ++i) p(i) = i;
 
@@ -1465,27 +1465,27 @@ FullMatrix<number>::householder(Vector<number2>& src)
 
   for (unsigned int j=0 ; j<n() ; ++j)
   {
-    double sigma = 0;
+    number2 sigma = 0;
     unsigned int i;
     for (i=j ; i<m() ; ++i) sigma += el(i,j)*el(i,j);
     if (fabs(sigma) < 1.e-15) return;
-    double s = el(j,j);
+    number2 s = el(j,j);
     s = (s<0) ? sqrt(sigma) : -sqrt(sigma);
-    double dj = s;
+    number2 dj = s;
 
-    double beta = 1./(s*el(j,j)-sigma);
+    number2 beta = 1./(s*el(j,j)-sigma);
     el(j,j) -= s;
 
     for (unsigned int k=j+1 ; k<n() ; ++k)
     {
-      double sum = 0.;
+      number2 sum = 0.;
       for (i=j ; i<m() ; ++i) sum += el(i,j)*el(i,k);
       sum *= beta;
 
       for (i=j ; i<m() ; ++i) el(i,k) += sum*el(i,j);
     }
 
-    double sum = 0.;
+    number2 sum = 0.;
     for (i=j ; i<m() ; ++i) sum += el(i,j)*src(i);
     sum *= beta;
 
@@ -1506,7 +1506,7 @@ FullMatrix<number>::least_squares(Vector<number2>& dst, Vector<number2>& src)
   householder(src);
   backward(dst, src);
 
-  double sum = 0.;
+  number2 sum = 0.;
   for (unsigned int i=n() ; i<m() ; ++i) sum += src(i) * src(i);
   return sqrt(sum);
 }
index bfef8d7f7dea29524903b45314263faa083b1765..a9bb565977ec420cf826c68f89ba17f10be321d1 100644 (file)
@@ -190,7 +190,7 @@ SparseMatrix<number>::vmult (Vector<somenumber>& dst, const Vector<somenumber>&
   somenumber   *dst_ptr = &dst(0);
   for (unsigned int row=0; row<n_rows; ++row)
     {
-      double s = 0.;
+      somenumber s = 0.;
       const number *const val_end_of_row = &val[cols->rowstart[row+1]];
       while (val_ptr != val_end_of_row)
        s += *val_ptr++ * src(*colnum_ptr++);
@@ -233,13 +233,13 @@ SparseMatrix<number>::matrix_norm (const Vector<somenumber>& v) const
   Assert(m() == v.size(), ExcDimensionsDontMatch(m(),v.size()));
   Assert(n() == v.size(), ExcDimensionsDontMatch(n(),v.size()));
 
-  double sum = 0.;
+  somenumber sum = 0.;
   const unsigned int n_rows = m();
   const number *val_ptr = &val[cols->rowstart[0]];
   const int *colnum_ptr = &cols->colnums[cols->rowstart[0]];
   for (unsigned int row=0; row<n_rows; ++row)
     {
-      double s = 0.;
+      somenumber s = 0.;
       const number *val_end_of_row = &val[cols->rowstart[row+1]];
       while (val_ptr != val_end_of_row)
        s += *val_ptr++ * v(*colnum_ptr++);
@@ -303,7 +303,7 @@ SparseMatrix<number>::residual (Vector<somenumber>& dst,
   Assert(m() == b.size(), ExcDimensionsDontMatch(m(),b.size()));
   Assert(n() == u.size(), ExcDimensionsDontMatch(n(),u.size()));
 
-  double s,norm=0.;   
+  somenumber s,norm=0.;   
   
   for (unsigned int i=0;i<m();i++)
     {
@@ -461,7 +461,7 @@ SparseMatrix<number>::SSOR (Vector<somenumber>& dst, const number om) const
   int p;
   const unsigned int  n = dst.size();
   unsigned int  j;
-  double s;
+  somenumber s;
   
   for (unsigned int i=0; i<n; i++)
     {
index 4df57c18c7faee664b4abdaa12f02cadd322e431..1bfd73193a06ac2bbffe50c075b840db35474359 100644 (file)
@@ -7,7 +7,7 @@
 
 
 #include <lac/sparsematrix.h>
-#include <lac/sparsematrix.templates.h>
+//#include <lac/sparsematrix.templates.h>
 #include <lac/ivector.h>
 
 #include <iostream>
@@ -487,122 +487,3 @@ SparseMatrixStruct::n_nonzero_elements () const {
 };
 
 
-
-
-
-
-
-
-// explicit instantiations for "float"
-template class SparseMatrix<float>;
-
-template SparseMatrix<float> & SparseMatrix<float>::copy_from (const SparseMatrix<float> &);
-template SparseMatrix<float> & SparseMatrix<float>::copy_from (const SparseMatrix<double> &);
-
-template void SparseMatrix<float>::add_scaled (const float, const SparseMatrix<float> &);
-template void SparseMatrix<float>::add_scaled (const float, const SparseMatrix<double> &);
-
-template void SparseMatrix<float>::vmult (Vector<float> &, const Vector<float> &) const;
-template void SparseMatrix<float>::vmult (Vector<double> &, const Vector<double> &) const;
-
-template void SparseMatrix<float>::Tvmult (Vector<float> &, const Vector<float> &) const;
-template void SparseMatrix<float>::Tvmult (Vector<double> &, const Vector<double> &) const;
-
-template double SparseMatrix<float>::matrix_norm (const Vector<float> &) const;
-template double SparseMatrix<float>::matrix_norm (const Vector<double> &) const;
-
-template double SparseMatrix<float>::residual (Vector<float> &,
-                                              const Vector<float> &,
-                                              const Vector<float> &) const;
-template double SparseMatrix<float>::residual (Vector<double> &,
-                                              const Vector<double> &,
-                                              const Vector<double> &) const;
-
-template void SparseMatrix<float>::precondition_SSOR (Vector<float> &,
-                                                     const Vector<float> &,
-                                                     float) const;
-template void SparseMatrix<float>::precondition_SSOR (Vector<double> &,
-                                                     const Vector<double> &,
-                                                     float) const;
-
-template void SparseMatrix<float>::precondition_SOR (Vector<float> &,
-                                                    const Vector<float> &,
-                                                    float) const;
-template void SparseMatrix<float>::precondition_SOR (Vector<double> &,
-                                                    const Vector<double> &,
-                                                    float) const;
-
-template void SparseMatrix<float>::precondition_Jacobi (Vector<float> &,
-                                                       const Vector<float> &,
-                                                       float) const;
-template void SparseMatrix<float>::precondition_Jacobi (Vector<double> &,
-                                                       const Vector<double> &,
-                                                       float) const;
-
-template void SparseMatrix<float>::SOR (Vector<float> &,
-                                       float) const;
-template void SparseMatrix<float>::SOR (Vector<double> &,
-                                       float) const;
-
-template void SparseMatrix<float>::SSOR (Vector<float> &,
-                                        float) const;
-template void SparseMatrix<float>::SSOR (Vector<double> &,
-                                        float) const;
-
-
-
-// explicit instantiations for "double"
-template class SparseMatrix<double>;
-
-template SparseMatrix<double> & SparseMatrix<double>::copy_from (const SparseMatrix<float> &);
-template SparseMatrix<double> & SparseMatrix<double>::copy_from (const SparseMatrix<double> &);
-
-template void SparseMatrix<double>::add_scaled (const double, const SparseMatrix<float> &);
-template void SparseMatrix<double>::add_scaled (const double, const SparseMatrix<double> &);
-
-template void SparseMatrix<double>::vmult (Vector<float> &, const Vector<float> &) const;
-template void SparseMatrix<double>::vmult (Vector<double> &, const Vector<double> &) const;
-
-template void SparseMatrix<double>::Tvmult (Vector<float> &, const Vector<float> &) const;
-template void SparseMatrix<double>::Tvmult (Vector<double> &, const Vector<double> &) const;
-
-template double SparseMatrix<double>::matrix_norm (const Vector<float> &) const;
-template double SparseMatrix<double>::matrix_norm (const Vector<double> &) const;
-
-template double SparseMatrix<double>::residual (Vector<float> &,
-                                               const Vector<float> &,
-                                               const Vector<float> &) const;
-template double SparseMatrix<double>::residual (Vector<double> &,
-                                               const Vector<double> &,
-                                               const Vector<double> &) const;
-
-template void SparseMatrix<double>::precondition_SSOR (Vector<float> &,
-                                                      const Vector<float> &,
-                                                      double) const;
-template void SparseMatrix<double>::precondition_SSOR (Vector<double> &,
-                                                      const Vector<double> &,
-                                                      double) const;
-
-template void SparseMatrix<double>::precondition_SOR (Vector<float> &,
-                                                     const Vector<float> &,
-                                                     double) const;
-template void SparseMatrix<double>::precondition_SOR (Vector<double> &,
-                                                     const Vector<double> &,
-                                                     double) const;
-
-template void SparseMatrix<double>::precondition_Jacobi (Vector<float> &,
-                                                        const Vector<float> &,
-                                                        double) const;
-template void SparseMatrix<double>::precondition_Jacobi (Vector<double> &,
-                                                        const Vector<double> &,
-                                                        double) const;
-
-template void SparseMatrix<double>::SOR (Vector<float> &,
-                                        double) const;
-template void SparseMatrix<double>::SOR (Vector<double> &,
-                                        double) const;
-
-template void SparseMatrix<double>::SSOR (Vector<float> &,
-                                         double) const;
-template void SparseMatrix<double>::SSOR (Vector<double> &,
-                                         double) const;
diff --git a/deal.II/lac/source/sparse_matrix.double.cc b/deal.II/lac/source/sparse_matrix.double.cc
new file mode 100644 (file)
index 0000000..d3adea7
--- /dev/null
@@ -0,0 +1,88 @@
+// $Id$
+
+// SparseMatrix template instantiation
+
+
+/* Instantiation is controlled by preprocessor symbols:
+ *
+ * 1. TYPEMAT : numerical type used in the matrix
+ * 2. TYPE2 : numerical type for function arguments
+ */
+
+#include <cmath>
+#include <lac/sparsematrix.templates.h>
+
+#define TYPEMAT double
+
+template class SparseMatrix<TYPEMAT>;
+
+#define TYPE2 float
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
+
+#undef TYPE2
+#define TYPE2 double
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
diff --git a/deal.II/lac/source/sparse_matrix.float.cc b/deal.II/lac/source/sparse_matrix.float.cc
new file mode 100644 (file)
index 0000000..5be1b54
--- /dev/null
@@ -0,0 +1,88 @@
+// $Id$
+
+// SparseMatrix template instantiation
+
+
+/* Instantiation is controlled by preprocessor symbols:
+ *
+ * 1. TYPEMAT : numerical type used in the matrix
+ * 2. TYPE2 : numerical type for function arguments
+ */
+
+#include <cmath>
+#include <lac/sparsematrix.templates.h>
+
+#define TYPEMAT float
+
+template class SparseMatrix<TYPEMAT>;
+
+#define TYPE2 float
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
+
+#undef TYPE2
+#define TYPE2 double
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
diff --git a/deal.II/lac/source/sparse_matrix.long_double.cc b/deal.II/lac/source/sparse_matrix.long_double.cc
new file mode 100644 (file)
index 0000000..6678c95
--- /dev/null
@@ -0,0 +1,196 @@
+// $Id$
+
+// SparseMatrix template instantiation
+
+
+/* Instantiation is controlled by preprocessor symbols:
+ *
+ * 1. TYPEMAT : numerical type used in the matrix
+ * 2. TYPE2 : numerical type for function arguments
+ */
+
+#include <cmath>
+#include <lac/sparsematrix.templates.h>
+
+#define TYPEMAT long double
+
+template class SparseMatrix<TYPEMAT>;
+
+#define TYPE2 float
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
+
+#undef TYPE2
+#define TYPE2 double
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
+
+#undef TYPE2
+#define TYPE2 long double
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
+
+#undef TYPEMAT
+#define TYPEMAT double
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
+
+#undef TYPEMAT
+#define TYPEMAT float
+
+template SparseMatrix<TYPEMAT> &
+SparseMatrix<TYPEMAT>::copy_from (const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::add_scaled (const TYPEMAT,
+                                                const SparseMatrix<TYPE2> &);
+
+template void SparseMatrix<TYPEMAT>::vmult (Vector<TYPE2> &,
+                                           const Vector<TYPE2> &) const;
+template void SparseMatrix<TYPEMAT>::Tvmult (Vector<TYPE2> &,
+                                            const Vector<TYPE2> &) const;
+
+template double
+SparseMatrix<TYPEMAT>::matrix_norm (const Vector<TYPE2> &) const;
+
+template double SparseMatrix<TYPEMAT>::residual (Vector<TYPE2> &,
+                                              const Vector<TYPE2> &,
+                                              const Vector<TYPE2> &) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SSOR (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &,
+                                                     TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &,
+                                                    TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
+                                                       const Vector<TYPE2> &,
+                                                       TYPEMAT) const;
+
+template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
diff --git a/deal.II/lac/source/vector.long_double.cc b/deal.II/lac/source/vector.long_double.cc
new file mode 100644 (file)
index 0000000..2fc8e16
--- /dev/null
@@ -0,0 +1,13 @@
+// $Id$
+
+#include <lac/vector.templates.h>
+
+// explicit instantiations
+template class Vector<long double>;
+template Vector<long double>& Vector<long double>::operator=(const Vector<float>&);
+template Vector<long double>& Vector<long double>::operator=(const Vector<double>&);
+template Vector<double>& Vector<double>::operator=(const Vector<long double>&);
+template Vector<float>& Vector<float>::operator=(const Vector<long double>&);
+// see the .h file for why these functions are disabled.
+// template Vector<float>::Vector (const Vector<double>& v);
+// template Vector<double>::Vector (const Vector<float>& v);

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.