quadraticEqnI.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29
31{}
32
33
35:
36 VectorSpace<quadraticEqn, scalar, 3>(Foam::zero{})
37{}
38
39
41(
42 const scalar a,
43 const scalar b,
44 const scalar c
45)
46{
47 this->v_[A] = a;
48 this->v_[B] = b;
49 this->v_[C] = c;
50}
51
52
53// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
54
55inline Foam::scalar Foam::quadraticEqn::a() const
56{
57 return this->v_[A];
58}
59
60
61inline Foam::scalar Foam::quadraticEqn::b() const
62{
63 return this->v_[B];
64}
65
66
67inline Foam::scalar Foam::quadraticEqn::c() const
68{
69 return this->v_[C];
70}
71
72
73inline Foam::scalar& Foam::quadraticEqn::a()
74{
75 return this->v_[A];
76}
77
78
79inline Foam::scalar& Foam::quadraticEqn::b()
80{
81 return this->v_[B];
82}
83
84
85inline Foam::scalar& Foam::quadraticEqn::c()
86{
87 return this->v_[C];
88}
89
90
91inline Foam::scalar Foam::quadraticEqn::value(const scalar x) const
92{
93 return x*(x*a() + b()) + c();
94}
95
96
97inline Foam::scalar Foam::quadraticEqn::derivative(const scalar x) const
98{
99 return x*2*a() + b();
100}
101
102
103inline Foam::scalar Foam::quadraticEqn::error(const scalar x) const
104{
105 return
106 SMALL*mag(x)*(mag(x*a()) + mag(b()))
107 + SMALL*(mag(x*(x*a() + b())) + mag(c()));
108}
109
110
111// ************************************************************************* //
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
Graphite solid properties.
Definition: C.H:53
Templated vector space.
Definition: VectorSpace.H:79
Class to handle errors and exceptions in a simple, consistent stream-based manner.
Definition: error.H:77
Container to encapsulate various operations for quadratic equation of the forms with real coefficient...
Definition: quadraticEqn.H:91
scalar value(const scalar x) const
Evaluate the quadratic equation at x.
Definition: quadraticEqnI.H:91
scalar c() const
Definition: quadraticEqnI.H:67
scalar derivative(const scalar x) const
Evaluate the derivative of the quadratic equation at x.
Definition: quadraticEqnI.H:97
scalar a() const
Definition: quadraticEqnI.H:55
scalar b() const
Definition: quadraticEqnI.H:61
quadraticEqn()
Construct null.
Definition: quadraticEqnI.H:30
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Definition: zero.H:63
Namespace for OpenFOAM.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
volScalarField & b
Definition: createFields.H:27