Changeset 999 in Sophya for trunk/SophyaLib

Timestamp:

May 3, 2000, 6:39:26 PM (25 years ago)

Author:

ansari

Message:

mise en template des inline + new GausPiv+Determinant cmv 3/5/00

Location:

trunk/SophyaLib/TArray

Files:

• sopemtx.cc (modified) (1 diff)
• sopemtx.h (modified) (4 diffs)

trunk/SophyaLib/TArray/sopemtx.cc

@@ -447 +447 @@
 //! Gaussian pivoting
 /*!
-  Diagonalize matrix \b a, doing the same operations on matrix \b b
-  \return determinat of \b a.
-
-  If \b b is identity matrix, return inverse of \b a
+  Do Gauss pivoting of \b a, doing the same operations on matrix \b b
+  \return determinant of \b a.
+  \verbatim
+  Solve linear system A(n,n) * X(n,m) = B(n,m)
+  and put solution X in B for return.
+  \endverbatim
+  \warning If \b b is identity matrix, return inverse of \b a
+  \warning matrix \b a and \b b are modified.
  */  
 template <class T>

trunk/SophyaLib/TArray/sopemtx.h

@@ -23 +23 @@
   \sa TMatrix TArray
 */
-
 //! Class for simple operation on TMatrix
 template <class T>
@@ -47 +46 @@
 // Resolution du systeme A*C = B
 //------------------------------------------------------------
-
-/*! \ingroup TArray \fn LinSolveInPlace(TMatrix<r_4>&,TVector<r_4>&)
+/*! \ingroup TArray \fn LinSolveInPlace(TMatrix<T>&,TVector<T>&)
     \brief  Solve A*C = B for C in place and return determinant
 */
-inline r_4 LinSolveInPlace(TMatrix<r_4>& a, TVector<r_4>& b)
+template <class T>
+inline T LinSolveInPlace(TMatrix<T>& a, TVector<T>& b)
 {
 if(a.NCols() != b.NRows() || a.NCols() != a.NRows())
-  throw(SzMismatchError("LinSolveInPlace(TMatrix<r_4>,TVector<r_4>) size mismatch"));
-return SimpleMatrixOperation<r_4>::GausPiv(a,b);
-}
-
-/*! \ingroup TArray \fn LinSolveInPlace(TMatrix<r_8>&,TVector<r_8>&)
-    \brief  Solve A*X = B in place and return determinant 
-*/
-inline r_8 LinSolveInPlace(TMatrix<r_8>& a, TVector<r_8>& b)
-{
-if(a.NCols() != b.NRows() || a.NCols() != a.NRows())
-  throw(SzMismatchError("LinSolveInPlace(TMatrix<r_8>,TVector<r_8>) size mismatch"));
-return SimpleMatrixOperation<r_8>::GausPiv(a,b);
-}
-
-/*! \ingroup TArray 
-    \fn LinSolveInPlace(TMatrix< complex<r_4> >&, TVector< complex<r_4> >&)
-    \brief Solve A*X = B in place and return determinant */
-inline complex<r_4> LinSolveInPlace(TMatrix< complex<r_4> >& a, TVector< complex<r_4> >& b)
-{
-if(a.NCols() != b.NRows() || a.NCols() != a.NRows())
-  throw(SzMismatchError("LinSolveInPlace(TMatrix< complex<r_4> >,TVector< complex<r_4> >) size mismatch"));
-return SimpleMatrixOperation< complex<r_4> >::GausPiv(a,b);
-}
-
-/*! \ingroup TArray 
-    \fn LinSolveInPlace(TMatrix< complex<r_8> >&, TVector< complex<r_8> >&)
-    \brief  Solve A*X = B in place and return determinant
-*/
-inline complex<r_8> LinSolveInPlace(TMatrix< complex<r_8> >& a, TVector< complex<r_8> >& b)
-{
-if(a.NCols() != b.NRows() || a.NCols() != a.NRows())
-  throw(SzMismatchError("LinSolveInPlace(TMatrix< complex<r_8> >,TVector< complex<r_8> >) size mismatch"));
-return SimpleMatrixOperation< complex<r_8> >::GausPiv(a,b);
+  throw(SzMismatchError("LinSolveInPlace(TMatrix<T>,TVector<T>) size mismatch"));
+return SimpleMatrixOperation<T>::GausPiv(a,b);
 }
 
@@ -92 +60 @@
 // Resolution du systeme A*C = B, avec C retourne dans B
 //------------------------------------------------------------
-
 /*! \ingroup TArray 
-    \fn LinSolve(const TMatrix<r_4>&, const TVector<r_4>&, TVector<r_4>&)
+    \fn LinSolve(const TMatrix<T>&, const TVector<T>&, TVector<T>&)
     \brief Solve A*C = B and return C and determinant
 */
-
-inline r_4 LinSolve(const TMatrix<r_4>& a, const TVector<r_4>& b, TVector<r_4>& c) {
+template <class T>
+inline T LinSolve(const TMatrix<T>& a, const TVector<T>& b, TVector<T>& c) {
   if(a.NCols()!=b.NRows() || a.NCols()!=a.NRows())
-    throw(SzMismatchError("LinSolve(TMatrix<r_4>,TVector<r_4>) size mismatch"));
-    c = b; TMatrix<r_4> a1(a);
-  return SimpleMatrixOperation<r_4>::GausPiv(a1,c);
-}
-
-/*! \ingroup TArray 
-    \fn LinSolve(const TMatrix<r_8>&, const TVector<r_8>&, TVector<r_8>&)
-    \brief Solve A*C = B and return C and determinant
-*/
-inline r_8 LinSolve(const TMatrix<r_8>& a, const TVector<r_8>& b, TVector<r_8>& c) {
-  if(a.NCols()!=b.NRows() || a.NCols()!=a.NRows())
-    throw(SzMismatchError("LinSolve(TMatrix<r_8>,TVector<r_8>) size mismatch"));
-    c = b; TMatrix<r_8> a1(a);
-  return SimpleMatrixOperation<r_8>::GausPiv(a1,c);
-}
-
-/*! \ingroup TArray 
-    \fn LinSolve(const TMatrix< complex<r_4> >&, const TVector< complex<r_4> >&, TVector< complex<r_4> >&)
-    \brief Solve A*C = B and return C and determinant
-*/
-inline complex<r_4> LinSolve(const TMatrix< complex<r_4> >& a, const TVector< complex<r_4> >& b, TVector< complex<r_4> >& c) {
-  if(a.NCols()!=b.NRows() || a.NCols()!=a.NRows())
-    throw(SzMismatchError("LinSolve(TMatrix< complex<r_4> >,TVector< complex<r_4> >) size mismatch"));
-    c = b; TMatrix< complex<r_4> > a1(a);
-  return SimpleMatrixOperation< complex<r_4> >::GausPiv(a1,c);
-}
-
-/*! \ingroup TArray 
-    \fn LinSolve(const TMatrix< complex<r_8> >&, const TVector< complex<r_8> >&, TVector< complex<r_8> >&)
-    \brief  Solve A*C = B and return C and determinant
-*/
-inline complex<r_8> LinSolve(const TMatrix< complex<r_8> >& a, const TVector< complex<r_8> >& b, TVector< complex<r_8> >& c) {
-  if(a.NCols()!=b.NRows() || a.NCols()!=a.NRows())
-    throw(SzMismatchError("LinSolve(TMatrix< complex<r_8> >,TVector< complex<r_8> >) size mismatch"));
-    c = b; TMatrix< complex<r_8> > a1(a);
-  return SimpleMatrixOperation< complex<r_8> >::GausPiv(a1,c);
+    throw(SzMismatchError("LinSolve(TMatrix<T>,TVector<T>) size mismatch"));
+    c = b; TMatrix<T> a1(a);
+  return SimpleMatrixOperation<T>::GausPiv(a1,c);
 }
 
@@ -151 +85 @@
 
 /*! \ingroup TArray 
-    \fn Inverse(TMatrix<r_4> const &)
+    \fn Inverse(TMatrix<T> const &)
     \brief To inverse a TMatrix
 */
-inline TMatrix<r_4> Inverse(TMatrix<r_4> const & A)
-  {return SimpleMatrixOperation<r_4>::Inverse(A);}
+template <class T>
+inline TMatrix<T> Inverse(TMatrix<T> const & A)
+  {return SimpleMatrixOperation<T>::Inverse(A);}
+
 /*! \ingroup TArray 
-    \fn Inverse(TMatrix<r_8> const &)
-    \brief To inverse a TMatrix
+    \fn Determinant(TMatrix<T> const &)
+    \brief Give the TMatrix determinant
 */
-inline TMatrix<r_8> Inverse(TMatrix<r_8> const & A)
-  {return SimpleMatrixOperation<r_8>::Inverse(A);}
+template <class T>
+inline T Determinant(TMatrix<T> const & A) {
+  TMatrix<T> a(A,false);
+  TMatrix<T> b(a.NCols(),a.NRows());  b = IdentityMatrix(1.);
+  return SimpleMatrixOperation<T>::GausPiv(a,b);
+}
+
 /*! \ingroup TArray 
-    \fn Inverse(TMatrix< complex<r_4> > const &)
-    \brief To inverse a TMatrix (complex numbers)
+    \fn GausPiv(TMatrix<T> const &,TMatrix<T> &)
+    \brief Gauss pivoting on matrix \b A, doing the same operations
+           on matrix \b B and return determinant of \b A.
+    \sa SOPHYA::SimpleMatrixOperation::GausPiv(TMatrix<T>&,TMatrix<T>&)
 */
-inline TMatrix< complex<r_4> > Inverse(TMatrix< complex<r_4> > const & A)
-  {return SimpleMatrixOperation< complex<r_4> >::Inverse(A);}
-/*! \ingroup TArray 
-    \fn Inverse(TMatrix< complex<r_8> > const &)
-    \brief To inverse a TMatrix (complex numbers)
-*/
-inline TMatrix< complex<r_8> > Inverse(TMatrix< complex<r_8> > const & A)
-  {return SimpleMatrixOperation< complex<r_8> >::Inverse(A);}
+template <class T>
+inline T GausPiv(TMatrix<T> const & A,TMatrix<T> & B) {
+  TMatrix<T> a(A,false);
+  return SimpleMatrixOperation<T>::GausPiv(a,B);
+}
+
 
 } // Fin du namespace