FreeFOAM The Cross-Platform CFD Toolkit
Co.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  Co
26 
27 Description
28  Calculates and writes the Co number as a surfaceScalarField obtained
29  from field phi.
30 
31  The -noWrite option just outputs the max values without writing the
32  field.
33 
34 Usage
35 
36  - Co [OPTIONS]
37 
38  @param -noWrite \n
39  Suppress output to files.
40 
41  @param -dict <dictionary name>\n
42  Use named dictionary instead of system/controlDict.
43 
44  @param -noZero \n
45  Ignore timestep 0.
46 
47  @param -constant \n
48  Include the constant directory.
49 
50  @param -time <time>\n
51  Apply only to specific time.
52 
53  @param -latestTime \n
54  Only apply to latest time step.
55 
56  @param -case <dir>\n
57  Case directory.
58 
59  @param -parallel \n
60  Run in parallel.
61 
62  @param -help \n
63  Display help message.
64 
65  @param -doc \n
66  Display Doxygen API documentation page for this application.
67 
68  @param -srcDoc \n
69  Display Doxygen source documentation page for this application.
70 
71 \*---------------------------------------------------------------------------*/
72 
73 #include <postCalc/calc.H>
74 #include <finiteVolume/fvc.H>
75 
76 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77 
78 namespace Foam
79 {
80  tmp<volScalarField> Co(const surfaceScalarField& Cof)
81  {
82  const fvMesh& mesh = Cof.mesh();
83 
84  tmp<volScalarField> tCo
85  (
86  new volScalarField
87  (
88  IOobject
89  (
90  "Co",
91  mesh.time().timeName(),
92  mesh
93  ),
94  mesh,
95  dimensionedScalar("0", Cof.dimensions(), 0)
96  )
97  );
98 
99  volScalarField& Co = tCo();
100 
101  // Set local references to mesh data
102  const unallocLabelList& owner = mesh.owner();
103  const unallocLabelList& neighbour = mesh.neighbour();
104 
105  forAll(owner, facei)
106  {
107  label own = owner[facei];
108  label nei = neighbour[facei];
109 
110  Co[own] = max(Co[own], Cof[facei]);
111  Co[nei] = max(Co[nei], Cof[facei]);
112  }
113 
114  forAll(Co.boundaryField(), patchi)
115  {
116  Co.boundaryField()[patchi] = Cof.boundaryField()[patchi];
117  }
118 
119  return tCo;
120  }
121 }
122 
123 
124 void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh)
125 {
126  bool writeResults = !args.optionFound("noWrite");
127 
128  IOobject phiHeader
129  (
130  "phi",
131  runTime.timeName(),
132  mesh,
133  IOobject::MUST_READ
134  );
135 
136  if (phiHeader.headerOk())
137  {
138  autoPtr<surfaceScalarField> CoPtr;
139 
140  Info<< " Reading phi" << endl;
141  surfaceScalarField phi(phiHeader, mesh);
142  Info<< " Calculating Co" << endl;
143 
144  if (phi.dimensions() == dimensionSet(1, 0, -1, 0, 0))
145  {
146  // compressible
148  (
149  IOobject
150  (
151  "rho",
152  runTime.timeName(),
153  mesh,
154  IOobject::MUST_READ
155  ),
156  mesh
157  );
158 
159  CoPtr.set
160  (
162  (
163  IOobject
164  (
165  "Cof",
166  runTime.timeName(),
167  mesh,
168  IOobject::NO_READ
169  ),
170  (
171  mesh.surfaceInterpolation::deltaCoeffs()
172  * (mag(phi)/(fvc::interpolate(rho)*mesh.magSf()))
173  * runTime.deltaT()
174  )
175  )
176  );
177  }
178  else if (phi.dimensions() == dimensionSet(0, 3, -1, 0, 0))
179  {
180  // incompressible
181  CoPtr.set
182  (
184  (
185  IOobject
186  (
187  "Cof",
188  runTime.timeName(),
189  mesh,
190  IOobject::NO_READ
191  ),
192  (
193  mesh.surfaceInterpolation::deltaCoeffs()
194  * (mag(phi)/mesh.magSf())
195  * runTime.deltaT()
196  )
197  )
198  );
199  }
200  else
201  {
202  FatalErrorIn(args.executable())
203  << "Incorrect dimensions of phi: " << phi.dimensions()
204  << abort(FatalError);
205  }
206 
207  Info<< "Co max : " << max(CoPtr()).value() << endl;
208 
209  if (writeResults)
210  {
211  CoPtr().write();
212  Co(CoPtr())().write();
213  }
214  }
215  else
216  {
217  Info<< " No phi" << endl;
218  }
219 
220  Info<< "\nEnd\n" << endl;
221 }
222 
223 // ************************ vim: set sw=4 sts=4 et: ************************ //