OpenFOAM logo
Open Source CFD Toolkit

fvMatrix.H

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2005 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software; you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM; if not, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     fvMatrix<Type>
00027 
00028 Description
00029     Finite-Volume matrix.
00030 
00031 SourceFiles
00032     fvMatrix.C
00033     fvMatrixSolve.C
00034 
00035 \*---------------------------------------------------------------------------*/
00036 
00037 #ifndef fvMatrix_H
00038 #define fvMatrix_H
00039 
00040 #include "volFields.H"
00041 #include "surfaceFields.H"
00042 #include "lduMatrix.H"
00043 #include "tmp.H"
00044 #include "autoPtr.H"
00045 #include "dimensionedTypes.H"
00046 #include "className.H"
00047 
00048 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00049 
00050 namespace Foam
00051 {
00052 
00053 // * * * * * * Forward declaration of template friend fuctions * * * * * * * //
00054 
00055 template<class Type>
00056 class fvMatrix;
00057 
00058 template<class Type>
00059 tmp<GeometricField<Type,fvPatchField,volMesh> > operator&
00060 (
00061     const fvMatrix<Type>&,
00062     const GeometricField<Type,fvPatchField,volMesh>&
00063 );
00064 
00065 template<class Type>
00066 tmp<GeometricField<Type,fvPatchField,volMesh> > operator&
00067 (
00068     const fvMatrix<Type>&,
00069     const tmp<GeometricField<Type,fvPatchField,volMesh> >&
00070 );
00071 
00072 template<class Type>
00073 tmp<GeometricField<Type,fvPatchField,volMesh> > operator&
00074 (
00075     const tmp<fvMatrix<Type> >&,
00076     const GeometricField<Type,fvPatchField,volMesh>&
00077 );
00078 
00079 template<class Type>
00080 tmp<GeometricField<Type,fvPatchField,volMesh> > operator&
00081 (
00082     const tmp<fvMatrix<Type> >&,
00083     const tmp<GeometricField<Type,fvPatchField,volMesh> >&
00084 );
00085 
00086 template<class Type>
00087 Ostream& operator<<(Ostream&, const fvMatrix<Type>&);
00088 
00089 
00090 /*---------------------------------------------------------------------------*\
00091                            Class fvMatrix Declaration
00092 \*---------------------------------------------------------------------------*/
00093 
00094 template<class Type>
00095 class fvMatrix
00096 :
00097     public refCount,
00098     public lduMatrix
00099 {
00100     // Private data
00101 
00102         // Reference to GeometricField<Type, fvPatchField, volMesh>
00103         GeometricField<Type, fvPatchField, volMesh>& psi_;
00104 
00105         //- Dimension set
00106         dimensionSet dimensions_;
00107 
00108         //- Source term
00109         Field<Type> source_;
00110 
00111         //- Boundary scalar field containing pseudo-matrix coeffs
00112         //  for internal cells
00113         FieldField<Field, Type> internalCoeffs_;
00114 
00115         //- Boundary scalar field containing pseudo-matrix coeffs
00116         //  for boundary cells
00117         FieldField<Field, Type> boundaryCoeffs_;
00118 
00119 
00120         //- Face flux field for non-orthogonal correction
00121         mutable GeometricField<Type, fvPatchField, surfaceMesh>
00122             *faceFluxCorrectionPtr_;
00123 
00124 
00125     // Private member functions
00126 
00127         //- Add patch contribution to internal field
00128         template<class Type2>
00129         void addToInternalField
00130         (
00131             const unallocLabelList& addr,
00132             const Field<Type2>& pf,
00133             Field<Type2>& intf
00134         ) const;
00135 
00136         template<class Type2>
00137         void addToInternalField
00138         (
00139             const unallocLabelList& addr,
00140             const tmp<Field<Type2> >& tpf,
00141             Field<Type2>& intf
00142         ) const;
00143 
00144         //- Subtract patch contribution from internal field
00145         template<class Type2>
00146         void subtractFromInternalField
00147         (
00148             const unallocLabelList& addr,
00149             const Field<Type2>& pf,
00150             Field<Type2>& intf
00151         ) const;
00152 
00153         template<class Type2>
00154         void subtractFromInternalField
00155         (
00156             const unallocLabelList& addr,
00157             const tmp<Field<Type2> >& tpf,
00158             Field<Type2>& intf
00159         ) const;
00160 
00161 
00162         // Matrix completion functionality
00163 
00164             void addBoundaryDiag
00165             (
00166                 scalarField& diag,
00167                 const direction cmpt
00168             ) const;
00169 
00170             void addCmptAvBoundaryDiag(scalarField& diag) const;
00171 
00172             void addBoundarySource(Field<Type>& source) const;
00173 
00174             void addBoundarySource
00175             (
00176                 scalarField& source,
00177                 const direction cmpt
00178             ) const;
00179 
00180 
00181 public:
00182 
00183     ClassName("fvMatrix");
00184 
00185 
00186     // Constructors
00187 
00188         //- Construct given a field to solve for
00189         fvMatrix
00190         (
00191             GeometricField<Type, fvPatchField, volMesh>&,
00192             const dimensionSet&
00193         );
00194 
00195         //- Construct as copy
00196         fvMatrix(const fvMatrix<Type>&);
00197 
00198         //- Construct as copy of tmp<fvMatrix<Type> > deleting argument
00199 #       ifdef ConstructFromTmp
00200         fvMatrix(const tmp<fvMatrix<Type> >&);
00201 #       endif
00202 
00203         //- Construct from Istream given field to solve for
00204         fvMatrix(GeometricField<Type, fvPatchField, volMesh>&, Istream&);
00205 
00206 
00207     // Destructor
00208 
00209         virtual ~fvMatrix();
00210 
00211 
00212     // Member functions
00213 
00214         // Access
00215 
00216             const GeometricField<Type, fvPatchField, volMesh>& psi() const
00217             {
00218                 return psi_;
00219             }
00220 
00221             GeometricField<Type, fvPatchField, volMesh>& psi()
00222             {
00223                 return psi_;
00224             }
00225 
00226             const dimensionSet& dimensions() const
00227             {
00228                 return dimensions_;
00229             }
00230 
00231             Field<Type>& source()
00232             {
00233                 return source_;
00234             }
00235 
00236             const Field<Type>& source() const
00237             {
00238                 return source_;
00239             }
00240 
00241             //- fvBoundary scalar field containing pseudo-matrix coeffs
00242             //  for internal cells
00243             FieldField<Field, Type>& internalCoeffs()
00244             {
00245                 return internalCoeffs_;
00246             }
00247 
00248             //- fvBoundary scalar field containing pseudo-matrix coeffs
00249             //  for boundary cells
00250             FieldField<Field, Type>& boundaryCoeffs()
00251             {
00252                 return boundaryCoeffs_;
00253             }
00254 
00255 
00256             //- Declare return type of the faceFluxCorrectionPtr() function
00257             typedef GeometricField<Type, fvPatchField, surfaceMesh>
00258                 *surfaceTypeFieldPtr;
00259 
00260             //- Return pointer to face-flux non-orthogonal correction field
00261             surfaceTypeFieldPtr& faceFluxCorrectionPtr()
00262             {
00263                 return faceFluxCorrectionPtr_;
00264             }
00265 
00266 
00267         // Operations
00268 
00269             //- Set solution in given cells and eliminate corresponding
00270             //  equations from the matrix
00271             void setValues
00272             (
00273                 const labelList& cells,
00274                 const Field<Type>& values
00275             );
00276 
00277             //- Set reference level for solution
00278             void setReference
00279             (
00280                 const label cell,
00281                 const Type& value
00282             );
00283 
00284             //- Set reference level for a component of the solution
00285             //  on a given patch face
00286             void setComponentReference
00287             (
00288                 const label patchi,
00289                 const label facei,
00290                 const direction cmpt,
00291                 const scalar value
00292             );
00293 
00294             //- Relax matrix (for steady-state solution).
00295             //  alpha = 1 : diagonally equal
00296             //  alpha < 1 :    ,,      dominant
00297             //  alpha = 0 : do nothing
00298             void relax(const scalar alpha);
00299 
00300             //- Relax matrix (for steadty-state solution).
00301             //  alpha is read from controlDict
00302             void relax();
00303 
00304             //- Solve returning the solution statistics.
00305             //  Solver controls read Istream
00306             lduMatrix::solverPerformance solve(Istream&);
00307 
00308             //- Solve returning the solution statistics.
00309             //  Solver controls read from fvSolution
00310             lduMatrix::solverPerformance solve();
00311 
00312             //- Return the matrix residual
00313             tmp<Field<Type> > residual() const;
00314 
00315             //- Return the matrix diagonal
00316             tmp<scalarField> D() const;
00317 
00318             //- Return the central coefficient
00319             tmp<volScalarField> A() const;
00320 
00321             //- Return the H operation source
00322             tmp<GeometricField<Type, fvPatchField, volMesh> > H() const;
00323 
00324             //- Return the face-flux field from the matrix
00325             tmp<GeometricField<Type, fvPatchField, surfaceMesh> > flux() const;
00326 
00327 
00328     // Member operators
00329 
00330         void operator=(const fvMatrix<Type>&);
00331         void operator=(const tmp<fvMatrix<Type> >&);
00332 
00333         void negate();
00334 
00335         void operator+=(const fvMatrix<Type>&);
00336         void operator+=(const tmp<fvMatrix<Type> >&);
00337 
00338         void operator-=(const fvMatrix<Type>&);
00339         void operator-=(const tmp<fvMatrix<Type> >&);
00340 
00341         void operator+=(const GeometricField<Type,fvPatchField,volMesh>&);
00342         void operator+=(const tmp<GeometricField<Type,fvPatchField,volMesh> >&);
00343 
00344         void operator-=(const GeometricField<Type,fvPatchField,volMesh>&);
00345         void operator-=(const tmp<GeometricField<Type,fvPatchField,volMesh> >&);
00346 
00347         void operator+=(const dimensioned<Type>&);
00348         void operator-=(const dimensioned<Type>&);
00349 
00350         void operator*=(const volScalarField&);
00351         void operator*=(const tmp<volScalarField>&);
00352 
00353         void operator*=(const dimensioned<scalar>&);
00354 
00355 
00356     // Friend operators
00357 
00358         friend tmp<GeometricField<Type,fvPatchField,volMesh> > operator& <Type>
00359         (
00360             const fvMatrix<Type>&,
00361             const GeometricField<Type,fvPatchField,volMesh>&
00362         );
00363 
00364         friend tmp<GeometricField<Type,fvPatchField,volMesh> > operator& <Type>
00365         (
00366             const fvMatrix<Type>&,
00367             const tmp<GeometricField<Type,fvPatchField,volMesh> >&
00368         );
00369 
00370         friend tmp<GeometricField<Type,fvPatchField,volMesh> > operator& <Type>
00371         (
00372             const tmp<fvMatrix<Type> >&,
00373             const GeometricField<Type,fvPatchField,volMesh>&
00374         );
00375 
00376         friend tmp<GeometricField<Type,fvPatchField,volMesh> > operator& <Type>
00377         (
00378             const tmp<fvMatrix<Type> >&,
00379             const tmp<GeometricField<Type,fvPatchField,volMesh> >&
00380         );
00381 
00382 
00383     // Ostream operator
00384 
00385         friend Ostream& operator<< <Type>
00386         (
00387             Ostream&,
00388             const fvMatrix<Type>&
00389         );
00390 };
00391 
00392 
00393 // * * * * * * * * * * * * * * * Global functions  * * * * * * * * * * * * * //
00394 
00395 template<class Type>
00396 void checkMethod
00397 (
00398     const fvMatrix<Type>&,
00399     const fvMatrix<Type>&,
00400     const char*
00401 );
00402 
00403 template<class Type>
00404 void checkMethod
00405 (
00406     const fvMatrix<Type>&,
00407     const GeometricField<Type, fvPatchField, volMesh>&,
00408     const char*
00409 );
00410 
00411 template<class Type>
00412 void checkMethod
00413 (
00414     const fvMatrix<Type>&,
00415     const dimensioned<Type>&,
00416     const char*
00417 );
00418 
00419 
00420 //- Solve returning the solution statistics given convergence tolerance
00421 //  Solver controls read Istream
00422 template<class Type>
00423 lduMatrix::solverPerformance solve(fvMatrix<Type>&, Istream&);
00424 
00425 
00426 //- Solve returning the solution statistics given convergence tolerance,
00427 //  deleting temporary matrix after solution.
00428 //  Solver controls read Istream
00429 template<class Type>
00430 lduMatrix::solverPerformance solve(const tmp<fvMatrix<Type> >&, Istream&);
00431 
00432 
00433 //- Solve returning the solution statistics given convergence tolerance
00434 //  Solver controls read fvSolution
00435 template<class Type>
00436 lduMatrix::solverPerformance solve(fvMatrix<Type>&);
00437 
00438 
00439 //- Solve returning the solution statistics given convergence tolerance,
00440 //  deleting temporary matrix after solution.
00441 //  Solver controls read fvSolution
00442 template<class Type>
00443 lduMatrix::solverPerformance solve(const tmp<fvMatrix<Type> >&);
00444 
00445 
00446 // * * * * * * * * * * * * * * * Global operators  * * * * * * * * * * * * * //
00447 
00448 template<class Type>
00449 tmp<fvMatrix<Type> > operator-
00450 (
00451     const fvMatrix<Type>&
00452 );
00453 
00454 template<class Type>
00455 tmp<fvMatrix<Type> > operator-
00456 (
00457     const tmp<fvMatrix<Type> >&
00458 );
00459 
00460 template<class Type>
00461 tmp<fvMatrix<Type> > operator+
00462 (
00463     const fvMatrix<Type>&,
00464     const fvMatrix<Type>&
00465 );
00466 
00467 template<class Type>
00468 tmp<fvMatrix<Type> > operator+
00469 (
00470     const tmp<fvMatrix<Type> >&,
00471     const fvMatrix<Type>&
00472 );
00473 
00474 template<class Type>
00475 tmp<fvMatrix<Type> > operator+
00476 (
00477     const fvMatrix<Type>&,
00478     const tmp<fvMatrix<Type> >&
00479 );
00480 
00481 template<class Type>
00482 tmp<fvMatrix<Type> > operator+
00483 (
00484     const tmp<fvMatrix<Type> >&,
00485     const tmp<fvMatrix<Type> >&
00486 );
00487 
00488 template<class Type>
00489 tmp<fvMatrix<Type> > operator-
00490 (
00491     const fvMatrix<Type>&,
00492     const fvMatrix<Type>&
00493 );
00494 
00495 template<class Type>
00496 tmp<fvMatrix<Type> > operator-
00497 (
00498     const tmp<fvMatrix<Type> >&,
00499     const fvMatrix<Type>&
00500 );
00501 
00502 template<class Type>
00503 tmp<fvMatrix<Type> > operator-
00504 (
00505     const fvMatrix<Type>&,
00506     const tmp<fvMatrix<Type> >&
00507 );
00508 
00509 template<class Type>
00510 tmp<fvMatrix<Type> > operator-
00511 (
00512     const tmp<fvMatrix<Type> >&,
00513     const tmp<fvMatrix<Type> >&
00514 );
00515 
00516 template<class Type>
00517 tmp<fvMatrix<Type> > operator==
00518 (
00519     const fvMatrix<Type>&,
00520     const fvMatrix<Type>&
00521 );
00522 
00523 template<class Type>
00524 tmp<fvMatrix<Type> > operator==
00525 (
00526     const tmp<fvMatrix<Type> >&,
00527     const fvMatrix<Type>&
00528 );
00529 
00530 template<class Type>
00531 tmp<fvMatrix<Type> > operator==
00532 (
00533     const fvMatrix<Type>&,
00534     const tmp<fvMatrix<Type> >&
00535 );
00536 
00537 template<class Type>
00538 tmp<fvMatrix<Type> > operator==
00539 (
00540     const tmp<fvMatrix<Type> >&,
00541     const tmp<fvMatrix<Type> >&
00542 );
00543 
00544 template<class Type>
00545 tmp<fvMatrix<Type> > operator+
00546 (
00547     const fvMatrix<Type>&,
00548     const GeometricField<Type, fvPatchField, volMesh>&
00549 );
00550 
00551 template<class Type>
00552 tmp<fvMatrix<Type> > operator+
00553 (
00554     const tmp<fvMatrix<Type> >&,
00555     const GeometricField<Type, fvPatchField, volMesh>&
00556 );
00557 
00558 template<class Type>
00559 tmp<fvMatrix<Type> > operator+
00560 (
00561     const fvMatrix<Type>&,
00562     const tmp<GeometricField<Type, fvPatchField, volMesh> >&
00563 );
00564 
00565 template<class Type>
00566 tmp<fvMatrix<Type> > operator+
00567 (
00568     const tmp<fvMatrix<Type> >&,
00569     const tmp<GeometricField<Type, fvPatchField, volMesh> >&
00570 );
00571 
00572 template<class Type>
00573 tmp<fvMatrix<Type> > operator+
00574 (
00575     const GeometricField<Type, fvPatchField, volMesh>&,
00576     const fvMatrix<Type>&
00577 );
00578 
00579 template<class Type>
00580 tmp<fvMatrix<Type> > operator+
00581 (
00582     const GeometricField<Type, fvPatchField, volMesh>&,
00583     const tmp<fvMatrix<Type> >&
00584 );
00585 
00586 template<class Type>
00587 tmp<fvMatrix<Type> > operator+
00588 (
00589     const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
00590     const fvMatrix<Type>&
00591 );
00592 
00593 template<class Type>
00594 tmp<fvMatrix<Type> > operator+
00595 (
00596     const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
00597     const tmp<fvMatrix<Type> >&
00598 );
00599 
00600 template<class Type>
00601 tmp<fvMatrix<Type> > operator-
00602 (
00603     const fvMatrix<Type>&,
00604     const GeometricField<Type, fvPatchField, volMesh>&
00605 );
00606 
00607 template<class Type>
00608 tmp<fvMatrix<Type> > operator-
00609 (
00610     const tmp<fvMatrix<Type> >&,
00611     const GeometricField<Type, fvPatchField, volMesh>&
00612 );
00613 
00614 template<class Type>
00615 tmp<fvMatrix<Type> > operator-
00616 (
00617     const fvMatrix<Type>&,
00618     const tmp<GeometricField<Type, fvPatchField, volMesh> >&
00619 );
00620 
00621 template<class Type>
00622 tmp<fvMatrix<Type> > operator-
00623 (
00624     const tmp<fvMatrix<Type> >&,
00625     const tmp<GeometricField<Type, fvPatchField, volMesh> >&
00626 );
00627 
00628 template<class Type>
00629 tmp<fvMatrix<Type> > operator-
00630 (
00631     const GeometricField<Type, fvPatchField, volMesh>&,
00632     const fvMatrix<Type>&
00633 );
00634 
00635 template<class Type>
00636 tmp<fvMatrix<Type> > operator-
00637 (
00638     const GeometricField<Type, fvPatchField, volMesh>&,
00639     const tmp<fvMatrix<Type> >&
00640 );
00641 
00642 template<class Type>
00643 tmp<fvMatrix<Type> > operator-
00644 (
00645     const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
00646     const fvMatrix<Type>&
00647 );
00648 
00649 template<class Type>
00650 tmp<fvMatrix<Type> > operator-
00651 (
00652     const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
00653     const tmp<fvMatrix<Type> >&
00654 );
00655 
00656 template<class Type>
00657 tmp<fvMatrix<Type> > operator+
00658 (
00659     const tmp<fvMatrix<Type> >&,
00660     const dimensioned<Type>&
00661 );
00662 
00663 template<class Type>
00664 tmp<fvMatrix<Type> > operator+
00665 (
00666     const dimensioned<Type>&,
00667     const tmp<fvMatrix<Type> >&
00668 );
00669 
00670 template<class Type>
00671 tmp<fvMatrix<Type> > operator-
00672 (
00673     const tmp<fvMatrix<Type> >&,
00674     const dimensioned<Type>&
00675 );
00676 
00677 template<class Type>
00678 tmp<fvMatrix<Type> > operator-
00679 (
00680     const dimensioned<Type>&,
00681     const tmp<fvMatrix<Type> >&
00682 );
00683 
00684 template<class Type>
00685 tmp<fvMatrix<Type> > operator==
00686 (
00687     const fvMatrix<Type>&,
00688     const GeometricField<Type, fvPatchField, volMesh>&
00689 );
00690 
00691 template<class Type>
00692 tmp<fvMatrix<Type> > operator==
00693 (
00694     const tmp<fvMatrix<Type> >&,
00695     const GeometricField<Type, fvPatchField, volMesh>&
00696 );
00697 
00698 template<class Type>
00699 tmp<fvMatrix<Type> > operator==
00700 (
00701     const fvMatrix<Type>&,
00702     const tmp<GeometricField<Type, fvPatchField, volMesh> >&
00703 );
00704 
00705 template<class Type>
00706 tmp<fvMatrix<Type> > operator==
00707 (
00708     const tmp<fvMatrix<Type> >&,
00709     const tmp<GeometricField<Type, fvPatchField, volMesh> >&
00710 );
00711 
00712 template<class Type>
00713 tmp<fvMatrix<Type> > operator==
00714 (
00715     const fvMatrix<Type>&,
00716     const dimensioned<Type>&
00717 );
00718 
00719 template<class Type>
00720 tmp<fvMatrix<Type> > operator==
00721 (
00722     const tmp<fvMatrix<Type> >&,
00723     const dimensioned<Type>&
00724 );
00725 
00726 
00727 template<class Type>
00728 tmp<fvMatrix<Type> > operator*
00729 (
00730     const volScalarField&,
00731     const fvMatrix<Type>&
00732 );
00733 
00734 template<class Type>
00735 tmp<fvMatrix<Type> > operator*
00736 (
00737     const volScalarField&,
00738     const tmp<fvMatrix<Type> >&
00739 );
00740 
00741 template<class Type>
00742 tmp<fvMatrix<Type> > operator*
00743 (
00744     const tmp<volScalarField>&,
00745     const fvMatrix<Type>&
00746 );
00747 
00748 template<class Type>
00749 tmp<fvMatrix<Type> > operator*
00750 (
00751     const tmp<volScalarField>&,
00752     const tmp<fvMatrix<Type> >&
00753 );
00754 
00755 
00756 template<class Type>
00757 tmp<fvMatrix<Type> > operator*
00758 (
00759     const dimensioned<scalar>&,
00760     const fvMatrix<Type>&
00761 );
00762 
00763 template<class Type>
00764 tmp<fvMatrix<Type> > operator*
00765 (
00766     const dimensioned<scalar>&,
00767     const tmp<fvMatrix<Type> >&
00768 );
00769 
00770 
00771 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00772 
00773 } // End namespace Foam
00774 
00775 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00776 
00777 #ifdef NoRepository
00778 #   include "fvMatrix.C"
00779 #endif
00780 
00781 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00782 
00783 #endif
00784 
00785 // ************************************************************************* //
For further information go to www.openfoam.org