effectivenessHeatExchangerSource.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) 2013-2021 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::fv::effectivenessHeatExchangerSource
26 
27 Description
28  Heat exchanger source model, in which the heat exchanger is defined as a
29  selection of cells.
30 
31  The total heat exchange source is given by:
32  \f[
33  Q_t = e(\phi, \dot{m}_2) (T_2 - T_1) \phi c_p
34  \f]
35 
36  where:
37  \vartable
38  Q_t | total heat source
39  e(\phi,\dot{m}_2) | effectiveness table
40  \phi | net mass flux entering heat exchanger [kg/s]
41  \dot{m}_2 | secondary mass flow rate [kg/s]
42  T_1 | primary inlet temperature [K]
43  T_2 | secondary inlet temperature [K]
44  c_p | specific heat capacity [J/kg/K]
45  \endvartable
46 
47 
48  The distribution inside the hear exchanger is given by:
49  \f[
50  Q_c = \frac{V_c |U_c| (T_c - T_{ref})}{\sum(V_c |U_c| (T_c - T_{ref}))}
51  \f]
52 
53  where:
54  \vartable
55  Q_c | source for cell
56  V_c | volume of the cell [m^3]
57  U_c | local cell velocity [m/s]
58  T_c | local call temperature [K]
59  T_{ref} | min or max(T) in cell zone depending on the sign of Q_t [K]
60  \endvartable
61 
62 Usage
63  Example usage:
64  \verbatim
65  effectivenessHeatExchangerSource1
66  {
67  type effectivenessHeatExchangerSource;
68 
69  selectionMode cellZone;
70  cellZone porosity;
71 
72  secondaryMassFlowRate 1.0;
73  secondaryInletT 336;
74  primaryInletT 293;
75 
76  faceZone facesZoneInletOriented;
77 
78  effectiveness <function2>;
79  }
80  \endverbatim
81 
82 Note
83  - The effectiveness Function2 is described in terms of the primary and
84  secondary mass flow rates and has the same units as the secondary mass
85  flow rate and kg/s for phi.
86  - faceZone is the faces at the inlet of the cellzone, it needs to be
87  created with flip map flags. It is used to integrate the net mass flow
88  rate into the heat exchanger
89 
90 SourceFiles
91  effectivenessHeatExchangerSource.C
92 
93 \*---------------------------------------------------------------------------*/
94 
95 #ifndef effectivenessHeatExchangerSource_H
96 #define effectivenessHeatExchangerSource_H
97 
98 #include "fvModel.H"
99 #include "fvCellSet.H"
100 #include "Function2.H"
101 
102 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
103 
104 namespace Foam
105 {
106 namespace fv
107 {
108 
109 /*---------------------------------------------------------------------------*\
110  Class effectivenessHeatExchangerSource Declaration
111 \*---------------------------------------------------------------------------*/
112 
113 class effectivenessHeatExchangerSource
114 :
115  public fvModel
116 {
117  // Private data
118 
119  //- The set of cells the fvConstraint applies to
120  fvCellSet set_;
121 
122  //- Secondary flow mass rate [kg/s]
123  scalar secondaryMassFlowRate_;
124 
125  //- Inlet secondary temperature [K]
126  scalar secondaryInletT_;
127 
128  //- Primary air temperature at the heat exchanger inlet [K]
129  scalar primaryInletT_;
130 
131  //- 2D function for effectiveness [kg/s]
132  // function of primary and secondary mass flow rates [kg/s]
133  autoPtr<Function2<scalar>> eTable_;
134 
135  //- Name of velocity field; default = U
136  word UName_;
137 
138  //- Name of temperature field; default = T
139  word TName_;
140 
141  //- Name of the flux
142  word phiName_;
143 
144  //- Name of the faceZone at the heat exchange inlet
145  word faceZoneName_;
146 
147  //- Id for the face zone
148  label zoneID_;
149 
150  //- Local list of face IDs
151  labelList faceId_;
152 
153  //- Local list of patch ID per face
154  labelList facePatchId_;
155 
156  //- List of +1/-1 representing face flip map (1 use as is, -1 negate)
157  labelList faceSign_;
158 
159  //- Area of the face zone
160  scalar faceZoneArea_;
161 
162 
163 private:
164 
165  // Private Member Functions
166 
167  //- Non-virtual read
168  void readCoeffs();
169 
170  //- Set the zone information
171  void setZone();
172 
173  //- Calculate total area of faceZone across processors
174  void calculateTotalArea(scalar& area) const;
175 
176 
177 public:
178 
179  //- Runtime type information
180  TypeName("effectivenessHeatExchangerSource");
181 
182 
183  // Constructors
184 
185  //- Construct from components
187  (
188  const word& name,
189  const word& modelType,
190  const dictionary& dict,
191  const fvMesh& mesh
192  );
193 
194  //- Disallow default bitwise copy construction
196  (
198  );
199 
200 
201  //- Destructor
203  {}
204 
205 
206  // Member Functions
207 
208  //- Return the list of fields for which the fvModel adds source term
209  // to the transport equation
210  virtual wordList addSupFields() const;
211 
212  //- Explicit and implicit source for compressible equation
213  virtual void addSup
214  (
215  const volScalarField& rho,
216  fvMatrix<scalar>& eqn,
217  const word& fieldName
218  ) const;
219 
220  //- Update for mesh changes
221  virtual void updateMesh(const mapPolyMesh&);
222 
223  //- Read dictionary
224  virtual bool read(const dictionary& dict);
225 
226 
227  // Member Operators
228 
229  //- Disallow default bitwise assignment
230  void operator=(const effectivenessHeatExchangerSource&) = delete;
231 };
232 
233 
234 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
235 
236 } // End namespace fv
237 } // End namespace Foam
238 
239 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240 
241 #endif
242 
243 // ************************************************************************* //
effectivenessHeatExchangerSource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
dictionary dict
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
const word & name() const
Return const access to the source name.
Definition: fvModelI.H:28
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
void operator=(const effectivenessHeatExchangerSource &)=delete
Disallow default bitwise assignment.
virtual bool read(const dictionary &dict)
Read dictionary.
const fvMesh & mesh() const
Return const access to the mesh database.
Definition: fvModelI.H:34
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
List< label > labelList
A List of labels.
Definition: labelList.H:56
TypeName("effectivenessHeatExchangerSource")
Runtime type information.
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
Heat exchanger source model, in which the heat exchanger is defined as a selection of cells...
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual void updateMesh(const mapPolyMesh &)
Update for mesh changes.
virtual void addSup(const volScalarField &rho, fvMatrix< scalar > &eqn, const word &fieldName) const
Explicit and implicit source for compressible equation.
Namespace for OpenFOAM.