forces.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) 2011-2016 OpenFOAM Foundation
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 \*---------------------------------------------------------------------------*/
25 
26 #include "forces.H"
27 #include "fvcGrad.H"
28 #include "porosityModel.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace functionObjects
38 {
39  defineTypeNameAndDebug(forces, 0);
40 
41  addToRunTimeSelectionTable(functionObject, forces, dictionary);
42 }
43 }
44 
45 
46 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
47 
49 (
50  const dictionary& dict
51 ) const
52 {
53  DynamicList<word> names(1);
54 
55  const word forceType(dict.lookup("type"));
56  names.append(forceType);
57 
58  if (dict.found("binData"))
59  {
60  const dictionary& binDict(dict.subDict("binData"));
61  label nb = readLabel(binDict.lookup("nBin"));
62  if (nb > 0)
63  {
64  names.append(forceType + "_bins");
65  }
66  }
67 
68  return names;
69 }
70 
71 
73 {
74  if (i == 0)
75  {
76  // force data
77 
78  writeHeader(file(i), "Forces");
79  writeHeaderValue(file(i), "CofR", coordSys_.origin());
80  writeCommented(file(i), "Time");
81 
82  file(i)
83  << "forces(pressure viscous porous) "
84  << "moment(pressure viscous porous)";
85 
86  if (localSystem_)
87  {
88  file(i)
89  << tab
90  << "localForces(pressure,viscous,porous) "
91  << "localMoments(pressure,viscous,porous)";
92  }
93  }
94  else if (i == 1)
95  {
96  // bin data
97 
98  writeHeader(file(i), "Force bins");
99  writeHeaderValue(file(i), "bins", nBin_);
100  writeHeaderValue(file(i), "start", binMin_);
101  writeHeaderValue(file(i), "delta", binDx_);
102  writeHeaderValue(file(i), "direction", binDir_);
103 
104  vectorField binPoints(nBin_);
105  writeCommented(file(i), "x co-ords :");
106  forAll(binPoints, pointi)
107  {
108  binPoints[pointi] = (binMin_ + (pointi + 1)*binDx_)*binDir_;
109  file(i) << tab << binPoints[pointi].x();
110  }
111  file(i) << nl;
112 
113  writeCommented(file(i), "y co-ords :");
114  forAll(binPoints, pointi)
115  {
116  file(i) << tab << binPoints[pointi].y();
117  }
118  file(i) << nl;
119 
120  writeCommented(file(i), "z co-ords :");
121  forAll(binPoints, pointi)
122  {
123  file(i) << tab << binPoints[pointi].z();
124  }
125  file(i) << nl;
126 
127  writeCommented(file(i), "Time");
128 
129  for (label j = 0; j < nBin_; j++)
130  {
131  const word jn('(' + Foam::name(j) + ')');
132  const word f("forces" + jn + "[pressure,viscous,porous]");
133  const word m("moments" + jn + "[pressure,viscous,porous]");
134 
135  file(i)<< tab << f << tab << m;
136  }
137  if (localSystem_)
138  {
139  for (label j = 0; j < nBin_; j++)
140  {
141  const word jn('(' + Foam::name(j) + ')');
142  const word f("localForces" + jn + "[pressure,viscous,porous]");
143  const word m("localMoments" + jn + "[pressure,viscous,porous]");
144 
145  file(i)<< tab << f << tab << m;
146  }
147  }
148  }
149  else
150  {
152  << "Unhandled file index: " << i
153  << abort(FatalError);
154  }
155 
156  file(i)<< endl;
157 }
158 
159 
161 {
162  if (initialised_)
163  {
164  return;
165  }
166 
167  if (directForceDensity_)
168  {
169  if (!obr_.foundObject<volVectorField>(fDName_))
170  {
172  << "Could not find " << fDName_ << " in database."
173  << exit(FatalError);
174  }
175  }
176  else
177  {
178  if
179  (
180  !obr_.foundObject<volVectorField>(UName_)
181  || !obr_.foundObject<volScalarField>(pName_)
182 
183  )
184  {
186  << "Could not find " << UName_ << ", " << pName_
187  << exit(FatalError);
188  }
189 
190  if
191  (
192  rhoName_ != "rhoInf"
193  && !obr_.foundObject<volScalarField>(rhoName_)
194  )
195  {
197  << "Could not find " << rhoName_
198  << exit(FatalError);
199  }
200  }
201 
202  initialised_ = true;
203 }
204 
205 
208 {
209  typedef compressible::turbulenceModel cmpTurbModel;
210  typedef incompressible::turbulenceModel icoTurbModel;
211 
212  if (obr_.foundObject<cmpTurbModel>(cmpTurbModel::propertiesName))
213  {
214  const cmpTurbModel& turb =
215  obr_.lookupObject<cmpTurbModel>(cmpTurbModel::propertiesName);
216 
217  return turb.devRhoReff();
218  }
219  else if (obr_.foundObject<icoTurbModel>(icoTurbModel::propertiesName))
220  {
222  obr_.lookupObject<icoTurbModel>(icoTurbModel::propertiesName);
223 
224  return rho()*turb.devReff();
225  }
226  else if (obr_.foundObject<fluidThermo>(fluidThermo::dictName))
227  {
228  const fluidThermo& thermo =
229  obr_.lookupObject<fluidThermo>(fluidThermo::dictName);
230 
231  const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
232 
233  return -thermo.mu()*dev(twoSymm(fvc::grad(U)));
234  }
235  else if
236  (
237  obr_.foundObject<transportModel>("transportProperties")
238  )
239  {
240  const transportModel& laminarT =
241  obr_.lookupObject<transportModel>("transportProperties");
242 
243  const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
244 
245  return -rho()*laminarT.nu()*dev(twoSymm(fvc::grad(U)));
246  }
247  else if (obr_.foundObject<dictionary>("transportProperties"))
248  {
249  const dictionary& transportProperties =
250  obr_.lookupObject<dictionary>("transportProperties");
251 
252  dimensionedScalar nu(transportProperties.lookup("nu"));
253 
254  const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
255 
256  return -rho()*nu*dev(twoSymm(fvc::grad(U)));
257  }
258  else
259  {
261  << "No valid model for viscous stress calculation"
262  << exit(FatalError);
263 
264  return volSymmTensorField::null();
265  }
266 }
267 
268 
270 {
271  if (obr_.foundObject<fluidThermo>(basicThermo::dictName))
272  {
273  const fluidThermo& thermo =
274  obr_.lookupObject<fluidThermo>(basicThermo::dictName);
275 
276  return thermo.mu();
277  }
278  else if
279  (
280  obr_.foundObject<transportModel>("transportProperties")
281  )
282  {
283  const transportModel& laminarT =
284  obr_.lookupObject<transportModel>("transportProperties");
285 
286  return rho()*laminarT.nu();
287  }
288  else if (obr_.foundObject<dictionary>("transportProperties"))
289  {
290  const dictionary& transportProperties =
291  obr_.lookupObject<dictionary>("transportProperties");
292 
294  (
295  "nu",
296  dimViscosity,
297  transportProperties.lookup("nu")
298  );
299 
300  return rho()*nu;
301  }
302  else
303  {
305  << "No valid model for dynamic viscosity calculation"
306  << exit(FatalError);
307 
308  return volScalarField::null();
309  }
310 }
311 
312 
314 {
315  if (rhoName_ == "rhoInf")
316  {
317  const fvMesh& mesh = refCast<const fvMesh>(obr_);
318 
319  return tmp<volScalarField>
320  (
321  new volScalarField
322  (
323  IOobject
324  (
325  "rho",
326  mesh.time().timeName(),
327  mesh
328  ),
329  mesh,
330  dimensionedScalar("rho", dimDensity, rhoRef_)
331  )
332  );
333  }
334  else
335  {
336  return(obr_.lookupObject<volScalarField>(rhoName_));
337  }
338 }
339 
340 
342 {
343  if (p.dimensions() == dimPressure)
344  {
345  return 1.0;
346  }
347  else
348  {
349  if (rhoName_ != "rhoInf")
350  {
352  << "Dynamic pressure is expected but kinematic is provided."
353  << exit(FatalError);
354  }
355 
356  return rhoRef_;
357  }
358 }
359 
360 
362 (
363  const vectorField& Md,
364  const vectorField& fN,
365  const vectorField& fT,
366  const vectorField& fP,
367  const vectorField& d
368 )
369 {
370  if (nBin_ == 1)
371  {
372  force_[0][0] += sum(fN);
373  force_[1][0] += sum(fT);
374  force_[2][0] += sum(fP);
375  moment_[0][0] += sum(Md^fN);
376  moment_[1][0] += sum(Md^fT);
377  moment_[2][0] += sum(Md^fP);
378  }
379  else
380  {
381  scalarField dd((d & binDir_) - binMin_);
382 
383  forAll(dd, i)
384  {
385  label bini = min(max(floor(dd[i]/binDx_), 0), force_[0].size() - 1);
386 
387  force_[0][bini] += fN[i];
388  force_[1][bini] += fT[i];
389  force_[2][bini] += fP[i];
390  moment_[0][bini] += Md[i]^fN[i];
391  moment_[1][bini] += Md[i]^fT[i];
392  moment_[2][bini] += Md[i]^fP[i];
393  }
394  }
395 }
396 
397 
399 {
400  Log << type() << " " << name() << " write:" << nl
401  << " sum of forces:" << nl
402  << " pressure : " << sum(force_[0]) << nl
403  << " viscous : " << sum(force_[1]) << nl
404  << " porous : " << sum(force_[2]) << nl
405  << " sum of moments:" << nl
406  << " pressure : " << sum(moment_[0]) << nl
407  << " viscous : " << sum(moment_[1]) << nl
408  << " porous : " << sum(moment_[2])
409  << endl;
410 
411  writeTime(file(0));
412  file(0) << tab << setw(1) << '('
413  << sum(force_[0]) << setw(1) << ' '
414  << sum(force_[1]) << setw(1) << ' '
415  << sum(force_[2]) << setw(3) << ") ("
416  << sum(moment_[0]) << setw(1) << ' '
417  << sum(moment_[1]) << setw(1) << ' '
418  << sum(moment_[2]) << setw(1) << ')'
419  << endl;
420 
421  if (localSystem_)
422  {
423  vectorField localForceN(coordSys_.localVector(force_[0]));
424  vectorField localForceT(coordSys_.localVector(force_[1]));
425  vectorField localForceP(coordSys_.localVector(force_[2]));
426  vectorField localMomentN(coordSys_.localVector(moment_[0]));
427  vectorField localMomentT(coordSys_.localVector(moment_[1]));
428  vectorField localMomentP(coordSys_.localVector(moment_[2]));
429 
430  writeTime(file(0));
431  file(0) << tab << setw(1) << '('
432  << sum(localForceN) << setw(1) << ' '
433  << sum(localForceT) << setw(1) << ' '
434  << sum(localForceP) << setw(3) << ") ("
435  << sum(localMomentN) << setw(1) << ' '
436  << sum(localMomentT) << setw(1) << ' '
437  << sum(localMomentP) << setw(1) << ')'
438  << endl;
439  }
440 }
441 
442 
444 {
445  if (nBin_ == 1)
446  {
447  return;
448  }
449 
450  List<Field<vector>> f(force_);
451  List<Field<vector>> m(moment_);
452 
453  if (binCumulative_)
454  {
455  for (label i = 1; i < f[0].size(); i++)
456  {
457  f[0][i] += f[0][i-1];
458  f[1][i] += f[1][i-1];
459  f[2][i] += f[2][i-1];
460 
461  m[0][i] += m[0][i-1];
462  m[1][i] += m[1][i-1];
463  m[2][i] += m[2][i-1];
464  }
465  }
466 
467  writeTime(file(1));
468 
469  forAll(f[0], i)
470  {
471  file(1)
472  << tab << setw(1) << '('
473  << f[0][i] << setw(1) << ' '
474  << f[1][i] << setw(1) << ' '
475  << f[2][i] << setw(3) << ") ("
476  << m[0][i] << setw(1) << ' '
477  << m[1][i] << setw(1) << ' '
478  << m[2][i] << setw(1) << ')';
479  }
480 
481  if (localSystem_)
482  {
483  List<Field<vector>> lf(3);
484  List<Field<vector>> lm(3);
485  lf[0] = coordSys_.localVector(force_[0]);
486  lf[1] = coordSys_.localVector(force_[1]);
487  lf[2] = coordSys_.localVector(force_[2]);
488  lm[0] = coordSys_.localVector(moment_[0]);
489  lm[1] = coordSys_.localVector(moment_[1]);
490  lm[2] = coordSys_.localVector(moment_[2]);
491 
492  if (binCumulative_)
493  {
494  for (label i = 1; i < lf[0].size(); i++)
495  {
496  lf[0][i] += lf[0][i-1];
497  lf[1][i] += lf[1][i-1];
498  lf[2][i] += lf[2][i-1];
499  lm[0][i] += lm[0][i-1];
500  lm[1][i] += lm[1][i-1];
501  lm[2][i] += lm[2][i-1];
502  }
503  }
504 
505  forAll(lf[0], i)
506  {
507  file(1)
508  << tab << setw(1) << '('
509  << lf[0][i] << setw(1) << ' '
510  << lf[1][i] << setw(1) << ' '
511  << lf[2][i] << setw(3) << ") ("
512  << lm[0][i] << setw(1) << ' '
513  << lm[1][i] << setw(1) << ' '
514  << lm[2][i] << setw(1) << ')';
515  }
516  }
517 
518  file(1) << endl;
519 }
520 
521 
522 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
523 
525 (
526  const word& name,
527  const Time& runTime,
528  const dictionary& dict
529 )
530 :
531  writeFiles(name, runTime, dict, name),
532  force_(3),
533  moment_(3),
534  patchSet_(),
535  pName_(word::null),
536  UName_(word::null),
537  rhoName_(word::null),
538  directForceDensity_(false),
539  fDName_(""),
540  rhoRef_(VGREAT),
541  pRef_(0),
542  coordSys_(),
543  localSystem_(false),
544  porosity_(false),
545  nBin_(1),
546  binDir_(Zero),
547  binDx_(0.0),
548  binMin_(GREAT),
549  binPoints_(),
550  binCumulative_(true),
551  initialised_(false)
552 {
553  if (!isA<fvMesh>(obr_))
554  {
556  << "objectRegistry is not an fvMesh" << exit(FatalError);
557  }
558 
559  read(dict);
560  resetNames(createFileNames(dict));
561 }
562 
563 
565 (
566  const word& name,
567  const objectRegistry& obr,
568  const dictionary& dict
569 )
570 :
571  writeFiles(name, obr, dict, name),
572  force_(3),
573  moment_(3),
574  patchSet_(),
575  pName_(word::null),
576  UName_(word::null),
577  rhoName_(word::null),
578  directForceDensity_(false),
579  fDName_(""),
580  rhoRef_(VGREAT),
581  pRef_(0),
582  coordSys_(),
583  localSystem_(false),
584  porosity_(false),
585  nBin_(1),
586  binDir_(Zero),
587  binDx_(0.0),
588  binMin_(GREAT),
589  binPoints_(),
590  binCumulative_(true),
591  initialised_(false)
592 {
593  if (!isA<fvMesh>(obr_))
594  {
596  << "objectRegistry is not an fvMesh" << exit(FatalError);
597  }
598 
599  read(dict);
600  resetNames(createFileNames(dict));
601 }
602 
603 
604 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
605 
607 {}
608 
609 
610 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
611 
613 {
614  writeFiles::read(dict);
615 
616  initialised_ = false;
617 
618  Log << type() << " " << name() << ":" << nl;
619 
620  directForceDensity_ = dict.lookupOrDefault("directForceDensity", false);
621 
622  const fvMesh& mesh = refCast<const fvMesh>(obr_);
623  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
624 
625  patchSet_ = pbm.patchSet(wordReList(dict.lookup("patches")));
626 
627  if (directForceDensity_)
628  {
629  // Optional entry for fDName
630  fDName_ = dict.lookupOrDefault<word>("fD", "fD");
631  }
632  else
633  {
634  // Optional entries U and p
635  pName_ = dict.lookupOrDefault<word>("p", "p");
636  UName_ = dict.lookupOrDefault<word>("U", "U");
637  rhoName_ = dict.lookupOrDefault<word>("rho", "rho");
638 
639  // Reference density needed for incompressible calculations
640  if (rhoName_ == "rhoInf")
641  {
642  rhoRef_ = readScalar(dict.lookup("rhoInf"));
643  }
644 
645  // Reference pressure, 0 by default
646  pRef_ = dict.lookupOrDefault<scalar>("pRef", 0.0);
647  }
648 
649  coordSys_.clear();
650 
651  // Centre of rotation for moment calculations
652  // specified directly, from coordinate system, or implicitly (0 0 0)
653  if (!dict.readIfPresent<point>("CofR", coordSys_.origin()))
654  {
655  coordSys_ = coordinateSystem(obr_, dict);
656  localSystem_ = true;
657  }
658 
659  dict.readIfPresent("porosity", porosity_);
660  if (porosity_)
661  {
662  Log << " Including porosity effects" << endl;
663  }
664  else
665  {
666  Log << " Not including porosity effects" << endl;
667  }
668 
669  if (dict.found("binData"))
670  {
671  const dictionary& binDict(dict.subDict("binData"));
672  binDict.lookup("nBin") >> nBin_;
673 
674  if (nBin_ < 0)
675  {
677  << "Number of bins (nBin) must be zero or greater"
678  << exit(FatalIOError);
679  }
680  else if ((nBin_ == 0) || (nBin_ == 1))
681  {
682  nBin_ = 1;
683  forAll(force_, i)
684  {
685  force_[i].setSize(1);
686  moment_[i].setSize(1);
687  }
688  }
689 
690  if (nBin_ > 1)
691  {
692  binDict.lookup("direction") >> binDir_;
693  binDir_ /= mag(binDir_);
694 
695  binMin_ = GREAT;
696  scalar binMax = -GREAT;
697  forAllConstIter(labelHashSet, patchSet_, iter)
698  {
699  label patchi = iter.key();
700  const polyPatch& pp = pbm[patchi];
701  scalarField d(pp.faceCentres() & binDir_);
702  binMin_ = min(min(d), binMin_);
703  binMax = max(max(d), binMax);
704  }
705  reduce(binMin_, minOp<scalar>());
706  reduce(binMax, maxOp<scalar>());
707 
708  // slightly boost binMax so that region of interest is fully
709  // within bounds
710  binMax = 1.0001*(binMax - binMin_) + binMin_;
711 
712  binDx_ = (binMax - binMin_)/scalar(nBin_);
713 
714  // create the bin points used for writing
715  binPoints_.setSize(nBin_);
716  forAll(binPoints_, i)
717  {
718  binPoints_[i] = (i + 0.5)*binDir_*binDx_;
719  }
720 
721  binDict.lookup("cumulative") >> binCumulative_;
722 
723  // allocate storage for forces and moments
724  forAll(force_, i)
725  {
726  force_[i].setSize(nBin_);
727  moment_[i].setSize(nBin_);
728  }
729  }
730  }
731 
732  if (nBin_ == 1)
733  {
734  // allocate storage for forces and moments
735  force_[0].setSize(1);
736  force_[1].setSize(1);
737  force_[2].setSize(1);
738  moment_[0].setSize(1);
739  moment_[1].setSize(1);
740  moment_[2].setSize(1);
741  }
742 
743  return true;
744 }
745 
746 
748 {
749  initialise();
750 
751  force_[0] = Zero;
752  force_[1] = Zero;
753  force_[2] = Zero;
754 
755  moment_[0] = Zero;
756  moment_[1] = Zero;
757  moment_[2] = Zero;
758 
759  if (directForceDensity_)
760  {
761  const volVectorField& fD = obr_.lookupObject<volVectorField>(fDName_);
762 
763  const fvMesh& mesh = fD.mesh();
764 
765  const surfaceVectorField::Boundary& Sfb =
766  mesh.Sf().boundaryField();
767 
768  forAllConstIter(labelHashSet, patchSet_, iter)
769  {
770  label patchi = iter.key();
771 
772  vectorField Md
773  (
774  mesh.C().boundaryField()[patchi] - coordSys_.origin()
775  );
776 
777  scalarField sA(mag(Sfb[patchi]));
778 
779  // Normal force = surfaceUnitNormal*(surfaceNormal & forceDensity)
780  vectorField fN
781  (
782  Sfb[patchi]/sA
783  *(
784  Sfb[patchi] & fD.boundaryField()[patchi]
785  )
786  );
787 
788  // Tangential force (total force minus normal fN)
789  vectorField fT(sA*fD.boundaryField()[patchi] - fN);
790 
791  //- Porous force
792  vectorField fP(Md.size(), Zero);
793 
794  applyBins(Md, fN, fT, fP, mesh.C().boundaryField()[patchi]);
795  }
796  }
797  else
798  {
799  const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
800  const volScalarField& p = obr_.lookupObject<volScalarField>(pName_);
801 
802  const fvMesh& mesh = U.mesh();
803 
804  const surfaceVectorField::Boundary& Sfb =
805  mesh.Sf().boundaryField();
806 
807  tmp<volSymmTensorField> tdevRhoReff = devRhoReff();
808  const volSymmTensorField::Boundary& devRhoReffb
809  = tdevRhoReff().boundaryField();
810 
811  // Scale pRef by density for incompressible simulations
812  scalar pRef = pRef_/rho(p);
813 
814  forAllConstIter(labelHashSet, patchSet_, iter)
815  {
816  label patchi = iter.key();
817 
818  vectorField Md
819  (
820  mesh.C().boundaryField()[patchi] - coordSys_.origin()
821  );
822 
823  vectorField fN
824  (
825  rho(p)*Sfb[patchi]*(p.boundaryField()[patchi] - pRef)
826  );
827 
828  vectorField fT(Sfb[patchi] & devRhoReffb[patchi]);
829 
830  vectorField fP(Md.size(), Zero);
831 
832  applyBins(Md, fN, fT, fP, mesh.C().boundaryField()[patchi]);
833  }
834  }
835 
836  if (porosity_)
837  {
838  const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
839  const volScalarField rho(this->rho());
840  const volScalarField mu(this->mu());
841 
842  const fvMesh& mesh = U.mesh();
843 
844  const HashTable<const porosityModel*> models =
845  obr_.lookupClass<porosityModel>();
846 
847  if (models.empty())
848  {
850  << "Porosity effects requested, but no porosity models found "
851  << "in the database"
852  << endl;
853  }
854 
856  {
857  // non-const access required if mesh is changing
858  porosityModel& pm = const_cast<porosityModel&>(*iter());
859 
860  vectorField fPTot(pm.force(U, rho, mu));
861 
862  const labelList& cellZoneIDs = pm.cellZoneIDs();
863 
864  forAll(cellZoneIDs, i)
865  {
866  label zoneI = cellZoneIDs[i];
867  const cellZone& cZone = mesh.cellZones()[zoneI];
868 
869  const vectorField d(mesh.C(), cZone);
870  const vectorField fP(fPTot, cZone);
871  const vectorField Md(d - coordSys_.origin());
872 
873  const vectorField fDummy(Md.size(), Zero);
874 
875  applyBins(Md, fDummy, fDummy, fP, d);
876  }
877  }
878  }
879 
884 }
885 
886 
888 {
889  return sum(force_[0]) + sum(force_[1]) + sum(force_[2]);
890 }
891 
892 
894 {
895  return sum(moment_[0]) + sum(moment_[1]) + sum(moment_[2]);
896 }
897 
898 
900 {
901  return true;
902 }
903 
904 
906 {
907  calcForcesMoment();
908 
909  if (Pstream::master())
910  {
912 
913  writeForces();
914 
915  writeBins();
916 
917  Log << endl;
918  }
919 
920  return true;
921 }
922 
923 
924 // ************************************************************************* //
Base class for other coordinate system specifications.
virtual void writeFileHeader(const label i)
Output file header information.
Definition: forces.C:72
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual tmp< volScalarField > mu() const =0
Dynamic viscosity of mixture [kg/m/s].
U
Definition: pEqn.H:83
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const char tab
Definition: Ostream.H:261
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
functionObject base class for writing files
Definition: writeFiles.H:56
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
virtual vector forceEff() const
Return the total force.
Definition: forces.C:887
virtual tmp< vectorField > force(const volVectorField &U, const volScalarField &rho, const volScalarField &mu)
Return the force over the cell zone(s)
const dimensionSet dimViscosity
compressible::turbulenceModel & turb
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
const volVectorField & C() const
Return cell centres as volVectorField.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
virtual bool read(const dictionary &)
Read the forces data.
Definition: forces.C:612
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
Macros for easy insertion into run-time selection tables.
void writeBins()
Helper function to write bin data.
Definition: forces.C:443
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
tmp< volSymmTensorField > devRhoReff() const
Return the effective viscous stress (laminar + turbulent).
Definition: forces.C:207
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
Templated abstract base class for single-phase incompressible turbulence models.
tmp< volScalarField > mu() const
Dynamic viscosity field.
Definition: forces.C:269
virtual bool execute()
Execute, currently does nothing.
Definition: forces.C:899
virtual bool read(const dictionary &)
Read optional controls.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dynamicFvMesh & mesh
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Calculate the gradient of the given field.
A class for handling words, derived from string.
Definition: word.H:59
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:49
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:292
static const word null
An empty word.
Definition: word.H:77
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
virtual void calcForcesMoment()
Calculate the forces and moments.
Definition: forces.C:747
static const zero Zero
Definition: zero.H:91
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
void applyBins(const vectorField &Md, const vectorField &fN, const vectorField &fT, const vectorField &fP, const vectorField &d)
Accumulate bin data.
Definition: forces.C:362
An STL-conforming hash table.
Definition: HashTable.H:61
const dimensionSet dimPressure
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
label readLabel(Istream &is)
Definition: label.H:64
const dimensionSet & dimensions() const
Return dimensions.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Foam::polyBoundaryMesh.
forces(const forces &)
Disallow default bitwise copy construct.
static const char nl
Definition: Ostream.H:262
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
void initialise()
Initialise the fields.
Definition: forces.C:160
scalar pRef
Definition: createFields.H:19
labelList f(nPoints)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
A subset of mesh cells.
Definition: cellZone.H:61
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual vector momentEff() const
Return the total moment.
Definition: forces.C:893
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimDensity
static const GeometricField< symmTensor, fvPatchField, volMesh > & null()
Return a null geometric field.
void writeForces()
Helper function to write force data.
Definition: forces.C:398
const dimensionedScalar mu
Atomic mass unit.
virtual tmp< volSymmTensorField > devReff() const
Return the effective stress tensor.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:469
label patchi
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
#define WarningInFunction
Report a warning using Foam::Warning.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const Mesh & mesh() const
Return mesh.
Base-class for all transport models used by the incompressible turbulence models. ...
addToRunTimeSelectionTable(functionObject, blendingFactor, dictionary)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual bool write()
Write function.
Definition: writeFiles.C:197
#define Log
Report write to Foam::Info if the local log switch is true.
dimensioned< scalar > mag(const dimensioned< Type > &)
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:284
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:54
Registry of regIOobjects.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const labelList & cellZoneIDs() const
Return const access to the cell zone IDs.
virtual ~forces()
Destructor.
Definition: forces.C:606
tmp< volScalarField > rho() const
Return rho if specified otherwise rhoRef.
Definition: forces.C:313
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
wordList createFileNames(const dictionary &dict) const
Create file names for forces and bins.
Definition: forces.C:49
volScalarField & nu
Top level model for porosity models.
Definition: porosityModel.H:55
virtual bool write()
Write the forces.
Definition: forces.C:905
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
IOerror FatalIOError