createFields.H
Go to the documentation of this file.
3 
4 Info<< "Reading field D\n" << endl;
6 (
7  IOobject
8  (
9  "D",
10  runTime.timeName(),
11  mesh,
12  IOobject::MUST_READ,
13  IOobject::AUTO_WRITE
14  ),
15  mesh
16 );
17 
18 
19 autoPtr<volScalarField> Tptr(NULL);
20 
21 if (thermalStress)
22 {
23  Info<< "Reading field T\n" << endl;
24  Tptr.reset
25  (
26  new volScalarField
27  (
28  IOobject
29  (
30  "T",
31  runTime.timeName(),
32  mesh,
33  IOobject::MUST_READ,
34  IOobject::AUTO_WRITE
35  ),
36  mesh
37  )
38  );
39 }
40 
41 
42 Info<< "Calculating stress field sigmaD\n" << endl;
43 volSymmTensorField sigmaD
44 (
45  IOobject
46  (
47  "sigmaD",
48  runTime.timeName(),
49  mesh,
50  IOobject::NO_READ,
51  IOobject::NO_WRITE
52  ),
53  mu*twoSymm(fvc::grad(D)) + lambda*(I*tr(fvc::grad(D)))
54 );
55 
56 Info<< "Calculating explicit part of div(sigma) divSigmaExp\n" << endl;
57 volVectorField divSigmaExp
58 (
59  IOobject
60  (
61  "divSigmaExp",
62  runTime.timeName(),
63  mesh,
64  IOobject::NO_READ,
65  IOobject::NO_WRITE
66  ),
67  fvc::div(sigmaD)
68 );
69 
71 {
72  divSigmaExp -= fvc::laplacian(2*mu + lambda, D, "laplacian(DD,D)");
73 }
74 else
75 {
76  divSigmaExp -= fvc::div((2*mu + lambda)*fvc::grad(D), "div(sigmaD)");
77 }
78 
79 mesh.setFluxRequired(D.name());
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Switch compactNormalStress(stressControl.lookup("compactNormalStress"))
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
volScalarField mu(IOobject("mu", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), laminarTransport.nu()*rhoInfValue)
dynamicFvMesh & mesh
static const Identity< scalar > I
Definition: Identity.H:93
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
Info<< "Reading field D\n"<< endl;volVectorField D(IOobject("D", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);autoPtr< volScalarField > Tptr(NULL)
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
messageStream Info