advectionDiffusionPatchDistMethod.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) 2015-2017 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"
28 #include "fvcGrad.H"
29 #include "fvcDiv.H"
30 #include "fvmDiv.H"
31 #include "fvmLaplacian.H"
32 #include "fvmSup.H"
34 
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
42 namespace patchDistMethods
43 {
46 }
47 }
48 
49 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
50 
51 Foam::patchDistMethods::advectionDiffusion::advectionDiffusion
52 (
53  const dictionary& dict,
54  const fvMesh& mesh,
55  const labelHashSet& patchIDs
56 )
57 :
58  patchDistMethod(mesh, patchIDs),
59  coeffs_(dict.optionalSubDict(type() + "Coeffs")),
60  pdmPredictor_
61  (
63  (
64  coeffs_,
65  mesh,
66  patchIDs
67  )
68  ),
69  epsilon_(coeffs_.lookupOrDefault<scalar>("epsilon", 0.1)),
70  tolerance_(coeffs_.lookupOrDefault<scalar>("tolerance", 1e-3)),
71  maxIter_(coeffs_.lookupOrDefault<int>("maxIter", 10)),
72  predicted_(false)
73 {}
74 
75 
76 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
77 
79 {
80  return correct(y, const_cast<volVectorField&>(volVectorField::null()));
81 }
82 
83 
85 (
86  volScalarField& y,
88 )
89 {
90  if (!predicted_)
91  {
92  pdmPredictor_->correct(y);
93  predicted_ = true;
94  }
95 
97  (
98  IOobject
99  (
100  "ny",
101  mesh_.time().timeName(),
102  mesh_,
105  false
106  ),
107  mesh_,
109  patchTypes<vector>(mesh_, patchIDs_)
110  );
111 
112  const fvPatchList& patches = mesh_.boundary();
113  volVectorField::Boundary& nybf = ny.boundaryFieldRef();
114 
115  forAllConstIter(labelHashSet, patchIDs_, iter)
116  {
117  label patchi = iter.key();
118  nybf[patchi] == -patches[patchi].nf();
119  }
120 
121  int iter = 0;
122  scalar initialResidual = 0;
123 
124  do
125  {
126  ny = fvc::grad(y);
127  ny /= (mag(ny) + SMALL);
128 
130  nf /= (mag(nf) + SMALL);
131 
132  surfaceScalarField yPhi("yPhi", mesh_.Sf() & nf);
133 
134  fvScalarMatrix yEqn
135  (
136  fvm::div(yPhi, y)
137  - fvm::Sp(fvc::div(yPhi), y)
138  - epsilon_*y*fvm::laplacian(y)
139  ==
140  dimensionedScalar("1", dimless, 1.0)
141  );
142 
143  yEqn.relax();
144  initialResidual = yEqn.solve().initialResidual();
145 
146  } while (initialResidual > tolerance_ && ++iter < maxIter_);
147 
148  // Only calculate n if the field is defined
149  if (notNull(n))
150  {
151  n = -ny;
152  }
153 
154  return true;
155 }
156 
157 
158 // ************************************************************************* //
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
static autoPtr< patchDistMethod > New(const dictionary &dict, const fvMesh &mesh, const labelHashSet &patchIDs)
const double e
Elementary charge.
Definition: doubleFloat.H:78
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
addToRunTimeSelectionTable(patchDistMethod, advectionDiffusion, dictionary)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Calculate the matrix for the laplacian of the field.
patches[0]
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Macros for easy insertion into run-time selection tables.
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:759
Calculate the gradient of the given field.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
static const zero Zero
Definition: zero.H:91
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Calculate the divergence of the given field.
virtual bool correct(volScalarField &y)
Correct the given distance-to-patch field.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
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
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
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.
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:46
static const GeometricField< vector, fvPatchField, volMesh > & null()
Return a null geometric field.
Calculate the matrix for the divergence of the given field and flux.
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
defineTypeNameAndDebug(advectionDiffusion, 0)
dimensioned< scalar > mag(const dimensioned< Type > &)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Specialisation of patchDist for wall distance calculation.
Calculation of approximate distance to nearest patch for all cells and boundary by solving the Eikona...
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.