Loading...
Searching...
No Matches
Gosman.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2014-2018 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#include "Gosman.H"
29#include "phasePair.H"
30#include "phaseCompressibleTurbulenceModel.H"
32#include "dragModel.H"
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
39{
42 (
44 Gosman,
46 );
47}
48}
49
50
51// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
52
54(
55 const dictionary& dict,
56 const phasePair& pair
57)
58:
60 sigma_("sigma", dimless, dict)
61{}
62
63
64// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
65
67{}
68
69
70// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
71
74{
75 const fvMesh& mesh(pair_.phase1().mesh());
76 const dragModel&
77 drag
78 (
79 mesh.lookupObject<dragModel>
80 (
81 IOobject::groupName(dragModel::typeName, pair_.name())
82 )
83 );
84
85 return
86 0.75
87 *drag.CdRe()
88 *pair_.dispersed()
89 *pair_.continuous().nu()
90 *continuousTurbulence().nut()
91 /(
92 sigma_
93 *sqr(pair_.dispersed().d())
94 )
95 *pair_.continuous().rho();
96}
97
98
99// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition dictionary.H:133
Mesh data needed to do the Finite Volume discretisation.
Definition fvMesh.H:85
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition phasePair.H:52
A class for managing temporary objects.
Definition tmp.H:75
const phasePair & pair_
Phase pair.
const phaseCompressibleTurbulenceModel & continuousTurbulence() const
Return a reference to the turbulence model for the continuous phase.
turbulentDispersionModel(const dictionary &dict, const phasePair &pair)
Construct from a dictionary and a phase pair.
Turbulent dispersion model of Gosman et al.
Definition Gosman.H:63
virtual tmp< volScalarField > D() const
Turbulent diffusivity.
Definition Gosman.C:66
Gosman(const dictionary &dict, const phasePair &pair)
Construct from a dictionary and a phase pair.
Definition Gosman.C:47
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition className.H:142
dynamicFvMesh & mesh
Namespace for OpenFOAM.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dictionary dict
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag