binaryNode.C
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) 2016-2022 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 \*---------------------------------------------------------------------------*/
25 
26 #include "binaryNode.H"
27 #include "ISAT.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
32 :
33  leafLeft_(nullptr),
34  leafRight_(nullptr),
35  nodeLeft_(nullptr),
36  nodeRight_(nullptr),
37  parent_(nullptr)
38 {}
39 
40 
42 (
43  chemPointISAT* elementLeft,
44  chemPointISAT* elementRight,
45  binaryNode* parent
46 )
47 :
48  leafLeft_(elementLeft),
49  leafRight_(elementRight),
50  nodeLeft_(nullptr),
51  nodeRight_(nullptr),
52  parent_(parent),
53  v_(elementLeft->completeSpaceSize(), 0)
54 {
55  calcV(*elementLeft, *elementRight, v_);
56  a_ = calcA(*elementLeft, *elementRight);
57 }
58 
59 
60 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
61 
63 (
64  const chemPointISAT& elementLeft,
65  const chemPointISAT& elementRight,
66  scalarField& v
67 )
68 {
69  // LT is the transpose of the L matrix
70  const scalarSquareMatrix& LT = elementLeft.LT();
71  const bool reduction = elementLeft.table().reduction();
72 
73  // Difference of composition in the full species domain
74  const scalarField& phil(elementLeft.phi());
75  const scalarField& phir(elementRight.phi());
76  const scalarField& scaleFactor(elementLeft.scaleFactor());
77  const scalar epsTol = elementLeft.tolerance();
78 
79  // v = LT.T()*LT*(phir - phil)
80  for (label i=0; i<elementLeft.completeSpaceSize(); i++)
81  {
82  label si = i;
83  bool outOfIndexI = true;
84  if (reduction)
85  {
86  if (i<elementLeft.completeSpaceSize() - 3)
87  {
88  si = elementLeft.completeToSimplifiedIndex()[i];
89  outOfIndexI = (si == -1);
90  }
91  else // temperature and pressure
92  {
93  outOfIndexI = false;
94  const label dif = i - (elementLeft.completeSpaceSize() - 3);
95  si = elementLeft.nActive() + dif;
96  }
97  }
98  if (!reduction || (reduction && !(outOfIndexI)))
99  {
100  v[i] = 0;
101  for (label j=0; j<elementLeft.completeSpaceSize(); j++)
102  {
103  label sj = j;
104  bool outOfIndexJ = true;
105 
106  if (reduction)
107  {
108  if (j < elementLeft.completeSpaceSize() - 3)
109  {
110  sj = elementLeft.completeToSimplifiedIndex()[j];
111  outOfIndexJ = (sj==-1);
112  }
113  else
114  {
115  outOfIndexJ = false;
116  const label dif =
117  j - (elementLeft.completeSpaceSize() - 3);
118  sj = elementLeft.nActive() + dif;
119  }
120  }
121 
122  if (!reduction ||(reduction && !(outOfIndexJ)))
123  {
124  // Since L is a lower triangular matrix k=0->min(i, j)
125  for (label k=0; k<=min(si, sj); k++)
126  {
127  v[i] += LT(k, si)*LT(k, sj)*(phir[j] - phil[j]);
128  }
129  }
130  }
131  }
132  else
133  {
134  // When it is an inactive species the diagonal element of LT is
135  // 1/(scaleFactor*epsTol)
136  v[i] = (phir[i] - phil[i])/sqr(scaleFactor[i]*epsTol);
137  }
138  }
139 }
140 
141 
142 // ************************************************************************* //
void calcV(const chemPointISAT &elementLeft, const chemPointISAT &elementRight, scalarField &v)
Compute vector v:
Definition: binaryNode.C:63
binaryNode * parent_
Parent node.
Definition: binaryNode.H:66
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Leaf of the binary tree. The chemPoint stores the composition &#39;phi&#39;, the mapping of this composition ...
chemPointISAT * leafRight_
Element on the right.
Definition: binaryNode.H:57
surfaceScalarField phir(fvc::flux(mixture.Udm()))
binaryNode * nodeLeft_
Node which follows on the left.
Definition: binaryNode.H:60
const scalar & tolerance() const
label nActive() const
label k
Boltzmann constant.
scalarField v_
Definition: binaryNode.H:68
binaryNode()
Construct null.
Definition: binaryNode.C:31
const List< label > & completeToSimplifiedIndex() const
const scalarField & phi() const
const chemistryTabulationMethods::ISAT & table() const
Access to the ISAT table.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
Node of the binary tree.
Definition: binaryNode.H:48
binaryNode * nodeRight_
Node which follows on the right.
Definition: binaryNode.H:63
scalar calcA(const chemPointISAT &elementLeft, const chemPointISAT &elementRight)
Compute a the product v^T.phih, with phih = (phi0 + phiq)/2.
Definition: binaryNodeI.H:30
const scalarSquareMatrix & LT() const
const scalarField & scaleFactor() const
bool reduction() const
Return true if reduction is applied to the state variables.
Definition: ISAT.H:241
label & completeSpaceSize()
chemPointISAT * leafLeft_
Element on the left.
Definition: binaryNode.H:54