FreeFOAM The Cross-Platform CFD Toolkit
stressComponents.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) 1991-2010 OpenCFD Ltd.
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  stressComponents
26 
27 Description
28  Calculates and writes the scalar fields of the six components of the stress
29  tensor sigma for each time.
30 
31 Usage
32 
33  - stressComponents [OPTIONS]
34 
35  @param -noZero \n
36  Ignore timestep 0.
37 
38  @param -constant \n
39  Include the constant directory.
40 
41  @param -time <time>\n
42  Apply only to specific time.
43 
44  @param -latestTime \n
45  Only apply to latest time step.
46 
47  @param -case <dir>\n
48  Case directory.
49 
50  @param -parallel \n
51  Run in parallel.
52 
53  @param -help \n
54  Display help message.
55 
56  @param -doc \n
57  Display Doxygen API documentation page for this application.
58 
59  @param -srcDoc \n
60  Display Doxygen source documentation page for this application.
61 
62 \*---------------------------------------------------------------------------*/
63 
64 #include <finiteVolume/fvCFD.H>
67 
68 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
69 
70 int main(int argc, char *argv[])
71 {
72  timeSelector::addOptions();
73 
74 # include <OpenFOAM/setRootCase.H>
75 # include <OpenFOAM/createTime.H>
76 
77  instantList timeDirs = timeSelector::select0(runTime, args);
78 
79 # include <OpenFOAM/createMesh.H>
80 
81  forAll(timeDirs, timeI)
82  {
83  runTime.setTime(timeDirs[timeI], timeI);
84 
85  Info<< "Time = " << runTime.timeName() << endl;
86 
87  IOobject Uheader
88  (
89  "U",
90  runTime.timeName(),
91  mesh,
92  IOobject::MUST_READ
93  );
94 
95  // Check U exists
96  if (Uheader.headerOk())
97  {
98  mesh.readUpdate();
99 
100  Info<< " Reading U" << endl;
101  volVectorField U(Uheader, mesh);
102 
103 # include <finiteVolume/createPhi.H>
104 
106 
107  volSymmTensorField sigma
108  (
109  IOobject
110  (
111  "sigma",
112  runTime.timeName(),
113  mesh,
114  IOobject::NO_READ,
115  IOobject::AUTO_WRITE
116  ),
118  );
119 
120 
121  volScalarField sigmaxx
122  (
123  IOobject
124  (
125  "sigmaxx",
126  runTime.timeName(),
127  mesh,
128  IOobject::NO_READ
129  ),
130  sigma.component(symmTensor::XX)
131  );
132  sigmaxx.write();
133 
134  volScalarField sigmayy
135  (
136  IOobject
137  (
138  "sigmayy",
139  runTime.timeName(),
140  mesh,
141  IOobject::NO_READ
142  ),
143  sigma.component(symmTensor::YY)
144  );
145  sigmayy.write();
146 
147  volScalarField sigmazz
148  (
149  IOobject
150  (
151  "sigmazz",
152  runTime.timeName(),
153  mesh,
154  IOobject::NO_READ
155  ),
156  sigma.component(symmTensor::ZZ)
157  );
158  sigmazz.write();
159 
160  volScalarField sigmaxy
161  (
162  IOobject
163  (
164  "sigmaxy",
165  runTime.timeName(),
166  mesh,
167  IOobject::NO_READ
168  ),
169  sigma.component(symmTensor::XY)
170  );
171  sigmaxy.write();
172 
173  volScalarField sigmaxz
174  (
175  IOobject
176  (
177  "sigmaxz",
178  runTime.timeName(),
179  mesh,
180  IOobject::NO_READ
181  ),
182  sigma.component(symmTensor::XZ)
183  );
184  sigmaxz.write();
185 
186  volScalarField sigmayz
187  (
188  IOobject
189  (
190  "sigmayz",
191  runTime.timeName(),
192  mesh,
193  IOobject::NO_READ
194  ),
195  sigma.component(symmTensor::YZ)
196  );
197  sigmayz.write();
198 
200  (
201  IOobject
202  (
203  "Ub",
204  runTime.timeName(),
205  mesh,
206  IOobject::NO_READ
207  ),
208  U,
209  zeroGradientFvPatchVectorField::typeName
210  );
212  Ub.write();
213 
214  volScalarField sigmaUn
215  (
216  IOobject
217  (
218  "sigmaUn",
219  runTime.timeName(),
220  mesh,
221  IOobject::NO_READ
222  ),
223  0.0*sigma.component(symmTensor::YZ)
224  );
225 
226  forAll(sigmaUn.boundaryField(), patchI)
227  {
228  sigmaUn.boundaryField()[patchI] =
229  (
230  mesh.boundary()[patchI].nf()
231  & sigma.boundaryField()[patchI]
232  )().component(vector::X);
233  }
234 
235  sigmaUn.write();
236  }
237  else
238  {
239  Info<< " No U" << endl;
240  }
241 
242  Info<< endl;
243  }
244 
245  Info<< "End" << endl;
246 
247  return 0;
248 }
249 
250 
251 // ************************ vim: set sw=4 sts=4 et: ************************ //