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-2016 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 Note
35  It might be better if this class were organised as a hierachy starting
36  from an empty matrix, then deriving diagonal, symmetric and asymmetric
37  matrices.
38 
39 SourceFiles
40  LduMatrixATmul.C
41  LduMatrix.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 "Field.H"
55 #include "FieldField.H"
57 #include "SolverPerformance.H"
58 #include "typeInfo.H"
59 #include "autoPtr.H"
60 #include "runTimeSelectionTables.H"
61 
62 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
63 
64 namespace Foam
65 {
66 
67 // Forward declaration of friend functions and operators
68 
69 template<class Type, class DType, class LUType>
70 class LduMatrix;
71 
72 template<class Type, class DType, class LUType>
73 Ostream& operator<<
74 (
75  Ostream&,
77 );
78 
79 
80 /*---------------------------------------------------------------------------*\
81  Class LduMatrix Declaration
82 \*---------------------------------------------------------------------------*/
83 
84 template<class Type, class DType, class LUType>
85 class LduMatrix
86 {
87  // private data
88 
89  //- LDU mesh reference
90  const lduMesh& lduMesh_;
91 
92  //- Diagonal coefficients
93  Field<DType> *diagPtr_;
94 
95  //- Off-diagonal coefficients
96  Field<LUType> *upperPtr_, *lowerPtr_;
97 
98  //- Source
99  Field<Type> *sourcePtr_;
100 
101  //- Field interfaces (processor patches etc.)
103 
104  //- Off-diagonal coefficients for interfaces
105  FieldField<Field, LUType> interfacesUpper_, interfacesLower_;
106 
107 
108 public:
110  friend class SolverPerformance<Type>;
111 
112  //- Abstract base-class for LduMatrix solvers
113  class solver
114  {
115  protected:
116 
117  // Protected data
121 
122  //- Dictionary of controls
124 
125  //- Maximum number of iterations in the solver
126  label maxIter_;
127 
128  //- Minimum number of iterations in the solver
129  label minIter_;
130 
131  //- Final convergence tolerance
132  Type tolerance_;
133 
134  //- Convergence tolerance relative to the initial
135  Type relTol_;
136 
137 
138  // Protected Member Functions
139 
140  //- Read a control parameter from controlDict
141  template<class T>
142  inline void readControl
143  (
144  const dictionary& controlDict,
145  T& control,
146  const word& controlName
147  );
148 
149 
150  //- Read the control parameters from the controlDict_
151  virtual void readControls();
152 
153 
154  public:
155 
156  //- Runtime type information
157  virtual const word& type() const = 0;
158 
159 
160  // Declare run-time constructor selection tables
161 
163  (
164  autoPtr,
165  solver,
166  symMatrix,
167  (
168  const word& fieldName,
169  const LduMatrix<Type, DType, LUType>& matrix,
170  const dictionary& solverDict
171  ),
172  (
173  fieldName,
174  matrix,
175  solverDict
176  )
177  );
178 
180  (
181  autoPtr,
182  solver,
183  asymMatrix,
184  (
185  const word& fieldName,
186  const LduMatrix<Type, DType, LUType>& matrix,
187  const dictionary& solverDict
188  ),
189  (
190  fieldName,
191  matrix,
192  solverDict
193  )
194  );
195 
196 
197  // Constructors
198 
199  solver
200  (
201  const word& fieldName,
202  const LduMatrix<Type, DType, LUType>& matrix,
203  const dictionary& solverDict
204  );
205 
206 
207  // Selectors
208 
209  //- Return a new solver
210  static autoPtr<solver> New
211  (
212  const word& fieldName,
213  const LduMatrix<Type, DType, LUType>& matrix,
214  const dictionary& solverDict
215  );
216 
217 
218  // Destructor
220  virtual ~solver()
221  {}
222 
223 
224  // Member functions
225 
226  // Access
228  const word& fieldName() const
229  {
230  return fieldName_;
231  }
234  {
235  return matrix_;
236  }
237 
238 
239  //- Read and reset the solver parameters from the given dictionary
240  virtual void read(const dictionary& solverDict);
241 
243  (
245  ) const = 0;
246 
247  //- Return the matrix norm used to normalise the residual for the
248  // stopping criterion
249  Type normFactor
250  (
251  const Field<Type>& psi,
252  const Field<Type>& Apsi,
253  Field<Type>& tmpField
254  ) const;
255  };
256 
257 
258  //- Abstract base-class for LduMatrix smoothers
259  class smoother
260  {
261  protected:
262 
263  // Protected data
267 
268 
269  public:
270 
271  //- Runtime type information
272  virtual const word& type() const = 0;
273 
274 
275  // Declare run-time constructor selection tables
276 
278  (
279  autoPtr,
280  smoother,
281  symMatrix,
282  (
283  const word& fieldName,
284  const LduMatrix<Type, DType, LUType>& matrix
285  ),
286  (
287  fieldName,
288  matrix
289  )
290  );
291 
293  (
294  autoPtr,
295  smoother,
296  asymMatrix,
297  (
298  const word& fieldName,
299  const LduMatrix<Type, DType, LUType>& matrix
300  ),
301  (
302  fieldName,
303  matrix
304  )
305  );
306 
307 
308  // Constructors
309 
310  smoother
311  (
312  const word& fieldName,
313  const LduMatrix<Type, DType, LUType>& matrix
314  );
315 
316 
317  // Selectors
318 
319  //- Return a new smoother
320  static autoPtr<smoother> New
321  (
322  const word& fieldName,
323  const LduMatrix<Type, DType, LUType>& matrix,
324  const dictionary& smootherDict
325  );
326 
327 
328  // Destructor
330  virtual ~smoother()
331  {}
332 
333 
334  // Member functions
335 
336  // Access
338  const word& fieldName() const
339  {
340  return fieldName_;
341  }
344  {
345  return matrix_;
346  }
347 
348 
349  //- Smooth the solution for a given number of sweeps
350  virtual void smooth
351  (
352  Field<Type>& psi,
353  const label nSweeps
354  ) const = 0;
355  };
356 
357 
358  //- Abstract base-class for LduMatrix preconditioners
359  class preconditioner
360  {
361  protected:
362 
363  // Protected data
364 
365  //- Reference to the base-solver this preconditioner is used with
366  const solver& solver_;
367 
368 
369  public:
370 
371  //- Runtime type information
372  virtual const word& type() const = 0;
373 
374 
375  // Declare run-time constructor selection tables
376 
378  (
379  autoPtr,
381  symMatrix,
382  (
383  const solver& sol,
384  const dictionary& preconditionerDict
385  ),
386  (sol, preconditionerDict)
387  );
388 
390  (
391  autoPtr,
393  asymMatrix,
394  (
395  const solver& sol,
396  const dictionary& preconditionerDict
397  ),
398  (sol, preconditionerDict)
399  );
400 
401 
402  // Constructors
403 
405  (
406  const solver& sol
407  )
408  :
409  solver_(sol)
410  {}
411 
412 
413  // Selectors
414 
415  //- Return a new preconditioner
417  (
418  const solver& sol,
419  const dictionary& preconditionerDict
420  );
421 
422 
423  // Destructor
425  virtual ~preconditioner()
426  {}
427 
428 
429  // Member functions
430 
431  //- Read and reset the preconditioner parameters
432  // from the given dictionary
433  virtual void read(const dictionary& preconditionerDict)
434  {}
435 
436  //- Return wA the preconditioned form of residual rA
437  virtual void precondition
438  (
439  Field<Type>& wA,
440  const Field<Type>& rA
441  ) const = 0;
442 
443  //- Return wT the transpose-matrix preconditioned form of
444  // residual rT.
445  // This is only required for preconditioning asymmetric matrices.
446  virtual void preconditionT
447  (
448  Field<Type>& wT,
449  const Field<Type>& rT
450  ) const
451  {
453  }
454  };
455 
456 
457  // Static data
458 
459  // Declare name of the class and its debug switch
460  ClassName("LduMatrix");
461 
462 
463  // Constructors
464 
465  //- Construct given an LDU addressed mesh.
466  // The coefficients are initially empty for subsequent setting.
467  LduMatrix(const lduMesh&);
468 
469  //- Construct as copy
471 
472  //- Construct as copy or re-use as specified.
474 
475  //- Construct given an LDU addressed mesh and an Istream
476  // from which the coefficients are read
477  LduMatrix(const lduMesh&, Istream&);
478 
479 
480  // Destructor
481 
482  ~LduMatrix();
483 
484 
485  // Member functions
486 
487  // Access to addressing
488 
489  //- Return the LDU mesh from which the addressing is obtained
490  const lduMesh& mesh() const
491  {
492  return lduMesh_;
493  }
494 
495  //- Return the LDU addressing
496  const lduAddressing& lduAddr() const
497  {
498  return lduMesh_.lduAddr();
499  }
500 
501  //- Return the patch evaluation schedule
502  const lduSchedule& patchSchedule() const
503  {
504  return lduAddr().patchSchedule();
505  }
506 
507  //- Return interfaces
509  {
510  return interfaces_;
511  }
512 
513  //- Return interfaces
515  {
516  return interfaces_;
517  }
518 
519 
520  // Access to coefficients
521 
522  Field<DType>& diag();
523  Field<LUType>& upper();
524  Field<LUType>& lower();
525  Field<Type>& source();
528  {
529  return interfacesUpper_;
530  }
533  {
534  return interfacesLower_;
535  }
536 
537 
538  const Field<DType>& diag() const;
539  const Field<LUType>& upper() const;
540  const Field<LUType>& lower() const;
541  const Field<Type>& source() const;
544  {
545  return interfacesUpper_;
546  }
549  {
550  return interfacesLower_;
551  }
552 
554  bool hasDiag() const
555  {
556  return (diagPtr_);
557  }
559  bool hasUpper() const
560  {
561  return (upperPtr_);
562  }
564  bool hasLower() const
565  {
566  return (lowerPtr_);
567  }
569  bool hasSource() const
570  {
571  return (sourcePtr_);
572  }
574  bool diagonal() const
575  {
576  return (diagPtr_ && !lowerPtr_ && !upperPtr_);
577  }
579  bool symmetric() const
580  {
581  return (diagPtr_ && (!lowerPtr_ && upperPtr_));
582  }
584  bool asymmetric() const
585  {
586  return (diagPtr_ && lowerPtr_ && upperPtr_);
587  }
588 
589 
590  // operations
591 
592  void sumDiag();
593  void negSumDiag();
594 
595  void sumMagOffDiag(Field<LUType>& sumOff) const;
596 
597  //- Matrix multiplication
598  void Amul(Field<Type>&, const tmp<Field<Type>>&) const;
599 
600  //- Matrix transpose multiplication
601  void Tmul(Field<Type>&, const tmp<Field<Type>>&) const;
602 
603 
604  //- Sum the coefficients on each row of the matrix
605  void sumA(Field<Type>&) const;
606 
607 
608  void residual(Field<Type>& rA, const Field<Type>& psi) const;
609 
610  tmp<Field<Type>> residual(const Field<Type>& psi) const;
611 
612 
613  //- Initialise the update of interfaced interfaces
614  // for matrix operations
616  (
617  const FieldField<Field, LUType>& interfaceCoeffs,
618  const Field<Type>& psiif,
619  Field<Type>& result
620  ) const;
621 
622  //- Update interfaced interfaces for matrix operations
624  (
625  const FieldField<Field, LUType>& interfaceCoeffs,
626  const Field<Type>& psiif,
627  Field<Type>& result
628  ) const;
629 
630 
631  tmp<Field<Type>> H(const Field<Type>&) const;
632  tmp<Field<Type>> H(const tmp<Field<Type>>&) const;
633 
634  tmp<Field<Type>> faceH(const Field<Type>&) const;
635  tmp<Field<Type>> faceH(const tmp<Field<Type>>&) const;
636 
637 
638  // Member operators
639 
641 
642  void negate();
643 
646 
647  void operator*=(const scalarField&);
648  void operator*=(scalar);
649 
650 
651  // Ostream operator
652 
653  friend Ostream& operator<< <Type, DType, LUType>
654  (
655  Ostream&,
657  );
658 };
659 
660 
661 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
662 
663 } // End namespace Foam
664 
665 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
667 #define makeLduMatrix(Type, DType, LUType) \
668  \
669 typedef Foam::LduMatrix<Type, DType, LUType> \
670  ldu##Type##DType##LUType##Matrix; \
671  \
672 defineNamedTemplateTypeNameAndDebug(ldu##Type##DType##LUType##Matrix, 0); \
673  \
674  \
675 typedef LduMatrix<Type, DType, LUType>::smoother \
676  ldu##Type##DType##LUType##Smoother; \
677  \
678 defineTemplateRunTimeSelectionTable \
679 ( \
680  ldu##Type##DType##LUType##Smoother, \
681  symMatrix \
682 ); \
683  \
684 defineTemplateRunTimeSelectionTable \
685 ( \
686  ldu##Type##DType##LUType##Smoother, \
687  asymMatrix \
688 ); \
689  \
690  \
691 typedef LduMatrix<Type, DType, LUType>::preconditioner \
692  ldu##Type##DType##LUType##Preconditioner; \
693  \
694 defineTemplateRunTimeSelectionTable \
695 ( \
696  ldu##Type##DType##LUType##Preconditioner, \
697  symMatrix \
698 ); \
699  \
700 defineTemplateRunTimeSelectionTable \
701 ( \
702  ldu##Type##DType##LUType##Preconditioner, \
703  asymMatrix \
704 ); \
705  \
706  \
707 typedef LduMatrix<Type, DType, LUType>::solver \
708  ldu##Type##DType##LUType##Solver; \
709  \
710 defineTemplateRunTimeSelectionTable \
711 ( \
712  ldu##Type##DType##LUType##Solver, \
713  symMatrix \
714 ); \
715  \
716 defineTemplateRunTimeSelectionTable \
717 ( \
718  ldu##Type##DType##LUType##Solver, \
719  asymMatrix \
720 );
721 
723 #define makeLduPreconditioner(Precon, Type, DType, LUType) \
724  \
725 typedef Precon<Type, DType, LUType> \
726  Precon##Type##DType##LUType##Preconditioner; \
727 defineNamedTemplateTypeNameAndDebug \
728 ( \
729  Precon##Type##DType##LUType##Preconditioner, \
730  0 \
731 );
733 #define makeLduSymPreconditioner(Precon, Type, DType, LUType) \
734  \
735 LduMatrix<Type, DType, LUType>::preconditioner:: \
736 addsymMatrixConstructorToTable<Precon##Type##DType##LUType##Preconditioner> \
737 add##Precon##Type##DType##LUType##PreconditionerSymMatrixConstructorToTable_;
739 #define makeLduAsymPreconditioner(Precon, Type, DType, LUType) \
740  \
741 LduMatrix<Type, DType, LUType>::preconditioner:: \
742 addasymMatrixConstructorToTable<Precon##Type##DType##LUType##Preconditioner> \
743 add##Precon##Type##DType##LUType##PreconditionerAsymMatrixConstructorToTable_;
744 
746 #define makeLduSmoother(Smoother, Type, DType, LUType) \
747  \
748 typedef Smoother<Type, DType, LUType> \
749  Smoother##Type##DType##LUType##Smoother; \
750  \
751 defineNamedTemplateTypeNameAndDebug \
752 ( \
753  Smoother##Type##DType##LUType##Smoother, \
754  0 \
755 );
757 #define makeLduSymSmoother(Smoother, Type, DType, LUType) \
758  \
759 LduMatrix<Type, DType, LUType>::smoother:: \
760  addsymMatrixConstructorToTable<Smoother##Type##DType##LUType##Smoother> \
761  add##Smoother##Type##DType##LUType##SymMatrixConstructorToTable_;
763 #define makeLduAsymSmoother(Smoother, Type, DType, LUType) \
764  \
765 LduMatrix<Type, DType, LUType>::smoother:: \
766  addasymMatrixConstructorToTable<Smoother##Type##DType##LUType##Smoother> \
767  add##Smoother##Type##DType##LUType##AsymMatrixConstructorToTable_;
768 
770 #define makeLduSolver(Solver, Type, DType, LUType) \
771  \
772 typedef Solver<Type, DType, LUType> \
773  Solver##Type##DType##LUType##Solver; \
774  \
775 defineNamedTemplateTypeNameAndDebug \
776 ( \
777  Solver##Type##DType##LUType##Solver, \
778  0 \
779 );
781 #define makeLduSymSolver(Solver, Type, DType, LUType) \
782  \
783 LduMatrix<Type, DType, LUType>::solver:: \
784  addsymMatrixConstructorToTable<Solver##Type##DType##LUType##Solver> \
785  add##Solver##Type##DType##LUType##SymMatrixConstructorToTable_;
787 #define makeLduAsymSolver(Solver, Type, DType, LUType) \
788  \
789 LduMatrix<Type, DType, LUType>::solver:: \
790  addasymMatrixConstructorToTable<Solver##Type##DType##LUType##Solver> \
791  add##Solver##Type##DType##LUType##AsymMatrixConstructorToTable_;
792 
793 
794 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
795 
796 #ifdef NoRepository
797  #include "LduMatrixI.H"
798  #include "LduMatrix.C"
799 #endif
800 
801 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
802 
803 #endif
804 
805 // ************************************************************************* //
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
virtual const word & type() const =0
Runtime type information.
runTime controlDict().lookup("adjustTimeStep") >> adjustTimeStep
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const lduMesh & mesh() const
Return the LDU mesh from which the addressing is obtained.
Definition: LduMatrix.H:489
Field< Type > & source()
Definition: LduMatrix.C:248
tmp< Field< Type > > H(const Field< Type > &) const
FieldField< Field, LUType > & interfacesLower()
Definition: LduMatrix.H:531
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
bool hasDiag() const
Definition: LduMatrix.H:553
LduMatrix(const lduMesh &)
Construct given an LDU addressed mesh.
Definition: LduMatrix.C:32
void Amul(Field< Type > &, const tmp< Field< Type >> &) const
Matrix multiplication.
void residual(Field< Type > &rA, const Field< Type > &psi) const
bool hasLower() const
Definition: LduMatrix.H:563
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:59
bool asymmetric() const
Definition: LduMatrix.H:583
Field< LUType > & lower()
Definition: LduMatrix.C:225
Type normFactor(const Field< Type > &psi, const Field< Type > &Apsi, Field< Type > &tmpField) const
Return the matrix norm used to normalise the residual for the.
label minIter_
Minimum number of iterations in the solver.
Definition: LduMatrix.H:128
bool hasSource() const
Definition: LduMatrix.H:568
Generic field type.
Definition: FieldField.H:51
void smooth(volScalarField &field, const scalar coeff)
Definition: fvcSmooth.C:35
dictionary controlDict_
Dictionary of controls.
Definition: LduMatrix.H:122
const word & fieldName() const
Definition: LduMatrix.H:227
ClassName("LduMatrix")
void Tmul(Field< Type > &, const tmp< Field< Type >> &) const
Matrix transpose multiplication.
const LduInterfaceFieldPtrsList< Type > & interfaces() const
Return interfaces.
Definition: LduMatrix.H:507
Abstract base-class for LduMatrix smoothers.
Definition: LduMatrix.H:258
A class for handling words, derived from string.
Definition: word.H:59
void readControl(const dictionary &controlDict, T &control, const word &controlName)
Read a control parameter from controlDict.
Definition: LduMatrixI.H:31
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics...
Abstract base-class for LduMatrix preconditioners.
Definition: LduMatrix.H:358
void operator*=(const scalarField &)
void initMatrixInterfaces(const FieldField< Field, LUType > &interfaceCoeffs, const Field< Type > &psiif, Field< Type > &result) const
Initialise the update of interfaced interfaces.
Field< LUType > & upper()
Definition: LduMatrix.C:202
const LduMatrix< Type, DType, LUType > & matrix() const
Definition: LduMatrix.H:232
void operator-=(const LduMatrix< Type, DType, LUType > &)
Abstract base-class for LduMatrix solvers.
Definition: LduMatrix.H:112
void sumA(Field< Type > &) const
Sum the coefficients on each row of the matrix.
solver(const word &fieldName, const LduMatrix< Type, DType, LUType > &matrix, const dictionary &solverDict)
const lduSchedule & patchSchedule() const
Return the patch evaluation schedule.
Definition: LduMatrix.H:501
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void updateMatrixInterfaces(const FieldField< Field, LUType > &interfaceCoeffs, const Field< Type > &psiif, Field< Type > &result) const
Update interfaced interfaces for matrix operations.
Field< DType > & diag()
Definition: LduMatrix.C:190
void sumMagOffDiag(Field< LUType > &sumOff) const
Type relTol_
Convergence tolerance relative to the initial.
Definition: LduMatrix.H:134
FieldField< Field, LUType > & interfacesUpper()
Definition: LduMatrix.H:526
LduMatrix is a general matrix class in which the coefficients are stored as three arrays...
Definition: LduMatrix.H:69
virtual void read(const dictionary &solverDict)
Read and reset the solver parameters from the given dictionary.
label maxIter_
Maximum number of iterations in the solver.
Definition: LduMatrix.H:125
const lduAddressing & lduAddr() const
Return the LDU addressing.
Definition: LduMatrix.H:495
void operator+=(const LduMatrix< Type, DType, LUType > &)
bool symmetric() const
Definition: LduMatrix.H:578
static autoPtr< solver > New(const word &fieldName, const LduMatrix< Type, DType, LUType > &matrix, const dictionary &solverDict)
Return a new solver.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
const volScalarField & psi
The class contains the addressing required by the lduMatrix: upper, lower and losort.
virtual SolverPerformance< Type > solve(Field< Type > &psi) const =0
const LduMatrix< Type, DType, LUType > & matrix_
Definition: LduMatrix.H:119
tmp< Field< Type > > faceH(const Field< Type > &) const
Macros to ease declaration of run-time selection tables.
virtual void readControls()
Read the control parameters from the controlDict_.
bool hasUpper() const
Definition: LduMatrix.H:558
A class for managing temporary objects.
Definition: PtrList.H:54
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:366
declareRunTimeSelectionTable(autoPtr, solver, symMatrix,(const word &fieldName, const LduMatrix< Type, DType, LUType > &matrix, const dictionary &solverDict),(fieldName, matrix, solverDict))
void operator=(const LduMatrix< Type, DType, LUType > &)
List of coupled interface fields to be used in coupling.
Type tolerance_
Final convergence tolerance.
Definition: LduMatrix.H:131
Namespace for OpenFOAM.
bool diagonal() const
Definition: LduMatrix.H:573
virtual const lduSchedule & patchSchedule() const =0