All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
createFields.H
Go to the documentation of this file.
2 
3 Info<< "Reading field D\n" << endl;
5 (
6  IOobject
7  (
8  "D",
9  runTime.timeName(),
10  mesh,
11  IOobject::MUST_READ,
12  IOobject::AUTO_WRITE
13  ),
14  mesh
15 );
16 
17 Info<< "Calculating stress field sigmaD\n" << endl;
18 volSymmTensorField sigmaD
19 (
20  IOobject
21  (
22  "sigmaD",
23  runTime.timeName(),
24  mesh,
25  IOobject::NO_READ,
26  IOobject::NO_WRITE
27  ),
28  mu*twoSymm(fvc::grad(D)) + lambda*(I*tr(fvc::grad(D)))
29 );
30 
31 Info<< "Calculating explicit part of div(sigma) divSigmaExp\n" << endl;
32 volVectorField divSigmaExp
33 (
34  IOobject
35  (
36  "divSigmaExp",
37  runTime.timeName(),
38  mesh,
39  IOobject::NO_READ,
40  IOobject::NO_WRITE
41  ),
42  fvc::div(sigmaD)
43 );
44 
46 {
47  divSigmaExp -= fvc::laplacian(2*mu + lambda, D, "laplacian(DD,D)");
48 }
49 else
50 {
51  divSigmaExp -= fvc::div((2*mu + lambda)*fvc::grad(D), "div(sigmaD)");
52 }
53 
54 mesh.setFluxRequired(D.name());
55 
56 #include "createFvModels.H"
57 #include "createFvConstraints.H"
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:61
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Switch compactNormalStress(stressControl.lookup("compactNormalStress"))
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
volScalarField mu(IOobject("mu", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), laminarTransport.nu() *rhoInfValue)
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
dynamicFvMesh & mesh
static const Identity< scalar > I
Definition: Identity.H:93
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
messageStream Info