From 67ffcd30d8b6f80db1c3787f6a1e77fbad6fc84d Mon Sep 17 00:00:00 2001 From: kanschat Date: Mon, 1 Nov 2010 03:30:38 +0000 Subject: [PATCH] undo commit that destroyed FullMatrix git-svn-id: https://svn.dealii.org/trunk@22562 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/base/table.h | 590 +++++++++++----------- deal.II/include/deal.II/lac/full_matrix.h | 98 ++-- deal.II/source/lac/lapack_full_matrix.cc | 64 +-- 3 files changed, 387 insertions(+), 365 deletions(-) diff --git a/deal.II/include/deal.II/base/table.h b/deal.II/include/deal.II/base/table.h index 3809abbb79..3daf142709 100644 --- a/deal.II/include/deal.II/base/table.h +++ b/deal.II/include/deal.II/base/table.h @@ -17,7 +17,6 @@ #include #include #include -#include #include #include @@ -77,16 +76,10 @@ namespace internal * template parameters. Specialization for accessors to constant * objects. */ - template struct Types + template struct Types { typedef const T value_type; typedef const TableBase TableType; - - typedef typename std::vector::const_iterator iterator; - typedef typename std::vector::const_iterator const_iterator; - - typedef typename std::vector::const_reference reference; - typedef typename std::vector::const_reference const_reference; }; /** @@ -95,18 +88,12 @@ namespace internal * template parameters. Specialization for accessors to non-constant * objects. */ - template struct Types + template struct Types { typedef T value_type; typedef TableBase TableType; - - typedef typename std::vector::iterator iterator; - typedef typename std::vector::const_iterator const_iterator; - - typedef typename std::vector::reference reference; - typedef typename std::vector::const_reference const_reference; }; - + /** * @internal @@ -144,17 +131,19 @@ namespace internal * generate objects of this class, as they are only thought as * temporaries for access to elements of the table class, not for * passing them around as arguments of functions, etc. - * + * * @author Wolfgang Bangerth, 2002 */ template class Accessor { public: - typedef typename Types::TableType TableType; - - typedef typename Types::iterator iterator; - typedef typename Types::const_iterator const_iterator; + /** + * Import two typedefs from the + * switch class above. + */ + typedef typename Types::value_type * pointer; + typedef typename Types::TableType TableType; private: /** @@ -166,7 +155,7 @@ namespace internal * data we may access. */ Accessor (const TableType &table, - const iterator data); + const pointer data); /** * Default constructor. Not @@ -174,9 +163,9 @@ namespace internal * private. */ Accessor (); - + public: - + /** * Copy constructor. This constructor * is public so that one can pass @@ -218,13 +207,13 @@ namespace internal * elements constant. */ const TableType &table; - const iterator data; + const pointer data; // declare some other classes // as friends. make sure to // work around bugs in some // compilers -#ifndef DEAL_II_NAMESP_TEMPL_FRIEND_BUG +#ifndef DEAL_II_NAMESP_TEMPL_FRIEND_BUG template friend class dealii::Table; template friend class Accessor; @@ -239,7 +228,7 @@ namespace internal }; - + /** * @internal * Accessor class for tables. This is the specialization for the last @@ -264,13 +253,12 @@ namespace internal * library algorithms. */ typedef typename Types::value_type value_type; - - typedef typename Types::iterator iterator; - typedef typename Types::const_iterator const_iterator; - - typedef typename Types::reference reference; - typedef typename Types::const_reference const_reference; - + typedef value_type* pointer; + typedef const value_type* const_pointer; + typedef value_type* iterator; + typedef const value_type* const_iterator; + typedef value_type& reference; + typedef const value_type& const_reference; typedef size_t size_type; typedef ptrdiff_t difference_type; @@ -281,13 +269,13 @@ namespace internal typedef typename Types::TableType TableType; private: - + /** - * Constructor. Take an iterator + * Constructor. Take a pointer * to the table object to know * about the sizes of the * various dimensions, and a - * iterator to the subset of + * pointer to the subset of * data we may access (which in * this particular case is only * one row). @@ -307,7 +295,7 @@ namespace internal * data consistency. */ Accessor (const TableType &table, - const iterator data); + const pointer data); /** * Default constructor. Not needed, @@ -315,7 +303,7 @@ namespace internal */ Accessor (); - public: + public: /** * Copy constructor. This constructor * is public so that one can pass @@ -330,13 +318,13 @@ namespace internal */ Accessor (const Accessor &a); - + /** * Index operator. Performs a * range check. */ reference operator [] (const unsigned int) const; - + /** * Return the length of one row, * i.e. the number of elements @@ -344,21 +332,21 @@ namespace internal * index of the table object. */ unsigned int size () const; - + /** * Return an iterator to the * first element of this * row. */ iterator begin () const; - + /** * Return an interator to the * element past the end of * this row. */ iterator end () const; - + private: /** * Store the data given to the @@ -369,7 +357,7 @@ namespace internal * elements constant. */ const TableType &table; - const iterator data; + const pointer data; // declare some other classes // as friends. make sure to @@ -389,7 +377,7 @@ namespace internal #endif }; } - + } // namespace internal @@ -478,7 +466,7 @@ class TableBase : public Subscriptor * dimensions to zero. */ TableBase (); - + /** * Constructor. Initialize the * array with the given @@ -500,12 +488,12 @@ class TableBase : public Subscriptor */ template TableBase (const TableBase &src); - + /** * Destructor. Free allocated memory. */ ~TableBase (); - + /** * Assignment operator. * Copy all elements of src @@ -519,7 +507,7 @@ class TableBase : public Subscriptor * operator which is not what we want. */ TableBase& operator = (const TableBase& src); - + /** * Copy operator. * Copy all elements of src @@ -532,7 +520,7 @@ class TableBase : public Subscriptor */ template TableBase& operator = (const TableBase &src); - + /** * Set all entries to their * default value (i.e. copy them @@ -541,7 +529,7 @@ class TableBase : public Subscriptor * size of the table, though. */ void reset_values (); - + /** * Set the dimensions of this object to * the sizes given in the argument, and @@ -559,13 +547,13 @@ class TableBase : public Subscriptor * i. */ unsigned int size (const unsigned int i) const; - + /** * Return the sizes of this * object in each direction. */ const TableIndices & size () const; - + /** * Return the number of elements * stored in this object, which @@ -582,7 +570,7 @@ class TableBase : public Subscriptor * n_elements()==0. */ bool empty () const; - + /** * Fill array with an array of * elements. The input array must @@ -615,14 +603,13 @@ class TableBase : public Subscriptor * the same value. */ void fill (const T& value); - + /** * Return a read-write reference * to the indicated element. */ - typename std::vector::reference - operator () (const TableIndices &indices); - + T & operator() (const TableIndices &indices); + /** * Return the value of the * indicated element as a @@ -636,8 +623,7 @@ class TableBase : public Subscriptor * don't know here whether * copying is expensive or not. */ - typename std::vector::const_reference - operator () (const TableIndices &indices) const; + const T & operator() (const TableIndices &indices) const; /** * Determine an estimate for the @@ -655,7 +641,7 @@ class TableBase : public Subscriptor * does no index checking. */ unsigned int position (const TableIndices &indices) const; - + /** * Return a read-write reference * to the indicated element. @@ -665,8 +651,8 @@ class TableBase : public Subscriptor * used internally and in * functions already checked. */ - typename std::vector::reference el (const TableIndices &indices); - + T & el (const TableIndices &indices); + /** * Return the value of the * indicated element as a @@ -685,8 +671,8 @@ class TableBase : public Subscriptor * don't know here whether * copying is expensive or not. */ - typename std::vector::const_reference el (const TableIndices &indices) const; - + const T & el (const TableIndices &indices) const; + /** * @deprecated This function * accesses data directly and @@ -699,13 +685,23 @@ class TableBase : public Subscriptor * cast from const), otherwise, * keep away! */ - typename std::vector::const_pointer data () const; - + const T* data () const; + protected: /** * Component-array. */ - std::vector values; + T* val; + + /** + * Size of array. This may be + * larger than the number of + * actually used elements, since + * we don't shrink the array upon + * calls to reinit unless the + * new size is zero. + */ + unsigned int val_size; /** * Size in each direction of the @@ -731,7 +727,7 @@ class TableBase : public Subscriptor * in partial specializations of this class, with fixed numbers of * dimensions. See there, and in the documentation of the base class * for more information. - * + * * @ingroup data * @author Wolfgang Bangerth, 2002 */ @@ -750,7 +746,7 @@ class Table : public TableBase * * For the rationale of this class, and a description of the * interface, see the base class. - * + * * @ingroup data * @author Wolfgang Bangerth, 2002 */ @@ -763,7 +759,7 @@ class Table<1,T> : public TableBase<1,T> * dimensions to zero. */ Table (); - + /** * Constructor. Pass down the * given dimension to the base @@ -778,7 +774,7 @@ class Table<1,T> : public TableBase<1,T> * data element. Returns a * read-only reference. */ - typename std::vector::const_reference + const T & operator [] (const unsigned int i) const; /** @@ -788,7 +784,7 @@ class Table<1,T> : public TableBase<1,T> * data element. Returns a * read-write reference. */ - typename std::vector::reference + T & operator [] (const unsigned int i); /** @@ -798,7 +794,7 @@ class Table<1,T> : public TableBase<1,T> * data element. Returns a * read-only reference. */ - typename std::vector::const_reference + const T & operator () (const unsigned int i) const; /** @@ -808,26 +804,26 @@ class Table<1,T> : public TableBase<1,T> * data element. Returns a * read-write reference. */ - typename std::vector::reference + T & operator () (const unsigned int i); /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference - operator () (const TableIndices<1> &indices); - + T & + operator() (const TableIndices<1> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference - operator () (const TableIndices<1> &indices) const; + const T & + operator() (const TableIndices<1> &indices) const; }; @@ -841,7 +837,7 @@ class Table<1,T> : public TableBase<1,T> * of the full matrix classes in this library, and to keep a minimal * compatibility with a predecessor class (vector2d), some * additional functions are provided. - * + * * @ingroup data * @author Wolfgang Bangerth, 2002 */ @@ -912,9 +908,9 @@ class Table<2,T> : public TableBase<2,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j) const; - + /** * Direct access to one element @@ -925,28 +921,28 @@ class Table<2,T> : public TableBase<2,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j); /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference - operator () (const TableIndices<2> &indices); - + T & + operator() (const TableIndices<2> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference - operator () (const TableIndices<2> &indices) const; - + const T & + operator() (const TableIndices<2> &indices) const; + /** * Number of rows. This function * really makes only sense since @@ -954,7 +950,7 @@ class Table<2,T> : public TableBase<2,T> * object here. */ unsigned int n_rows () const; - + /** * Number of columns. This function * really makes only sense since @@ -979,9 +975,9 @@ class Table<2,T> : public TableBase<2,T> * table classes for 2d arrays, * then called vector2d. */ - typename std::vector::reference el (const unsigned int i, + T & el (const unsigned int i, const unsigned int j); - + /** * Return the value of the * element (i,j) as a @@ -1006,7 +1002,7 @@ class Table<2,T> : public TableBase<2,T> * table classes for 2d arrays, * then called vector2d. */ - typename std::vector::const_reference el (const unsigned int i, + const T & el (const unsigned int i, const unsigned int j) const; }; @@ -1018,7 +1014,7 @@ class Table<2,T> : public TableBase<2,T> * * For the rationale of this class, and a description of the * interface, see the base class. - * + * * @ingroup data * @author Wolfgang Bangerth, 2002 */ @@ -1078,10 +1074,10 @@ class Table<3,T> : public TableBase<3,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j, const unsigned int k) const; - + /** * Direct access to one element @@ -1092,25 +1088,25 @@ class Table<3,T> : public TableBase<3,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j, const unsigned int k); /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference operator () (const TableIndices<3> &indices); - + T & operator() (const TableIndices<3> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference operator () (const TableIndices<3> &indices) const; + const T & operator() (const TableIndices<3> &indices) const; }; @@ -1121,7 +1117,7 @@ class Table<3,T> : public TableBase<3,T> * * For the rationale of this class, and a description of the * interface, see the base class. - * + * * @ingroup data * @author Wolfgang Bangerth, Ralf Hartmann 2002 */ @@ -1144,7 +1140,7 @@ class Table<4,T> : public TableBase<4,T> const unsigned int size2, const unsigned int size3, const unsigned int size4); - + /** * Access operator. Generate an * object that accesses the @@ -1182,11 +1178,11 @@ class Table<4,T> : public TableBase<4,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l) const; - + /** * Direct access to one element @@ -1197,28 +1193,28 @@ class Table<4,T> : public TableBase<4,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l); /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference - operator () (const TableIndices<4> &indices); - + T & + operator() (const TableIndices<4> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference - operator () (const TableIndices<4> &indices) const; + const T & + operator() (const TableIndices<4> &indices) const; }; @@ -1229,7 +1225,7 @@ class Table<4,T> : public TableBase<4,T> * * For the rationale of this class, and a description of the * interface, see the base class. - * + * * @ingroup data * @author Wolfgang Bangerth, Ralf Hartmann 2002 */ @@ -1291,12 +1287,12 @@ class Table<5,T> : public TableBase<5,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l, const unsigned int m) const; - + /** * Direct access to one element * of the table by specifying all @@ -1306,7 +1302,7 @@ class Table<5,T> : public TableBase<5,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l, @@ -1314,21 +1310,21 @@ class Table<5,T> : public TableBase<5,T> /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference - operator () (const TableIndices<5> &indices); - + T & + operator() (const TableIndices<5> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference - operator () (const TableIndices<5> &indices) const; + const T & + operator() (const TableIndices<5> &indices) const; }; @@ -1339,7 +1335,7 @@ class Table<5,T> : public TableBase<5,T> * * For the rationale of this class, and a description of the * interface, see the base class. - * + * * @ingroup data * @author Wolfgang Bangerth, Ralf Hartmann 2002 */ @@ -1402,12 +1398,12 @@ class Table<6,T> : public TableBase<6,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l, const unsigned int m, - const unsigned int n) const; + const unsigned int n) const; /** * Direct access to one element @@ -1418,7 +1414,7 @@ class Table<6,T> : public TableBase<6,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l, @@ -1427,21 +1423,21 @@ class Table<6,T> : public TableBase<6,T> /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference - operator () (const TableIndices<6> &indices); - + T & + operator() (const TableIndices<6> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference - operator () (const TableIndices<6> &indices) const; + const T & + operator() (const TableIndices<6> &indices) const; }; @@ -1451,7 +1447,7 @@ class Table<6,T> : public TableBase<6,T> * * For the rationale of this class, and a description of the * interface, see the base class. - * + * * @ingroup data * @author Wolfgang Bangerth, 2002, Ralf Hartmann 2004 */ @@ -1515,13 +1511,13 @@ class Table<7,T> : public TableBase<7,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l, const unsigned int m, const unsigned int n, - const unsigned int o) const; + const unsigned int o) const; /** * Direct access to one element @@ -1532,7 +1528,7 @@ class Table<7,T> : public TableBase<7,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j, const unsigned int k, const unsigned int l, @@ -1542,21 +1538,21 @@ class Table<7,T> : public TableBase<7,T> /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::reference - operator () (const TableIndices<7> &indices); - + T & + operator() (const TableIndices<7> &indices); + /** * Make the corresponding - * operator () from the TableBase + * operator() from the TableBase * base class available also in * this class. */ - typename std::vector::const_reference - operator () (const TableIndices<7> &indices) const; + const T & + operator() (const TableIndices<7> &indices) const; }; @@ -1615,7 +1611,7 @@ class TransposeTable : public TableBase<2,T> * This version of the function * only allows read access. */ - typename std::vector::const_reference operator () (const unsigned int i, + const T & operator () (const unsigned int i, const unsigned int j) const; /** @@ -1627,9 +1623,9 @@ class TransposeTable : public TableBase<2,T> * This version of the function * allows read-write access. */ - typename std::vector::reference operator () (const unsigned int i, + T & operator () (const unsigned int i, const unsigned int j); - + /** * Number of rows. This function * really makes only sense since @@ -1637,7 +1633,7 @@ class TransposeTable : public TableBase<2,T> * object here. */ unsigned int n_rows () const; - + /** * Number of columns. This function * really makes only sense since @@ -1662,9 +1658,9 @@ class TransposeTable : public TableBase<2,T> * table classes for 2d arrays, * then called vector2d. */ - typename std::vector::reference el (const unsigned int i, + T & el (const unsigned int i, const unsigned int j); - + /** * Return the value of the * element (i,j) as a @@ -1689,7 +1685,7 @@ class TransposeTable : public TableBase<2,T> * table classes for 2d arrays, * then called vector2d. */ - typename std::vector::const_reference el (const unsigned int i, + const T & el (const unsigned int i, const unsigned int j) const; }; @@ -1702,12 +1698,18 @@ class TransposeTable : public TableBase<2,T> template TableBase::TableBase () + : + val (0), + val_size (0) {} template TableBase::TableBase (const TableIndices &sizes) + : + val (0), + val_size (0) { reinit (sizes); } @@ -1717,11 +1719,13 @@ TableBase::TableBase (const TableIndices &sizes) template TableBase::TableBase (const TableBase &src) : - Subscriptor () + Subscriptor (), + val (0), + val_size (0) { reinit (src.table_size); if (src.n_elements() != 0) - std::copy (src.values.begin(), src.values.end(), values.begin()); + fill (src.data()); } @@ -1729,22 +1733,25 @@ TableBase::TableBase (const TableBase &src) template template TableBase::TableBase (const TableBase &src) + : + val (0), + val_size (0) { reinit (src.table_size); if (src.n_elements() != 0) - std::copy (src.values.begin(), src.values.end(), values.begin()); + fill (src.data()); } namespace internal { - namespace TableBaseAccessors + namespace TableBaseAccessors { template inline Accessor::Accessor (const TableType &table, - const iterator data) + const pointer data) : table (table), data (data) @@ -1759,9 +1766,9 @@ namespace internal table (a.table), data (a.data) {} + - - + template inline Accessor::Accessor () @@ -1774,9 +1781,9 @@ namespace internal // not need to be copied around Assert (false, ExcInternalError()); } + - - + template inline Accessor @@ -1790,7 +1797,7 @@ namespace internal // case i==0 if (i==0) return Accessor (table, data); - else + else { // note: P>1, otherwise the // specialization would have @@ -1798,7 +1805,7 @@ namespace internal unsigned int subobject_size = table.size()[N-1]; for (int p=P-1; p>1; --p) subobject_size *= table.size()[N-p]; - const iterator new_data = data + i*subobject_size; + const pointer new_data = data + i*subobject_size; return Accessor (table, new_data); } } @@ -1808,7 +1815,7 @@ namespace internal template inline Accessor::Accessor (const TableType &table, - const iterator data) + const pointer data) : table (table), data (data) @@ -1828,7 +1835,7 @@ namespace internal // not need to be copied around Assert (false, ExcInternalError()); } - + template @@ -1840,7 +1847,7 @@ namespace internal {} - + template inline typename Accessor::reference @@ -1852,7 +1859,7 @@ namespace internal } - + template inline unsigned int @@ -1888,7 +1895,10 @@ namespace internal template inline TableBase::~TableBase () -{} +{ + if (val != 0) + delete[] val; +} @@ -1899,8 +1909,8 @@ TableBase::operator = (const TableBase& m) { reinit (m.size()); if (!empty()) - std::copy (&m.values[0], &m.values[n_elements()], &values[0]); - + std::copy (&m.val[0], &m.val[n_elements()], &val[0]); + return *this; } @@ -1914,8 +1924,8 @@ TableBase::operator = (const TableBase& m) { reinit (m.size()); if (!empty()) - std::copy (&m.values[0], &m.values[n_elements()], &values[0]); - + std::copy (&m.val[0], &m.val[n_elements()], &val[0]); + return *this; } @@ -1927,7 +1937,7 @@ void TableBase::reset_values () { if (n_elements() != 0) - std::fill (values.begin(), values.end(), T()); + std::fill_n (val, n_elements(), T()); } @@ -1938,7 +1948,7 @@ void TableBase::fill (const T& value) { if (n_elements() != 0) - std::fill (values.begin(), values.end(), T()); + std::fill_n (val, n_elements(), value); } @@ -1951,14 +1961,19 @@ TableBase::reinit (const TableIndices &new_sizes, const bool fast) { table_size = new_sizes; - + const unsigned int new_size = n_elements(); - + // if zero size was given: free all // memory if (new_size == 0) { - values.resize (0); + if (val != 0) + delete[] val; + + val = 0; + val_size = 0; + // set all sizes to zero, even // if one was previously // nonzero. This simplifies @@ -1967,11 +1982,18 @@ TableBase::reinit (const TableIndices &new_sizes, return; } - + // if new size is nonzero: // if necessary allocate // additional memory - values.resize (new_size); + if (val_size::fill (const T2* entries) Assert (n_elements() != 0, ExcMessage("Trying to fill an empty matrix.")); - std::copy (entries, entries+n_elements(), values.begin()); + std::copy (entries, entries+n_elements(), val); } @@ -2056,7 +2078,7 @@ inline unsigned int TableBase::memory_consumption () const { - return sizeof(*this) + MemoryConsumption::memory_consumption(values); + return sizeof(*this) + val_size*sizeof(T); } @@ -2094,8 +2116,8 @@ TableBase::position (const TableIndices &indices) const template inline -typename std::vector::const_reference -TableBase::operator () (const TableIndices &indices) const +const T & +TableBase::operator() (const TableIndices &indices) const { for (unsigned int n=0; n::operator () (const TableIndices &indices) const template inline -typename std::vector::reference -TableBase::operator () (const TableIndices &indices) +T & +TableBase::operator() (const TableIndices &indices) { for (unsigned int n=0; n::operator () (const TableIndices &indices) template inline -typename std::vector::const_reference +const T & TableBase::el (const TableIndices &indices) const -{ - return values[position(indices)]; +{ + return val[position(indices)]; } template inline -typename std::vector::reference +T & TableBase::el (const TableIndices &indices) { - return values[position(indices)]; + return val[position(indices)]; } template inline -typename std::vector::const_pointer +const T * TableBase::data () const { - return &values[0]; + return val; } @@ -2167,68 +2189,68 @@ Table<1,T>::Table (const unsigned int size) template inline -typename std::vector::const_reference +const T & Table<1,T>::operator [] (const unsigned int i) const { Assert (i < this->table_size[0], ExcIndexRange (i, 0, this->table_size[0])); - return this->values[i]; + return this->val[i]; } template inline -typename std::vector::reference +T & Table<1,T>::operator [] (const unsigned int i) { Assert (i < this->table_size[0], ExcIndexRange (i, 0, this->table_size[0])); - return this->values[i]; + return this->val[i]; } template inline -typename std::vector::const_reference +const T & Table<1,T>::operator () (const unsigned int i) const { Assert (i < this->table_size[0], ExcIndexRange (i, 0, this->table_size[0])); - return this->values[i]; + return this->val[i]; } template inline -typename std::vector::reference +T & Table<1,T>::operator () (const unsigned int i) { Assert (i < this->table_size[0], ExcIndexRange (i, 0, this->table_size[0])); - return this->values[i]; + return this->val[i]; } template inline -typename std::vector::const_reference +const T & Table<1,T>::operator () (const TableIndices<1> &indices) const { - return TableBase<1,T>::operator () (indices); + return TableBase<1,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<1,T>::operator () (const TableIndices<1> &indices) { - return TableBase<1,T>::operator () (indices); + return TableBase<1,T>::operator() (indices); } @@ -2271,7 +2293,7 @@ Table<2,T>::operator [] (const unsigned int i) const Assert (i < this->table_size[0], ExcIndexRange (i, 0, this->table_size[0])); return internal::TableBaseAccessors::Accessor<2,T,true,1>(*this, - this->values.begin()+i*n_cols()); + this->val+i*n_cols()); } @@ -2284,14 +2306,14 @@ Table<2,T>::operator [] (const unsigned int i) Assert (i < this->table_size[0], ExcIndexRange (i, 0, this->table_size[0])); return internal::TableBaseAccessors::Accessor<2,T,false,1>(*this, - this->values.begin()+i*n_cols()); + this->val+i*n_cols()); } template inline -typename std::vector::const_reference +const T & Table<2,T>::operator () (const unsigned int i, const unsigned int j) const { @@ -2299,14 +2321,14 @@ Table<2,T>::operator () (const unsigned int i, ExcIndexRange (i, 0, this->table_size[0])); Assert (j < this->table_size[1], ExcIndexRange (j, 0, this->table_size[1])); - return this->values[i*this->table_size[1]+j]; + return this->val[i*this->table_size[1]+j]; } template inline -typename std::vector::reference +T & Table<2,T>::operator () (const unsigned int i, const unsigned int j) { @@ -2314,49 +2336,49 @@ Table<2,T>::operator () (const unsigned int i, ExcIndexRange (i, 0, this->table_size[0])); Assert (j < this->table_size[1], ExcIndexRange (j, 0, this->table_size[1])); - return this->values[i*this->table_size[1]+j]; + return this->val[i*this->table_size[1]+j]; } template inline -typename std::vector::const_reference +const T & Table<2,T>::operator () (const TableIndices<2> &indices) const { - return TableBase<2,T>::operator () (indices); + return TableBase<2,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<2,T>::operator () (const TableIndices<2> &indices) { - return TableBase<2,T>::operator () (indices); + return TableBase<2,T>::operator() (indices); } template inline -typename std::vector::const_reference +const T & Table<2,T>::el (const unsigned int i, const unsigned int j) const { - return this->values[i*this->table_size[1]+j]; + return this->val[i*this->table_size[1]+j]; } template inline -typename std::vector::reference +T & Table<2,T>::el (const unsigned int i, const unsigned int j) { - return this->values[i*this->table_size[1]+j]; + return this->val[i*this->table_size[1]+j]; } @@ -2414,7 +2436,7 @@ TransposeTable::reinit (const unsigned int size1, template inline -typename std::vector::const_reference +const T & TransposeTable::operator () (const unsigned int i, const unsigned int j) const { @@ -2422,14 +2444,14 @@ TransposeTable::operator () (const unsigned int i, ExcIndexRange (i, 0, this->table_size[1])); Assert (j < this->table_size[0], ExcIndexRange (j, 0, this->table_size[0])); - return this->values[j*this->table_size[1]+i]; + return this->val[j*this->table_size[1]+i]; } template inline -typename std::vector::reference +T & TransposeTable::operator () (const unsigned int i, const unsigned int j) { @@ -2437,29 +2459,29 @@ TransposeTable::operator () (const unsigned int i, ExcIndexRange (i, 0, this->table_size[1])); Assert (j < this->table_size[0], ExcIndexRange (j, 0, this->table_size[0])); - return this->values[j*this->table_size[1]+i]; + return this->val[j*this->table_size[1]+i]; } template inline -typename std::vector::const_reference +const T & TransposeTable::el (const unsigned int i, const unsigned int j) const { - return this->values[j*this->table_size[1]+i]; + return this->val[j*this->table_size[1]+i]; } template inline -typename std::vector::reference +T & TransposeTable::el (const unsigned int i, const unsigned int j) { - return this->values[j*this->table_size[1]+i]; + return this->val[j*this->table_size[1]+i]; } @@ -2516,7 +2538,7 @@ Table<3,T>::operator [] (const unsigned int i) const this->table_size[2]; return (internal::TableBaseAccessors::Accessor<3,T,true,2> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } @@ -2532,14 +2554,14 @@ Table<3,T>::operator [] (const unsigned int i) this->table_size[2]; return (internal::TableBaseAccessors::Accessor<3,T,false,2> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } template inline -typename std::vector::const_reference +const T & Table<3,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k) const @@ -2550,7 +2572,7 @@ Table<3,T>::operator () (const unsigned int i, ExcIndexRange (j, 0, this->table_size[1])); Assert (k < this->table_size[2], ExcIndexRange (k, 0, this->table_size[2])); - return this->values[(i*this->table_size[1]+j) + return this->val[(i*this->table_size[1]+j) *this->table_size[2] + k]; } @@ -2558,7 +2580,7 @@ Table<3,T>::operator () (const unsigned int i, template inline -typename std::vector::reference +T & Table<3,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k) @@ -2569,7 +2591,7 @@ Table<3,T>::operator () (const unsigned int i, ExcIndexRange (j, 0, this->table_size[1])); Assert (k < this->table_size[2], ExcIndexRange (k, 0, this->table_size[2])); - return this->values[(i*this->table_size[1]+j) + return this->val[(i*this->table_size[1]+j) *this->table_size[2] + k]; } @@ -2577,20 +2599,20 @@ Table<3,T>::operator () (const unsigned int i, template inline -typename std::vector::const_reference +const T & Table<3,T>::operator () (const TableIndices<3> &indices) const { - return TableBase<3,T>::operator () (indices); + return TableBase<3,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<3,T>::operator () (const TableIndices<3> &indices) { - return TableBase<3,T>::operator () (indices); + return TableBase<3,T>::operator() (indices); } @@ -2626,7 +2648,7 @@ Table<4,T>::operator [] (const unsigned int i) const this->table_size[3]; return (internal::TableBaseAccessors::Accessor<4,T,true,3> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } @@ -2643,14 +2665,14 @@ Table<4,T>::operator [] (const unsigned int i) this->table_size[3]; return (internal::TableBaseAccessors::Accessor<4,T,false,3> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } template inline -typename std::vector::const_reference +const T & Table<4,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -2664,7 +2686,7 @@ Table<4,T>::operator () (const unsigned int i, ExcIndexRange (k, 0, this->table_size[2])); Assert (l < this->table_size[3], ExcIndexRange (l, 0, this->table_size[3])); - return this->values[((i*this->table_size[1]+j) + return this->val[((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l]; } @@ -2673,7 +2695,7 @@ Table<4,T>::operator () (const unsigned int i, template inline -typename std::vector::reference +T & Table<4,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -2687,7 +2709,7 @@ Table<4,T>::operator () (const unsigned int i, ExcIndexRange (k, 0, this->table_size[2])); Assert (l < this->table_size[3], ExcIndexRange (l, 0, this->table_size[3])); - return this->values[((i*this->table_size[1]+j) + return this->val[((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l]; } @@ -2696,20 +2718,20 @@ Table<4,T>::operator () (const unsigned int i, template inline -typename std::vector::const_reference +const T & Table<4,T>::operator () (const TableIndices<4> &indices) const { - return TableBase<4,T>::operator () (indices); + return TableBase<4,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<4,T>::operator () (const TableIndices<4> &indices) { - return TableBase<4,T>::operator () (indices); + return TableBase<4,T>::operator() (indices); } @@ -2748,7 +2770,7 @@ Table<5,T>::operator [] (const unsigned int i) const this->table_size[4]; return (internal::TableBaseAccessors::Accessor<5,T,true,4> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } @@ -2766,14 +2788,14 @@ Table<5,T>::operator [] (const unsigned int i) this->table_size[4]; return (internal::TableBaseAccessors::Accessor<5,T,false,4> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } template inline -typename std::vector::const_reference +const T & Table<5,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -2790,7 +2812,7 @@ Table<5,T>::operator () (const unsigned int i, ExcIndexRange (l, 0, this->table_size[3])); Assert (m < this->table_size[4], ExcIndexRange (m, 0, this->table_size[4])); - return this->values[(((i*this->table_size[1]+j) + return this->val[(((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l) *this->table_size[4] + m]; @@ -2800,7 +2822,7 @@ Table<5,T>::operator () (const unsigned int i, template inline -typename std::vector::reference +T & Table<5,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -2817,7 +2839,7 @@ Table<5,T>::operator () (const unsigned int i, ExcIndexRange (l, 0, this->table_size[3])); Assert (m < this->table_size[4], ExcIndexRange (m, 0, this->table_size[4])); - return this->values[(((i*this->table_size[1]+j) + return this->val[(((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l) *this->table_size[4] + m]; @@ -2827,20 +2849,20 @@ Table<5,T>::operator () (const unsigned int i, template inline -typename std::vector::const_reference +const T & Table<5,T>::operator () (const TableIndices<5> &indices) const { - return TableBase<5,T>::operator () (indices); + return TableBase<5,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<5,T>::operator () (const TableIndices<5> &indices) { - return TableBase<5,T>::operator () (indices); + return TableBase<5,T>::operator() (indices); } @@ -2880,7 +2902,7 @@ Table<6,T>::operator [] (const unsigned int i) const this->table_size[5]; return (internal::TableBaseAccessors::Accessor<6,T,true,5> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } @@ -2899,14 +2921,14 @@ Table<6,T>::operator [] (const unsigned int i) this->table_size[5]; return (internal::TableBaseAccessors::Accessor<6,T,false,5> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } template inline -typename std::vector::const_reference +const T & Table<6,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -2926,7 +2948,7 @@ Table<6,T>::operator () (const unsigned int i, ExcIndexRange (m, 0, this->table_size[4])); Assert (n < this->table_size[5], ExcIndexRange (n, 0, this->table_size[5])); - return this->values[((((i*this->table_size[1]+j) + return this->val[((((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l) *this->table_size[4] + m) @@ -2937,7 +2959,7 @@ Table<6,T>::operator () (const unsigned int i, template inline -typename std::vector::reference +T & Table<6,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -2957,7 +2979,7 @@ Table<6,T>::operator () (const unsigned int i, ExcIndexRange (m, 0, this->table_size[4])); Assert (n < this->table_size[5], ExcIndexRange (n, 0, this->table_size[5])); - return this->values[((((i*this->table_size[1]+j) + return this->val[((((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l) *this->table_size[4] + m) @@ -2968,20 +2990,20 @@ Table<6,T>::operator () (const unsigned int i, template inline -typename std::vector::const_reference +const T & Table<6,T>::operator () (const TableIndices<6> &indices) const { - return TableBase<6,T>::operator () (indices); + return TableBase<6,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<6,T>::operator () (const TableIndices<6> &indices) { - return TableBase<6,T>::operator () (indices); + return TableBase<6,T>::operator() (indices); } @@ -3023,7 +3045,7 @@ Table<7,T>::operator [] (const unsigned int i) const this->table_size[6]; return (internal::TableBaseAccessors::Accessor<7,T,true,6> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } @@ -3043,14 +3065,14 @@ Table<7,T>::operator [] (const unsigned int i) this->table_size[6]; return (internal::TableBaseAccessors::Accessor<7,T,false,6> (*this, - this->values.begin() + i*subobject_size)); + this->val+i*subobject_size)); } template inline -typename std::vector::const_reference +const T & Table<7,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -3073,7 +3095,7 @@ Table<7,T>::operator () (const unsigned int i, ExcIndexRange (n, 0, this->table_size[5])); Assert (o < this->table_size[6], ExcIndexRange (o, 0, this->table_size[6])); - return this->values[(((((i*this->table_size[1]+j) + return this->val[(((((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l) *this->table_size[4] + m) @@ -3085,7 +3107,7 @@ Table<7,T>::operator () (const unsigned int i, template inline -typename std::vector::reference +T & Table<7,T>::operator () (const unsigned int i, const unsigned int j, const unsigned int k, @@ -3108,7 +3130,7 @@ Table<7,T>::operator () (const unsigned int i, ExcIndexRange (n, 0, this->table_size[5])); Assert (o < this->table_size[5], ExcIndexRange (o, 0, this->table_size[6])); - return this->values[(((((i*this->table_size[1]+j) + return this->val[(((((i*this->table_size[1]+j) *this->table_size[2] + k) *this->table_size[3] + l) *this->table_size[4] + m) @@ -3120,20 +3142,20 @@ Table<7,T>::operator () (const unsigned int i, template inline -typename std::vector::const_reference +const T & Table<7,T>::operator () (const TableIndices<7> &indices) const { - return TableBase<7,T>::operator () (indices); + return TableBase<7,T>::operator() (indices); } template inline -typename std::vector::reference +T & Table<7,T>::operator () (const TableIndices<7> &indices) { - return TableBase<7,T>::operator () (indices); + return TableBase<7,T>::operator() (indices); } diff --git a/deal.II/include/deal.II/lac/full_matrix.h b/deal.II/include/deal.II/lac/full_matrix.h index a62e83af85..42bcbe880b 100644 --- a/deal.II/include/deal.II/lac/full_matrix.h +++ b/deal.II/include/deal.II/lac/full_matrix.h @@ -90,9 +90,9 @@ class FullMatrix : public Table<2,number> */ typedef typename numbers::NumberTraits::real_type real_type; - + class const_iterator; - + /** * Accessor class for iterators */ @@ -127,7 +127,7 @@ class FullMatrix : public Table<2,number> * Value of this matrix entry. */ number value() const; - + protected: /** * The matrix accessed. @@ -159,11 +159,11 @@ class FullMatrix : public Table<2,number> public: /** * Constructor. - */ + */ const_iterator(const FullMatrix *matrix, const unsigned int row, const unsigned int col); - + /** * Prefix increment. */ @@ -211,7 +211,7 @@ class FullMatrix : public Table<2,number> * operator above. */ bool operator > (const const_iterator&) const; - + private: /** * Store an object of the @@ -224,7 +224,7 @@ class FullMatrix : public Table<2,number> * See also the base class Table. */ //@{ - + /** * Constructor. Initialize the * matrix as a square matrix with @@ -240,7 +240,7 @@ class FullMatrix : public Table<2,number> * allocated. */ explicit FullMatrix (const unsigned int n = 0); - + /** * Constructor. Initialize the * matrix as a rectangular @@ -248,8 +248,8 @@ class FullMatrix : public Table<2,number> */ FullMatrix (const unsigned int rows, const unsigned int cols); - - /** + + /** * Copy constructor. This * constructor does a deep copy * of the matrix. Therefore, it @@ -293,14 +293,14 @@ class FullMatrix : public Table<2,number> * @endverbatim */ explicit FullMatrix (const IdentityMatrix &id); - + /** * Assignment operator. */ FullMatrix & operator = (const FullMatrix&); - + /** * Variable assignment operator. */ @@ -334,7 +334,7 @@ class FullMatrix : public Table<2,number> */ FullMatrix & operator = (const IdentityMatrix &id); - + /** * Assignment operator for a * LapackFullMatrix. The calling matrix @@ -344,7 +344,7 @@ class FullMatrix : public Table<2,number> template FullMatrix & operator = (const LAPACKFullMatrix&); - + /** * Assignment from different * matrix classes. This @@ -355,7 +355,7 @@ class FullMatrix : public Table<2,number> */ template void copy_from (const MATRIX&); - + /** * Fill rectangular block. * @@ -380,7 +380,7 @@ class FullMatrix : public Table<2,number> const unsigned int dst_offset_j = 0, const unsigned int src_offset_i = 0, const unsigned int src_offset_j = 0); - + /** * Make function of base class @@ -388,7 +388,7 @@ class FullMatrix : public Table<2,number> */ template void fill (const number2*); - + /** * Fill with permutation of * another matrix. @@ -415,7 +415,7 @@ class FullMatrix : public Table<2,number> //@} ///@name Non-modifying operators //@{ - + /** * Comparison operator. Be * careful with this thing, it @@ -433,7 +433,7 @@ class FullMatrix : public Table<2,number> * m x n-matrix. */ unsigned int m () const; - + /** * Number of columns of this matrix. * To remember: this matrix is an @@ -520,7 +520,7 @@ class FullMatrix : public Table<2,number> * over rows). */ real_type linfty_norm () const; - + /** * Compute the Frobenius norm of * the matrix. Return value is @@ -534,7 +534,7 @@ class FullMatrix : public Table<2,number> * the vector space. */ real_type frobenius_norm () const; - + /** * Compute the relative norm of * the skew-symmetric part. The @@ -549,7 +549,7 @@ class FullMatrix : public Table<2,number> * accuracy, or not. */ real_type relative_symmetry_norm2 () const; - + /** * Computes the determinant of a * matrix. This is only @@ -571,7 +571,7 @@ class FullMatrix : public Table<2,number> * be quadratic for this function. */ number trace () const; - + /** * Output of the matrix in * user-defined format. @@ -629,7 +629,7 @@ class FullMatrix : public Table<2,number> const char *zero_string = " ", const double denominator = 1., const double threshold = 0.) const; - + /** * Determine an estimate for the * memory consumption (in bytes) @@ -640,7 +640,7 @@ class FullMatrix : public Table<2,number> //@} ///@name Iterator functions //@{ - + /** * STL-like iterator with the * first entry. @@ -651,7 +651,7 @@ class FullMatrix : public Table<2,number> * Final iterator. */ const_iterator end () const; - + /** * STL-like iterator with the * first entry of row r. @@ -666,7 +666,7 @@ class FullMatrix : public Table<2,number> //@} ///@name Modifying operators //@{ - + /** * Scale the entire matrix by a * fixed factor. @@ -888,7 +888,7 @@ class FullMatrix : public Table<2,number> const number b, const FullMatrix &B, const number c, - const FullMatrix &C); + const FullMatrix &C); /** * Symmetrize the matrix by @@ -903,15 +903,15 @@ class FullMatrix : public Table<2,number> void symmetrize (); /** - * A=Inverse(A). A must be a square matrix. - * Inversion of + * A=Inverse(A). A must be a square matrix. + * Inversion of * this matrix by Gauss-Jordan * algorithm with partial * pivoting. This process is * well-behaved for positive * definite matrices, but be * aware of round-off errors in - * the indefinite case. + * the indefinite case. * * In case deal.II was configured with * LAPACK, the functions Xgetrf and @@ -942,7 +942,7 @@ class FullMatrix : public Table<2,number> /** * Assign the Cholesky decomposition - * of the given matrix to *this. + * of the given matrix to *this. * The given matrix must be symmetric * positive definite. * @@ -954,33 +954,33 @@ class FullMatrix : public Table<2,number> void cholesky (const FullMatrix &A); /** - * *this(i,j) = $V(i) W(j)$ + * *this(i,j) = $V(i) W(j)$ * where $V,W$ * are vectors of the same length. */ template void outer_product (const Vector &V, const Vector &W); - + /** * Assign the left_inverse of the given matrix - * to *this. The calculation being + * to *this. The calculation being * performed is (AT*A)-1 * *AT. - */ + */ template void left_invert (const FullMatrix &M); /** * Assign the right_inverse of the given matrix - * to *this. The calculation being + * to *this. The calculation being * performed is AT*(A*AT) * -1. */ template void right_invert (const FullMatrix &M); - /** + /** * Fill matrix with elements * extracted from a tensor, * taking rows included between @@ -1017,7 +1017,7 @@ class FullMatrix : public Table<2,number> * that of the matrix coincide. */ template - void + void copy_to(Tensor<2,dim> &T, const unsigned int src_r_i=0, const unsigned int src_r_j=dim-1, @@ -1061,7 +1061,7 @@ class FullMatrix : public Table<2,number> void mmult (FullMatrix &C, const FullMatrix &B, const bool adding=false) const; - + /** * Matrix-matrix-multiplication using * transpose of this. @@ -1125,7 +1125,7 @@ class FullMatrix : public Table<2,number> void mTmult (FullMatrix &C, const FullMatrix &B, const bool adding=false) const; - + /** * Matrix-matrix-multiplication using * transpose of this and @@ -1158,7 +1158,7 @@ class FullMatrix : public Table<2,number> void TmTmult (FullMatrix &C, const FullMatrix &B, const bool adding=false) const; - + /** * Matrix-vector-multiplication. * @@ -1180,7 +1180,7 @@ class FullMatrix : public Table<2,number> void vmult (Vector &w, const Vector &v, const bool adding=false) const; - + /** * Adding Matrix-vector-multiplication. * w += A*v @@ -1191,7 +1191,7 @@ class FullMatrix : public Table<2,number> template void vmult_add (Vector &w, const Vector &v) const; - + /** * Transpose * matrix-vector-multiplication. @@ -1242,7 +1242,7 @@ class FullMatrix : public Table<2,number> void precondition_Jacobi (Vector &dst, const Vector &src, const number omega = 1.) const; - + /** * dst=b-A*x. Residual calculation, * returns the l2-norm @@ -1326,7 +1326,7 @@ class FullMatrix : public Table<2,number> */ DeclException0 (ExcMatrixNotPositiveDefinite); //@} - + friend class Accessor; }; @@ -1365,7 +1365,7 @@ FullMatrix::operator = (const number d) Assert (d==number(0), ExcScalarAssignmentOnlyForZeroValue()); if (this->n_elements() != 0) - memset (&this->values[0], 0, this->n_elements()*sizeof(number)); + memset (this->val, 0, this->n_elements()*sizeof(number)); return *this; } @@ -1476,7 +1476,7 @@ typename FullMatrix::const_iterator & FullMatrix::const_iterator::operator++ () { Assert (accessor.a_row < accessor.matrix->m(), ExcIteratorPastEnd()); - + ++accessor.a_col; if (accessor.a_col >= accessor.matrix->n()) { @@ -1597,7 +1597,7 @@ FullMatrix::print (STREAM &s, const unsigned int p) const { Assert (!this->empty(), ExcEmptyMatrix()); - + for (unsigned int i=0; im(); ++i) { for (unsigned int j=0; jn(); ++j) diff --git a/deal.II/source/lac/lapack_full_matrix.cc b/deal.II/source/lac/lapack_full_matrix.cc index f051a05374..74bede304e 100644 --- a/deal.II/source/lac/lapack_full_matrix.cc +++ b/deal.II/source/lac/lapack_full_matrix.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2005, 2006, 2008, 2010 by the deal.II authors +// Copyright (C) 2005, 2006, 2008 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -71,7 +71,7 @@ LAPACKFullMatrix::operator = (const FullMatrix& M) for (unsigned int i=0;in_rows();++i) for (unsigned int j=0;jn_cols();++j) (*this)(i,j) = M(i,j); - + state = LAPACKSupport::matrix; return *this; } @@ -82,10 +82,10 @@ LAPACKFullMatrix & LAPACKFullMatrix::operator = (const double d) { Assert (d==0, ExcScalarAssignmentOnlyForZeroValue()); - + if (this->n_elements() != 0) - std::fill (this->values.begin(), this->values.end(), number()); - + std::fill_n (this->val, this->n_elements(), number()); + state = LAPACKSupport::matrix; return *this; } @@ -101,12 +101,12 @@ LAPACKFullMatrix::vmult ( { Assert ((state == matrix) || (state == inverse_matrix), ExcState(state)); - + const int mm = this->n_rows(); const int nn = this->n_cols(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); - + gemv("N", &mm, &nn, &alpha, this->data(), &mm, v.val, &one, &beta, w.val, &one); } @@ -119,12 +119,12 @@ LAPACKFullMatrix::Tvmult ( const bool adding) const { Assert (state == matrix, ExcState(state)); - + const int mm = this->n_rows(); const int nn = this->n_cols(); const number alpha = 1.; const number beta = (adding ? 1. : 0.); - + gemv("T", &mm, &nn, &alpha, this->data(), &mm, v.val, &one, &beta, w.val, &one); } @@ -171,7 +171,7 @@ LAPACKFullMatrix::compute_lu_factorization() Assert(info >= 0, ExcInternalError()); Assert(info == 0, LACExceptions::ExcSingular()); - + state = lu; } @@ -180,7 +180,7 @@ template void LAPACKFullMatrix::invert() { - Assert(state == matrix || state == lu, + Assert(state == matrix || state == lu, ExcState(state)); const int mm = this->n_rows(); const int nn = this->n_cols(); @@ -203,7 +203,7 @@ LAPACKFullMatrix::invert() Assert(info >= 0, ExcInternalError()); Assert(info == 0, LACExceptions::ExcSingular()); - + state = inverse_matrix; } @@ -216,7 +216,7 @@ LAPACKFullMatrix::apply_lu_factorization(Vector& v, Assert(state == lu, ExcState(state)); Assert(this->n_rows() == this->n_cols(), LACExceptions::ExcNotQuadratic()); - + const char* trans = transposed ? &T : &N; const int nn = this->n_cols(); const number* values = this->data(); @@ -241,14 +241,14 @@ LAPACKFullMatrix::compute_eigenvalues( wi.resize(nn); if (right) vr.resize(nn*nn); if (left) vl.resize(nn*nn); - + number* values = const_cast (this->data()); - + int info = 0; int lwork = 1; const char * const jobvr = (right) ? (&V) : (&N); const char * const jobvl = (left) ? (&V) : (&N); - + // Optimal workspace query: // The LAPACK routine DGEEV requires @@ -257,7 +257,7 @@ LAPACKFullMatrix::compute_eigenvalues( // work.size>=4*nn. // However, to improve performance, a // somewhat larger workspace may be needed. - + // SOME implementations of the LAPACK routine // provide a workspace query call, // info:=0, lwork:=-1 @@ -271,7 +271,7 @@ LAPACKFullMatrix::compute_eigenvalues( #ifndef DEAL_II_LIBLAPACK_NOQUERYMODE lwork = -1; work.resize(1); - + geev(jobvl, jobvr, &nn, values, &nn, &wr[0], &wi[0], &vl[0], &nn, &vr[0], &nn, @@ -320,23 +320,23 @@ LAPACKFullMatrix::compute_generalized_eigenvalues_symmetric ( const int nn = this->n_cols(); Assert(static_cast(nn) == this->n_rows(), ExcNotQuadratic()); Assert(B.n_rows() == B.n_cols(), ExcNotQuadratic()); - Assert(static_cast(nn) == B.n_cols(), + Assert(static_cast(nn) == B.n_cols(), ExcDimensionMismatch (nn, B.n_cols())); - Assert(eigenvectors.size() <= static_cast(nn), + Assert(eigenvectors.size() <= static_cast(nn), ExcMessage ("eigenvectors.size() > matrix.n_cols()")); - + wr.resize(nn); - wi.resize(nn); //This is set purley for consistency reasons with the + wi.resize(nn); //This is set purley for consistency reasons with the //eigenvalues() function. - + number* values_A = const_cast (this->data()); number* values_B = const_cast (B.data()); - + int info = 0; int lwork = 1; const char * const jobz = (eigenvectors.size() > 0) ? (&V) : (&N); const char * const uplo = (&U); - + // Optimal workspace query: // The LAPACK routine DSYGV requires @@ -345,7 +345,7 @@ LAPACKFullMatrix::compute_generalized_eigenvalues_symmetric ( // work.size>=3*nn-1. // However, to improve performance, a // somewhat larger workspace may be needed. - + // SOME implementations of the LAPACK routine // provide a workspace query call, // info:=0, lwork:=-1 @@ -359,7 +359,7 @@ LAPACKFullMatrix::compute_generalized_eigenvalues_symmetric ( #ifndef DEAL_II_LIBLAPACK_NOQUERYMODE lwork = -1; work.resize(1); - + sygv (&itype, jobz, uplo, &nn, values_A, &nn, values_B, &nn, &wr[0], &work[0], &lwork, &info); @@ -410,7 +410,7 @@ void LAPACKFullMatrix::compute_lu_factorization() { Assert(false, ExcNeedsLAPACK()); - + } @@ -475,10 +475,10 @@ LAPACKFullMatrix::print_formatted ( const double threshold) const { unsigned int width = width_; - + Assert ((!this->empty()) || (this->n_cols()+this->n_rows()==0), ExcInternalError()); - + // set output format, but store old // state std::ios::fmtflags old_flags = out.flags(); @@ -494,8 +494,8 @@ LAPACKFullMatrix::print_formatted ( if (!width) width = precision+2; } - - for (unsigned int i=0; in_rows(); ++i) + + for (unsigned int i=0; in_rows(); ++i) { for (unsigned int j=0; jn_cols(); ++j) if (std::fabs(this->el(i,j)) > threshold) -- 2.39.5