FreeFOAM The Cross-Platform CFD Toolkit
fvScalarMatrix.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 "fvScalarMatrix.H"
28 
29 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
30 
31 template<>
33 (
34  const label patchi,
35  const label facei,
36  const direction,
37  const scalar value
38 )
39 {
40  if (psi_.needReference())
41  {
42  if (Pstream::master())
43  {
44  internalCoeffs_[patchi][facei] +=
45  diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]];
46 
47  boundaryCoeffs_[patchi][facei] +=
48  diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]]
49  *value;
50  }
51  }
52 }
53 
54 
55 template<>
58 (
59  const dictionary& solverControls
60 )
61 {
62  if (debug)
63  {
64  Info<< "fvMatrix<scalar>::solver(const dictionary& solverControls) : "
65  "solver for fvMatrix<scalar>"
66  << endl;
67  }
68 
69  scalarField saveDiag = diag();
70  addBoundaryDiag(diag(), 0);
71 
73  (
75  (
76  *this,
77  lduMatrix::solver::New
78  (
79  psi_.name(),
80  *this,
81  boundaryCoeffs_,
82  internalCoeffs_,
83  psi_.boundaryField().interfaces(),
84  solverControls
85  )
86  )
87  );
88 
89  diag() = saveDiag;
90 
91  return solverPtr;
92 }
93 
94 
95 template<>
97 (
98  const dictionary& solverControls
99 )
100 {
101  scalarField saveDiag = fvMat_.diag();
102  fvMat_.addBoundaryDiag(fvMat_.diag(), 0);
103 
104  scalarField totalSource = fvMat_.source();
105  fvMat_.addBoundarySource(totalSource, false);
106 
107  // assign new solver controls
108  solver_->read(solverControls);
109 
110  lduMatrix::solverPerformance solverPerf =
111  solver_->solve(fvMat_.psi().internalField(), totalSource);
112 
113  solverPerf.print();
114 
115  fvMat_.diag() = saveDiag;
116 
117  fvMat_.psi().correctBoundaryConditions();
118 
119  return solverPerf;
120 }
121 
122 
123 template<>
125 (
126  const dictionary& solverControls
127 )
128 {
129  if (debug)
130  {
131  Info<< "fvMatrix<scalar>::solve(const dictionary& solverControls) : "
132  "solving fvMatrix<scalar>"
133  << endl;
134  }
135 
136  scalarField saveDiag = diag();
137  addBoundaryDiag(diag(), 0);
138 
139  scalarField totalSource = source_;
140  addBoundarySource(totalSource, false);
141 
142  // Solver call
143  lduMatrix::solverPerformance solverPerf = lduMatrix::solver::New
144  (
145  psi_.name(),
146  *this,
147  boundaryCoeffs_,
148  internalCoeffs_,
149  psi_.boundaryField().interfaces(),
150  solverControls
151  )->solve(psi_.internalField(), totalSource);
152 
153  solverPerf.print();
154 
155  diag() = saveDiag;
156 
157  psi_.correctBoundaryConditions();
158 
159  return solverPerf;
160 }
161 
162 
163 template<>
165 {
166  scalarField boundaryDiag(psi_.size(), 0.0);
167  addBoundaryDiag(boundaryDiag, 0);
168 
169  tmp<scalarField> tres
170  (
171  lduMatrix::residual
172  (
173  psi_.internalField(),
174  source_ - boundaryDiag*psi_.internalField(),
175  boundaryCoeffs_,
176  psi_.boundaryField().interfaces(),
177  0
178  )
179  );
180 
181  addBoundarySource(tres());
182 
183  return tres;
184 }
185 
186 
187 template<>
189 {
190  tmp<volScalarField> tHphi
191  (
192  new volScalarField
193  (
194  IOobject
195  (
196  "H("+psi_.name()+')',
197  psi_.instance(),
198  psi_.mesh(),
199  IOobject::NO_READ,
200  IOobject::NO_WRITE
201  ),
202  psi_.mesh(),
203  dimensions_/dimVol,
204  zeroGradientFvPatchScalarField::typeName
205  )
206  );
207  volScalarField& Hphi = tHphi();
208 
209  Hphi.internalField() = (lduMatrix::H(psi_.internalField()) + source_);
210  addBoundarySource(Hphi.internalField());
211 
212  Hphi.internalField() /= psi_.mesh().V();
214 
215  return tHphi;
216 }
217 
218 
219 template<>
221 {
223  (
224  new volScalarField
225  (
226  IOobject
227  (
228  "H(1)",
229  psi_.instance(),
230  psi_.mesh(),
231  IOobject::NO_READ,
232  IOobject::NO_WRITE
233  ),
234  psi_.mesh(),
235  dimensions_/(dimVol*psi_.dimensions()),
236  zeroGradientFvPatchScalarField::typeName
237  )
238  );
239  volScalarField& H1_ = tH1();
240 
241  H1_.internalField() = lduMatrix::H1();
242  //addBoundarySource(Hphi.internalField());
243 
244  H1_.internalField() /= psi_.mesh().V();
245  H1_.correctBoundaryConditions();
246 
247  return tH1;
248 }
249 
250 
251 // ************************ vim: set sw=4 sts=4 et: ************************ //