readMechanicalProperties.H
Go to the documentation of this file.
1  Info<< "Reading mechanical properties\n" << endl;
2 
3  IOdictionary mechanicalProperties
4  (
5  IOobject
6  (
7  "mechanicalProperties",
8  runTime.constant(),
9  mesh,
10  IOobject::MUST_READ_IF_MODIFIED,
11  IOobject::NO_WRITE
12  )
13  );
14 
15  const dictionary& rhoDict(mechanicalProperties.subDict("rho"));
16  word rhoType(rhoDict.lookup("type"));
17 
18  autoPtr<volScalarField> rhoPtr;
19 
20  IOobject rhoIO
21  (
22  "rho",
23  runTime.timeName(0),
24  mesh,
25  IOobject::NO_READ,
26  IOobject::NO_WRITE
27  );
28 
29  if (rhoType == "uniform")
30  {
31  scalar rhoValue(readScalar(rhoDict.lookup("value")));
32 
33  rhoPtr.reset
34  (
35  new volScalarField
36  (
37  rhoIO,
38  mesh,
40  (
42  rhoValue
43  )
44  )
45  );
46  }
47  else if (rhoType == "field")
48  {
49  rhoIO.readOpt() = IOobject::MUST_READ;
50 
51  rhoPtr.reset
52  (
53  new volScalarField
54  (
55  rhoIO,
56  mesh
57  )
58  );
59  }
60  else
61  {
63  << "Valid type entries are uniform or field for rho"
64  << abort(FatalError);
65  }
66 
68 
69  const dictionary& EDict(mechanicalProperties.subDict("E"));
70  word EType(EDict.lookup("type"));
71 
72  autoPtr<volScalarField> EPtr;
73 
74  IOobject EHeader
75  (
76  "E",
77  runTime.timeName(0),
78  mesh,
79  IOobject::NO_READ,
80  IOobject::NO_WRITE
81  );
82 
83  if (EType == "uniform")
84  {
85  scalar rhoEValue(readScalar(EDict.lookup("value")));
86 
87  EPtr.reset
88  (
89  new volScalarField
90  (
91  EHeader,
92  mesh,
94  (
96  rhoEValue
97  )
98  )
99  );
100  }
101  else if (EType == "field")
102  {
103  EHeader.readOpt() = IOobject::MUST_READ;
104 
105  EPtr.reset
106  (
107  new volScalarField
108  (
109  EHeader,
110  mesh
111  )
112  );
113  }
114  else
115  {
117  << "Valid type entries are uniform or field for E"
118  << abort(FatalError);
119  }
120 
121  volScalarField& rhoE = EPtr();
122 
123  autoPtr<volScalarField> nuPtr;
124 
125  IOobject nuIO
126  (
127  "nu",
128  runTime.timeName(0),
129  mesh,
130  IOobject::NO_READ,
131  IOobject::NO_WRITE
132  );
133 
134  const dictionary& nuDict(mechanicalProperties.subDict("nu"));
135  word nuType(nuDict.lookup("type"));
136 
137  if (nuType == "uniform")
138  {
139  scalar nuValue(readScalar(nuDict.lookup("value")));
140  nuPtr.reset
141  (
142  new volScalarField
143  (
144  nuIO,
145  mesh,
147  (
148  dimless,
149  nuValue
150  )
151  )
152  );
153  }
154  else if (nuType == "field")
155  {
156  nuIO.readOpt() = IOobject::MUST_READ;
157  nuPtr.reset
158  (
159  new volScalarField
160  (
161  nuIO,
162  mesh
163  )
164  );
165  }
166  else
167  {
169  << "Valid type entries are uniform or field for nu"
170  << abort(FatalError);
171  }
172 
174 
175  Info<< "Normalising E : E/rho\n" << endl;
176  volScalarField E(rhoE/rho);
177 
178  Info<< "Calculating Lame's coefficients\n" << endl;
179 
180  volScalarField mu(E/(2.0*(1.0 + nu)));
181  volScalarField lambda(nu*E/((1.0 + nu)*(1.0 - 2.0*nu)));
182  volScalarField threeK(E/(1.0 - 2.0*nu));
183 
184  Switch planeStress(mechanicalProperties.lookup("planeStress"));
185 
186  if (planeStress)
187  {
188  Info<< "Plane Stress\n" << endl;
189 
190  lambda = nu*E/((1.0 + nu)*(1.0 - nu));
191  threeK = E/(1.0 - nu);
192  }
193  else
194  {
195  Info<< "Plane Strain\n" << endl;
196  }
#define readScalar
Definition: doubleScalar.C:38
IOobject EHeader("E", runTime.timeName(0), mesh, IOobject::NO_READ, IOobject::NO_WRITE)
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
autoPtr< volScalarField > EPtr
Info<< "Reading mechanical properties\"<< endl;IOdictionary mechanicalProperties(IOobject("mechanicalProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));const dictionary &rhoDict(mechanicalProperties.subDict("rho"));word rhoType(rhoDict.lookup("type"));autoPtr< volScalarField > rhoPtr
IOobject nuIO("nu", runTime.timeName(0), mesh, IOobject::NO_READ, IOobject::NO_WRITE)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
dynamicFvMesh & mesh
autoPtr< volScalarField > nuPtr
word nuType(nuDict.lookup("type"))
const dictionary & EDict(mechanicalProperties.subDict("E"))
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const dictionary & nuDict(mechanicalProperties.subDict("nu"))
const dimensionedScalar mu
Atomic mass unit.
word EType(EDict.lookup("type"))
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
IOobject rhoIO("rho", runTime.timeName(0), mesh, IOobject::NO_READ, IOobject::NO_WRITE)
volScalarField & nu