All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
linearEqnI.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) 2017-2018 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
29 {}
30 
31 
33 :
34  VectorSpace<linearEqn, scalar, 2>(Foam::zero())
35 {}
36 
37 
38 inline Foam::linearEqn::linearEqn(const scalar a, const scalar b)
39 {
40  this->v_[coefficient::a] = a;
41  this->v_[coefficient::b] = b;
42 }
43 
44 
45 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
46 
47 inline Foam::scalar Foam::linearEqn::a() const
48 {
49  return this->v_[coefficient::a];
50 }
51 
52 
53 inline Foam::scalar Foam::linearEqn::b() const
54 {
55  return this->v_[coefficient::b];
56 }
57 
58 
59 inline Foam::scalar& Foam::linearEqn::a()
60 {
61  return this->v_[coefficient::a];
62 }
63 
64 
65 inline Foam::scalar& Foam::linearEqn::b()
66 {
67  return this->v_[coefficient::b];
68 }
69 
70 
71 inline Foam::scalar Foam::linearEqn::value(const scalar x) const
72 {
73  return x*a() + b();
74 }
75 
76 
77 inline Foam::scalar Foam::linearEqn::derivative(const scalar x) const
78 {
79  return a();
80 }
81 
82 
83 inline Foam::scalar Foam::linearEqn::error(const scalar x) const
84 {
85  return small*(mag(x*a()) + mag(b()));
86 }
87 
88 
90 {
91  /*
92 
93  This function solves a linear equation of the following form:
94 
95  a*x + b = 0
96  x + b = 0
97 
98  */
99 
100  const scalar a = this->a();
101  const scalar b = this->b();
102 
103  if (a == 0)
104  {
105  return Roots<1>(rootType::nan, 0);
106  }
107 
108  if (mag(b/vGreat) >= mag(a))
109  {
110  return Roots<1>
111  (
112  sign(a) == sign(b)
115  0
116  );
117  }
118 
119  return Roots<1>(rootType::real, -b/a);
120 }
121 
122 
123 // ************************************************************************* //
dimensionedScalar sign(const dimensionedScalar &ds)
scalar b() const
Definition: linearEqnI.H:53
Templated vector space.
Definition: VectorSpace.H:53
Templated storage for the roots of polynomial equations, plus flags to indicate the nature of the roo...
Definition: Roots.H:64
scalar a() const
Definition: linearEqnI.H:47
linearEqn()
Construct null.
Definition: linearEqnI.H:28
scalar value(const scalar x) const
Evaluate at x.
Definition: linearEqnI.H:71
scalar error(const scalar x) const
Estimate the error of evaluation at x.
Definition: linearEqnI.H:83
Linear equation of the form a*x + b = 0.
Definition: linearEqn.H:48
Roots< 1 > roots() const
Get the roots.
Definition: linearEqnI.H:89
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:49
dimensioned< scalar > mag(const dimensioned< Type > &)
scalar derivative(const scalar x) const
Evaluate the derivative at x.
Definition: linearEqnI.H:77
scalar v_[Ncmpts]
The components of this vector space.
Definition: VectorSpace.H:84
Namespace for OpenFOAM.