FreeFOAM The Cross-Platform CFD Toolkit
simpleMatrix.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include <OpenFOAM/simpleMatrix.H>
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Type>
32 :
33  scalarSquareMatrix(mSize),
34  source_(mSize)
35 {}
36 
37 
38 template<class Type>
40 (
41  const label mSize,
42  const scalar coeffVal,
43  const Type& sourceVal
44 )
45 :
46  scalarSquareMatrix(mSize, mSize, coeffVal),
47  source_(mSize, sourceVal)
48 {}
49 
50 
51 template<class Type>
53 (
54  const scalarSquareMatrix& matrix,
55  const Field<Type>& source
56 )
57 :
58  scalarSquareMatrix(matrix),
59  source_(source)
60 {}
61 
62 
63 template<class Type>
65 :
67  source_(is)
68 {}
69 
70 
71 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
72 
73 template<class Type>
75 {
76  scalarSquareMatrix tmpMatrix = *this;
77  Field<Type> sourceSol = source_;
78 
79  Foam::solve(tmpMatrix, sourceSol);
80 
81  return sourceSol;
82 }
83 
84 
85 template<class Type>
87 {
88  scalarSquareMatrix luMatrix = *this;
89  Field<Type> sourceSol = source_;
90 
91  Foam::LUsolve(luMatrix, sourceSol);
92 
93  return sourceSol;
94 }
95 
96 
97 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
98 
99 template<class Type>
101 {
102  if (this == &m)
103  {
104  FatalErrorIn("simpleMatrix<Type>::operator=(const simpleMatrix<Type>&)")
105  << "Attempted assignment to self"
106  << abort(FatalError);
107  }
108 
109  if (n() != m.n())
110  {
111  FatalErrorIn("simpleMatrix<Type>::operator=(const simpleMatrix<Type>&)")
112  << "Different size matrices"
113  << abort(FatalError);
114  }
115 
116  if (source_.size() != m.source_.size())
117  {
118  FatalErrorIn("simpleMatrix<Type>::operator=(const simpleMatrix<Type>&)")
119  << "Different size source vectors"
120  << abort(FatalError);
121  }
122 
124  source_ = m.source_;
125 }
126 
127 
128 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
129 
130 template<class Type>
131 Foam::simpleMatrix<Type> Foam::operator+
132 (
133  const simpleMatrix<Type>& m1,
134  const simpleMatrix<Type>& m2
135 )
136 {
137  return simpleMatrix<Type>
138  (
139  static_cast<const scalarSquareMatrix&>(m1)
140  + static_cast<const scalarSquareMatrix&>(m2),
141  m1.source_ + m2.source_
142  );
143 }
144 
145 
146 template<class Type>
147 Foam::simpleMatrix<Type> Foam::operator-
148 (
149  const simpleMatrix<Type>& m1,
150  const simpleMatrix<Type>& m2
151 )
152 {
153  return simpleMatrix<Type>
154  (
155  static_cast<const scalarSquareMatrix&>(m1)
156  - static_cast<const scalarSquareMatrix&>(m2),
157  m1.source_ - m2.source_
158  );
159 }
160 
161 
162 template<class Type>
164 {
165  return simpleMatrix<Type>(s*m.matrix_, s*m.source_);
166 }
167 
168 
169 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
170 
171 template<class Type>
172 Foam::Ostream& Foam::operator<<(Ostream& os, const simpleMatrix<Type>& m)
173 {
174  os << static_cast<const scalarSquareMatrix&>(m) << nl << m.source_;
175  return os;
176 }
177 
178 
179 // ************************ vim: set sw=4 sts=4 et: ************************ //