motionDirectionalDiffusivity.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) 2011-2012 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 
27 #include "surfaceInterpolate.H"
30 #include "velocityMotionSolver.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(motionDirectionalDiffusivity, 0);
37 
39  (
40  motionDiffusivity,
41  motionDirectionalDiffusivity,
42  Istream
43  );
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
49 Foam::motionDirectionalDiffusivity::motionDirectionalDiffusivity
50 (
51  const fvMesh& mesh,
52  Istream& mdData
53 )
54 :
55  uniformDiffusivity(mesh, mdData),
56  diffusivityVector_(mdData)
57 {}
58 
59 
60 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
61 
63 {}
64 
65 
66 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
67 
69 {
70  static bool first = true;
71 
72  if (!first)
73  {
74  const volVectorField& cellMotionU =
75  mesh().lookupObject<volVectorField>("cellMotionU");
76 
78  (
79  IOobject
80  (
81  "D",
82  mesh().time().timeName(),
83  mesh()
84  ),
85  diffusivityVector_.y()*vector::one
86  + (diffusivityVector_.x() - diffusivityVector_.y())*cellMotionU
87  /(mag(cellMotionU) + dimensionedScalar("small", dimVelocity, SMALL)),
88  zeroGradientFvPatchVectorField::typeName
89  );
90  D.correctBoundaryConditions();
91 
92  const surfaceVectorField n(mesh().Sf()/mesh().magSf());
93  faceDiffusivity_ == (n & cmptMultiply(fvc::interpolate(D), n));
94  }
95  else
96  {
97  first = false;
98 
99  const velocityMotionSolver& mSolver =
100  mesh().lookupObject<velocityMotionSolver>("dynamicMeshDict");
101 
102  const_cast<velocityMotionSolver&>(mSolver).solve();
103  correct();
104  }
105 }
106 
107 
108 // ************************************************************************* //
Uniform uniform finite volume mesh motion diffusivity.
Virtual base class for velocity motion solver.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
virtual void correct()
Correct the motion diffusivity.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
dynamicFvMesh & mesh
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
word timeName
Definition: getTimeIndex.H:3
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/ubuntu/OpenFOAM-5.0/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:72
defineTypeNameAndDebug(combustionModel, 0)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
const dimensionSet dimVelocity