From 2bddc9e902ffcc023fd6d11036570a6e926a8361 Mon Sep 17 00:00:00 2001 From: Ralf Hartmann Date: Tue, 21 Feb 2006 17:08:42 +0000 Subject: [PATCH] Implementation of equ function for one, two and three matrices. git-svn-id: https://svn.dealii.org/trunk@12449 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/full_matrix.h | 30 ++ .../lac/include/lac/full_matrix.templates.h | 285 ++++++++++++++++++ deal.II/lac/source/full_matrix.double.cc | 11 + deal.II/lac/source/full_matrix.float.cc | 11 + 4 files changed, 337 insertions(+) diff --git a/deal.II/lac/include/lac/full_matrix.h b/deal.II/lac/include/lac/full_matrix.h index f0f55c63d1..1431160b60 100644 --- a/deal.II/lac/include/lac/full_matrix.h +++ b/deal.II/lac/include/lac/full_matrix.h @@ -793,6 +793,36 @@ class FullMatrix : public Table<2,number> */ void diagadd (const number s); + /** + * Assignment *this = + * a*A. + */ + template + void equ (const number a, + const FullMatrix &A); + + /** + * Assignment *this = a*A + + * b*B. + */ + template + void equ (const number a, + const FullMatrix &A, + const number b, + const FullMatrix &B); + + /** + * Assignment *this = a*A + + * b*B + c*C. + */ + template + void equ (const number a, + const FullMatrix &A, + const number b, + const FullMatrix &B, + const number c, + const FullMatrix &C); + /** * Symmetrize the matrix by * forming the mean value between diff --git a/deal.II/lac/include/lac/full_matrix.templates.h b/deal.II/lac/include/lac/full_matrix.templates.h index 7aa8b172cb..ec667beb4f 100644 --- a/deal.II/lac/include/lac/full_matrix.templates.h +++ b/deal.II/lac/include/lac/full_matrix.templates.h @@ -546,6 +546,291 @@ void FullMatrix::diagadd (const number src) } +template +template +void FullMatrix::equ (const number a, + const FullMatrix &A) +{ + Assert (!this->empty(), ExcEmptyMatrix()); + + Assert (m() == A.m(), ExcDimensionMismatch(m(), A.m())); + Assert (n() == A.n(), ExcDimensionMismatch(n(), A.n())); + + number* val = const_cast (this->data()); + const number2* Aval = A.data(); + + if ((n()==2) && (m()==2)) + { + val[0] = a * Aval[0]; + val[1] = a * Aval[1]; + val[2] = a * Aval[2]; + val[3] = a * Aval[3]; + } + else if ((n()==3) && (m()==3)) + { + val[0] = a * Aval[0]; + val[1] = a * Aval[1]; + val[2] = a * Aval[2]; + val[3] = a * Aval[3]; + val[4] = a * Aval[4]; + val[5] = a * Aval[5]; + val[6] = a * Aval[6]; + val[7] = a * Aval[7]; + val[8] = a * Aval[8]; + } + else if ((n()==4) && (m()==4)) + { + val[0] = a * Aval[0]; + val[1] = a * Aval[1]; + val[2] = a * Aval[2]; + val[3] = a * Aval[3]; + val[4] = a * Aval[4]; + val[5] = a * Aval[5]; + val[6] = a * Aval[6]; + val[7] = a * Aval[7]; + val[8] = a * Aval[8]; + val[9] = a * Aval[9]; + val[10] = a * Aval[10]; + val[11] = a * Aval[11]; + val[12] = a * Aval[12]; + val[13] = a * Aval[13]; + val[14] = a * Aval[14]; + val[15] = a * Aval[15]; + } + else if ((n()==5) && (m()==5)) + { + val[0] = a * Aval[0]; + val[1] = a * Aval[1]; + val[2] = a * Aval[2]; + val[3] = a * Aval[3]; + val[4] = a * Aval[4]; + val[5] = a * Aval[5]; + val[6] = a * Aval[6]; + val[7] = a * Aval[7]; + val[8] = a * Aval[8]; + val[9] = a * Aval[9]; + val[10] = a * Aval[10]; + val[11] = a * Aval[11]; + val[12] = a * Aval[12]; + val[13] = a * Aval[13]; + val[14] = a * Aval[14]; + val[15] = a * Aval[15]; + val[16] = a * Aval[16]; + val[17] = a * Aval[17]; + val[18] = a * Aval[18]; + val[19] = a * Aval[19]; + + val[20] = a * Aval[20]; + val[21] = a * Aval[21]; + val[22] = a * Aval[22]; + val[23] = a * Aval[23]; + val[24] = a * Aval[24]; + } + else + { + const unsigned int size = n()*m(); + for (unsigned int i=0; i +template +void +FullMatrix::equ (const number a, + const FullMatrix &A, + const number b, + const FullMatrix &B) +{ + Assert (!this->empty(), ExcEmptyMatrix()); + + Assert (m() == A.m(), ExcDimensionMismatch(m(), A.m())); + Assert (n() == A.n(), ExcDimensionMismatch(n(), A.n())); + Assert (m() == B.m(), ExcDimensionMismatch(m(), B.m())); + Assert (n() == B.n(), ExcDimensionMismatch(n(), B.n())); + + number* val = const_cast (this->data()); + const number2* Aval = A.data(); + const number2* Bval = B.data(); + + if ((n()==2) && (m()==2)) + { + val[0] = a * Aval[0] + b * Bval[0]; + val[1] = a * Aval[1] + b * Bval[1]; + val[2] = a * Aval[2] + b * Bval[2]; + val[3] = a * Aval[3] + b * Bval[3]; + } + else if ((n()==3) && (m()==3)) + { + val[0] = a * Aval[0] + b * Bval[0]; + val[1] = a * Aval[1] + b * Bval[1]; + val[2] = a * Aval[2] + b * Bval[2]; + val[3] = a * Aval[3] + b * Bval[3]; + val[4] = a * Aval[4] + b * Bval[4]; + val[5] = a * Aval[5] + b * Bval[5]; + val[6] = a * Aval[6] + b * Bval[6]; + val[7] = a * Aval[7] + b * Bval[7]; + val[8] = a * Aval[8] + b * Bval[8]; + } + else if ((n()==4) && (m()==4)) + { + val[0] = a * Aval[0] + b * Bval[0]; + val[1] = a * Aval[1] + b * Bval[1]; + val[2] = a * Aval[2] + b * Bval[2]; + val[3] = a * Aval[3] + b * Bval[3]; + val[4] = a * Aval[4] + b * Bval[4]; + val[5] = a * Aval[5] + b * Bval[5]; + val[6] = a * Aval[6] + b * Bval[6]; + val[7] = a * Aval[7] + b * Bval[7]; + val[8] = a * Aval[8] + b * Bval[8]; + val[9] = a * Aval[9] + b * Bval[9]; + val[10] = a * Aval[10] + b * Bval[10]; + val[11] = a * Aval[11] + b * Bval[11]; + val[12] = a * Aval[12] + b * Bval[12]; + val[13] = a * Aval[13] + b * Bval[13]; + val[14] = a * Aval[14] + b * Bval[14]; + val[15] = a * Aval[15] + b * Bval[15]; + } + else if ((n()==5) && (m()==5)) + { + val[0] = a * Aval[0] + b * Bval[0]; + val[1] = a * Aval[1] + b * Bval[1]; + val[2] = a * Aval[2] + b * Bval[2]; + val[3] = a * Aval[3] + b * Bval[3]; + val[4] = a * Aval[4] + b * Bval[4]; + val[5] = a * Aval[5] + b * Bval[5]; + val[6] = a * Aval[6] + b * Bval[6]; + val[7] = a * Aval[7] + b * Bval[7]; + val[8] = a * Aval[8] + b * Bval[8]; + val[9] = a * Aval[9] + b * Bval[9]; + val[10] = a * Aval[10] + b * Bval[10]; + val[11] = a * Aval[11] + b * Bval[11]; + val[12] = a * Aval[12] + b * Bval[12]; + val[13] = a * Aval[13] + b * Bval[13]; + val[14] = a * Aval[14] + b * Bval[14]; + val[15] = a * Aval[15] + b * Bval[15]; + val[16] = a * Aval[16] + b * Bval[16]; + val[17] = a * Aval[17] + b * Bval[17]; + val[18] = a * Aval[18] + b * Bval[18]; + val[19] = a * Aval[19] + b * Bval[19]; + + val[20] = a * Aval[20] + b * Bval[20]; + val[21] = a * Aval[21] + b * Bval[21]; + val[22] = a * Aval[22] + b * Bval[22]; + val[23] = a * Aval[23] + b * Bval[23]; + val[24] = a * Aval[24] + b * Bval[24]; + } + else + { + const unsigned int size = n()*m(); + for (unsigned int i=0; i +template +void +FullMatrix::equ (const number a, + const FullMatrix &A, + const number b, + const FullMatrix &B, + const number c, + const FullMatrix &C) +{ + Assert (!this->empty(), ExcEmptyMatrix()); + + Assert (m() == A.m(), ExcDimensionMismatch(m(), A.m())); + Assert (n() == A.n(), ExcDimensionMismatch(n(), A.n())); + Assert (m() == B.m(), ExcDimensionMismatch(m(), B.m())); + Assert (n() == B.n(), ExcDimensionMismatch(n(), B.n())); + Assert (m() == C.m(), ExcDimensionMismatch(m(), C.m())); + Assert (n() == C.n(), ExcDimensionMismatch(n(), C.n())); + + number* val = const_cast (this->data()); + const number2* Aval = A.data(); + const number2* Bval = B.data(); + const number2* Cval = C.data(); + + if ((n()==2) && (m()==2)) + { + val[0] = a * Aval[0] + b * Bval[0] + c * Cval[0]; + val[1] = a * Aval[1] + b * Bval[1] + c * Cval[1]; + val[2] = a * Aval[2] + b * Bval[2] + c * Cval[2]; + val[3] = a * Aval[3] + b * Bval[3] + c * Cval[3]; + } + else if ((n()==3) && (m()==3)) + { + val[0] = a * Aval[0] + b * Bval[0] + c * Cval[0]; + val[1] = a * Aval[1] + b * Bval[1] + c * Cval[1]; + val[2] = a * Aval[2] + b * Bval[2] + c * Cval[2]; + val[3] = a * Aval[3] + b * Bval[3] + c * Cval[3]; + val[4] = a * Aval[4] + b * Bval[4] + c * Cval[4]; + val[5] = a * Aval[5] + b * Bval[5] + c * Cval[5]; + val[6] = a * Aval[6] + b * Bval[6] + c * Cval[6]; + val[7] = a * Aval[7] + b * Bval[7] + c * Cval[7]; + val[8] = a * Aval[8] + b * Bval[8] + c * Cval[8]; + } + else if ((n()==4) && (m()==4)) + { + val[0] = a * Aval[0] + b * Bval[0] + c * Cval[0]; + val[1] = a * Aval[1] + b * Bval[1] + c * Cval[1]; + val[2] = a * Aval[2] + b * Bval[2] + c * Cval[2]; + val[3] = a * Aval[3] + b * Bval[3] + c * Cval[3]; + val[4] = a * Aval[4] + b * Bval[4] + c * Cval[4]; + val[5] = a * Aval[5] + b * Bval[5] + c * Cval[5]; + val[6] = a * Aval[6] + b * Bval[6] + c * Cval[6]; + val[7] = a * Aval[7] + b * Bval[7] + c * Cval[7]; + val[8] = a * Aval[8] + b * Bval[8] + c * Cval[8]; + val[9] = a * Aval[9] + b * Bval[9] + c * Cval[9]; + val[10] = a * Aval[10] + b * Bval[10] + c * Cval[10]; + val[11] = a * Aval[11] + b * Bval[11] + c * Cval[11]; + val[12] = a * Aval[12] + b * Bval[12] + c * Cval[12]; + val[13] = a * Aval[13] + b * Bval[13] + c * Cval[13]; + val[14] = a * Aval[14] + b * Bval[14] + c * Cval[14]; + val[15] = a * Aval[15] + b * Bval[15] + c * Cval[15]; + } + else if ((n()==5) && (m()==5)) + { + val[0] = a * Aval[0] + b * Bval[0] + c * Cval[0]; + val[1] = a * Aval[1] + b * Bval[1] + c * Cval[1]; + val[2] = a * Aval[2] + b * Bval[2] + c * Cval[2]; + val[3] = a * Aval[3] + b * Bval[3] + c * Cval[3]; + val[4] = a * Aval[4] + b * Bval[4] + c * Cval[4]; + val[5] = a * Aval[5] + b * Bval[5] + c * Cval[5]; + val[6] = a * Aval[6] + b * Bval[6] + c * Cval[6]; + val[7] = a * Aval[7] + b * Bval[7] + c * Cval[7]; + val[8] = a * Aval[8] + b * Bval[8] + c * Cval[8]; + val[9] = a * Aval[9] + b * Bval[9] + c * Cval[9]; + val[10] = a * Aval[10] + b * Bval[10] + c * Cval[10]; + val[11] = a * Aval[11] + b * Bval[11] + c * Cval[11]; + val[12] = a * Aval[12] + b * Bval[12] + c * Cval[12]; + val[13] = a * Aval[13] + b * Bval[13] + c * Cval[13]; + val[14] = a * Aval[14] + b * Bval[14] + c * Cval[14]; + val[15] = a * Aval[15] + b * Bval[15] + c * Cval[15]; + val[16] = a * Aval[16] + b * Bval[16] + c * Cval[16]; + val[17] = a * Aval[17] + b * Bval[17] + c * Cval[17]; + val[18] = a * Aval[18] + b * Bval[18] + c * Cval[18]; + val[19] = a * Aval[19] + b * Bval[19] + c * Cval[19]; + + val[20] = a * Aval[20] + b * Bval[20] + c * Cval[20]; + val[21] = a * Aval[21] + b * Bval[21] + c * Cval[21]; + val[22] = a * Aval[22] + b * Bval[22] + c * Cval[22]; + val[23] = a * Aval[23] + b * Bval[23] + c * Cval[23]; + val[24] = a * Aval[24] + b * Bval[24] + c * Cval[24]; + } + else + { + const unsigned int size = n()*m(); + for (unsigned int i=0; i template void FullMatrix::mmult (FullMatrix &dst, diff --git a/deal.II/lac/source/full_matrix.double.cc b/deal.II/lac/source/full_matrix.double.cc index 0e80525e22..0dd4504d8a 100644 --- a/deal.II/lac/source/full_matrix.double.cc +++ b/deal.II/lac/source/full_matrix.double.cc @@ -32,12 +32,23 @@ template FullMatrix& FullMatrix::operator =( template void FullMatrix::fill ( const FullMatrix&, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::add (const TYPEMAT, const FullMatrix&); +template void FullMatrix::add (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); +template void FullMatrix::add (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); template void FullMatrix::add_scaled (const TYPEMAT, const FullMatrix&); template void FullMatrix::add ( const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::Tadd (const TYPEMAT, const FullMatrix&); template void FullMatrix::Tadd ( const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); +template void FullMatrix::equ (const TYPEMAT, const FullMatrix&); +template void FullMatrix::equ (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); +template void FullMatrix::equ (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); template void FullMatrix::mmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::Tmmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::add_diag (const TYPEMAT, const FullMatrix&); diff --git a/deal.II/lac/source/full_matrix.float.cc b/deal.II/lac/source/full_matrix.float.cc index d77b4df006..4de71ff519 100644 --- a/deal.II/lac/source/full_matrix.float.cc +++ b/deal.II/lac/source/full_matrix.float.cc @@ -33,12 +33,23 @@ template FullMatrix& FullMatrix::operator =( template void FullMatrix::fill ( const FullMatrix&, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::add (const TYPEMAT, const FullMatrix&); +template void FullMatrix::add (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); +template void FullMatrix::add (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); template void FullMatrix::add_scaled (const TYPEMAT, const FullMatrix&); template void FullMatrix::add ( const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::Tadd (const TYPEMAT, const FullMatrix&); template void FullMatrix::Tadd ( const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); +template void FullMatrix::equ (const TYPEMAT, const FullMatrix&); +template void FullMatrix::equ (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); +template void FullMatrix::equ (const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&, + const TYPEMAT, const FullMatrix&); template void FullMatrix::mmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::Tmmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::add_diag (const TYPEMAT, const FullMatrix&); -- 2.39.5