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-2023 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<VolField<Type>> operator&
64 (
65  const fvMatrix<Type>&,
67 );
68 
69 template<class Type>
70 tmp<VolField<Type>> operator&
71 (
72  const fvMatrix<Type>&,
74 );
75 
76 template<class Type>
77 tmp<VolField<Type>> operator&
78 (
79  const fvMatrix<Type>&,
80  const tmp<VolField<Type>>&
81 );
82 
83 template<class Type>
84 tmp<VolField<Type>> operator&
85 (
86  const tmp<fvMatrix<Type>>&,
88 );
89 
90 template<class Type>
91 tmp<VolField<Type>> operator&
92 (
93  const tmp<fvMatrix<Type>>&,
95 );
96 
97 template<class Type>
98 tmp<VolField<Type>> operator&
99 (
100  const tmp<fvMatrix<Type>>&,
101  const tmp<VolField<Type>>&
102 );
103 
104 template<class Type>
106 
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 VolField<Type>
123  // Converted into a non-const reference at the point of solution.
124  const VolField<Type>& psi_;
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
142  mutable SurfaceField<Type>* faceFluxCorrectionPtr_;
143 
144 
145 protected:
146 
147  //- Declare friendship with the fvSolver class
148  friend class fvSolver;
149 
150  // Protected Member Functions
151 
152  //- Add patch contribution to internal field
153  template<class Type2>
154  void addToInternalField
155  (
156  const labelUList& addr,
157  const Field<Type2>& pf,
158  Field<Type2>& intf
159  ) const;
160 
161  template<class Type2>
162  void addToInternalField
163  (
164  const labelUList& addr,
165  const tmp<Field<Type2>>& tpf,
166  Field<Type2>& intf
167  ) const;
168 
169  //- Subtract patch contribution from internal field
170  template<class Type2>
172  (
173  const labelUList& addr,
174  const Field<Type2>& pf,
175  Field<Type2>& intf
176  ) const;
177 
178  template<class Type2>
180  (
181  const labelUList& addr,
182  const tmp<Field<Type2>>& tpf,
183  Field<Type2>& intf
184  ) const;
185 
186 
187  // Matrix completion functionality
188 
189  void addBoundaryDiag
190  (
191  scalarField& diag,
192  const direction cmpt
193  ) const;
194 
196 
197  void addBoundarySource
198  (
200  const bool couples=true
201  ) const;
202 
203  // Matrix manipulation functionality
204 
205  //- Set solution in the given cell to the specified value
206  void setValue
207  (
208  const label celli,
209  const Type& value
210  );
211 
212  //- Set solution in the given cell to the specified value
213  void setValue
214  (
215  const label celli,
216  const Type& value,
217  const scalar fraction,
218  const scalarField& ddtDiag
219  );
220 
221 
222 public:
223 
224  //- Solver class returned by the solver function
225  // used for systems in which it is useful to cache the solver for reuse
226  // e.g. if the solver is potentially expensive to construct (AMG) and can
227  // be used several times (PISO)
228  class fvSolver
229  {
230  fvMatrix<Type>& fvMat_;
231 
233 
234  public:
235 
236  // Constructors
237 
239  :
240  fvMat_(fvMat),
241  solver_(sol)
242  {}
243 
244 
245  // Member Functions
246 
247  //- Solve returning the solution statistics.
248  // Use the given solver controls
250 
251  //- Solve returning the solution statistics.
252  // Solver controls read from fvSolution
254  };
255 
256 
257  ClassName("fvMatrix");
258 
259 
260  // Constructors
261 
262  //- Construct given a field to solve for
263  fvMatrix
264  (
265  const VolField<Type>&,
266  const dimensionSet&
267  );
268 
269  //- Copy constructor
270  fvMatrix(const fvMatrix<Type>&);
271 
272  //- Copy constructor of tmp<fvMatrix<Type>> deleting argument
273  fvMatrix(const tmp<fvMatrix<Type>>&);
274 
275  //- Construct from Istream given field to solve for
276  fvMatrix(const VolField<Type>&, Istream&);
277 
278  //- Clone
279  tmp<fvMatrix<Type>> clone() const;
280 
281 
282  //- Destructor
283  virtual ~fvMatrix();
284 
285 
286  // Member Functions
287 
288  // Access
289 
291  {
292  return const_cast<VolField<Type>&>
293  (
294  psi_
295  );
296  }
297 
298  const VolField<Type>& psi() const
299  {
300  return psi_;
301  }
302 
303  const dimensionSet& dimensions() const
304  {
305  return dimensions_;
306  }
307 
309  {
310  return source_;
311  }
312 
313  const Field<Type>& source() const
314  {
315  return source_;
316  }
317 
318  //- fvBoundary scalar field containing pseudo-matrix coeffs
319  // for internal cells
321  {
322  return internalCoeffs_;
323  }
324 
325  //- fvBoundary scalar field containing pseudo-matrix coeffs
326  // for boundary cells
328  {
329  return boundaryCoeffs_;
330  }
331 
332  //- Return pointer to face-flux non-orthogonal correction field
334  {
335  return faceFluxCorrectionPtr_;
336  }
337 
338 
339  // Operations
340 
341  //- Set solution in given cells to the specified values
342  // and eliminate the corresponding equations from the matrix.
343  template<template<class> class ListType>
344  void setValues
345  (
346  const labelUList& cells,
347  const ListType<Type>& values
348  );
349 
350  //- Set solution in given cells to the specified values
351  // and eliminate the corresponding equations from the matrix.
352  template<template<class> class ListType>
353  void setValues
354  (
355  const labelUList& cells,
356  const ListType<Type>& values,
357  const scalarList& fractions,
358  const bool hasDdt = true
359  );
360 
361  //- Set reference level for solution
362  void setReference
363  (
364  const label celli,
365  const Type& value,
366  const bool forceReference = false
367  );
368 
369  //- Set reference level for a component of the solution
370  // on a given patch face
372  (
373  const label patchi,
374  const label facei,
375  const direction cmpt,
376  const scalar value
377  );
378 
379  //- Return the relaxation factor for this equation
380  scalar relaxationFactor() const;
381 
382  //- Relax matrix (for steady-state solution).
383  // alpha = 1 : diagonally equal
384  // alpha < 1 : diagonally dominant
385  // alpha = 0 : do nothing
386  // Note: Requires positive diagonal.
387  void relax(const scalar alpha);
388 
389  //- Relax matrix (for steady-state solution).
390  // alpha is read from controlDict
391  void relax();
392 
393  //- Manipulate based on a boundary field
394  void boundaryManipulate
395  (
396  typename VolField<Type>::
397  Boundary& values
398  );
399 
400  //- Construct and return the solver
401  // Use the given solver controls
403 
404  //- Construct and return the solver
405  // Solver controls read from fvSolution
407 
408  //- Solve segregated or coupled returning the solution statistics.
409  // Use the given solver controls
411 
412  //- Solve segregated returning the solution statistics.
413  // Use the given solver controls
415 
416  //- Solve coupled returning the solution statistics.
417  // Use the given solver controls
419 
420  //- Solve segregated or coupled returning the solution statistics.
421  // Solver controls read from fvSolution
423 
424  //- Solve returning the solution statistics.
425  // Solver controls read from fvSolution
427 
428  //- Return the matrix residual
429  tmp<Field<Type>> residual() const;
430 
431  //- Return the matrix scalar diagonal
432  tmp<scalarField> D() const;
433 
434  //- Return the matrix Type diagonal
435  tmp<Field<Type>> DD() const;
436 
437  //- Return the central coefficient
438  tmp<volScalarField> A() const;
439 
440  //- Return the explicit source
442 
443  //- Return the implicit source
445 
446  //- Return the H operation source
447  tmp<VolField<Type>> H() const;
448 
449  //- Return H(1)
450  tmp<volScalarField> H1() const;
451 
452  //- Return the face-flux field from the matrix
454  flux() const;
455 
456 
457  // Member Operators
458 
459  void operator=(const fvMatrix<Type>&);
460  void operator=(const tmp<fvMatrix<Type>>&);
461 
462  void negate();
463 
464  void operator+=(const fvMatrix<Type>&);
465  void operator+=(const tmp<fvMatrix<Type>>&);
466 
467  void operator-=(const fvMatrix<Type>&);
468  void operator-=(const tmp<fvMatrix<Type>>&);
469 
470  void operator+=
471  (
473  );
474  void operator+=
475  (
477  );
478  void operator+=
479  (
480  const tmp<VolField<Type>>&
481  );
482 
483  void operator-=
484  (
486  );
487  void operator-=
488  (
490  );
491  void operator-=
492  (
493  const tmp<VolField<Type>>&
494  );
495 
496  void operator+=(const dimensioned<Type>&);
497  void operator-=(const dimensioned<Type>&);
498 
499  void operator+=(const zero&);
500  void operator-=(const zero&);
501 
504  void operator*=(const tmp<volScalarField>&);
505  void operator*=(const dimensioned<scalar>&);
506 
509  void operator/=(const tmp<volScalarField>&);
510  void operator/=(const dimensioned<scalar>&);
511 
512 
513  // Friend operators
514 
516  operator& <Type>
517  (
518  const fvMatrix<Type>&,
520  );
521 
523  operator& <Type>
524  (
525  const fvMatrix<Type>&,
526  const tmp<VolField<Type>>&
527  );
528 
530  operator& <Type>
531  (
532  const tmp<fvMatrix<Type>>&,
534  );
535 
537  operator& <Type>
538  (
539  const tmp<fvMatrix<Type>>&,
540  const tmp<VolField<Type>>&
541  );
542 
543 
544  // Ostream operator
545 
546  friend Ostream& operator<< <Type>
547  (
548  Ostream&,
549  const fvMatrix<Type>&
550  );
551 };
552 
553 
554 // * * * * * * * * * * * * * * * Global functions * * * * * * * * * * * * * //
555 
556 template<class Type>
557 void checkMethod
558 (
559  const fvMatrix<Type>&,
560  const fvMatrix<Type>&,
561  const char*
562 );
563 
564 template<class Type>
565 void checkMethod
566 (
567  const fvMatrix<Type>&,
569  const char*
570 );
571 
572 template<class Type>
573 void checkMethod
574 (
575  const fvMatrix<Type>&,
576  const dimensioned<Type>&,
577  const char*
578 );
579 
580 
581 //- Solve returning the solution statistics given convergence tolerance
582 // Use the given solver controls
583 template<class Type>
585 
586 //- Solve returning the solution statistics given convergence tolerance,
587 // deleting temporary matrix after solution.
588 // Use the given solver controls
589 template<class Type>
591 (
592  const tmp<fvMatrix<Type>>&,
593  const word&
594 );
595 
596 //- Solve returning the solution statistics given convergence tolerance
597 // Solver controls read fvSolution
598 template<class Type>
600 
601 //- Solve returning the solution statistics given convergence tolerance,
602 // deleting temporary matrix after solution.
603 // Solver controls read fvSolution
604 template<class Type>
606 
607 //- Return the correction form of the given matrix
608 // by subtracting the matrix multiplied by the current field
609 template<class Type>
611 
612 //- Return the correction form of the given temporary matrix
613 // by subtracting the matrix multiplied by the current field
614 template<class Type>
616 
617 
618 // * * * * * * * * * * * * * * * Global operators * * * * * * * * * * * * * //
619 
620 template<class Type>
621 tmp<fvMatrix<Type>> operator==
622 (
623  const fvMatrix<Type>&,
624  const fvMatrix<Type>&
625 );
626 
627 template<class Type>
628 tmp<fvMatrix<Type>> operator==
629 (
630  const tmp<fvMatrix<Type>>&,
631  const fvMatrix<Type>&
632 );
633 
634 template<class Type>
635 tmp<fvMatrix<Type>> operator==
636 (
637  const fvMatrix<Type>&,
638  const tmp<fvMatrix<Type>>&
639 );
640 
641 template<class Type>
642 tmp<fvMatrix<Type>> operator==
643 (
644  const tmp<fvMatrix<Type>>&,
645  const tmp<fvMatrix<Type>>&
646 );
647 
648 
649 template<class Type>
650 tmp<fvMatrix<Type>> operator==
651 (
652  const fvMatrix<Type>&,
654 );
655 
656 template<class Type>
657 tmp<fvMatrix<Type>> operator==
658 (
659  const fvMatrix<Type>&,
661 );
662 
663 template<class Type>
664 tmp<fvMatrix<Type>> operator==
665 (
666  const fvMatrix<Type>&,
667  const tmp<VolField<Type>>&
668 );
669 
670 template<class Type>
671 tmp<fvMatrix<Type>> operator==
672 (
673  const tmp<fvMatrix<Type>>&,
675 );
676 
677 template<class Type>
678 tmp<fvMatrix<Type>> operator==
679 (
680  const tmp<fvMatrix<Type>>&,
682 );
683 
684 template<class Type>
685 tmp<fvMatrix<Type>> operator==
686 (
687  const tmp<fvMatrix<Type>>&,
688  const tmp<VolField<Type>>&
689 );
690 
691 template<class Type>
692 tmp<fvMatrix<Type>> operator==
693 (
694  const fvMatrix<Type>&,
695  const dimensioned<Type>&
696 );
697 
698 template<class Type>
699 tmp<fvMatrix<Type>> operator==
700 (
701  const tmp<fvMatrix<Type>>&,
702  const dimensioned<Type>&
703 );
704 
705 
706 template<class Type>
707 tmp<fvMatrix<Type>> operator==
708 (
709  const fvMatrix<Type>&,
710  const zero&
711 );
712 
713 template<class Type>
714 tmp<fvMatrix<Type>> operator==
715 (
716  const tmp<fvMatrix<Type>>&,
717  const zero&
718 );
719 
720 
721 template<class Type>
722 tmp<fvMatrix<Type>> operator-
723 (
724  const fvMatrix<Type>&
725 );
726 
727 template<class Type>
728 tmp<fvMatrix<Type>> operator-
729 (
730  const tmp<fvMatrix<Type>>&
731 );
732 
733 
734 template<class Type>
735 tmp<fvMatrix<Type>> operator+
736 (
737  const fvMatrix<Type>&,
738  const fvMatrix<Type>&
739 );
740 
741 template<class Type>
742 tmp<fvMatrix<Type>> operator+
743 (
744  const tmp<fvMatrix<Type>>&,
745  const fvMatrix<Type>&
746 );
747 
748 template<class Type>
749 tmp<fvMatrix<Type>> operator+
750 (
751  const fvMatrix<Type>&,
752  const tmp<fvMatrix<Type>>&
753 );
754 
755 template<class Type>
756 tmp<fvMatrix<Type>> operator+
757 (
758  const tmp<fvMatrix<Type>>&,
759  const tmp<fvMatrix<Type>>&
760 );
761 
762 
763 template<class Type>
764 tmp<fvMatrix<Type>> operator+
765 (
766  const fvMatrix<Type>&,
768 );
769 
770 template<class Type>
771 tmp<fvMatrix<Type>> operator+
772 (
773  const fvMatrix<Type>&,
775 );
776 
777 template<class Type>
778 tmp<fvMatrix<Type>> operator+
779 (
780  const fvMatrix<Type>&,
781  const tmp<VolField<Type>>&
782 );
783 
784 template<class Type>
785 tmp<fvMatrix<Type>> operator+
786 (
787  const tmp<fvMatrix<Type>>&,
789 );
790 
791 template<class Type>
792 tmp<fvMatrix<Type>> operator+
793 (
794  const tmp<fvMatrix<Type>>&,
796 );
797 
798 template<class Type>
799 tmp<fvMatrix<Type>> operator+
800 (
801  const tmp<fvMatrix<Type>>&,
802  const tmp<VolField<Type>>&
803 );
804 
805 template<class Type>
806 tmp<fvMatrix<Type>> operator+
807 (
809  const fvMatrix<Type>&
810 );
811 
812 template<class Type>
813 tmp<fvMatrix<Type>> operator+
814 (
816  const fvMatrix<Type>&
817 );
818 
819 template<class Type>
820 tmp<fvMatrix<Type>> operator+
821 (
822  const tmp<VolField<Type>>&,
823  const fvMatrix<Type>&
824 );
825 
826 template<class Type>
827 tmp<fvMatrix<Type>> operator+
828 (
830  const tmp<fvMatrix<Type>>&
831 );
832 
833 template<class Type>
834 tmp<fvMatrix<Type>> operator+
835 (
837  const tmp<fvMatrix<Type>>&
838 );
839 
840 template<class Type>
841 tmp<fvMatrix<Type>> operator+
842 (
843  const tmp<VolField<Type>>&,
844  const tmp<fvMatrix<Type>>&
845 );
846 
847 
848 template<class Type>
849 tmp<fvMatrix<Type>> operator+
850 (
851  const fvMatrix<Type>&,
852  const dimensioned<Type>&
853 );
854 
855 template<class Type>
856 tmp<fvMatrix<Type>> operator+
857 (
858  const tmp<fvMatrix<Type>>&,
859  const dimensioned<Type>&
860 );
861 
862 template<class Type>
863 tmp<fvMatrix<Type>> operator+
864 (
865  const dimensioned<Type>&,
866  const fvMatrix<Type>&
867 );
868 
869 template<class Type>
870 tmp<fvMatrix<Type>> operator+
871 (
872  const dimensioned<Type>&,
873  const tmp<fvMatrix<Type>>&
874 );
875 
876 
877 template<class Type>
878 tmp<fvMatrix<Type>> operator-
879 (
880  const fvMatrix<Type>&,
881  const fvMatrix<Type>&
882 );
883 
884 template<class Type>
885 tmp<fvMatrix<Type>> operator-
886 (
887  const tmp<fvMatrix<Type>>&,
888  const fvMatrix<Type>&
889 );
890 
891 template<class Type>
892 tmp<fvMatrix<Type>> operator-
893 (
894  const fvMatrix<Type>&,
895  const tmp<fvMatrix<Type>>&
896 );
897 
898 template<class Type>
899 tmp<fvMatrix<Type>> operator-
900 (
901  const tmp<fvMatrix<Type>>&,
902  const tmp<fvMatrix<Type>>&
903 );
904 
905 
906 template<class Type>
907 tmp<fvMatrix<Type>> operator-
908 (
909  const fvMatrix<Type>&,
911 );
912 
913 template<class Type>
914 tmp<fvMatrix<Type>> operator-
915 (
916  const fvMatrix<Type>&,
918 );
919 
920 template<class Type>
921 tmp<fvMatrix<Type>> operator-
922 (
923  const fvMatrix<Type>&,
924  const tmp<VolField<Type>>&
925 );
926 
927 template<class Type>
928 tmp<fvMatrix<Type>> operator-
929 (
930  const tmp<fvMatrix<Type>>&,
932 );
933 
934 template<class Type>
935 tmp<fvMatrix<Type>> operator-
936 (
937  const tmp<fvMatrix<Type>>&,
939 );
940 
941 template<class Type>
942 tmp<fvMatrix<Type>> operator-
943 (
944  const tmp<fvMatrix<Type>>&,
945  const tmp<VolField<Type>>&
946 );
947 
948 template<class Type>
949 tmp<fvMatrix<Type>> operator-
950 (
952  const fvMatrix<Type>&
953 );
954 
955 template<class Type>
956 tmp<fvMatrix<Type>> operator-
957 (
959  const fvMatrix<Type>&
960 );
961 
962 template<class Type>
963 tmp<fvMatrix<Type>> operator-
964 (
965  const tmp<VolField<Type>>&,
966  const fvMatrix<Type>&
967 );
968 
969 template<class Type>
970 tmp<fvMatrix<Type>> operator-
971 (
973  const tmp<fvMatrix<Type>>&
974 );
975 
976 template<class Type>
977 tmp<fvMatrix<Type>> operator-
978 (
980  const tmp<fvMatrix<Type>>&
981 );
982 
983 template<class Type>
984 tmp<fvMatrix<Type>> operator-
985 (
986  const tmp<VolField<Type>>&,
987  const tmp<fvMatrix<Type>>&
988 );
989 
990 
991 template<class Type>
992 tmp<fvMatrix<Type>> operator-
993 (
994  const fvMatrix<Type>&,
995  const dimensioned<Type>&
996 );
997 
998 template<class Type>
999 tmp<fvMatrix<Type>> operator-
1000 (
1001  const tmp<fvMatrix<Type>>&,
1002  const dimensioned<Type>&
1003 );
1004 
1005 template<class Type>
1006 tmp<fvMatrix<Type>> operator-
1007 (
1008  const dimensioned<Type>&,
1009  const fvMatrix<Type>&
1010 );
1011 
1012 template<class Type>
1013 tmp<fvMatrix<Type>> operator-
1014 (
1015  const dimensioned<Type>&,
1016  const tmp<fvMatrix<Type>>&
1017 );
1018 
1019 
1020 template<class Type>
1021 tmp<fvMatrix<Type>> operator*
1022 (
1023  const volScalarField::Internal&,
1024  const fvMatrix<Type>&
1025 );
1026 
1027 template<class Type>
1028 tmp<fvMatrix<Type>> operator*
1029 (
1031  const fvMatrix<Type>&
1032 );
1033 
1034 template<class Type>
1035 tmp<fvMatrix<Type>> operator*
1036 (
1037  const tmp<volScalarField>&,
1038  const fvMatrix<Type>&
1039 );
1040 
1041 template<class Type>
1042 tmp<fvMatrix<Type>> operator*
1043 (
1044  const volScalarField::Internal&,
1045  const tmp<fvMatrix<Type>>&
1046 );
1047 
1048 template<class Type>
1049 tmp<fvMatrix<Type>> operator*
1050 (
1052  const tmp<fvMatrix<Type>>&
1053 );
1054 
1055 template<class Type>
1056 tmp<fvMatrix<Type>> operator*
1057 (
1058  const tmp<volScalarField>&,
1059  const tmp<fvMatrix<Type>>&
1060 );
1061 
1062 
1063 template<class Type>
1064 tmp<fvMatrix<Type>> operator*
1065 (
1066  const dimensioned<scalar>&,
1067  const fvMatrix<Type>&
1068 );
1069 
1070 template<class Type>
1071 tmp<fvMatrix<Type>> operator*
1072 (
1073  const dimensioned<scalar>&,
1074  const tmp<fvMatrix<Type>>&
1075 );
1076 
1077 template<class Type>
1078 tmp<fvMatrix<Type>> operator*
1079 (
1080  const volScalarField::Internal&,
1081  const fvMatrix<Type>&
1082 );
1083 
1084 template<class Type>
1085 tmp<fvMatrix<Type>> operator*
1086 (
1088  const fvMatrix<Type>&
1089 );
1090 
1091 template<class Type>
1092 tmp<fvMatrix<Type>> operator*
1093 (
1094  const tmp<volScalarField>&,
1095  const fvMatrix<Type>&
1096 );
1097 
1098 template<class Type>
1099 tmp<fvMatrix<Type>> operator*
1100 (
1101  const volScalarField::Internal&,
1102  const tmp<fvMatrix<Type>>&
1103 );
1104 
1105 template<class Type>
1106 tmp<fvMatrix<Type>> operator*
1107 (
1109  const tmp<fvMatrix<Type>>&
1110 );
1111 
1112 template<class Type>
1113 tmp<fvMatrix<Type>> operator*
1114 (
1115  const tmp<volScalarField>&,
1116  const tmp<fvMatrix<Type>>&
1117 );
1118 
1119 
1120 template<class Type>
1121 tmp<fvMatrix<Type>> operator*
1122 (
1123  const dimensioned<scalar>&,
1124  const fvMatrix<Type>&
1125 );
1126 
1127 template<class Type>
1128 tmp<fvMatrix<Type>> operator*
1129 (
1130  const dimensioned<scalar>&,
1131  const tmp<fvMatrix<Type>>&
1132 );
1133 
1134 
1135 template<class Type>
1136 tmp<fvMatrix<Type>> operator/
1137 (
1138  const fvMatrix<Type>&,
1140 );
1141 
1142 template<class Type>
1143 tmp<fvMatrix<Type>> operator/
1144 (
1145  const fvMatrix<Type>&,
1147 );
1148 
1149 template<class Type>
1150 tmp<fvMatrix<Type>> operator/
1151 (
1152  const fvMatrix<Type>&,
1153  const tmp<volScalarField>&
1154 );
1155 
1156 template<class Type>
1157 tmp<fvMatrix<Type>> operator/
1158 (
1159  const tmp<fvMatrix<Type>>&,
1161 );
1162 
1163 template<class Type>
1164 tmp<fvMatrix<Type>> operator/
1165 (
1166  const tmp<fvMatrix<Type>>&,
1168 );
1169 
1170 template<class Type>
1171 tmp<fvMatrix<Type>> operator/
1172 (
1173  const tmp<fvMatrix<Type>>&,
1174  const tmp<volScalarField>&
1175 );
1176 
1177 
1178 template<class Type>
1179 tmp<fvMatrix<Type>> operator/
1180 (
1181  const fvMatrix<Type>&,
1182  const dimensioned<scalar>&
1183 );
1184 
1185 template<class Type>
1186 tmp<fvMatrix<Type>> operator/
1187 (
1188  const tmp<fvMatrix<Type>>&,
1189  const dimensioned<scalar>&
1190 );
1191 
1192 
1193 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1194 
1195 } // End namespace Foam
1196 
1197 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1198 
1199 #ifdef NoRepository
1200  #include "fvMatrix.C"
1201 #endif
1202 
1203 // Specialisation for scalars
1204 #include "fvScalarMatrix.H"
1205 
1206 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1207 
1208 #endif
1209 
1210 // ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic field type.
Definition: FieldField.H:77
Pre-declare SubField and related Field type.
Definition: Field.H:83
Generic GeometricField class.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics.
A List with indirect addressing.
Definition: UIndirectList.H:60
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Dimension set for the base types.
Definition: dimensionSet.H:125
Generic dimensioned Type class.
Solver class returned by the solver function.
Definition: fvMatrix.H:228
SolverPerformance< Type > solve()
Solve returning the solution statistics.
fvSolver(fvMatrix< Type > &fvMat, autoPtr< lduMatrix::solver > sol)
Definition: fvMatrix.H:237
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
void operator-=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1137
autoPtr< fvSolver > solver()
Construct and return the solver.
tmp< Field< Type > > DD() const
Return the matrix Type diagonal.
Definition: fvMatrix.C:784
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:808
void addToInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Add patch contribution to internal field.
Definition: fvMatrix.C:40
void setValue(const label celli, const Type &value)
Set solution in the given cell to the specified value.
Definition: fvMatrix.C:179
void boundaryManipulate(typename VolField< Type >::Boundary &values)
Manipulate based on a boundary field.
Definition: fvMatrix.C:762
SurfaceField< Type > *& faceFluxCorrectionPtr()
Return pointer to face-flux non-orthogonal correction field.
Definition: fvMatrix.H:332
void relax()
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:754
SolverPerformance< Type > solve()
Solve returning the solution statistics.
tmp< volScalarField > H1() const
Return H(1)
Definition: fvMatrix.C:923
tmp< VolField< Type > > H() const
Return the H operation source.
Definition: fvMatrix.C:868
void operator+=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1103
Field< Type > & source()
Definition: fvMatrix.H:307
virtual ~fvMatrix()
Destructor.
Definition: fvMatrix.C:482
VolField< Type > & psi()
Definition: fvMatrix.H:289
tmp< SurfaceField< Type > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:964
FieldField< Field, Type > & internalCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:319
FieldField< Field, Type > & boundaryCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:326
void operator*=(const volScalarField::Internal &)
Definition: fvMatrix.C:1272
void addCmptAvBoundaryDiag(scalarField &diag) const
Definition: fvMatrix.C:131
void addBoundarySource(Field< Type > &source, const bool couples=true) const
Definition: fvMatrix.C:147
fvMatrix(const VolField< Type > &, const dimensionSet &)
Construct given a field to solve for.
Definition: fvMatrix.C:285
tmp< volScalarField::Internal > Sp() const
Return the implicit source.
Definition: fvMatrix.C:847
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:563
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
void operator/=(const volScalarField::Internal &)
Definition: fvMatrix.C:1346
ClassName("fvMatrix")
SolverPerformance< Type > solveCoupled(const dictionary &)
Solve coupled returning the solution statistics.
tmp< Field< Type > > residual() const
Return the matrix residual.
scalar relaxationFactor() const
Return the relaxation factor for this equation.
Definition: fvMatrix.C:578
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
Definition: fvMatrix.C:113
void negate()
Definition: fvMatrix.C:1088
void setValues(const labelUList &cells, const ListType< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:502
tmp< fvMatrix< Type > > clone() const
Clone.
Definition: fvMatrix.C:470
SolverPerformance< Type > solveSegregated(const dictionary &)
Solve segregated returning the solution statistics.
tmp< VolInternalField< Type > > Su() const
Return the explicit source.
Definition: fvMatrix.C:829
tmp< scalarField > D() const
Return the matrix scalar diagonal.
Definition: fvMatrix.C:775
void subtractFromInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Subtract patch contribution from internal field.
Definition: fvMatrix.C:77
const dimensionSet & dimensions() const
Definition: fvMatrix.H:302
void operator=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1044
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
Definition: lduMatrix.H:80
scalarField & diag()
Definition: lduMatrix.C:186
Reference counter for various OpenFOAM components.
Definition: refCount.H:50
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:50
Macro definitions for declaring ClassName(), NamespaceName(), etc.
A scalar instance of fvMatrix.
label patchi
const cellShapeList & cells
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Namespace for OpenFOAM.
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
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
void checkMethod(const fvMatrix< Type > &, const fvMatrix< Type > &, const char *)
Definition: fvMatrix.C:1422
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
Ostream & operator<<(Ostream &os, const fvConstraints &constraints)
uint8_t direction
Definition: direction.H:45
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
dictionary fractions(initialConditions.subDict("fractions"))
Foam::surfaceFields.