lduMatrix.H
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) 2011-2015 OpenFOAM Foundation
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 Class
25  Foam::lduMatrix
26 
27 Description
28  lduMatrix is a general matrix class in which the coefficients are
29  stored as three arrays, one for the upper triangle, one for the
30  lower triangle and a third for the diagonal.
31 
32  Addressing arrays must be supplied for the upper and lower triangles.
33 
34  It might be better if this class were organised as a hierachy starting
35  from an empty matrix, then deriving diagonal, symmetric and asymmetric
36  matrices.
37 
38 SourceFiles
39  lduMatrixATmul.C
40  lduMatrix.C
41  lduMatrixTemplates.C
42  lduMatrixOperations.C
43  lduMatrixSolver.C
44  lduMatrixPreconditioner.C
45  lduMatrixTests.C
46  lduMatrixUpdateMatrixInterfaces.C
47 
48 \*---------------------------------------------------------------------------*/
49 
50 #ifndef lduMatrix_H
51 #define lduMatrix_H
52 
53 #include "lduMesh.H"
54 #include "primitiveFieldsFwd.H"
55 #include "FieldField.H"
57 #include "typeInfo.H"
58 #include "autoPtr.H"
59 #include "runTimeSelectionTables.H"
60 #include "solverPerformance.H"
61 #include "InfoProxy.H"
62 
63 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
64 
65 namespace Foam
66 {
67 
68 // Forward declaration of friend functions and operators
69 
70 class lduMatrix;
71 Ostream& operator<<(Ostream&, const lduMatrix&);
72 
73 
74 /*---------------------------------------------------------------------------*\
75  Class lduMatrix Declaration
76 \*---------------------------------------------------------------------------*/
77 
78 class lduMatrix
79 {
80  // private data
81 
82  //- LDU mesh reference
83  const lduMesh& lduMesh_;
84 
85  //- Coefficients (not including interfaces)
86  scalarField *lowerPtr_, *diagPtr_, *upperPtr_;
87 
88 
89 public:
90 
91  //- Abstract base-class for lduMatrix solvers
92  class solver
93  {
94  protected:
95 
96  // Protected data
97 
103 
104  //- Dictionary of controls
106 
107  //- Maximum number of iterations in the solver
108  label maxIter_;
109 
110  //- Minimum number of iterations in the solver
111  label minIter_;
112 
113  //- Final convergence tolerance
114  scalar tolerance_;
115 
116  //- Convergence tolerance relative to the initial
117  scalar relTol_;
118 
119 
120  // Protected Member Functions
121 
122  //- Read the control parameters from the controlDict_
123  virtual void readControls();
124 
125 
126  public:
127 
128  //- Runtime type information
129  virtual const word& type() const = 0;
130 
131 
132  // Declare run-time constructor selection tables
133 
135  (
136  autoPtr,
137  solver,
138  symMatrix,
139  (
140  const word& fieldName,
141  const lduMatrix& matrix,
142  const FieldField<Field, scalar>& interfaceBouCoeffs,
143  const FieldField<Field, scalar>& interfaceIntCoeffs,
144  const lduInterfaceFieldPtrsList& interfaces,
145  const dictionary& solverControls
146  ),
147  (
148  fieldName,
149  matrix,
150  interfaceBouCoeffs,
151  interfaceIntCoeffs,
152  interfaces,
153  solverControls
154  )
155  );
156 
158  (
159  autoPtr,
160  solver,
161  asymMatrix,
162  (
163  const word& fieldName,
164  const lduMatrix& matrix,
165  const FieldField<Field, scalar>& interfaceBouCoeffs,
166  const FieldField<Field, scalar>& interfaceIntCoeffs,
167  const lduInterfaceFieldPtrsList& interfaces,
168  const dictionary& solverControls
169  ),
170  (
171  fieldName,
172  matrix,
173  interfaceBouCoeffs,
174  interfaceIntCoeffs,
175  interfaces,
176  solverControls
177  )
178  );
179 
180 
181  // Constructors
182 
183  solver
184  (
185  const word& fieldName,
186  const lduMatrix& matrix,
187  const FieldField<Field, scalar>& interfaceBouCoeffs,
188  const FieldField<Field, scalar>& interfaceIntCoeffs,
189  const lduInterfaceFieldPtrsList& interfaces,
190  const dictionary& solverControls
191  );
192 
193  // Selectors
194 
195  //- Return a new solver
196  static autoPtr<solver> New
197  (
198  const word& fieldName,
199  const lduMatrix& matrix,
200  const FieldField<Field, scalar>& interfaceBouCoeffs,
201  const FieldField<Field, scalar>& interfaceIntCoeffs,
202  const lduInterfaceFieldPtrsList& interfaces,
203  const dictionary& solverControls
204  );
205 
206 
207 
208  //- Destructor
209  virtual ~solver()
210  {}
211 
212 
213  // Member functions
214 
215  // Access
217  const word& fieldName() const
218  {
219  return fieldName_;
220  }
222  const lduMatrix& matrix() const
223  {
224  return matrix_;
225  }
228  {
229  return interfaceBouCoeffs_;
230  }
233  {
234  return interfaceIntCoeffs_;
235  }
238  {
239  return interfaces_;
240  }
241 
242 
243  //- Read and reset the solver parameters from the given stream
244  virtual void read(const dictionary&);
245 
246  virtual solverPerformance solve
247  (
248  scalarField& psi,
249  const scalarField& source,
250  const direction cmpt=0
251  ) const = 0;
252 
253  //- Return the matrix norm used to normalise the residual for the
254  // stopping criterion
255  scalar normFactor
256  (
257  const scalarField& psi,
258  const scalarField& source,
259  const scalarField& Apsi,
260  scalarField& tmpField
261  ) const;
262  };
263 
264 
265  //- Abstract base-class for lduMatrix smoothers
266  class smoother
267  {
268  protected:
269 
270  // Protected data
277 
278 
279  public:
280 
281  //- Find the smoother name (directly or from a sub-dictionary)
282  static word getName(const dictionary&);
283 
284  //- Runtime type information
285  virtual const word& type() const = 0;
286 
287 
288  // Declare run-time constructor selection tables
289 
291  (
292  autoPtr,
293  smoother,
294  symMatrix,
295  (
296  const word& fieldName,
297  const lduMatrix& matrix,
298  const FieldField<Field, scalar>& interfaceBouCoeffs,
299  const FieldField<Field, scalar>& interfaceIntCoeffs,
300  const lduInterfaceFieldPtrsList& interfaces
301  ),
302  (
303  fieldName,
304  matrix,
305  interfaceBouCoeffs,
306  interfaceIntCoeffs,
307  interfaces
308  )
309  );
310 
312  (
313  autoPtr,
314  smoother,
315  asymMatrix,
316  (
317  const word& fieldName,
318  const lduMatrix& matrix,
319  const FieldField<Field, scalar>& interfaceBouCoeffs,
320  const FieldField<Field, scalar>& interfaceIntCoeffs,
321  const lduInterfaceFieldPtrsList& interfaces
322  ),
323  (
324  fieldName,
325  matrix,
326  interfaceBouCoeffs,
327  interfaceIntCoeffs,
328  interfaces
329  )
330  );
331 
332 
333  // Constructors
334 
335  smoother
336  (
337  const word& fieldName,
338  const lduMatrix& matrix,
339  const FieldField<Field, scalar>& interfaceBouCoeffs,
340  const FieldField<Field, scalar>& interfaceIntCoeffs,
341  const lduInterfaceFieldPtrsList& interfaces
342  );
343 
344 
345  // Selectors
346 
347  //- Return a new smoother
348  static autoPtr<smoother> New
349  (
350  const word& fieldName,
351  const lduMatrix& matrix,
352  const FieldField<Field, scalar>& interfaceBouCoeffs,
353  const FieldField<Field, scalar>& interfaceIntCoeffs,
354  const lduInterfaceFieldPtrsList& interfaces,
355  const dictionary& solverControls
356  );
357 
358 
359  //- Destructor
360  virtual ~smoother()
361  {}
362 
363 
364  // Member functions
365 
366  // Access
368  const word& fieldName() const
369  {
370  return fieldName_;
371  }
373  const lduMatrix& matrix() const
374  {
375  return matrix_;
376  }
379  {
380  return interfaceBouCoeffs_;
381  }
384  {
385  return interfaceIntCoeffs_;
386  }
389  {
390  return interfaces_;
391  }
392 
393 
394  //- Smooth the solution for a given number of sweeps
395  virtual void smooth
396  (
397  scalarField& psi,
398  const scalarField& source,
399  const direction cmpt,
400  const label nSweeps
401  ) const = 0;
402  };
403 
404 
405  //- Abstract base-class for lduMatrix preconditioners
406  class preconditioner
407  {
408  protected:
409 
410  // Protected data
411 
412  //- Reference to the base-solver this preconditioner is used with
413  const solver& solver_;
414 
415 
416  public:
417 
418  //- Find the preconditioner name (directly or from a sub-dictionary)
419  static word getName(const dictionary&);
420 
421  //- Runtime type information
422  virtual const word& type() const = 0;
423 
424 
425  // Declare run-time constructor selection tables
426 
428  (
429  autoPtr,
431  symMatrix,
432  (
433  const solver& sol,
434  const dictionary& solverControls
435  ),
436  (sol, solverControls)
437  );
438 
440  (
441  autoPtr,
443  asymMatrix,
444  (
445  const solver& sol,
446  const dictionary& solverControls
447  ),
448  (sol, solverControls)
449  );
450 
451 
452  // Constructors
453 
455  (
456  const solver& sol
457  )
458  :
459  solver_(sol)
460  {}
461 
462 
463  // Selectors
464 
465  //- Return a new preconditioner
467  (
468  const solver& sol,
469  const dictionary& solverControls
470  );
471 
472 
473  //- Destructor
474  virtual ~preconditioner()
475  {}
476 
477 
478  // Member functions
479 
480  //- Read and reset the preconditioner parameters
481  // from the given stream
482  virtual void read(const dictionary&)
483  {}
484 
485  //- Return wA the preconditioned form of residual rA
486  virtual void precondition
487  (
488  scalarField& wA,
489  const scalarField& rA,
490  const direction cmpt=0
491  ) const = 0;
492 
493  //- Return wT the transpose-matrix preconditioned form of
494  // residual rT.
495  // This is only required for preconditioning asymmetric matrices.
496  virtual void preconditionT
497  (
498  scalarField& wT,
499  const scalarField& rT,
500  const direction cmpt=0
501  ) const
502  {
504  (
505  type() +"::preconditionT"
506  "(scalarField& wT, const scalarField& rT, "
507  "const direction cmpt)"
508  );
509  }
510  };
511 
512 
513  // Static data
514 
515  // Declare name of the class and its debug switch
516  ClassName("lduMatrix");
517 
518 
519  // Constructors
520 
521  //- Construct given an LDU addressed mesh.
522  // The coefficients are initially empty for subsequent setting.
523  lduMatrix(const lduMesh&);
524 
525  //- Construct as copy
526  lduMatrix(const lduMatrix&);
527 
528  //- Construct as copy or re-use as specified.
529  lduMatrix(lduMatrix&, bool reUse);
530 
531  //- Construct given an LDU addressed mesh and an Istream
532  // from which the coefficients are read
533  lduMatrix(const lduMesh&, Istream&);
534 
535 
536  //- Destructor
537  ~lduMatrix();
538 
539 
540  // Member functions
541 
542  // Access to addressing
543 
544  //- Return the LDU mesh from which the addressing is obtained
545  const lduMesh& mesh() const
546  {
547  return lduMesh_;
548  }
549 
550  //- Return the LDU addressing
551  const lduAddressing& lduAddr() const
552  {
553  return lduMesh_.lduAddr();
554  }
555 
556  //- Return the patch evaluation schedule
557  const lduSchedule& patchSchedule() const
558  {
559  return lduAddr().patchSchedule();
560  }
561 
562 
563  // Access to coefficients
564 
565  scalarField& lower();
566  scalarField& diag();
567  scalarField& upper();
568 
569  // Size with externally provided sizes (for constructing with 'fake'
570  // mesh in GAMG)
571 
572  scalarField& lower(const label size);
573  scalarField& diag(const label nCoeffs);
574  scalarField& upper(const label nCoeffs);
575 
576 
577  const scalarField& lower() const;
578  const scalarField& diag() const;
579  const scalarField& upper() const;
581  bool hasDiag() const
582  {
583  return (diagPtr_);
584  }
586  bool hasUpper() const
587  {
588  return (upperPtr_);
589  }
591  bool hasLower() const
592  {
593  return (lowerPtr_);
594  }
596  bool diagonal() const
597  {
598  return (diagPtr_ && !lowerPtr_ && !upperPtr_);
599  }
601  bool symmetric() const
602  {
603  return (diagPtr_ && (!lowerPtr_ && upperPtr_));
604  }
606  bool asymmetric() const
607  {
608  return (diagPtr_ && lowerPtr_ && upperPtr_);
609  }
610 
611 
612  // operations
613 
614  void sumDiag();
615  void negSumDiag();
616 
617  void sumMagOffDiag(scalarField& sumOff) const;
618 
619  //- Matrix multiplication with updated interfaces.
620  void Amul
621  (
622  scalarField&,
623  const tmp<scalarField>&,
626  const direction cmpt
627  ) const;
628 
629  //- Matrix transpose multiplication with updated interfaces.
630  void Tmul
631  (
632  scalarField&,
633  const tmp<scalarField>&,
636  const direction cmpt
637  ) const;
638 
639 
640  //- Sum the coefficients on each row of the matrix
641  void sumA
642  (
643  scalarField&,
646  ) const;
647 
648 
649  void residual
650  (
651  scalarField& rA,
652  const scalarField& psi,
653  const scalarField& source,
654  const FieldField<Field, scalar>& interfaceBouCoeffs,
655  const lduInterfaceFieldPtrsList& interfaces,
656  const direction cmpt
657  ) const;
658 
660  (
661  const scalarField& psi,
662  const scalarField& source,
663  const FieldField<Field, scalar>& interfaceBouCoeffs,
664  const lduInterfaceFieldPtrsList& interfaces,
665  const direction cmpt
666  ) const;
667 
668 
669  //- Initialise the update of interfaced interfaces
670  // for matrix operations
672  (
673  const FieldField<Field, scalar>& interfaceCoeffs,
674  const lduInterfaceFieldPtrsList& interfaces,
675  const scalarField& psiif,
676  scalarField& result,
677  const direction cmpt
678  ) const;
679 
680  //- Update interfaced interfaces for matrix operations
682  (
683  const FieldField<Field, scalar>& interfaceCoeffs,
684  const lduInterfaceFieldPtrsList& interfaces,
685  const scalarField& psiif,
686  scalarField& result,
687  const direction cmpt
688  ) const;
689 
690 
691  template<class Type>
692  tmp<Field<Type> > H(const Field<Type>&) const;
693 
694  template<class Type>
695  tmp<Field<Type> > H(const tmp<Field<Type> >&) const;
696 
697  tmp<scalarField> H1() const;
698 
699  template<class Type>
700  tmp<Field<Type> > faceH(const Field<Type>&) const;
701 
702  template<class Type>
703  tmp<Field<Type> > faceH(const tmp<Field<Type> >&) const;
704 
705 
706  // Info
707 
708  //- Return info proxy.
709  // Used to print matrix information to a stream
710  InfoProxy<lduMatrix> info() const
711  {
712  return *this;
713  }
714 
715 
716  // Member operators
717 
718  void operator=(const lduMatrix&);
719 
720  void negate();
721 
722  void operator+=(const lduMatrix&);
723  void operator-=(const lduMatrix&);
724 
725  void operator*=(const scalarField&);
726  void operator*=(scalar);
727 
728 
729  // Ostream operator
730 
731  friend Ostream& operator<<(Ostream&, const lduMatrix&);
732  friend Ostream& operator<<(Ostream&, const InfoProxy<lduMatrix>&);
733 };
734 
735 
736 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
737 
738 } // End namespace Foam
739 
740 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
741 
742 #ifdef NoRepository
743 # include "lduMatrixTemplates.C"
744 #endif
745 
746 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
747 
748 #endif
749 
750 // ************************************************************************* //
virtual void read(const dictionary &)
Read and reset the solver parameters from the given stream.
const lduMatrix & matrix_
Definition: lduMatrix.H:98
lduInterfaceFieldPtrsList interfaces_
Definition: lduMatrix.H:101
virtual void readControls()
Read the control parameters from the controlDict_.
unsigned char direction
Definition: direction.H:43
void sumMagOffDiag(scalarField &sumOff) const
tmp< Field< Type > > H(const Field< Type > &) const
virtual const word & type() const =0
Runtime type information.
void sumA(scalarField &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &) const
Sum the coefficients on each row of the matrix.
Generic field type.
Definition: FieldField.H:51
virtual ~solver()
Destructor.
Definition: lduMatrix.H:208
void operator*=(const scalarField &)
scalarField & lower()
Definition: lduMatrix.C:165
Macros to ease declaration of run-time selection tables.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void Amul(scalarField &, const tmp< scalarField > &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &, const direction cmpt) const
Matrix multiplication with updated interfaces.
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics...
label maxIter_
Maximum number of iterations in the solver.
Definition: lduMatrix.H:107
scalar normFactor(const scalarField &psi, const scalarField &source, const scalarField &Apsi, scalarField &tmpField) const
Return the matrix norm used to normalise the residual for the.
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
label minIter_
Minimum number of iterations in the solver.
Definition: lduMatrix.H:110
ClassName("lduMatrix")
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:51
declareRunTimeSelectionTable(autoPtr, solver, symMatrix,( const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls ),( fieldName, matrix, interfaceBouCoeffs, interfaceIntCoeffs, interfaces, solverControls ))
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Ostream & operator<<(Ostream &, const edgeMesh &)
Definition: edgeMeshIO.C:133
Namespace for OpenFOAM.
const lduAddressing & lduAddr() const
Return the LDU addressing.
Definition: lduMatrix.H:550
virtual const lduAddressing & lduAddr() const =0
Return ldu addressing.
scalar relTol_
Convergence tolerance relative to the initial.
Definition: lduMatrix.H:116
const FieldField< Field, scalar > & interfaceIntCoeffs_
Definition: lduMatrix.H:100
void Tmul(scalarField &, const tmp< scalarField > &, const FieldField< Field, scalar > &, const lduInterfaceFieldPtrsList &, const direction cmpt) const
Matrix transpose multiplication with updated interfaces.
tmp< scalarField > H1() const
Forward declarations of the specialisations of Field<T> for scalar, vector and tensor.
void updateMatrixInterfaces(const FieldField< Field, scalar > &interfaceCoeffs, const lduInterfaceFieldPtrsList &interfaces, const scalarField &psiif, scalarField &result, const direction cmpt) const
Update interfaced interfaces for matrix operations.
Abstract base-class for lduMatrix smoothers.
Definition: lduMatrix.H:265
void residual(scalarField &rA, const scalarField &psi, const scalarField &source, const FieldField< Field, scalar > &interfaceBouCoeffs, const lduInterfaceFieldPtrsList &interfaces, const direction cmpt) const
A helper class for outputting values to Ostream.
Definition: InfoProxy.H:45
InfoProxy< lduMatrix > info() const
Return info proxy.
Definition: lduMatrix.H:709
bool asymmetric() const
Definition: lduMatrix.H:605
bool hasUpper() const
Definition: lduMatrix.H:585
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
Definition: lduMatrix.H:77
const FieldField< Field, scalar > & interfaceBouCoeffs() const
Definition: lduMatrix.H:226
bool hasLower() const
Definition: lduMatrix.H:590
static autoPtr< solver > New(const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls)
Return a new solver.
Abstract base-class for lduMatrix preconditioners.
Definition: lduMatrix.H:405
const volScalarField & psi
Definition: createFields.H:24
const lduInterfaceFieldPtrsList & interfaces() const
Definition: lduMatrix.H:236
lduMatrix(const lduMesh &)
Construct given an LDU addressed mesh.
Definition: lduMatrix.C:40
scalarField & diag()
Definition: lduMatrix.C:183
const lduMatrix & matrix() const
Definition: lduMatrix.H:221
friend Ostream & operator<<(Ostream &, const lduMatrix &)
void operator=(const lduMatrix &)
bool diagonal() const
Definition: lduMatrix.H:595
scalar tolerance_
Final convergence tolerance.
Definition: lduMatrix.H:113
dictionary controlDict_
Dictionary of controls.
Definition: lduMatrix.H:104
virtual const lduSchedule & patchSchedule() const =0
lduMatrix member H operations.
bool symmetric() const
Definition: lduMatrix.H:600
bool hasDiag() const
Definition: lduMatrix.H:580
Abstract base-class for lduMatrix solvers.
Definition: lduMatrix.H:91
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
scalarField & upper()
Definition: lduMatrix.C:194
const FieldField< Field, scalar > & interfaceBouCoeffs_
Definition: lduMatrix.H:99
const lduSchedule & patchSchedule() const
Return the patch evaluation schedule.
Definition: lduMatrix.H:556
const lduMesh & mesh() const
Return the LDU mesh from which the addressing is obtained.
Definition: lduMatrix.H:544
The class contains the addressing required by the lduMatrix: upper, lower and losort.
void operator-=(const lduMatrix &)
const FieldField< Field, scalar > & interfaceIntCoeffs() const
Definition: lduMatrix.H:231
#define notImplemented(functionName)
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:356
tmp< Field< Type > > faceH(const Field< Type > &) const
virtual solverPerformance solve(scalarField &psi, const scalarField &source, const direction cmpt=0) const =0
void operator+=(const lduMatrix &)
~lduMatrix()
Destructor.
Definition: lduMatrix.C:146
void initMatrixInterfaces(const FieldField< Field, scalar > &interfaceCoeffs, const lduInterfaceFieldPtrsList &interfaces, const scalarField &psiif, scalarField &result, const direction cmpt) const
Initialise the update of interfaced interfaces.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
A class for managing temporary objects.
Definition: PtrList.H:118
void smooth(volScalarField &field, const scalar coeff)
Definition: fvcSmooth.C:35
solver(const word &fieldName, const lduMatrix &matrix, const FieldField< Field, scalar > &interfaceBouCoeffs, const FieldField< Field, scalar > &interfaceIntCoeffs, const lduInterfaceFieldPtrsList &interfaces, const dictionary &solverControls)
const word & fieldName() const
Definition: lduMatrix.H:216