Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
transportModels
twoPhaseInterfaceProperties
alphaContactAngle
dynamicAlphaContactAngle
dynamicAlphaContactAngleFvPatchScalarField.C
Go to the documentation of this file.
1
/*---------------------------------------------------------------------------*\
2
========= |
3
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4
\\ / O peration |
5
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
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 "
dynamicAlphaContactAngleFvPatchScalarField.H
"
27
#include <
OpenFOAM/addToRunTimeSelectionTable.H
>
28
#include <
finiteVolume/fvPatchFieldMapper.H
>
29
#include <
finiteVolume/volMesh.H
>
30
31
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32
33
Foam::dynamicAlphaContactAngleFvPatchScalarField::
34
dynamicAlphaContactAngleFvPatchScalarField
35
(
36
const
fvPatch
&
p
,
37
const
DimensionedField<scalar, volMesh>
& iF
38
)
39
:
40
alphaContactAngleFvPatchScalarField
(p, iF),
41
theta0_(0.0),
42
uTheta_(0.0),
43
thetaA_(0.0),
44
thetaR_(0.0)
45
{}
46
47
48
Foam::dynamicAlphaContactAngleFvPatchScalarField::
49
dynamicAlphaContactAngleFvPatchScalarField
50
(
51
const
dynamicAlphaContactAngleFvPatchScalarField
& gcpsf,
52
const
fvPatch
& p,
53
const
DimensionedField<scalar, volMesh>
& iF,
54
const
fvPatchFieldMapper
& mapper
55
)
56
:
57
alphaContactAngleFvPatchScalarField
(gcpsf, p, iF, mapper),
58
theta0_(gcpsf.theta0_),
59
uTheta_(gcpsf.uTheta_),
60
thetaA_(gcpsf.thetaA_),
61
thetaR_(gcpsf.thetaR_)
62
{}
63
64
65
Foam::dynamicAlphaContactAngleFvPatchScalarField::
66
dynamicAlphaContactAngleFvPatchScalarField
67
(
68
const
fvPatch
& p,
69
const
DimensionedField<scalar, volMesh>
& iF,
70
const
dictionary
& dict
71
)
72
:
73
alphaContactAngleFvPatchScalarField
(p, iF, dict),
74
theta0_(
readScalar
(dict.
lookup
(
"theta0"
))),
75
uTheta_(
readScalar
(dict.
lookup
(
"uTheta"
))),
76
thetaA_(
readScalar
(dict.
lookup
(
"thetaA"
))),
77
thetaR_(
readScalar
(dict.
lookup
(
"thetaR"
)))
78
{
79
evaluate();
80
}
81
82
83
Foam::dynamicAlphaContactAngleFvPatchScalarField::
84
dynamicAlphaContactAngleFvPatchScalarField
85
(
86
const
dynamicAlphaContactAngleFvPatchScalarField
& gcpsf
87
)
88
:
89
alphaContactAngleFvPatchScalarField
(gcpsf),
90
theta0_(gcpsf.theta0_),
91
uTheta_(gcpsf.uTheta_),
92
thetaA_(gcpsf.thetaA_),
93
thetaR_(gcpsf.thetaR_)
94
{}
95
96
97
Foam::dynamicAlphaContactAngleFvPatchScalarField::
98
dynamicAlphaContactAngleFvPatchScalarField
99
(
100
const
dynamicAlphaContactAngleFvPatchScalarField
& gcpsf,
101
const
DimensionedField<scalar, volMesh>
& iF
102
)
103
:
104
alphaContactAngleFvPatchScalarField
(gcpsf, iF),
105
theta0_(gcpsf.theta0_),
106
uTheta_(gcpsf.uTheta_),
107
thetaA_(gcpsf.thetaA_),
108
thetaR_(gcpsf.thetaR_)
109
{}
110
111
112
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113
114
Foam::tmp<Foam::scalarField>
115
Foam::dynamicAlphaContactAngleFvPatchScalarField::theta
116
(
117
const
fvPatchVectorField
& Up,
118
const
fvsPatchVectorField
& nHat
119
)
const
120
{
121
if
(uTheta_ < SMALL)
122
{
123
return
tmp<scalarField>
(
new
scalarField
(size(), theta0_));
124
}
125
126
vectorField
nf = patch().nf();
127
128
// Calculated the component of the velocity parallel to the wall
129
vectorField
Uwall = Up.
patchInternalField
() - Up;
130
Uwall -= (nf & Uwall)*nf;
131
132
// Find the direction of the interface parallel to the wall
133
vectorField
nWall = nHat - (nf & nHat)*nf;
134
135
// Normalise nWall
136
nWall /= (
mag
(nWall) + SMALL);
137
138
// Calculate Uwall resolved normal to the interface parallel to
139
// the interface
140
scalarField
uwall = nWall & Uwall;
141
142
return
theta0_ + (thetaA_ - thetaR_)*
tanh
(uwall/uTheta_);
143
}
144
145
146
void
Foam::dynamicAlphaContactAngleFvPatchScalarField::write
(
Ostream
& os)
const
147
{
148
alphaContactAngleFvPatchScalarField::write
(os);
149
os.
writeKeyword
(
"theta0"
) << theta0_ <<
token::END_STATEMENT
<<
nl
;
150
os.
writeKeyword
(
"uTheta"
) << uTheta_ <<
token::END_STATEMENT
<<
nl
;
151
os.
writeKeyword
(
"thetaA"
) << thetaA_ <<
token::END_STATEMENT
<<
nl
;
152
os.
writeKeyword
(
"thetaR"
) << thetaR_ <<
token::END_STATEMENT
<<
nl
;
153
writeEntry(
"value"
, os);
154
}
155
156
157
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
158
159
namespace
Foam
160
{
161
makePatchTypeField
162
(
163
fvPatchScalarField
,
164
dynamicAlphaContactAngleFvPatchScalarField
165
);
166
}
167
168
169
// ************************ vim: set sw=4 sts=4 et: ************************ //