fvMatrix.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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::fvMatrix
26 
27 Description
28  A special matrix type and solver, designed for finite volume
29  solutions of scalar equations.
30  Face addressing is used to make all matrix assembly
31  and solution loops vectorise.
32 
33 SourceFiles
34  fvMatrix.C
35  fvMatrixSolve.C
36  fvScalarMatrix.C
37 
38 \*---------------------------------------------------------------------------*/
39 
40 #ifndef fvMatrix_H
41 #define fvMatrix_H
42 
43 #include "volFields.H"
44 #include "surfaceFields.H"
45 #include "lduMatrix.H"
46 #include "tmp.H"
47 #include "autoPtr.H"
48 #include "dimensionedTypes.H"
49 #include "zero.H"
50 #include "className.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 // Forward declaration of friend functions and operators
58 
59 template<class Type>
60 class fvMatrix;
61 
62 template<class Type>
63 tmp<GeometricField<Type, fvPatchField, volMesh>> operator&
64 (
65  const fvMatrix<Type>&,
66  const DimensionedField<Type, volMesh>&
67 );
68 
69 template<class Type>
70 tmp<GeometricField<Type, fvPatchField, volMesh>> operator&
71 (
72  const fvMatrix<Type>&,
73  const tmp<DimensionedField<Type, volMesh>>&
74 );
75 
76 template<class Type>
77 tmp<GeometricField<Type, fvPatchField, volMesh>> operator&
78 (
79  const fvMatrix<Type>&,
80  const tmp<GeometricField<Type, fvPatchField, volMesh>>&
81 );
82 
83 template<class Type>
84 tmp<GeometricField<Type, fvPatchField, volMesh>> operator&
85 (
86  const tmp<fvMatrix<Type>>&,
87  const DimensionedField<Type, volMesh>&
88 );
89 
90 template<class Type>
91 tmp<GeometricField<Type, fvPatchField, volMesh>> operator&
92 (
93  const tmp<fvMatrix<Type>>&,
94  const tmp<DimensionedField<Type, volMesh>>&
95 );
96 
97 template<class Type>
98 tmp<GeometricField<Type, fvPatchField, volMesh>> operator&
99 (
100  const tmp<fvMatrix<Type>>&,
101  const tmp<GeometricField<Type, fvPatchField, volMesh>>&
102 );
103 
104 template<class Type>
105 Ostream& operator<<(Ostream&, const fvMatrix<Type>&);
107 template<class T> class UIndirectList;
108 
109 
110 /*---------------------------------------------------------------------------*\
111  Class fvMatrix Declaration
112 \*---------------------------------------------------------------------------*/
113 
114 template<class Type>
115 class fvMatrix
116 :
117  public tmp<fvMatrix<Type>>::refCount,
118  public lduMatrix
119 {
120  // Private Data
121 
122  //- Const reference to GeometricField<Type, fvPatchField, volMesh>
123  // Converted into a non-const reference at the point of solution.
125 
126  //- Dimension set
127  dimensionSet dimensions_;
128 
129  //- Source term
130  Field<Type> source_;
131 
132  //- Boundary scalar field containing pseudo-matrix coeffs
133  // for internal cells
134  FieldField<Field, Type> internalCoeffs_;
135 
136  //- Boundary scalar field containing pseudo-matrix coeffs
137  // for boundary cells
138  FieldField<Field, Type> boundaryCoeffs_;
139 
140 
141  //- Face flux field for non-orthogonal correction
143  *faceFluxCorrectionPtr_;
144 
145 
146 protected:
147 
148  //- Declare friendship with the fvSolver class
149  friend class fvSolver;
150 
151  // Protected Member Functions
152 
153  //- Add patch contribution to internal field
154  template<class Type2>
155  void addToInternalField
156  (
157  const labelUList& addr,
158  const Field<Type2>& pf,
159  Field<Type2>& intf
160  ) const;
161 
162  template<class Type2>
163  void addToInternalField
164  (
165  const labelUList& addr,
166  const tmp<Field<Type2>>& tpf,
167  Field<Type2>& intf
168  ) const;
169 
170  //- Subtract patch contribution from internal field
171  template<class Type2>
173  (
174  const labelUList& addr,
175  const Field<Type2>& pf,
176  Field<Type2>& intf
177  ) const;
178 
179  template<class Type2>
181  (
182  const labelUList& addr,
183  const tmp<Field<Type2>>& tpf,
184  Field<Type2>& intf
185  ) const;
186 
187 
188  // Matrix completion functionality
189 
190  void addBoundaryDiag
191  (
192  scalarField& diag,
193  const direction cmpt
194  ) const;
195 
197 
198  void addBoundarySource
199  (
201  const bool couples=true
202  ) const;
203 
204  // Matrix manipulation functionality
205 
206  //- Set solution in given cells to the specified values
207  template<template<class> class ListType>
208  void setValuesFromList
209  (
210  const labelUList& cells,
211  const ListType<Type>& values
212  );
213 
214 
215 public:
216 
217  //- Solver class returned by the solver function
218  // used for systems in which it is useful to cache the solver for reuse
219  // e.g. if the solver is potentially expensive to construct (AMG) and can
220  // be used several times (PISO)
221  class fvSolver
222  {
223  fvMatrix<Type>& fvMat_;
224 
226 
227  public:
228 
229  // Constructors
231  fvSolver(fvMatrix<Type>& fvMat, autoPtr<lduMatrix::solver> sol)
232  :
233  fvMat_(fvMat),
234  solver_(sol)
235  {}
236 
237 
238  // Member Functions
239 
240  //- Solve returning the solution statistics.
241  // Use the given solver controls
243 
244  //- Solve returning the solution statistics.
245  // Solver controls read from fvSolution
247  };
248 
249 
250  ClassName("fvMatrix");
251 
252 
253  // Constructors
254 
255  //- Construct given a field to solve for
256  fvMatrix
257  (
259  const dimensionSet&
260  );
261 
262  //- Copy constructor
263  fvMatrix(const fvMatrix<Type>&);
264 
265  //- Copy constructor of tmp<fvMatrix<Type>> deleting argument
266  fvMatrix(const tmp<fvMatrix<Type>>&);
267 
268  //- Construct from Istream given field to solve for
270 
271  //- Clone
272  tmp<fvMatrix<Type>> clone() const;
273 
274 
275  //- Destructor
276  virtual ~fvMatrix();
277 
278 
279  // Member Functions
280 
281  // Access
284  {
285  return psi_;
286  }
288  const dimensionSet& dimensions() const
289  {
290  return dimensions_;
291  }
294  {
295  return source_;
296  }
298  const Field<Type>& source() const
299  {
300  return source_;
301  }
302 
303  //- fvBoundary scalar field containing pseudo-matrix coeffs
304  // for internal cells
306  {
307  return internalCoeffs_;
308  }
309 
310  //- fvBoundary scalar field containing pseudo-matrix coeffs
311  // for boundary cells
313  {
314  return boundaryCoeffs_;
315  }
316 
317 
318  //- Declare return type of the faceFluxCorrectionPtr() function
321 
322  //- Return pointer to face-flux non-orthogonal correction field
323  surfaceTypeFieldPtr& faceFluxCorrectionPtr()
324  {
325  return faceFluxCorrectionPtr_;
326  }
327 
328 
329  // Operations
330 
331  //- Set solution in given cells to the specified values
332  // and eliminate the corresponding equations from the matrix.
333  void setValues
334  (
335  const labelUList& cells,
336  const UList<Type>& values
337  );
338 
339  //- Set solution in given cells to the specified values
340  // and eliminate the corresponding equations from the matrix.
341  void setValues
342  (
343  const labelUList& cells,
344  const UIndirectList<Type>& values
345  );
346 
347  //- Set reference level for solution
348  void setReference
349  (
350  const label celli,
351  const Type& value,
352  const bool forceReference = false
353  );
354 
355  //- Set reference level for a component of the solution
356  // on a given patch face
358  (
359  const label patchi,
360  const label facei,
361  const direction cmpt,
362  const scalar value
363  );
364 
365  //- Relax matrix (for steady-state solution).
366  // alpha = 1 : diagonally equal
367  // alpha < 1 : diagonally dominant
368  // alpha = 0 : do nothing
369  // Note: Requires positive diagonal.
370  void relax(const scalar alpha);
371 
372  //- Relax matrix (for steady-state solution).
373  // alpha is read from controlDict
374  void relax();
375 
376  //- Manipulate based on a boundary field
377  void boundaryManipulate
378  (
380  Boundary& values
381  );
382 
383  //- Construct and return the solver
384  // Use the given solver controls
386 
387  //- Construct and return the solver
388  // Solver controls read from fvSolution
390 
391  //- Solve segregated or coupled returning the solution statistics.
392  // Use the given solver controls
394 
395  //- Solve segregated returning the solution statistics.
396  // Use the given solver controls
398 
399  //- Solve coupled returning the solution statistics.
400  // Use the given solver controls
402 
403  //- Solve segregated or coupled returning the solution statistics.
404  // Solver controls read from fvSolution
406 
407  //- Solve returning the solution statistics.
408  // Solver controls read from fvSolution
410 
411  //- Return the matrix residual
412  tmp<Field<Type>> residual() const;
413 
414  //- Return the matrix scalar diagonal
415  tmp<scalarField> D() const;
416 
417  //- Return the matrix Type diagonal
418  tmp<Field<Type>> DD() const;
419 
420  //- Return the central coefficient
421  tmp<volScalarField> A() const;
422 
423  //- Return the H operation source
425 
426  //- Return H(1)
427  tmp<volScalarField> H1() const;
428 
429  //- Return the face-flux field from the matrix
431  flux() const;
432 
433 
434  // Member Operators
435 
436  void operator=(const fvMatrix<Type>&);
437  void operator=(const tmp<fvMatrix<Type>>&);
438 
439  void negate();
440 
441  void operator+=(const fvMatrix<Type>&);
442  void operator+=(const tmp<fvMatrix<Type>>&);
443 
444  void operator-=(const fvMatrix<Type>&);
445  void operator-=(const tmp<fvMatrix<Type>>&);
446 
447  void operator+=
448  (
450  );
451  void operator+=
452  (
454  );
455  void operator+=
456  (
458  );
459 
460  void operator-=
461  (
463  );
464  void operator-=
465  (
467  );
468  void operator-=
469  (
471  );
472 
473  void operator+=(const dimensioned<Type>&);
474  void operator-=(const dimensioned<Type>&);
475 
476  void operator+=(const zero&);
477  void operator-=(const zero&);
478 
481  void operator*=(const tmp<volScalarField>&);
482 
483  void operator*=(const dimensioned<scalar>&);
484 
485 
486  // Friend operators
487 
489  operator& <Type>
490  (
491  const fvMatrix<Type>&,
493  );
494 
496  operator& <Type>
497  (
498  const fvMatrix<Type>&,
500  );
501 
503  operator& <Type>
504  (
505  const tmp<fvMatrix<Type>>&,
507  );
508 
510  operator& <Type>
511  (
512  const tmp<fvMatrix<Type>>&,
514  );
515 
516 
517  // Ostream operator
518 
519  friend Ostream& operator<< <Type>
520  (
521  Ostream&,
522  const fvMatrix<Type>&
523  );
524 };
525 
526 
527 // * * * * * * * * * * * * * * * Global functions * * * * * * * * * * * * * //
528 
529 template<class Type>
530 void checkMethod
531 (
532  const fvMatrix<Type>&,
533  const fvMatrix<Type>&,
534  const char*
535 );
536 
537 template<class Type>
538 void checkMethod
539 (
540  const fvMatrix<Type>&,
542  const char*
543 );
544 
545 template<class Type>
546 void checkMethod
547 (
548  const fvMatrix<Type>&,
549  const dimensioned<Type>&,
550  const char*
551 );
552 
553 
554 //- Solve returning the solution statistics given convergence tolerance
555 // Use the given solver controls
556 template<class Type>
557 SolverPerformance<Type> solve(fvMatrix<Type>&, const word&);
558 
559 //- Solve returning the solution statistics given convergence tolerance,
560 // deleting temporary matrix after solution.
561 // Use the given solver controls
562 template<class Type>
564 (
565  const tmp<fvMatrix<Type>>&,
566  const word&
567 );
568 
569 //- Solve returning the solution statistics given convergence tolerance
570 // Solver controls read fvSolution
571 template<class Type>
572 SolverPerformance<Type> solve(fvMatrix<Type>&);
573 
574 //- Solve returning the solution statistics given convergence tolerance,
575 // deleting temporary matrix after solution.
576 // Solver controls read fvSolution
577 template<class Type>
578 SolverPerformance<Type> solve(const tmp<fvMatrix<Type>>&);
579 
580 //- Return the correction form of the given matrix
581 // by subtracting the matrix multiplied by the current field
582 template<class Type>
583 tmp<fvMatrix<Type>> correction(const fvMatrix<Type>&);
584 
585 //- Return the correction form of the given temporary matrix
586 // by subtracting the matrix multiplied by the current field
587 template<class Type>
588 tmp<fvMatrix<Type>> correction(const tmp<fvMatrix<Type>>&);
589 
590 
591 // * * * * * * * * * * * * * * * Global operators * * * * * * * * * * * * * //
592 
593 template<class Type>
594 tmp<fvMatrix<Type>> operator==
595 (
596  const fvMatrix<Type>&,
597  const fvMatrix<Type>&
598 );
599 
600 template<class Type>
601 tmp<fvMatrix<Type>> operator==
602 (
603  const tmp<fvMatrix<Type>>&,
604  const fvMatrix<Type>&
605 );
606 
607 template<class Type>
608 tmp<fvMatrix<Type>> operator==
609 (
610  const fvMatrix<Type>&,
611  const tmp<fvMatrix<Type>>&
612 );
613 
614 template<class Type>
615 tmp<fvMatrix<Type>> operator==
616 (
617  const tmp<fvMatrix<Type>>&,
618  const tmp<fvMatrix<Type>>&
619 );
620 
621 
622 template<class Type>
623 tmp<fvMatrix<Type>> operator==
624 (
625  const fvMatrix<Type>&,
627 );
628 
629 template<class Type>
630 tmp<fvMatrix<Type>> operator==
631 (
632  const fvMatrix<Type>&,
634 );
635 
636 template<class Type>
637 tmp<fvMatrix<Type>> operator==
638 (
639  const fvMatrix<Type>&,
641 );
642 
643 template<class Type>
644 tmp<fvMatrix<Type>> operator==
645 (
646  const tmp<fvMatrix<Type>>&,
648 );
649 
650 template<class Type>
651 tmp<fvMatrix<Type>> operator==
652 (
653  const tmp<fvMatrix<Type>>&,
655 );
656 
657 template<class Type>
658 tmp<fvMatrix<Type>> operator==
659 (
660  const tmp<fvMatrix<Type>>&,
662 );
663 
664 template<class Type>
665 tmp<fvMatrix<Type>> operator==
666 (
667  const fvMatrix<Type>&,
668  const dimensioned<Type>&
669 );
670 
671 template<class Type>
672 tmp<fvMatrix<Type>> operator==
673 (
674  const tmp<fvMatrix<Type>>&,
675  const dimensioned<Type>&
676 );
677 
678 
679 template<class Type>
680 tmp<fvMatrix<Type>> operator==
681 (
682  const fvMatrix<Type>&,
683  const zero&
684 );
685 
686 template<class Type>
687 tmp<fvMatrix<Type>> operator==
688 (
689  const tmp<fvMatrix<Type>>&,
690  const zero&
691 );
692 
693 
694 template<class Type>
695 tmp<fvMatrix<Type>> operator-
696 (
697  const fvMatrix<Type>&
698 );
699 
700 template<class Type>
701 tmp<fvMatrix<Type>> operator-
702 (
703  const tmp<fvMatrix<Type>>&
704 );
705 
706 
707 template<class Type>
708 tmp<fvMatrix<Type>> operator+
709 (
710  const fvMatrix<Type>&,
711  const fvMatrix<Type>&
712 );
713 
714 template<class Type>
715 tmp<fvMatrix<Type>> operator+
716 (
717  const tmp<fvMatrix<Type>>&,
718  const fvMatrix<Type>&
719 );
720 
721 template<class Type>
722 tmp<fvMatrix<Type>> operator+
723 (
724  const fvMatrix<Type>&,
725  const tmp<fvMatrix<Type>>&
726 );
727 
728 template<class Type>
729 tmp<fvMatrix<Type>> operator+
730 (
731  const tmp<fvMatrix<Type>>&,
732  const tmp<fvMatrix<Type>>&
733 );
734 
735 
736 template<class Type>
737 tmp<fvMatrix<Type>> operator+
738 (
739  const fvMatrix<Type>&,
741 );
742 
743 template<class Type>
744 tmp<fvMatrix<Type>> operator+
745 (
746  const fvMatrix<Type>&,
748 );
749 
750 template<class Type>
751 tmp<fvMatrix<Type>> operator+
752 (
753  const fvMatrix<Type>&,
755 );
756 
757 template<class Type>
758 tmp<fvMatrix<Type>> operator+
759 (
760  const tmp<fvMatrix<Type>>&,
762 );
763 
764 template<class Type>
765 tmp<fvMatrix<Type>> operator+
766 (
767  const tmp<fvMatrix<Type>>&,
769 );
770 
771 template<class Type>
772 tmp<fvMatrix<Type>> operator+
773 (
774  const tmp<fvMatrix<Type>>&,
776 );
777 
778 template<class Type>
779 tmp<fvMatrix<Type>> operator+
780 (
782  const fvMatrix<Type>&
783 );
784 
785 template<class Type>
786 tmp<fvMatrix<Type>> operator+
787 (
789  const fvMatrix<Type>&
790 );
791 
792 template<class Type>
793 tmp<fvMatrix<Type>> operator+
794 (
796  const fvMatrix<Type>&
797 );
798 
799 template<class Type>
800 tmp<fvMatrix<Type>> operator+
801 (
802  const DimensionedField<Type, volMesh>&,
803  const tmp<fvMatrix<Type>>&
804 );
805 
806 template<class Type>
807 tmp<fvMatrix<Type>> operator+
808 (
809  const tmp<DimensionedField<Type, volMesh>>&,
810  const tmp<fvMatrix<Type>>&
811 );
812 
813 template<class Type>
814 tmp<fvMatrix<Type>> operator+
815 (
816  const tmp<GeometricField<Type, fvPatchField, volMesh>>&,
817  const tmp<fvMatrix<Type>>&
818 );
819 
820 
821 template<class Type>
822 tmp<fvMatrix<Type>> operator+
823 (
824  const fvMatrix<Type>&,
825  const dimensioned<Type>&
826 );
827 
828 template<class Type>
829 tmp<fvMatrix<Type>> operator+
830 (
831  const tmp<fvMatrix<Type>>&,
832  const dimensioned<Type>&
833 );
834 
835 template<class Type>
836 tmp<fvMatrix<Type>> operator+
837 (
838  const dimensioned<Type>&,
839  const fvMatrix<Type>&
840 );
841 
842 template<class Type>
843 tmp<fvMatrix<Type>> operator+
844 (
845  const dimensioned<Type>&,
846  const tmp<fvMatrix<Type>>&
847 );
848 
849 
850 template<class Type>
851 tmp<fvMatrix<Type>> operator-
852 (
853  const fvMatrix<Type>&,
854  const fvMatrix<Type>&
855 );
856 
857 template<class Type>
858 tmp<fvMatrix<Type>> operator-
859 (
860  const tmp<fvMatrix<Type>>&,
861  const fvMatrix<Type>&
862 );
863 
864 template<class Type>
865 tmp<fvMatrix<Type>> operator-
866 (
867  const fvMatrix<Type>&,
868  const tmp<fvMatrix<Type>>&
869 );
870 
871 template<class Type>
872 tmp<fvMatrix<Type>> operator-
873 (
874  const tmp<fvMatrix<Type>>&,
875  const tmp<fvMatrix<Type>>&
876 );
877 
878 
879 template<class Type>
880 tmp<fvMatrix<Type>> operator-
881 (
882  const fvMatrix<Type>&,
883  const DimensionedField<Type, volMesh>&
884 );
885 
886 template<class Type>
887 tmp<fvMatrix<Type>> operator-
888 (
889  const fvMatrix<Type>&,
890  const tmp<DimensionedField<Type, volMesh>>&
891 );
892 
893 template<class Type>
894 tmp<fvMatrix<Type>> operator-
895 (
896  const fvMatrix<Type>&,
897  const tmp<GeometricField<Type, fvPatchField, volMesh>>&
898 );
899 
900 template<class Type>
901 tmp<fvMatrix<Type>> operator-
902 (
903  const tmp<fvMatrix<Type>>&,
904  const DimensionedField<Type, volMesh>&
905 );
906 
907 template<class Type>
908 tmp<fvMatrix<Type>> operator-
909 (
910  const tmp<fvMatrix<Type>>&,
911  const tmp<DimensionedField<Type, volMesh>>&
912 );
913 
914 template<class Type>
915 tmp<fvMatrix<Type>> operator-
916 (
917  const tmp<fvMatrix<Type>>&,
918  const tmp<GeometricField<Type, fvPatchField, volMesh>>&
919 );
920 
921 template<class Type>
922 tmp<fvMatrix<Type>> operator-
923 (
924  const DimensionedField<Type, volMesh>&,
925  const fvMatrix<Type>&
926 );
927 
928 template<class Type>
929 tmp<fvMatrix<Type>> operator-
930 (
931  const tmp<DimensionedField<Type, volMesh>>&,
932  const fvMatrix<Type>&
933 );
934 
935 template<class Type>
936 tmp<fvMatrix<Type>> operator-
937 (
938  const tmp<GeometricField<Type, fvPatchField, volMesh>>&,
939  const fvMatrix<Type>&
940 );
941 
942 template<class Type>
943 tmp<fvMatrix<Type>> operator-
944 (
945  const DimensionedField<Type, volMesh>&,
946  const tmp<fvMatrix<Type>>&
947 );
948 
949 template<class Type>
950 tmp<fvMatrix<Type>> operator-
951 (
952  const tmp<DimensionedField<Type, volMesh>>&,
953  const tmp<fvMatrix<Type>>&
954 );
955 
956 template<class Type>
957 tmp<fvMatrix<Type>> operator-
958 (
959  const tmp<GeometricField<Type, fvPatchField, volMesh>>&,
960  const tmp<fvMatrix<Type>>&
961 );
962 
963 
964 template<class Type>
965 tmp<fvMatrix<Type>> operator-
966 (
967  const fvMatrix<Type>&,
968  const dimensioned<Type>&
969 );
970 
971 template<class Type>
972 tmp<fvMatrix<Type>> operator-
973 (
974  const tmp<fvMatrix<Type>>&,
975  const dimensioned<Type>&
976 );
977 
978 template<class Type>
979 tmp<fvMatrix<Type>> operator-
980 (
981  const dimensioned<Type>&,
982  const fvMatrix<Type>&
983 );
984 
985 template<class Type>
986 tmp<fvMatrix<Type>> operator-
987 (
988  const dimensioned<Type>&,
989  const tmp<fvMatrix<Type>>&
990 );
991 
992 
993 template<class Type>
994 tmp<fvMatrix<Type>> operator*
995 (
997  const fvMatrix<Type>&
998 );
999 
1000 template<class Type>
1001 tmp<fvMatrix<Type>> operator*
1002 (
1004  const fvMatrix<Type>&
1005 );
1006 
1007 template<class Type>
1008 tmp<fvMatrix<Type>> operator*
1009 (
1010  const tmp<volScalarField>&,
1011  const fvMatrix<Type>&
1012 );
1013 
1014 template<class Type>
1015 tmp<fvMatrix<Type>> operator*
1016 (
1017  const volScalarField::Internal&,
1018  const tmp<fvMatrix<Type>>&
1019 );
1020 
1021 template<class Type>
1022 tmp<fvMatrix<Type>> operator*
1023 (
1024  const tmp<volScalarField::Internal>&,
1025  const tmp<fvMatrix<Type>>&
1026 );
1027 
1028 template<class Type>
1029 tmp<fvMatrix<Type>> operator*
1030 (
1031  const tmp<volScalarField>&,
1032  const tmp<fvMatrix<Type>>&
1033 );
1034 
1035 
1036 template<class Type>
1037 tmp<fvMatrix<Type>> operator*
1038 (
1039  const dimensioned<scalar>&,
1040  const fvMatrix<Type>&
1041 );
1042 
1043 template<class Type>
1044 tmp<fvMatrix<Type>> operator*
1045 (
1046  const dimensioned<scalar>&,
1047  const tmp<fvMatrix<Type>>&
1048 );
1049 
1050 
1051 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1052 
1053 } // End namespace Foam
1054 
1055 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1056 
1057 #ifdef NoRepository
1058  #include "fvMatrix.C"
1059 #endif
1060 
1061 // Specialisation for scalars
1062 #include "fvScalarMatrix.H"
1063 
1064 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1065 
1066 #endif
1067 
1068 // ************************************************************************* //
Foam::surfaceFields.
void addToInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Add patch contribution to internal field.
Definition: fvMatrix.C:38
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
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
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
Definition: fvMatrix.C:111
Reference counter for various OpenFOAM components.
Definition: refCount.H:49
void relax()
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:672
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
void negate()
Definition: fvMatrix.C:994
FieldField< Field, Type > & internalCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:304
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:282
uint8_t direction
Definition: direction.H:45
void addCmptAvBoundaryDiag(scalarField &diag) const
Definition: fvMatrix.C:129
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:506
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
Definition: fvMatrix.C:760
void operator-=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1043
Solver class returned by the solver function.
Definition: fvMatrix.H:220
Generic GeometricField class.
Generic dimensioned Type class.
tmp< volScalarField > H1() const
Return H(1)
Definition: fvMatrix.C:822
void subtractFromInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Subtract patch contribution from internal field.
Definition: fvMatrix.C:75
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
Generic field type.
Definition: FieldField.H:51
tmp< scalarField > D() const
Return the matrix scalar diagonal.
Definition: fvMatrix.C:705
void operator=(const fvMatrix< Type > &)
Definition: fvMatrix.C:950
SolverPerformance< Type > solveCoupled(const dictionary &)
Solve coupled returning the solution statistics.
surfaceTypeFieldPtr & faceFluxCorrectionPtr()
Return pointer to face-flux non-orthogonal correction field.
Definition: fvMatrix.H:322
fvSolver(fvMatrix< Type > &fvMat, autoPtr< lduMatrix::solver > sol)
Definition: fvMatrix.H:230
Dimension set for the base types.
Definition: dimensionSet.H:120
SolverPerformance< Type > solve()
Solve returning the solution statistics.
const cellShapeList & cells
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics...
void setValues(const labelUList &cells, const UList< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:484
FieldField< Field, Type > & boundaryCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:311
virtual ~fvMatrix()
Destructor.
Definition: fvMatrix.C:465
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
tmp< fvMatrix< Type > > clone() const
Clone.
Definition: fvMatrix.C:453
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Field< Type > & source()
Definition: fvMatrix.H:292
void operator*=(const volScalarField::Internal &)
Definition: fvMatrix.C:1178
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void operator+=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1009
ClassName("fvMatrix")
tmp< Field< Type > > DD() const
Return the matrix Type diagonal.
Definition: fvMatrix.C:714
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
Definition: lduMatrix.H:79
GeometricField< Type, fvsPatchField, surfaceMesh > * surfaceTypeFieldPtr
Declare return type of the faceFluxCorrectionPtr() function.
Definition: fvMatrix.H:319
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:738
label patchi
A scalar instance of fvMatrix.
Macro definitions for declaring ClassName(), NamespaceName(), etc.
A List with indirect addressing.
Definition: fvMatrix.H:106
autoPtr< fvSolver > solver()
Construct and return the solver.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void setValuesFromList(const labelUList &cells, const ListType< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:178
void checkMethod(const fvMatrix< Type > &, const fvMatrix< Type > &, const char *)
Definition: fvMatrix.C:1251
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:49
void addBoundarySource(Field< Type > &source, const bool couples=true) const
Definition: fvMatrix.C:145
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
void setComponentReference(const label patchi, const label facei, const direction cmpt, const scalar value)
Set reference level for a component of the solution.
Definition: fvMatrixSolve.C:34
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:863
SolverPerformance< Type > solveSegregated(const dictionary &)
Solve segregated returning the solution statistics.
void boundaryManipulate(typename GeometricField< Type, fvPatchField, volMesh >::Boundary &values)
Manipulate based on a boundary field.
Definition: fvMatrix.C:692
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
scalarField & diag()
Definition: lduMatrix.C:186
tmp< Field< Type > > residual() const
Return the matrix residual.
fvMatrix(const GeometricField< Type, fvPatchField, volMesh > &, const dimensionSet &)
Construct given a field to solve for.
Definition: fvMatrix.C:268
Namespace for OpenFOAM.
const dimensionSet & dimensions() const
Definition: fvMatrix.H:287