FreeFOAM The Cross-Platform CFD Toolkit
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 | Copyright (C) 1991-2010 OpenCFD Ltd.
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 <finiteVolume/volFields.H>
45 #include <OpenFOAM/lduMatrix.H>
46 #include <OpenFOAM/tmp.H>
47 #include <OpenFOAM/autoPtr.H>
49 #include <OpenFOAM/zeroField.H>
50 #include <OpenFOAM/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>&);
106 
107 
108 /*---------------------------------------------------------------------------*\
109  Class fvMatrix Declaration
110 \*---------------------------------------------------------------------------*/
111 
112 template<class Type>
113 class fvMatrix
114 :
115  public refCount,
116  public lduMatrix
117 {
118  // Private data
119 
120  // Reference to GeometricField<Type, fvPatchField, volMesh>
122 
123  //- Dimension set
124  dimensionSet dimensions_;
125 
126  //- Source term
127  Field<Type> source_;
128 
129  //- Boundary scalar field containing pseudo-matrix coeffs
130  // for internal cells
131  FieldField<Field, Type> internalCoeffs_;
132 
133  //- Boundary scalar field containing pseudo-matrix coeffs
134  // for boundary cells
135  FieldField<Field, Type> boundaryCoeffs_;
136 
137 
138  //- Face flux field for non-orthogonal correction
140  *faceFluxCorrectionPtr_;
141 
142 
143  // Private member functions
144 
145  //- Add patch contribution to internal field
146  template<class Type2>
147  void addToInternalField
148  (
149  const unallocLabelList& addr,
150  const Field<Type2>& pf,
151  Field<Type2>& intf
152  ) const;
153 
154  template<class Type2>
155  void addToInternalField
156  (
157  const unallocLabelList& addr,
158  const tmp<Field<Type2> >& tpf,
159  Field<Type2>& intf
160  ) const;
161 
162  //- Subtract patch contribution from internal field
163  template<class Type2>
164  void subtractFromInternalField
165  (
166  const unallocLabelList& addr,
167  const Field<Type2>& pf,
168  Field<Type2>& intf
169  ) const;
170 
171  template<class Type2>
172  void subtractFromInternalField
173  (
174  const unallocLabelList& addr,
175  const tmp<Field<Type2> >& tpf,
176  Field<Type2>& intf
177  ) const;
178 
179 
180  // Matrix completion functionality
181 
182  void addBoundaryDiag
183  (
184  scalarField& diag,
185  const direction cmpt
186  ) const;
187 
188  void addCmptAvBoundaryDiag(scalarField& diag) const;
189 
190  void addBoundarySource
191  (
193  const bool couples=true
194  ) const;
195 
196 
197 public:
198 
199  //- Solver class returned by the solver function
200  // used for systems in which it is useful to cache the solver for reuse
201  // e.g. if the solver is potentialy expensive to construct (AMG) and can
202  // be used several times (PISO)
203  class fvSolver
204  {
205  fvMatrix<Type>& fvMat_;
206 
208 
209  public:
210 
211  // Constructors
212 
214  :
215  fvMat_(fvMat),
216  solver_(sol)
217  {}
218 
219 
220  // Member functions
221 
222  //- Solve returning the solution statistics.
223  // Use the given solver controls
225 
226  //- Solve returning the solution statistics.
227  // Solver controls read from fvSolution
229  };
230 
231 
232  ClassName("fvMatrix");
233 
234 
235  // Constructors
236 
237  //- Construct given a field to solve for
238  fvMatrix
239  (
241  const dimensionSet&
242  );
243 
244  //- Construct as copy
245  fvMatrix(const fvMatrix<Type>&);
246 
247  //- Construct as copy of tmp<fvMatrix<Type> > deleting argument
248 # ifdef ConstructFromTmp
249  fvMatrix(const tmp<fvMatrix<Type> >&);
250 # endif
251 
252  //- Construct from Istream given field to solve for
254 
255 
256  // Destructor
257 
258  virtual ~fvMatrix();
259 
260 
261  // Member functions
262 
263  // Access
264 
266  {
267  return psi_;
268  }
269 
271  {
272  return psi_;
273  }
274 
275  const dimensionSet& dimensions() const
276  {
277  return dimensions_;
278  }
279 
281  {
282  return source_;
283  }
284 
285  const Field<Type>& source() const
286  {
287  return source_;
288  }
289 
290  //- fvBoundary scalar field containing pseudo-matrix coeffs
291  // for internal cells
293  {
294  return internalCoeffs_;
295  }
296 
297  //- fvBoundary scalar field containing pseudo-matrix coeffs
298  // for boundary cells
300  {
301  return boundaryCoeffs_;
302  }
303 
304 
305  //- Declare return type of the faceFluxCorrectionPtr() function
308 
309  //- Return pointer to face-flux non-orthogonal correction field
311  {
312  return faceFluxCorrectionPtr_;
313  }
314 
315 
316  // Operations
317 
318  //- Set solution in given cells and eliminate corresponding
319  // equations from the matrix
320  void setValues
321  (
322  const labelList& cells,
323  const Field<Type>& values
324  );
325 
326  //- Set reference level for solution
327  void setReference
328  (
329  const label celli,
330  const Type& value,
331  const bool forceReference = false
332  );
333 
334  //- Set reference level for a component of the solution
335  // on a given patch face
337  (
338  const label patchi,
339  const label facei,
340  const direction cmpt,
341  const scalar value
342  );
343 
344  //- Relax matrix (for steady-state solution).
345  // alpha = 1 : diagonally equal
346  // alpha < 1 : diagonally dominant
347  // alpha = 0 : do nothing
348  void relax(const scalar alpha);
349 
350  //- Relax matrix (for steady-state solution).
351  // alpha is read from controlDict
352  void relax();
353 
354  //- Manipulate based on a boundary field
355  void boundaryManipulate
356  (
358  GeometricBoundaryField& values
359  );
360 
361  //- Construct and return the solver
362  // Use the given solver controls
364 
365  //- Construct and return the solver
366  // Solver controls read from fvSolution
368 
369  //- Solve returning the solution statistics.
370  // Use the given solver controls
372 
373  //- Solve returning the solution statistics.
374  // Solver controls read from fvSolution
376 
377  //- Return the matrix residual
378  tmp<Field<Type> > residual() const;
379 
380  //- Return the matrix scalar diagonal
381  tmp<scalarField> D() const;
382 
383  //- Return the matrix Type diagonal
384  tmp<Field<Type> > DD() const;
385 
386  //- Return the central coefficient
387  tmp<volScalarField> A() const;
388 
389  //- Return the H operation source
391 
392  //- Return H(1)
393  tmp<volScalarField> H1() const;
394 
395  //- Return the face-flux field from the matrix
397  flux() const;
398 
399 
400  // Member operators
401 
402  void operator=(const fvMatrix<Type>&);
403  void operator=(const tmp<fvMatrix<Type> >&);
404 
405  void negate();
406 
407  void operator+=(const fvMatrix<Type>&);
408  void operator+=(const tmp<fvMatrix<Type> >&);
409 
410  void operator-=(const fvMatrix<Type>&);
411  void operator-=(const tmp<fvMatrix<Type> >&);
412 
413  void operator+=
414  (
416  );
417  void operator+=
418  (
420  );
421  void operator+=
422  (
424  );
425 
426  void operator-=
427  (
429  );
430  void operator-=
431  (
433  );
434  void operator-=
435  (
437  );
438 
439  void operator+=(const dimensioned<Type>&);
440  void operator-=(const dimensioned<Type>&);
441 
442  void operator+=(const zeroField&);
443  void operator-=(const zeroField&);
444 
447  void operator*=(const tmp<volScalarField>&);
448 
449  void operator*=(const dimensioned<scalar>&);
450 
451 
452  // Friend operators
453 
455  operator& <Type>
456  (
457  const fvMatrix<Type>&,
459  );
460 
462  operator& <Type>
463  (
464  const fvMatrix<Type>&,
466  );
467 
469  operator& <Type>
470  (
471  const tmp<fvMatrix<Type> >&,
473  );
474 
476  operator& <Type>
477  (
478  const tmp<fvMatrix<Type> >&,
480  );
481 
482 
483  // Ostream operator
484 
485  friend Ostream& operator<< <Type>
486  (
487  Ostream&,
488  const fvMatrix<Type>&
489  );
490 };
491 
492 
493 // * * * * * * * * * * * * * * * Global functions * * * * * * * * * * * * * //
494 
495 template<class Type>
496 void checkMethod
497 (
498  const fvMatrix<Type>&,
499  const fvMatrix<Type>&,
500  const char*
501 );
502 
503 template<class Type>
504 void checkMethod
505 (
506  const fvMatrix<Type>&,
507  const DimensionedField<Type, volMesh>&,
508  const char*
509 );
510 
511 template<class Type>
512 void checkMethod
513 (
514  const fvMatrix<Type>&,
515  const dimensioned<Type>&,
516  const char*
517 );
518 
519 
520 //- Solve returning the solution statistics given convergence tolerance
521 // Use the given solver controls
522 template<class Type>
523 lduMatrix::solverPerformance solve(fvMatrix<Type>&, const dictionary&);
524 
525 
526 //- Solve returning the solution statistics given convergence tolerance,
527 // deleting temporary matrix after solution.
528 // Use the given solver controls
529 template<class Type>
530 lduMatrix::solverPerformance solve
531 (
532  const tmp<fvMatrix<Type> >&,
533  const dictionary&
534 );
535 
536 
537 //- Solve returning the solution statistics given convergence tolerance
538 // Solver controls read fvSolution
539 template<class Type>
540 lduMatrix::solverPerformance solve(fvMatrix<Type>&);
541 
542 
543 //- Solve returning the solution statistics given convergence tolerance,
544 // deleting temporary matrix after solution.
545 // Solver controls read fvSolution
546 template<class Type>
547 lduMatrix::solverPerformance solve(const tmp<fvMatrix<Type> >&);
548 
549 
550 //- Return the correction form of the given matrix
551 // by subtracting the matrix multiplied by the current field
552 template<class Type>
553 tmp<fvMatrix<Type> > correction(const fvMatrix<Type>&);
554 
555 
556 //- Return the correction form of the given temporary matrix
557 // by subtracting the matrix multiplied by the current field
558 template<class Type>
559 tmp<fvMatrix<Type> > correction(const tmp<fvMatrix<Type> >&);
560 
561 
562 // * * * * * * * * * * * * * * * Global operators * * * * * * * * * * * * * //
563 
564 template<class Type>
565 tmp<fvMatrix<Type> > operator==
566 (
567  const fvMatrix<Type>&,
568  const fvMatrix<Type>&
569 );
570 
571 template<class Type>
572 tmp<fvMatrix<Type> > operator==
573 (
574  const tmp<fvMatrix<Type> >&,
575  const fvMatrix<Type>&
576 );
577 
578 template<class Type>
579 tmp<fvMatrix<Type> > operator==
580 (
581  const fvMatrix<Type>&,
582  const tmp<fvMatrix<Type> >&
583 );
584 
585 template<class Type>
586 tmp<fvMatrix<Type> > operator==
587 (
588  const tmp<fvMatrix<Type> >&,
589  const tmp<fvMatrix<Type> >&
590 );
591 
592 
593 template<class Type>
594 tmp<fvMatrix<Type> > operator==
595 (
596  const fvMatrix<Type>&,
597  const DimensionedField<Type, volMesh>&
598 );
599 
600 template<class Type>
601 tmp<fvMatrix<Type> > operator==
602 (
603  const fvMatrix<Type>&,
604  const tmp<DimensionedField<Type, volMesh> >&
605 );
606 
607 template<class Type>
608 tmp<fvMatrix<Type> > operator==
609 (
610  const fvMatrix<Type>&,
611  const tmp<GeometricField<Type, fvPatchField, volMesh> >&
612 );
613 
614 template<class Type>
615 tmp<fvMatrix<Type> > operator==
616 (
617  const tmp<fvMatrix<Type> >&,
618  const DimensionedField<Type, volMesh>&
619 );
620 
621 template<class Type>
622 tmp<fvMatrix<Type> > operator==
623 (
624  const tmp<fvMatrix<Type> >&,
625  const tmp<DimensionedField<Type, volMesh> >&
626 );
627 
628 template<class Type>
629 tmp<fvMatrix<Type> > operator==
630 (
631  const tmp<fvMatrix<Type> >&,
632  const tmp<GeometricField<Type, fvPatchField, volMesh> >&
633 );
634 
635 template<class Type>
636 tmp<fvMatrix<Type> > operator==
637 (
638  const fvMatrix<Type>&,
639  const dimensioned<Type>&
640 );
641 
642 template<class Type>
643 tmp<fvMatrix<Type> > operator==
644 (
645  const tmp<fvMatrix<Type> >&,
646  const dimensioned<Type>&
647 );
648 
649 
650 template<class Type>
651 tmp<fvMatrix<Type> > operator==
652 (
653  const fvMatrix<Type>&,
654  const zeroField&
655 );
656 
657 template<class Type>
658 tmp<fvMatrix<Type> > operator==
659 (
660  const tmp<fvMatrix<Type> >&,
661  const zeroField&
662 );
663 
664 
665 template<class Type>
666 tmp<fvMatrix<Type> > operator-
667 (
668  const fvMatrix<Type>&
669 );
670 
671 template<class Type>
672 tmp<fvMatrix<Type> > operator-
673 (
674  const tmp<fvMatrix<Type> >&
675 );
676 
677 
678 template<class Type>
679 tmp<fvMatrix<Type> > operator+
680 (
681  const fvMatrix<Type>&,
682  const fvMatrix<Type>&
683 );
684 
685 template<class Type>
686 tmp<fvMatrix<Type> > operator+
687 (
688  const tmp<fvMatrix<Type> >&,
689  const fvMatrix<Type>&
690 );
691 
692 template<class Type>
693 tmp<fvMatrix<Type> > operator+
694 (
695  const fvMatrix<Type>&,
696  const tmp<fvMatrix<Type> >&
697 );
698 
699 template<class Type>
700 tmp<fvMatrix<Type> > operator+
701 (
702  const tmp<fvMatrix<Type> >&,
703  const tmp<fvMatrix<Type> >&
704 );
705 
706 
707 template<class Type>
708 tmp<fvMatrix<Type> > operator+
709 (
710  const fvMatrix<Type>&,
711  const DimensionedField<Type, volMesh>&
712 );
713 
714 template<class Type>
715 tmp<fvMatrix<Type> > operator+
716 (
717  const fvMatrix<Type>&,
718  const tmp<DimensionedField<Type, volMesh> >&
719 );
720 
721 template<class Type>
722 tmp<fvMatrix<Type> > operator+
723 (
724  const fvMatrix<Type>&,
725  const tmp<GeometricField<Type, fvPatchField, volMesh> >&
726 );
727 
728 template<class Type>
729 tmp<fvMatrix<Type> > operator+
730 (
731  const tmp<fvMatrix<Type> >&,
732  const DimensionedField<Type, volMesh>&
733 );
734 
735 template<class Type>
736 tmp<fvMatrix<Type> > operator+
737 (
738  const tmp<fvMatrix<Type> >&,
739  const tmp<DimensionedField<Type, volMesh> >&
740 );
741 
742 template<class Type>
743 tmp<fvMatrix<Type> > operator+
744 (
745  const tmp<fvMatrix<Type> >&,
746  const tmp<GeometricField<Type, fvPatchField, volMesh> >&
747 );
748 
749 template<class Type>
750 tmp<fvMatrix<Type> > operator+
751 (
752  const DimensionedField<Type, volMesh>&,
753  const fvMatrix<Type>&
754 );
755 
756 template<class Type>
757 tmp<fvMatrix<Type> > operator+
758 (
759  const tmp<DimensionedField<Type, volMesh> >&,
760  const fvMatrix<Type>&
761 );
762 
763 template<class Type>
764 tmp<fvMatrix<Type> > operator+
765 (
766  const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
767  const fvMatrix<Type>&
768 );
769 
770 template<class Type>
771 tmp<fvMatrix<Type> > operator+
772 (
773  const DimensionedField<Type, volMesh>&,
774  const tmp<fvMatrix<Type> >&
775 );
776 
777 template<class Type>
778 tmp<fvMatrix<Type> > operator+
779 (
780  const tmp<DimensionedField<Type, volMesh> >&,
781  const tmp<fvMatrix<Type> >&
782 );
783 
784 template<class Type>
785 tmp<fvMatrix<Type> > operator+
786 (
787  const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
788  const tmp<fvMatrix<Type> >&
789 );
790 
791 
792 template<class Type>
793 tmp<fvMatrix<Type> > operator+
794 (
795  const fvMatrix<Type>&,
796  const dimensioned<Type>&
797 );
798 
799 template<class Type>
800 tmp<fvMatrix<Type> > operator+
801 (
802  const tmp<fvMatrix<Type> >&,
803  const dimensioned<Type>&
804 );
805 
806 template<class Type>
807 tmp<fvMatrix<Type> > operator+
808 (
809  const dimensioned<Type>&,
810  const fvMatrix<Type>&
811 );
812 
813 template<class Type>
814 tmp<fvMatrix<Type> > operator+
815 (
816  const dimensioned<Type>&,
817  const tmp<fvMatrix<Type> >&
818 );
819 
820 
821 template<class Type>
822 tmp<fvMatrix<Type> > operator-
823 (
824  const fvMatrix<Type>&,
825  const fvMatrix<Type>&
826 );
827 
828 template<class Type>
829 tmp<fvMatrix<Type> > operator-
830 (
831  const tmp<fvMatrix<Type> >&,
832  const fvMatrix<Type>&
833 );
834 
835 template<class Type>
836 tmp<fvMatrix<Type> > operator-
837 (
838  const fvMatrix<Type>&,
839  const tmp<fvMatrix<Type> >&
840 );
841 
842 template<class Type>
843 tmp<fvMatrix<Type> > operator-
844 (
845  const tmp<fvMatrix<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 DimensionedField<Type, volMesh>&
855 );
856 
857 template<class Type>
858 tmp<fvMatrix<Type> > operator-
859 (
860  const fvMatrix<Type>&,
861  const tmp<DimensionedField<Type, volMesh> >&
862 );
863 
864 template<class Type>
865 tmp<fvMatrix<Type> > operator-
866 (
867  const fvMatrix<Type>&,
868  const tmp<GeometricField<Type, fvPatchField, volMesh> >&
869 );
870 
871 template<class Type>
872 tmp<fvMatrix<Type> > operator-
873 (
874  const tmp<fvMatrix<Type> >&,
875  const DimensionedField<Type, volMesh>&
876 );
877 
878 template<class Type>
879 tmp<fvMatrix<Type> > operator-
880 (
881  const tmp<fvMatrix<Type> >&,
882  const tmp<DimensionedField<Type, volMesh> >&
883 );
884 
885 template<class Type>
886 tmp<fvMatrix<Type> > operator-
887 (
888  const tmp<fvMatrix<Type> >&,
889  const tmp<GeometricField<Type, fvPatchField, volMesh> >&
890 );
891 
892 template<class Type>
893 tmp<fvMatrix<Type> > operator-
894 (
895  const DimensionedField<Type, volMesh>&,
896  const fvMatrix<Type>&
897 );
898 
899 template<class Type>
900 tmp<fvMatrix<Type> > operator-
901 (
902  const tmp<DimensionedField<Type, volMesh> >&,
903  const fvMatrix<Type>&
904 );
905 
906 template<class Type>
907 tmp<fvMatrix<Type> > operator-
908 (
909  const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
910  const fvMatrix<Type>&
911 );
912 
913 template<class Type>
914 tmp<fvMatrix<Type> > operator-
915 (
916  const DimensionedField<Type, volMesh>&,
917  const tmp<fvMatrix<Type> >&
918 );
919 
920 template<class Type>
921 tmp<fvMatrix<Type> > operator-
922 (
923  const tmp<DimensionedField<Type, volMesh> >&,
924  const tmp<fvMatrix<Type> >&
925 );
926 
927 template<class Type>
928 tmp<fvMatrix<Type> > operator-
929 (
930  const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
931  const tmp<fvMatrix<Type> >&
932 );
933 
934 
935 template<class Type>
936 tmp<fvMatrix<Type> > operator-
937 (
938  const fvMatrix<Type>&,
939  const dimensioned<Type>&
940 );
941 
942 template<class Type>
943 tmp<fvMatrix<Type> > operator-
944 (
945  const tmp<fvMatrix<Type> >&,
946  const dimensioned<Type>&
947 );
948 
949 template<class Type>
950 tmp<fvMatrix<Type> > operator-
951 (
952  const dimensioned<Type>&,
953  const fvMatrix<Type>&
954 );
955 
956 template<class Type>
957 tmp<fvMatrix<Type> > operator-
958 (
959  const dimensioned<Type>&,
960  const tmp<fvMatrix<Type> >&
961 );
962 
963 
964 template<class Type>
965 tmp<fvMatrix<Type> > operator*
966 (
967  const DimensionedField<scalar, volMesh>&,
968  const fvMatrix<Type>&
969 );
970 
971 template<class Type>
972 tmp<fvMatrix<Type> > operator*
973 (
974  const tmp<DimensionedField<scalar, volMesh> >&,
975  const fvMatrix<Type>&
976 );
977 
978 template<class Type>
979 tmp<fvMatrix<Type> > operator*
980 (
981  const tmp<volScalarField>&,
982  const fvMatrix<Type>&
983 );
984 
985 template<class Type>
986 tmp<fvMatrix<Type> > operator*
987 (
988  const DimensionedField<scalar, volMesh>&,
989  const tmp<fvMatrix<Type> >&
990 );
991 
992 template<class Type>
993 tmp<fvMatrix<Type> > operator*
994 (
995  const tmp<DimensionedField<scalar, volMesh> >&,
996  const tmp<fvMatrix<Type> >&
997 );
998 
999 template<class Type>
1000 tmp<fvMatrix<Type> > operator*
1001 (
1002  const tmp<volScalarField>&,
1003  const tmp<fvMatrix<Type> >&
1004 );
1005 
1006 
1007 template<class Type>
1008 tmp<fvMatrix<Type> > operator*
1009 (
1010  const dimensioned<scalar>&,
1011  const fvMatrix<Type>&
1012 );
1013 
1014 template<class Type>
1015 tmp<fvMatrix<Type> > operator*
1016 (
1017  const dimensioned<scalar>&,
1018  const tmp<fvMatrix<Type> >&
1019 );
1020 
1021 
1022 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1023 
1024 } // End namespace Foam
1025 
1026 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1027 
1028 #ifdef NoRepository
1029 # include <finiteVolume/fvMatrix.C>
1030 #endif
1031 
1032 // Specialisation for scalars
1034 
1035 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1036 
1037 #endif
1038 
1039 // ************************ vim: set sw=4 sts=4 et: ************************ //