dnsFoam.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-2016 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 Application
25  dnsFoam
26 
27 Description
28  Direct numerical simulation solver for boxes of isotropic turbulence.
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "fvCFD.H"
33 #include "Kmesh.H"
34 #include "UOprocess.H"
35 #include "fft.H"
36 #include "calcEk.H"
37 #include "graph.H"
38 #include "pisoControl.H"
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 int main(int argc, char *argv[])
43 {
44  #include "postProcess.H"
45 
46  #include "setRootCase.H"
47  #include "createTime.H"
48  #include "createMeshNoClear.H"
49  #include "createControl.H"
50  #include "createFields.H"
51  #include "initContinuityErrs.H"
52 
53  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
54 
55  Info<< nl << "Starting time loop" << endl;
56 
57  while (runTime.loop())
58  {
59  Info<< "Time = " << runTime.timeName() << nl << endl;
60 
61  force.primitiveFieldRef() = ReImSum
62  (
63  fft::reverseTransform
64  (
65  K/(mag(K) + 1.0e-6) ^ forceGen.newField(), K.nn()
66  )
67  );
68 
69  #include "globalProperties.H"
70 
72  (
73  fvm::ddt(U)
74  + fvm::div(phi, U)
75  - fvm::laplacian(nu, U)
76  ==
77  force
78  );
79 
80  solve(UEqn == -fvc::grad(p));
81 
82 
83  // --- PISO loop
84  while (piso.correct())
85  {
86  volScalarField rAU(1.0/UEqn.A());
90  (
91  "phiHbyA",
94  );
95 
96  // Update the pressure BCs to ensure flux consistency
98 
99  fvScalarMatrix pEqn
100  (
102  );
103 
104  pEqn.solve(mesh.solver(p.select(piso.finalInnerIter())));
105 
106  phi = phiHbyA - pEqn.flux();
107 
108  #include "continuityErrs.H"
109 
110  U = HbyA - rAU*fvc::grad(p);
111  U.correctBoundaryConditions();
112  }
113 
114  runTime.write();
115 
116  if (runTime.writeTime())
117  {
118  calcEk(U, K).write
119  (
120  runTime.path()/"graphs"/runTime.timeName(),
121  "Ek",
122  runTime.graphFormat()
123  );
124  }
125 
126  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
127  << " ClockTime = " << runTime.elapsedClockTime() << " s"
128  << nl << endl;
129  }
130 
131  Info<< "End\n" << endl;
132 
133  return 0;
134 }
135 
136 
137 // ************************************************************************* //
surfaceScalarField & phi
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
U
Definition: pEqn.H:83
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
void constrainPressure(volScalarField &p, const RhoType &rho, const volVectorField &U, const surfaceScalarField &phiHbyA, const RAUType &rhorAU, const MRFType &MRF)
tmp< GeometricField< typename flux< Type >::type, fvsPatchField, surfaceMesh > > ddtCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: fvcDdt.C:155
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
pisoControl piso(mesh)
tmp< surfaceScalarField > interpolate(const RhoType &rho)
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
graph calcEk(const volVectorField &U, const Kmesh &K)
Definition: calcEk.C:27
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
CGAL::Exact_predicates_exact_constructions_kernel K
HbyA
Definition: pcEqn.H:74
scalarField ReImSum(const UList< complex > &cf)
Definition: complexFields.C:84
phiHbyA
Definition: pcEqn.H:73
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const dimensionedScalar e
Elementary charge.
Definition: doubleFloat.H:78
dynamicFvMesh & mesh
rhoEqn solve()
tmp< volVectorField > constrainHbyA(const tmp< volVectorField > &tHbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:33
static const char nl
Definition: Ostream.H:262
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU))
fvVectorMatrix & UEqn
Definition: UEqn.H:13
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
void write(Ostream &, const word &format) const
Write graph to stream in given format.
Definition: graph.C:221
volScalarField & p
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Execute application functionObjects to post-process existing results.
volScalarField & nu
volScalarField rAU(1.0/UEqn.A())