FreeFOAM The Cross-Platform CFD Toolkit
fvmSup.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 <finiteVolume/volFields.H>
28 #include <finiteVolume/fvMatrix.H>
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 template<class Type>
35 (
38 )
39 {
40  const fvMesh& mesh = vf.mesh();
41 
42  tmp<fvMatrix<Type> > tfvm
43  (
44  new fvMatrix<Type>
45  (
46  vf,
47  dimVol*su.dimensions()
48  )
49  );
50  fvMatrix<Type>& fvm = tfvm();
51 
52  fvm.source() -= mesh.V()*su.field();
53 
54  return tfvm;
55 }
56 
57 template<class Type>
60 (
63 )
64 {
65  tmp<fvMatrix<Type> > tfvm = fvm::Su(tsu(), vf);
66  tsu.clear();
67  return tfvm;
68 }
69 
70 template<class Type>
73 (
76 )
77 {
78  tmp<fvMatrix<Type> > tfvm = fvm::Su(tsu(), vf);
79  tsu.clear();
80  return tfvm;
81 }
82 
83 template<class Type>
86 (
87  const zeroField&,
89 )
90 {
91  return zeroField();
92 }
93 
94 
95 template<class Type>
98 (
101 )
102 {
103  const fvMesh& mesh = vf.mesh();
104 
105  tmp<fvMatrix<Type> > tfvm
106  (
107  new fvMatrix<Type>
108  (
109  vf,
110  dimVol*sp.dimensions()*vf.dimensions()
111  )
112  );
113  fvMatrix<Type>& fvm = tfvm();
114 
115  fvm.diag() += mesh.V()*sp.field();
116 
117  return tfvm;
118 }
119 
120 template<class Type>
123 (
126 )
127 {
128  tmp<fvMatrix<Type> > tfvm = fvm::Sp(tsp(), vf);
129  tsp.clear();
130  return tfvm;
131 }
132 
133 template<class Type>
136 (
137  const tmp<volScalarField>& tsp,
139 )
140 {
141  tmp<fvMatrix<Type> > tfvm = fvm::Sp(tsp(), vf);
142  tsp.clear();
143  return tfvm;
144 }
145 
146 
147 template<class Type>
150 (
151  const dimensionedScalar& sp,
153 )
154 {
155  const fvMesh& mesh = vf.mesh();
156 
157  tmp<fvMatrix<Type> > tfvm
158  (
159  new fvMatrix<Type>
160  (
161  vf,
162  dimVol*sp.dimensions()*vf.dimensions()
163  )
164  );
165  fvMatrix<Type>& fvm = tfvm();
166 
167  fvm.diag() += mesh.V()*sp.value();
168 
169  return tfvm;
170 }
171 
172 template<class Type>
175 (
176  const zeroField&,
178 )
179 {
180  return zeroField();
181 }
182 
183 
184 template<class Type>
187 (
190 )
191 {
192  const fvMesh& mesh = vf.mesh();
193 
194  tmp<fvMatrix<Type> > tfvm
195  (
196  new fvMatrix<Type>
197  (
198  vf,
199  dimVol*susp.dimensions()*vf.dimensions()
200  )
201  );
202  fvMatrix<Type>& fvm = tfvm();
203 
204  fvm.diag() += mesh.V()*max(susp.field(), scalar(0));
205 
206  fvm.source() -= mesh.V()*min(susp.field(), scalar(0))
207  *vf.internalField();
208 
209  return tfvm;
210 }
211 
212 template<class Type>
215 (
216  const tmp<DimensionedField<scalar, volMesh> >& tsusp,
218 )
219 {
220  tmp<fvMatrix<Type> > tfvm = fvm::SuSp(tsusp(), vf);
221  tsusp.clear();
222  return tfvm;
223 }
224 
225 template<class Type>
228 (
229  const tmp<volScalarField>& tsusp,
231 )
232 {
233  tmp<fvMatrix<Type> > tfvm = fvm::SuSp(tsusp(), vf);
234  tsusp.clear();
235  return tfvm;
236 }
237 
238 template<class Type>
241 (
242  const zeroField&,
244 )
245 {
246  return zeroField();
247 }
248 
249 
250 // ************************ vim: set sw=4 sts=4 et: ************************ //