fvMatrix.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2016-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::fvMatrix
29
30Description
31 A special matrix type and solver, designed for finite volume
32 solutions of scalar equations.
33 Face addressing is used to make all matrix assembly
34 and solution loops vectorise.
35
36SourceFiles
37 fvMatrix.C
38 fvMatrixSolve.C
39 fvScalarMatrix.C
40
41\*---------------------------------------------------------------------------*/
42
43#ifndef Foam_fvMatrix_H
44#define Foam_fvMatrix_H
45
46#include "volFields.H"
47#include "surfaceFields.H"
48#include "lduMatrix.H"
49#include "tmp.H"
50#include "autoPtr.H"
51#include "dimensionedTypes.H"
52#include "zero.H"
53#include "className.H"
55#include "lduMesh.H"
56
57// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
58
59namespace Foam
60{
61
62// Forward Declarations
63template<class Type> class fvMatrix;
64template<class T> class UIndirectList;
65
66template<class Type>
68
69
70template<class Type>
72(
73 const fvMatrix<Type>&,
75);
76
77template<class Type>
79(
80 const fvMatrix<Type>&,
82);
83
84template<class Type>
86(
87 const fvMatrix<Type>&,
89);
90
91template<class Type>
93(
94 const tmp<fvMatrix<Type>>&,
96);
97
98template<class Type>
100(
101 const tmp<fvMatrix<Type>>&,
103);
104
105template<class Type>
107(
108 const tmp<fvMatrix<Type>>&,
110);
111
112
113/*---------------------------------------------------------------------------*\
114 Class fvMatrix Declaration
115\*---------------------------------------------------------------------------*/
116
117template<class Type>
118class fvMatrix
119:
120 public refCount,
121 public lduMatrix
122{
123public:
124
125 // Public Types
126
127 //- The geometric field type for psi
128 typedef
131
132 //- Field type for face flux (for non-orthogonal correction)
133 typedef
136
139
140
141private:
142
143 // Private Data
144
145 //- Const reference to field
146 // Converted into a non-const reference at the point of solution.
147 const psiFieldType& psi_;
148
149 //- Originating fvMatrices when assembling matrices. Empty if not used.
150 PtrList<fvMatrix<Type>> subMatrices_;
151
152 //- Is the fvMatrix using implicit formulation
153 bool useImplicit_;
154
155 //- Name of the lduAssembly
156 word lduAssemblyName_;
157
158 //- Number of fvMatrices added to this
159 label nMatrix_;
160
161 //- Dimension set
162 dimensionSet dimensions_;
163
164 //- Source term
165 Field<Type> source_;
166
167 //- Boundary scalar field containing pseudo-matrix coeffs
168 //- for internal cells
169 FieldField<Field, Type> internalCoeffs_;
170
171 //- Boundary scalar field containing pseudo-matrix coeffs
172 //- for boundary cells
173 FieldField<Field, Type> boundaryCoeffs_;
174
175 //- Face flux field for non-orthogonal correction
176 mutable faceFluxFieldType* faceFluxCorrectionPtr_;
177
178
179protected:
180
181 //- Declare friendship with the fvSolver class
182 friend class fvSolver;
183
184
185 // Protected Member Functions
186
187 //- Add patch contribution to internal field
188 template<class Type2>
190 (
191 const labelUList& addr,
192 const Field<Type2>& pf,
193 Field<Type2>& intf
194 ) const;
195
196 template<class Type2>
198 (
199 const labelUList& addr,
200 const tmp<Field<Type2>>& tpf,
201 Field<Type2>& intf
202 ) const;
203
204 //- Subtract patch contribution from internal field
205 template<class Type2>
207 (
208 const labelUList& addr,
209 const Field<Type2>& pf,
210 Field<Type2>& intf
211 ) const;
212
213 template<class Type2>
215 (
216 const labelUList& addr,
217 const tmp<Field<Type2>>& tpf,
218 Field<Type2>& intf
219 ) const;
220
221
222 // Implicit helper functions
223
224 //- Name the implicit assembly addressing
225 // \return true if assembly is used
226 bool checkImplicit(const label fieldi = 0);
227
228
229 // Matrix completion functionality
230
231 void addBoundaryDiag
232 (
234 const direction cmpt
235 ) const;
236
238
240 (
242 const bool couples=true
243 ) const;
244
245
246 // Matrix manipulation functionality
247
248 //- Set solution in given cells to the specified values
249 template<template<class> class ListType>
251 (
252 const labelUList& cellLabels,
253 const ListType<Type>& values
254 );
255
256
257public:
258
259 //- Solver class returned by the solver function
260 //- used for systems in which it is useful to cache the solver for reuse.
261 // E.g. if the solver is potentially expensive to construct (AMG) and can
262 // be used several times (PISO)
263 class fvSolver
264 {
265 fvMatrix<Type>& fvMat_;
266
268
269 public:
270
271 // Constructors
274 :
275 fvMat_(fvMat),
276 solver_(std::move(sol))
277 {}
278
279
280 // Member Functions
281
282 //- Solve returning the solution statistics.
283 // Solver controls read from dictionary
284 SolverPerformance<Type> solve(const dictionary& solverControls);
285
286 //- Solve returning the solution statistics.
287 // Solver controls read from fvSolution
289 };
290
291
292 // Runtime information
293 ClassName("fvMatrix");
294
295
296 // Constructors
297
298 //- Construct given a field to solve for
300 (
302 const dimensionSet& ds
303 );
304
305 //- Copy construct
306 fvMatrix(const fvMatrix<Type>&);
307
308 //- Copy/move construct from tmp<fvMatrix<Type>>
309 fvMatrix(const tmp<fvMatrix<Type>>&);
310
311 //- Deprecated(2022-05) - construct with dimensionSet instead
312 // \deprecated(2022-05) - construct with dimensionSet instead
313 FOAM_DEPRECATED_FOR(2022-05, "Construct with dimensionSet")
315 (
316 const GeometricField<Type, fvPatchField, volMesh>& psi,
317 Istream& is
318 )
319 :
320 fvMatrix<Type>(psi, dimensionSet(is))
321 {}
322
323 //- Construct and return a clone
325 {
326 return tmp<fvMatrix<Type>>::New(*this);
327 }
328
329
330 //- Destructor
331 virtual ~fvMatrix();
332
333
334 // Member Functions
335
336 // Access
337
338 // Coupling
340 label nMatrices() const
341 {
342 return (nMatrix_ == 0 ? 1 : nMatrix_);
343 }
345 const fvMatrix<Type>& matrix(const label i) const
346 {
347 return (nMatrix_ == 0 ? *this : subMatrices_[i]);
348 }
350 fvMatrix<Type>& matrix(const label i)
351 {
352 return (nMatrix_ == 0 ? *this : subMatrices_[i]);
353 }
355 label globalPatchID
356 (
357 const label fieldi,
358 const label patchi
359 ) const
360 {
361 if (!lduMeshPtr())
362 {
363 return patchi;
364 }
365 else
366 {
367 return lduMeshPtr()->patchMap()[fieldi][patchi];
368 }
369 }
370
371 //- Transfer lower, upper, diag and source to this fvMatrix
373
374 //- Create or update ldu assembly
376
377 //- Access to lduPrimitiveMeshAssembly
379
380 //- Const Access to lduPrimitiveMeshAssembly
382
383 //- Manipulate matrix
384 void manipulateMatrix(direction cmp);
385
386 //- Manipulate boundary/internal coeffs for coupling
388
389 //- Set interfaces
390 void setInterfaces
391 (
393 PtrDynList<lduInterfaceField>& newInterfaces
394 );
395
396 //- Add internal and boundary contribution to local patches
398 (
399 label fieldi,
400 const FieldField<Field, Type>& fluxContrib,
402 bool internal
403 ) const;
404
405 //- Return optional lduAdressing
407 {
408 return *lduMeshPtr();
409 }
410
411 //- Return psi
413 (
414 const label i = 0
415 ) const
416 {
417 return
418 (
419 (i == 0 && nMatrix_ == 0) ? psi_ : matrix(i).psi()
420 );
421 }
422
425 (
426 const label i = 0
427 )
428 {
429 return
430 (
431 (i == 0 && nMatrix_ == 0)
432 ? const_cast
433 <
435 >(psi_)
436 : const_cast
437 <
439 >
440 (
441 matrix(i).psi()
442 )
443 );
444 }
445
446 //- Clear multiple fvMatrices
447 void clear()
448 {
449 subMatrices_.clear();
450 nMatrix_ = 0;
451 }
452
454 const dimensionSet& dimensions() const noexcept
455 {
456 return dimensions_;
457 }
460 {
461 return source_;
462 }
464 const Field<Type>& source() const noexcept
465 {
466 return source_;
467 }
468
469 //- fvBoundary scalar field containing pseudo-matrix coeffs
470 //- for internal cells
472 {
473 return internalCoeffs_;
474 }
475
476 //- fvBoundary scalar field containing pseudo-matrix coeffs
477 //- for internal cells
479 {
480 return internalCoeffs_;
481 }
482
483 //- fvBoundary scalar field containing pseudo-matrix coeffs
484 //- for boundary cells
486 {
487 return boundaryCoeffs_;
488 }
489
490 //- fvBoundary scalar field containing pseudo-matrix coeffs
491 //- for boundary cells
493 {
494 return boundaryCoeffs_;
495 }
496
497 //- Declare return type of the faceFluxCorrectionPtr() function
500
501 //- Return pointer to face-flux non-orthogonal correction field
503 {
504 return faceFluxCorrectionPtr_;
505 }
506
507 //- True if face-flux non-orthogonal correction field exists
509 {
510 return bool(faceFluxCorrectionPtr_);
511 }
512
513
514 // Operations
515
516 //- Set solution in given cells to the specified value
517 //- and eliminate the corresponding equations from the matrix.
518 void setValues
519 (
520 const labelUList& cellLabels,
521 const Type& value
522 );
523
524 //- Set solution in given cells to the specified values
525 //- and eliminate the corresponding equations from the matrix.
526 void setValues
527 (
528 const labelUList& cellLabels,
529 const UList<Type>& values
530 );
531
532 //- Set solution in given cells to the specified values
533 //- and eliminate the corresponding equations from the matrix.
534 void setValues
535 (
536 const labelUList& cellLabels,
537 const UIndirectList<Type>& values
538 );
539
540 //- Set reference level for solution
541 void setReference
542 (
543 const label celli,
544 const Type& value,
545 const bool forceReference = false
546 );
547
548 //- Set reference level for solution
549 void setReferences
550 (
551 const labelUList& cellLabels,
552 const Type& value,
553 const bool forceReference = false
554 );
555
556 //- Set reference level for solution
557 void setReferences
558 (
559 const labelUList& cellLabels,
560 const UList<Type>& values,
561 const bool forceReference = false
562 );
563
564 //- Set reference level for a component of the solution
565 //- on a given patch face
567 (
568 const label patchi,
569 const label facei,
570 const direction cmpt,
571 const scalar value
572 );
573
574 //- Add fvMatrix
576
577 //- Relax matrix (for steady-state solution).
578 // alpha = 1 : diagonally equal
579 // alpha < 1 : diagonally dominant
580 // alpha = 0 : do nothing
581 // Note: Requires positive diagonal.
582 void relax(const scalar alpha);
583
584 //- Relax matrix (for steady-state solution).
585 // alpha is read from controlDict
586 void relax();
587
588 //- Manipulate based on a boundary field
590 (
592 Boundary& values
593 );
594
595 //- Construct and return the solver
596 // Use the given solver controls
598
599 //- Construct and return the solver
600 // Solver controls read from fvSolution
602
603 //- Solve segregated or coupled returning the solution statistics.
604 // Use the given solver controls
606
607 //- Solve segregated returning the solution statistics.
608 // Use the given solver controls
610
611 //- Solve coupled returning the solution statistics.
612 // Use the given solver controls
614
615 //- Solve returning the solution statistics.
616 // Use the given solver controls
618
619 //- Solve returning the solution statistics.
620 // Solver controls read from fvSolution
622
623 //- Return the matrix residual
624 tmp<Field<Type>> residual() const;
625
626 //- Return the matrix scalar diagonal
627 tmp<scalarField> D() const;
628
629 //- Return the matrix Type diagonal
630 tmp<Field<Type>> DD() const;
631
632 //- Return the central coefficient
633 tmp<volScalarField> A() const;
634
635 //- Return the H operation source
637
638 //- Return H(1)
639 tmp<volScalarField> H1() const;
640
641 //- Return the face-flux field from the matrix
643 flux() const;
644
645 //- Return the solver dictionary taking into account finalIteration
646 const dictionary& solverDict() const;
647
648
649 // Member Operators
650
651 void operator=(const fvMatrix<Type>&);
652 void operator=(const tmp<fvMatrix<Type>>&);
653
654 //- Inplace negate
655 void negate();
656
657 void operator+=(const fvMatrix<Type>&);
658 void operator+=(const tmp<fvMatrix<Type>>&);
659
660 void operator-=(const fvMatrix<Type>&);
661 void operator-=(const tmp<fvMatrix<Type>>&);
662
665 void operator+=
666 (
668 );
669
672 void operator-=
673 (
675 );
676
677 void operator+=(const dimensioned<Type>&);
678 void operator-=(const dimensioned<Type>&);
679
680 void operator+=(const Foam::zero);
681 void operator-=(const Foam::zero);
682
685 void operator*=(const tmp<volScalarField>&);
686
687 void operator*=(const dimensioned<scalar>&);
688
689
690 // Friend Operators
691
693 operator& <Type>
694 (
695 const fvMatrix<Type>&,
697 );
698
700 operator& <Type>
701 (
702 const fvMatrix<Type>&,
704 );
705
707 operator& <Type>
708 (
709 const tmp<fvMatrix<Type>>&,
711 );
712
714 operator& <Type>
715 (
716 const tmp<fvMatrix<Type>>&,
718 );
719
720
721 // Ostream Operator
723 friend Ostream& operator<< <Type>
724 (
725 Ostream&,
726 const fvMatrix<Type>&
727 );
728};
729
730
731// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
732
733template<class Type>
734void checkMethod
735(
736 const fvMatrix<Type>&,
737 const fvMatrix<Type>&,
738 const char*
739);
740
741template<class Type>
742void checkMethod
743(
744 const fvMatrix<Type>&,
746 const char*
747);
748
749template<class Type>
750void checkMethod
751(
752 const fvMatrix<Type>&,
753 const dimensioned<Type>&,
754 const char*
755);
756
757
758//- Solve returning the solution statistics given convergence tolerance
759// Use the given solver controls
760template<class Type>
762
763
764//- Solve returning the solution statistics given convergence tolerance,
765//- deleting temporary matrix after solution.
766// Use the given solver controls
767template<class Type>
769(
770 const tmp<fvMatrix<Type>>&,
771 const dictionary&
772);
773
774
775//- Solve returning the solution statistics given convergence tolerance
776// Solver controls read fvSolution
777template<class Type>
779
780
781//- Solve returning the solution statistics given convergence tolerance,
782//- deleting temporary matrix after solution.
783// Solver controls read fvSolution
784template<class Type>
786
787
788//- Return the correction form of the given matrix
789//- by subtracting the matrix multiplied by the current field
790template<class Type>
792
793
794//- Return the correction form of the given temporary matrix
795//- by subtracting the matrix multiplied by the current field
796template<class Type>
798
799
800// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
801
802template<class Type>
803tmp<fvMatrix<Type>> operator==
804(
805 const fvMatrix<Type>&,
806 const fvMatrix<Type>&
807);
808
809template<class Type>
810tmp<fvMatrix<Type>> operator==
811(
812 const tmp<fvMatrix<Type>>&,
813 const fvMatrix<Type>&
814);
815
816template<class Type>
817tmp<fvMatrix<Type>> operator==
818(
819 const fvMatrix<Type>&,
820 const tmp<fvMatrix<Type>>&
821);
822
823template<class Type>
824tmp<fvMatrix<Type>> operator==
825(
826 const tmp<fvMatrix<Type>>&,
827 const tmp<fvMatrix<Type>>&
828);
829
830
831template<class Type>
832tmp<fvMatrix<Type>> operator==
833(
834 const fvMatrix<Type>&,
836);
837
838template<class Type>
839tmp<fvMatrix<Type>> operator==
840(
841 const fvMatrix<Type>&,
843);
844
845template<class Type>
846tmp<fvMatrix<Type>> operator==
847(
848 const fvMatrix<Type>&,
850);
851
852template<class Type>
853tmp<fvMatrix<Type>> operator==
854(
855 const tmp<fvMatrix<Type>>&,
857);
858
859template<class Type>
860tmp<fvMatrix<Type>> operator==
861(
862 const tmp<fvMatrix<Type>>&,
864);
865
866template<class Type>
867tmp<fvMatrix<Type>> operator==
868(
869 const tmp<fvMatrix<Type>>&,
871);
872
873template<class Type>
874tmp<fvMatrix<Type>> operator==
875(
876 const fvMatrix<Type>&,
877 const dimensioned<Type>&
878);
879
880template<class Type>
881tmp<fvMatrix<Type>> operator==
882(
883 const tmp<fvMatrix<Type>>&,
884 const dimensioned<Type>&
885);
886
887
888template<class Type>
889tmp<fvMatrix<Type>> operator==
890(
891 const fvMatrix<Type>&,
892 const Foam::zero
893);
894
895template<class Type>
896tmp<fvMatrix<Type>> operator==
897(
898 const tmp<fvMatrix<Type>>&,
899 const Foam::zero
900);
901
902
903//- Unary negation
904template<class Type>
905tmp<fvMatrix<Type>> operator-
906(
907 const fvMatrix<Type>&
908);
909
910//- Unary negation
911template<class Type>
912tmp<fvMatrix<Type>> operator-
913(
914 const tmp<fvMatrix<Type>>&
915);
916
917
918template<class Type>
919tmp<fvMatrix<Type>> operator+
920(
921 const fvMatrix<Type>&,
922 const fvMatrix<Type>&
923);
924
925template<class Type>
926tmp<fvMatrix<Type>> operator+
927(
928 const tmp<fvMatrix<Type>>&,
929 const fvMatrix<Type>&
930);
931
932template<class Type>
933tmp<fvMatrix<Type>> operator+
934(
935 const fvMatrix<Type>&,
936 const tmp<fvMatrix<Type>>&
937);
938
939template<class Type>
940tmp<fvMatrix<Type>> operator+
941(
942 const tmp<fvMatrix<Type>>&,
943 const tmp<fvMatrix<Type>>&
944);
945
946
947template<class Type>
948tmp<fvMatrix<Type>> operator+
949(
950 const fvMatrix<Type>&,
952);
953
954template<class Type>
955tmp<fvMatrix<Type>> operator+
956(
957 const fvMatrix<Type>&,
959);
960
961template<class Type>
962tmp<fvMatrix<Type>> operator+
963(
964 const fvMatrix<Type>&,
966);
967
968template<class Type>
969tmp<fvMatrix<Type>> operator+
970(
971 const tmp<fvMatrix<Type>>&,
973);
974
975template<class Type>
976tmp<fvMatrix<Type>> operator+
977(
978 const tmp<fvMatrix<Type>>&,
980);
981
982template<class Type>
983tmp<fvMatrix<Type>> operator+
984(
985 const tmp<fvMatrix<Type>>&,
987);
988
989template<class Type>
990tmp<fvMatrix<Type>> operator+
991(
993 const fvMatrix<Type>&
994);
995
996template<class Type>
997tmp<fvMatrix<Type>> operator+
998(
1000 const fvMatrix<Type>&
1001);
1002
1003template<class Type>
1004tmp<fvMatrix<Type>> operator+
1005(
1007 const fvMatrix<Type>&
1008);
1009
1010template<class Type>
1011tmp<fvMatrix<Type>> operator+
1012(
1014 const tmp<fvMatrix<Type>>&
1015);
1016
1017template<class Type>
1018tmp<fvMatrix<Type>> operator+
1019(
1021 const tmp<fvMatrix<Type>>&
1022);
1023
1024template<class Type>
1025tmp<fvMatrix<Type>> operator+
1026(
1028 const tmp<fvMatrix<Type>>&
1029);
1030
1031
1032template<class Type>
1033tmp<fvMatrix<Type>> operator+
1034(
1035 const fvMatrix<Type>&,
1036 const dimensioned<Type>&
1037);
1038
1039template<class Type>
1040tmp<fvMatrix<Type>> operator+
1041(
1042 const tmp<fvMatrix<Type>>&,
1043 const dimensioned<Type>&
1044);
1045
1046template<class Type>
1047tmp<fvMatrix<Type>> operator+
1048(
1049 const dimensioned<Type>&,
1050 const fvMatrix<Type>&
1051);
1052
1053template<class Type>
1054tmp<fvMatrix<Type>> operator+
1055(
1056 const dimensioned<Type>&,
1057 const tmp<fvMatrix<Type>>&
1058);
1059
1060
1061template<class Type>
1062tmp<fvMatrix<Type>> operator-
1063(
1064 const fvMatrix<Type>&,
1065 const fvMatrix<Type>&
1066);
1067
1068template<class Type>
1069tmp<fvMatrix<Type>> operator-
1070(
1071 const tmp<fvMatrix<Type>>&,
1072 const fvMatrix<Type>&
1073);
1074
1075template<class Type>
1076tmp<fvMatrix<Type>> operator-
1077(
1078 const fvMatrix<Type>&,
1079 const tmp<fvMatrix<Type>>&
1080);
1081
1082template<class Type>
1083tmp<fvMatrix<Type>> operator-
1084(
1085 const tmp<fvMatrix<Type>>&,
1086 const tmp<fvMatrix<Type>>&
1087);
1088
1089
1090template<class Type>
1091tmp<fvMatrix<Type>> operator-
1092(
1093 const fvMatrix<Type>&,
1095);
1096
1097template<class Type>
1098tmp<fvMatrix<Type>> operator-
1099(
1100 const fvMatrix<Type>&,
1102);
1103
1104template<class Type>
1105tmp<fvMatrix<Type>> operator-
1106(
1107 const fvMatrix<Type>&,
1109);
1110
1111template<class Type>
1112tmp<fvMatrix<Type>> operator-
1113(
1114 const tmp<fvMatrix<Type>>&,
1116);
1117
1118template<class Type>
1119tmp<fvMatrix<Type>> operator-
1120(
1121 const tmp<fvMatrix<Type>>&,
1123);
1124
1125template<class Type>
1126tmp<fvMatrix<Type>> operator-
1127(
1128 const tmp<fvMatrix<Type>>&,
1130);
1131
1132template<class Type>
1133tmp<fvMatrix<Type>> operator-
1134(
1136 const fvMatrix<Type>&
1137);
1138
1139template<class Type>
1140tmp<fvMatrix<Type>> operator-
1141(
1143 const fvMatrix<Type>&
1144);
1145
1146template<class Type>
1147tmp<fvMatrix<Type>> operator-
1148(
1150 const fvMatrix<Type>&
1151);
1152
1153template<class Type>
1154tmp<fvMatrix<Type>> operator-
1155(
1157 const tmp<fvMatrix<Type>>&
1158);
1159
1160template<class Type>
1161tmp<fvMatrix<Type>> operator-
1162(
1164 const tmp<fvMatrix<Type>>&
1165);
1166
1167template<class Type>
1168tmp<fvMatrix<Type>> operator-
1169(
1171 const tmp<fvMatrix<Type>>&
1172);
1173
1174
1175template<class Type>
1176tmp<fvMatrix<Type>> operator-
1177(
1178 const fvMatrix<Type>&,
1179 const dimensioned<Type>&
1180);
1181
1182template<class Type>
1183tmp<fvMatrix<Type>> operator-
1184(
1185 const tmp<fvMatrix<Type>>&,
1186 const dimensioned<Type>&
1187);
1188
1189template<class Type>
1190tmp<fvMatrix<Type>> operator-
1191(
1192 const dimensioned<Type>&,
1193 const fvMatrix<Type>&
1194);
1195
1196template<class Type>
1197tmp<fvMatrix<Type>> operator-
1198(
1199 const dimensioned<Type>&,
1200 const tmp<fvMatrix<Type>>&
1201);
1202
1203
1204template<class Type>
1205tmp<fvMatrix<Type>> operator*
1206(
1208 const fvMatrix<Type>&
1209);
1210
1211template<class Type>
1212tmp<fvMatrix<Type>> operator*
1213(
1215 const fvMatrix<Type>&
1216);
1217
1218template<class Type>
1219tmp<fvMatrix<Type>> operator*
1220(
1221 const tmp<volScalarField>&,
1222 const fvMatrix<Type>&
1223);
1224
1225template<class Type>
1226tmp<fvMatrix<Type>> operator*
1227(
1229 const tmp<fvMatrix<Type>>&
1230);
1231
1232template<class Type>
1233tmp<fvMatrix<Type>> operator*
1234(
1236 const tmp<fvMatrix<Type>>&
1237);
1238
1239template<class Type>
1240tmp<fvMatrix<Type>> operator*
1241(
1242 const tmp<volScalarField>&,
1243 const tmp<fvMatrix<Type>>&
1244);
1245
1246
1247template<class Type>
1248tmp<fvMatrix<Type>> operator*
1249(
1250 const dimensioned<scalar>&,
1251 const fvMatrix<Type>&
1252);
1253
1254template<class Type>
1255tmp<fvMatrix<Type>> operator*
1256(
1257 const dimensioned<scalar>&,
1258 const tmp<fvMatrix<Type>>&
1259);
1260
1261
1262// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1263
1264} // End namespace Foam
1265
1266// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1267
1268#ifdef NoRepository
1269 #include "fvMatrix.C"
1270#endif
1271
1272// Specialisation for scalars
1273#include "fvScalarMatrix.H"
1274
1275// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1276
1277#endif
1278
1279// ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A field of fields is a PtrList of fields with reference counting.
Definition: FieldField.H:80
Generic templated field type.
Definition: Field.H:82
Generic GeometricField class.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
A dynamically resizable PtrList with allocation management.
Definition: PtrDynList.H:64
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
void clear()
Clear the PtrList. Delete allocated entries and set size to zero.
Definition: PtrListI.H:97
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics.
A List with indirect addressing. Like IndirectList but does not store addressing.
Definition: IndirectList.H:79
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Definition: dimensionSet.H:109
Generic dimensioned Type class.
SolverPerformance< Type > solve()
Solve returning the solution statistics.
fvSolver(fvMatrix< Type > &fvMat, autoPtr< lduMatrix::solver > &&sol)
Definition: fvMatrix.H:272
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
void operator-=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1672
autoPtr< fvSolver > solver()
Construct and return the solver.
tmp< Field< Type > > DD() const
Return the matrix Type diagonal.
Definition: fvMatrix.C:1279
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:1443
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:1303
void manipulateMatrix(direction cmp)
Manipulate matrix.
Definition: fvMatrix.C:790
void addToInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Add patch contribution to internal field.
Definition: fvMatrix.C:49
label globalPatchID(const label fieldi, const label patchi) const
Definition: fvMatrix.H:355
void mapContributions(label fieldi, const FieldField< Field, Type > &fluxContrib, FieldField< Field, Type > &contrib, bool internal) const
Add internal and boundary contribution to local patches.
Definition: fvMatrix.C:548
void setReferences(const labelUList &cellLabels, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:1018
const Field< Type > & source() const noexcept
Definition: fvMatrix.H:463
const FieldField< Field, Type > & internalCoeffs() const noexcept
Definition: fvMatrix.H:470
void transferFvMatrixCoeffs()
Transfer lower, upper, diag and source to this fvMatrix.
Definition: fvMatrix.C:813
GeometricField< Type, fvsPatchField, surfaceMesh > * faceFluxFieldPtrType
Declare return type of the faceFluxCorrectionPtr() function.
Definition: fvMatrix.H:498
void relax()
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:1240
SolverPerformance< Type > solve()
Solve returning the solution statistics.
tmp< volScalarField > H1() const
Return H(1)
Definition: fvMatrix.C:1394
void operator+=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1634
const fvMatrix< Type > & matrix(const label i) const
Definition: fvMatrix.H:344
void createOrUpdateLduPrimitiveAssembly()
Create or update ldu assembly.
Definition: fvMatrix.C:906
faceFluxFieldPtrType & faceFluxCorrectionPtr()
Return pointer to face-flux non-orthogonal correction field.
Definition: fvMatrix.H:501
FieldField< Field, Type > & internalCoeffs() noexcept
Definition: fvMatrix.H:477
const FieldField< Field, Type > & boundaryCoeffs() const noexcept
Definition: fvMatrix.H:484
GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0)
Definition: fvMatrix.H:424
virtual ~fvMatrix()
Destructor.
Definition: fvMatrix.C:456
GeometricField< Type, fvPatchField, volMesh > psiFieldType
The geometric field type for psi.
Definition: fvMatrix.H:129
GeometricField< Type, fvsPatchField, surfaceMesh > faceFluxFieldType
Field type for face flux (for non-orthogonal correction)
Definition: fvMatrix.H:134
void setValues(const labelUList &cellLabels, const Type &value)
Definition: fvMatrix.C:969
label nMatrices() const
Definition: fvMatrix.H:339
const dimensionSet & dimensions() const noexcept
Definition: fvMatrix.H:453
const lduPrimitiveMeshAssembly & lduMeshAssembly()
Return optional lduAdressing.
Definition: fvMatrix.H:405
void operator*=(const volScalarField::Internal &)
Definition: fvMatrix.C:1805
void addCmptAvBoundaryDiag(scalarField &diag) const
Definition: fvMatrix.C:152
void addBoundarySource(Field< Type > &source, const bool couples=true) const
Definition: fvMatrix.C:177
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:1002
void setComponentReference(const label patchi, const label facei, const direction cmpt, const scalar value)
Definition: fvMatrixSolve.C:38
void boundaryManipulate(typename GeometricField< Type, fvPatchField, volMesh >::Boundary &values)
Manipulate based on a boundary field.
Definition: fvMatrix.C:1257
void setInterfaces(lduInterfaceFieldPtrsList &, PtrDynList< lduInterfaceField > &newInterfaces)
Set interfaces.
Definition: fvMatrix.C:470
SolverPerformance< Type > solveSegregatedOrCoupled(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:62
ClassName("fvMatrix")
SolverPerformance< Type > solveCoupled(const dictionary &)
Solve coupled returning the solution statistics.
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
Definition: fvMatrix.C:1332
tmp< Field< Type > > residual() const
Return the matrix residual.
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
Definition: fvMatrix.C:124
void negate()
Inplace negate.
Definition: fvMatrix.C:1619
tmp< fvMatrix< Type > > clone() const
Construct and return a clone.
Definition: fvMatrix.H:323
SolverPerformance< Type > solveSegregated(const dictionary &)
Solve segregated returning the solution statistics.
void clear()
Clear multiple fvMatrices.
Definition: fvMatrix.H:446
void addFvMatrix(fvMatrix< Type > &matrix)
Add fvMatrix.
Definition: fvMatrix.C:1063
tmp< scalarField > D() const
Return the matrix scalar diagonal.
Definition: fvMatrix.C:1270
void subtractFromInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Subtract patch contribution from internal field.
Definition: fvMatrix.C:87
const GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0) const
Return psi.
Definition: fvMatrix.H:412
bool checkImplicit(const label fieldi=0)
Name the implicit assembly addressing.
Definition: fvMatrix.C:307
lduPrimitiveMeshAssembly * lduMeshPtr()
Access to lduPrimitiveMeshAssembly.
Definition: fvMatrix.C:880
void operator=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1574
bool hasFaceFluxCorrection() const noexcept
True if face-flux non-orthogonal correction field exists.
Definition: fvMatrix.H:507
void setValuesFromList(const labelUList &cellLabels, const ListType< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:227
FieldField< Field, Type > & boundaryCoeffs() noexcept
Definition: fvMatrix.H:491
Field< Type > & source() noexcept
Definition: fvMatrix.H:458
fvMatrix< Type > & matrix(const label i)
Definition: fvMatrix.H:349
const dictionary & solverDict() const
Return the solver dictionary taking into account finalIteration.
Definition: fvMatrix.C:1558
void setBounAndInterCoeffs()
Manipulate boundary/internal coeffs for coupling.
Definition: fvMatrix.C:660
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:82
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
Definition: lduMatrix.H:84
scalarField & diag()
Definition: lduMatrix.C:192
An assembly of lduMatrix that is specific inter-region coupling through mapped patches.
const labelListList & patchMap() const
Return patchMap.
Reference counter for various OpenFOAM components.
Definition: refCount.H:51
A class for managing temporary objects.
Definition: tmp.H:65
Mesh data needed to do the Finite Volume discretisation.
Definition: volMesh.H:54
A class for handling words, derived from Foam::string.
Definition: word.H:68
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Definition: zero.H:63
Macro definitions for declaring ClassName(), NamespaceName(), etc.
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
Definition: className.H:67
const volScalarField & psi
bool
Definition: EEqn.H:20
A scalar instance of fvMatrix.
Namespace for OpenFOAM.
void checkMethod(const faMatrix< Type > &, const faMatrix< Type > &, const char *)
Definition: faMatrix.C:1050
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
uint8_t direction
Definition: direction.H:56
const direction noexcept
Definition: Scalar.H:223
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
volScalarField & alpha
CEqn solve()
#define FOAM_DEPRECATED_FOR(since, replacement)
Definition: stdFoam.H:52
Foam::surfaceFields.