STDMD.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2020 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Class
27  Foam::functionObjects::STDMD
28 
29 Group
30  grpFieldFunctionObjects
31 
32 Description
33  (Beta Release) STDMD (i.e. Streaming Total Dynamic Mode Decomposition) is
34  a variant of a data-driven dimensionality reduction method.
35 
36  STDMD is being used as a mathematical post-processing tool to compute
37  a set of dominant modes out of a given flow (or dataset) each of which is
38  associated with a constant frequency and decay rate, so that dynamic
39  features of a given flow may become interpretable, and tractable.
40  Among other Dynamic Mode Decomposition (DMD) variants, STDMD is presumed
41  to provide the general DMD method capabilities alongside economised and
42  feasible memory and CPU usage.
43 
44  The code implementation corresponds to Figs. 15-16 of the first citation
45  below, more broadly to Section 2.4.
46 
47  References:
48  \verbatim
49  STDMD and mode-sorting algorithms (tags:K, HRDC, KZ, HWR):
50  Kiewat, M. (2019).
51  Streaming modal decomposition approaches for vehicle aerodynamics.
52  PhD thesis. Munich: Technical University of Munich.
53  URL:mediatum.ub.tum.de/doc/1482652/1482652.pdf
54 
55  Hemati, M. S., Rowley, C. W.,
56  Deem, E. A., & Cattafesta, L. N. (2017).
57  De-biasing the dynamic mode decomposition
58  for applied Koopman spectral analysis of noisy datasets.
59  Theoretical and Computational Fluid Dynamics, 31(4), 349-368.
60  DOI:10.1007/s00162-017-0432-2
61 
62  Kou, J., & Zhang, W. (2017).
63  An improved criterion to select
64  dominant modes from dynamic mode decomposition.
65  European Journal of Mechanics-B/Fluids, 62, 109-129.
66  DOI:10.1016/j.euromechflu.2016.11.015
67 
68  Hemati, M. S., Williams, M. O., & Rowley, C. W. (2014).
69  Dynamic mode decomposition for large and streaming datasets.
70  Physics of Fluids, 26(11), 111701.
71  DOI:10.1063/1.4901016
72 
73  Parallel classical Gram-Schmidt process (tag:Ka):
74  Katagiri, T. (2003).
75  Performance evaluation of parallel
76  Gram-Schmidt re-orthogonalization methods.
77  In: Palma J. M. L. M., Sousa A. A., Dongarra J., Hernández V. (eds)
78  High Performance Computing for Computational Science — VECPAR 2002.
79  Lecture Notes in Computer Science, vol 2565, p. 302-314.
80  Berlin, Heidelberg: Springer.
81  DOI:10.1007/3-540-36569-9_19
82 
83  Parallel direct tall-skinny QR decomposition (tags:BGD, DGHL):
84  Benson, A. R., Gleich, D. F., & Demmel, J. (2013).
85  Direct QR factorizations for
86  tall-and-skinny matrices in MapReduce architectures.
87  2013 IEEE International Conference on Big Data.
88  DOI:10.1109/bigdata.2013.6691583
89 
90  Demmel, J., Grigori, L., Hoemmen, M., & Langou, J. (2012).
91  Communication-optimal parallel
92  and sequential QR and LU factorizations.
93  SIAM Journal on Scientific Computing, 34(1), A206-A239.
94  DOI:10.1137/080731992
95 
96  DMD properties:
97  Brunton S. L. (2018).
98  Dynamic mode decomposition overview.
99  Seattle, Washington: University of Washington.
100  youtu.be/sQvrK8AGCAo (Retrieved:24-04-20)
101  \endverbatim
102 
103  Operands:
104  \table
105  Operand | Type | Location
106  input | {vol,surface}<Type>Field <!--
107  --> | $FOAM_CASE/<time>/<inpField>
108  output file | dat <!--
109  --> | $FOAM_CASE/postProcessing/<FO>/<time>/<file>(s)
110  output field | volScalarField(s) | $FOAM_CASE/<time>/<outField>(s)
111  \endtable
112 
113  where \c <Type>=Scalar/Vector/SphericalTensor/SymmTensor/Tensor.
114 
115  Output fields:
116  \verbatim
117  modeReal<modeIndex><field> | Real part of a mode field
118  modeImag<modeIndex><field> | Imaginary part of a mode field
119  \endverbatim
120 
121  Output files:
122  \verbatim
123  uSTDMD.dat | Unfiltered STDMD output
124  STDMD.dat | Filtered STDMD output
125  \endverbatim
126 
127  wherein for each mode, the following quantities are output into files:
128  \verbatim
129  Frequency
130  Magnitude
131  Amplitude (real part)
132  Amplitude (imaginary part)
133  Eigenvalue (real part)
134  Eigenvalue (imaginary part)
135  \endverbatim
136 
137 Note
138  Operations on boundary fields, e.g. \c wallShearStress, are currently not
139  available.
140 
141 Usage
142  Minimal example by using \c system/controlDict.functions:
143  \verbatim
144  STDMD1
145  {
146  // Mandatory entries (unmodifiable)
147  type STDMD;
148  libs (fieldFunctionObjects);
149  field <inpField>;
150 
151  // Conditional mandatory entries (unmodifiable)
152  // either of the below
153  stdmdInterval 5.5;
154  executeInterval 10;
155 
156  // Optional entries (unmodifiable)
157  modeSorter kiewat;
158  nModes 50;
159  maxRank 50;
160  nGramSchmidt 5;
161  fMin 0;
162  fMax 1000000000;
163 
164  // Optional entries (run-time modifiable, yet not recommended)
165  testEigen false;
166  dumpEigen false;
167  minBasis 0.00000001;
168  minEVal 0.00000001;
169  absTol 0.001;
170  relTol 0.0001;
171  sortLimiter 500;
172 
173  // Optional (inherited) entries
174  ...
175  }
176  \endverbatim
177 
178  where the entries mean:
179  \table
180  Property | Description | Type | Req'd | Dflt
181  type | Type name: STDMD | word | yes | -
182  libs | Library name: fieldFunctionObjects | word | yes | -
183  field | Name of the operand field | word | yes | -
184  stdmdInterval | STDMD time-step size [s] | scalar| conditional | <!--
185  --> executeInterval*(current time-step of the simulation)
186  modeSorter | Mode-sorting algorithm variant | word | no | firstSnap
187  nModes | Number of output modes in input freq range | label | no | GREAT
188  maxRank | Max columns in accumulated matrices | label | no | GREAT
189  nGramSchmidt | Number of Gram-Schmidt iterations | label | no | 5
190  fMin | Min (non-negative) output frequency | label | no | 0
191  fMax | Max output frequency | label | no | GREAT
192  testEigen | Flag to verify eigendecompositions | bool | no | false
193  dumpEigen | Flag to log operands of eigendecomps | bool | no | false
194  minBasis | Orthogonal basis expansion threshold | scalar| no | 1e-8
195  minEVal | Min EVal for below EVals are omitted | scalar| no | 1e-8
196  absTol | Min abs tol in eigendecomposition tests | scalar| no | 1e-4
197  relTol | Relative tol in eigendecomposition tests | scalar| no | 1e-6
198  sortLimiter | Maximum allowable magnitude for <!--
199  --> mag(eigenvalue)*(number of STDMD steps) to be used in <!--
200  --> modeSorter={kiewat,kouZhang} to avoid overflow errors <!--
201  --> | scalar | no | 500.0
202  \endtable
203 
204  Options for the \c modeSorter entry:
205  \verbatim
206  kiewat | Modified weighted amplitude scaling method
207  kouZhang | Weighted amplitude scaling method
208  firstSnap | Method of first snapshot amplitude magnitude
209  \endverbatim
210 
211  The inherited entries are elaborated in:
212  - \link functionObject.H \endlink
213  - \link writeFile.H \endlink
214 
215  Usage by \c postProcess utility is not available.
216 
217 Note
218  - To specify the STDMD time-step size (not necessarily equal to the
219  time step of the simulation), entries of either \c stdmdInterval or
220  \c executeInterval must be available (highest to lowest precedence). While
221  \c stdmdInterval allows users to directly specify the STDMD time-step size
222  in seconds, in absence of \c stdmdInterval, for convenience,
223  \c executeInterval allows users to compute the STDMD time-step internally
224  by multiplying itself with the current time-step size of the simulation.
225  - Limitation: Restart is currently not available since intermediate writing
226  of STDMD matrices are not supported.
227  - Limitation: Non-physical input (e.g. full-zero fields) may upset STDMD.
228  - Warning: DMD is an active research area at the time of writing; therefore,
229  there would be cases whereat oddities might be encountered.
230  - Warning: This STDMD release is the \b beta release; therefore,
231  small-to-medium changes in input/output interfaces and internal structures
232  should be expected in the next versions.
233  - Warning: In \c modeSorter algorithms of \c kiewat and \c kouZhang, there
234  exists a function of \f$power(x,y)\f$ where \c x is the magnitude of an
235  eigenvalue, and \c y is the number of STDMD snapshots. This function poses
236  a risk of overflow errors since, for example, if \c x ends up above 1.5 with
237  500 snapshots or more, this function automatically throws the floating point
238  error meaning overflow. Therefore, the domain-range of this function is
239  heuristically constrained by the optional entry \c sortLimiter which skips
240  the evaluation of this function for a given
241  mag(eigenvalue)*(no. of STDMD steps), i.e. x*y, whose magnitude is larger
242  than \c sortLimiter.
243 
244 See also
245  - Foam::functionObjects::fvMeshFunctionObject
246  - Foam::functionObjects::writeFile
247  - ExtendedCodeGuide::functionObjects::field::STDMD
248 
249 SourceFiles
250  STDMD.C
251  STDMDTemplates.C
252 
253 \*---------------------------------------------------------------------------*/
254 
255 #ifndef functionObjects_STDMD_H
256 #define functionObjects_STDMD_H
257 
258 #include "fvMeshFunctionObject.H"
259 #include "writeFile.H"
260 #include "EigenMatrix.H"
261 #include "QRMatrix.H"
262 
263 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
264 
265 namespace Foam
266 {
267 namespace functionObjects
268 {
269 
270 /*---------------------------------------------------------------------------*\
271  Class STDMD Declaration
272 \*---------------------------------------------------------------------------*/
273 
274 class STDMD
275 :
276  public fvMeshFunctionObject,
277  public writeFile
278 {
279 
280  typedef SquareMatrix<scalar> SMatrix;
281  typedef SquareMatrix<complex> SCMatrix;
282  typedef RectangularMatrix<scalar> RMatrix;
283  typedef RectangularMatrix<complex> RCMatrix;
284 
285 
286  // Private Enumerations
287 
288  //- Options for the mode-sorting algorithm
289  enum modeSorterType
290  {
291  KIEWAT,
292  KOU_ZHANG,
293  FIRST_SNAPSHOT
294  };
295 
296  //- Names for modeSorterType
297  static const Enum<modeSorterType> modeSorterTypeNames;
298 
299 
300  // Private Data
301 
302  //- Mode-sorting algorithm (default=modeSorterType::FIRST_SNAPSHOT)
303  const enum modeSorterType modeSorter_;
304 
305  //- Name of the operand volume or surface field
306  const word fieldName_;
307 
308  //- Flag: First execution-step initialisation
309  bool initialised_;
310 
311  //- Flag: Verify eigendecompositions by using theoretical relations
312  bool testEigen_;
313 
314  //- Flag: Write operands of eigendecompositions to log
315  // To activate it, testEigen=true
316  bool dumpEigen_;
317 
318  //- Number of output modes within input frequency range
319  //- starting from the most energetic mode
320  const label nModes_;
321 
322  //- Maximum allowable matrix column for Qz_ and Gz_
323  // Qz_ is assumed to always have full-rank, thus Qz_.n() = rank
324  const label maxRank_;
325 
326  //- Number of maximum iterations of the classical Gram-Schmidt procedure
327  const label nGramSchmidt_;
328 
329  //- Min frequency: Output only entries corresponding to freqs >= fMin
330  const label fMin_;
331 
332  //- Max frequency: Output only entries corresponding to freqs <= fMax
333  const label fMax_;
334 
335  //- Number of base components of input field, e.g. 3 for vector
336  label nComps_;
337 
338  //- Number of elements in a snapshot
339  // A snapshot is an input dataset to be processed per execution-step
340  label nSnap_;
341 
342  //- Current execution-step index of STDMD,
343  //- not necessarily that of the simulation
344  label currIndex_;
345 
346  //- Maximum allowable magnitude for mag(eigenvalue)*(number of
347  //- STDMD steps) to be used in modeSorter={kiewat,kouZhang} to avoid
348  //- overflow errors
349  scalar sortLimiter_;
350 
351  //- Min absolute tolerance being used in eigendecomposition tests
352  scalar absTol_;
353 
354  //- Relative tolerance being used in eigendecomposition tests
355  scalar relTol_;
356 
357  //- Min value to execute orthogonal basis expansion of Qz_ and Gz_
358  scalar minBasis_;
359 
360  //- STDMD time-step size that equals to
361  //- (executeInterval of STDMD)*(deltaT of simulation) [s]
362  scalar dt_;
363 
364  //- Min EVal magnitude threshold below where EVals are omitted
365  scalar minMagEVal_;
366 
367  //- L2-norm of column vector z_
368  scalar zNorm_;
369 
370  //- L2-norm of column vector ez_
371  scalar ezNorm_;
372 
373  //- Augmented snapshot matrix (effectively a column vector) (K:Eq. 60)
374  // Upper half z_ = current-time snapshot slot
375  // Lower half z_ = previous-time snapshot slot
376  RMatrix z_;
377 
378  //- Working copy of the augmented snapshot matrix z_
379  //- being used in the classical Gram-Schmidt process
380  RMatrix ez_;
381 
382  //- First-processed snapshot required by the mode-sorting
383  //- algorithms at the final output computations (K:p. 43)
384  RMatrix X1_;
385 
386  //- Accumulated-in-time unitary similarity matrix produced by the
387  //- orthonormal decomposition of the augmented snapshot matrix z_
388  // (K:Eq. 60)
389  RMatrix Qz_;
390 
391  //- Accumulated-in-time (squared) upper triangular matrix produced by
392  //- the orthonormal decomposition of the augmented snapshot matrix z_
393  // (K:Eq. 64)
394  SMatrix Gz_;
395 
396  //- Upper half of Qz_
397  RMatrix QzUH_;
398 
399  //- Lower half of Qz_
400  RMatrix QzLH_;
401 
402  //- Moore-Penrose pseudo-inverse of R produced by
403  //- the QR decomposition of the last time-step QzUH_
404  RMatrix RxInv_;
405 
406  //- Projected STDMD operator (K:Eq. 78)
407  SMatrix Ap_;
408 
409  //- Output eigenvectors
410  RCMatrix oEVecs_;
411 
412  //- Output eigenvalues
413  List<complex> oEVals_;
414 
415  //- Output amplitudes
416  List<complex> oAmps_;
417 
418  //- Output (non-negative) frequencies
419  List<scalar> oFreqs_;
420 
421  //- Indices of oFreqs_ where freqs are
422  //- non-negative and within [fMin_, fMax_]
423  DynamicList<label> iFreqs_;
424 
425  //- Output (descending) magnitudes of (complex) amplitudes
426  List<scalar> oMags_;
427 
428  //- Indices of oMags_
429  List<label> iMags_;
430 
431 
432  // Private Member Functions
433 
434  // Process
435 
436  //- Move previous-time snapshot into previous-time slot in z_
437  //- and copy new current-time snapshot into current-time slot in z_
438  template<class GeoFieldType>
439  bool getSnapshot();
440 
441  //- Get the input field type to be processed by snapshot()
442  template<class Type>
443  bool getSnapshotType();
444 
445  //- Get the number of base components of input field
446  template<class GeoFieldType>
447  bool getComps();
448 
449  //- Return (parallel) L2-norm of a given column vector
450  scalar parnorm(const RMatrix& colVector) const;
451 
452  //- Move the current-time snapshot into the previous-time snapshot in z_
453  //- and copy the new field into the current-time snapshot
454  void snapshot();
455 
456  //- Initialise all working matrices at the first execution-step
457  void init();
458 
459  //- Initialise Qz_, Gz_ (both require the first two snapshots) and X1_
460  void initBasis();
461 
462  //- Execute (parallel) classical Gram-Schmidt
463  //- process to orthonormalise ez_ (Ka:Fig. 5)
464  void GramSchmidt();
465 
466  //- Expand orthonormal bases Qz_ and Gz_ by stacking a column
467  //- (ez_/ezNorm_) to Qz_, and a row (Zero) and column (Zero)
468  //- to Gz_ if (minBasis_ < ezNorm_/zNorm_)
469  void expandBasis();
470 
471  //- Update Gz_ before the potential orthonormal basis compression
472  void updateGz();
473 
474  //- Compress orthonormal basis for Qz_ and Gz_ if (maxRank_ < Qz_.n())
475  void compressBasis();
476 
477 
478  // Postprocess
479 
480  //- Return a new List containing elems of List at indices
481  template<class Type>
482  void filterIndexed
483  (
484  List<Type>& lst,
485  const UList<label>& indices
486  );
487 
488  //- Return a new Matrix containing columns of Matrix at indices
489  template<class MatrixType>
490  void filterIndexed
491  (
492  MatrixType& lst,
493  const UList<label>& indices
494  );
495 
496  //- Compute global Ap (K:Eq. 78)
497  void calcAp();
498 
499  //- Compute eigenvalues and eigenvectors
500  void calcEigen();
501 
502  //- Weak-type floating-point comparison
503  // bit.ly/2Trdbgs (Eq. 2), and PEP-485
504  bool close
505  (
506  const scalar s1,
507  const scalar s2,
508  const scalar absTol = 0, //<! comparisons near zero
509  const scalar relTol = 1e-8 //<! e.g. vals the same within 8 decimals
510  ) const;
511 
512  //- Test real/complex eigenvalues by using
513  //- the theoretical relation: (sum(eigenvalues) - trace(A) ~ 0)
514  void testEigenvalues
515  (
516  const SquareMatrix<scalar>& A,
517  const DiagonalMatrix<scalar>& EValsRe
518  ) const;
519 
520  //- Test real eigenvectors by using the theoretical relation:
521  //- ((A & EVec - EVal*EVec).norm() ~ 0)
522  void testEigenvectors
523  (
524  const SquareMatrix<scalar>& A,
525  const DiagonalMatrix<scalar>& EValsRe,
526  const SquareMatrix<scalar>& EVecs
527  ) const;
528 
529  //- Test complex eigenvectors by using the theoretical relation:
530  //- ((A & EVec - EVal*EVec).norm() ~ 0)
531  void testEigenvectors
532  (
533  const SquareMatrix<scalar>& A,
534  const List<complex>& EVals,
535  const RectangularMatrix<complex>& EVecs
536  ) const;
537 
538  //- Remove any eigenvalues whose magnitude is smaller than
539  //- minMagEVal_ while keeping the order of elements the same
540  void filterEVals();
541 
542  //- Compute and filter frequencies (K:Eq. 81)
543  void calcFreqs();
544 
545  //- Compute frequency indices
546  // Locate indices where oFreqs_ are
547  // in [fMin_, fMax_], and store them in iFreqs_ indices
548  void calcFreqI();
549 
550  //- Compute amplitudes
551  void calcAmps();
552 
553  //- Compute magnitudes
554  // Includes advanced sorting methods using
555  // the chosen weighted amplitude scaling method
556  void calcMags();
557 
558  //- Compute the ordered magnitude indices
559  void calcMagI();
560 
561  //- Compute modes
562  void calcModes();
563 
564  //- Eigenvalue weighted amplitude scaling (KZ:Eq. 33)
565  //- Modified eigenvalue weighted amplitude scaling (K)
566  scalar sorter
567  (
568  const List<scalar>& weight,
569  const complex& amplitude,
570  const complex& eval,
571  const scalar modeNorm
572  ) const;
573 
574  //- Output file header information
575  virtual void writeFileHeader(Ostream& os) const;
576 
577  //- Filter objects according to iFreqs_ and iMags_
578  void filterOutput();
579 
580  //- Write unfiltered/filtered data
581  void writeOutput(OFstream& os) const;
582 
583  //- Compute STDMD output
584  void calcOutput();
585 
586 
587 public:
588 
589  //- Runtime type information
590  TypeName("STDMD");
591 
592 
593  // Constructors
594 
595  //- No default construct
596  STDMD() = delete;
597 
598  //- Construct from Time and dictionary
599  STDMD
600  (
601  const word& name,
602  const Time& runTime,
603  const dictionary& dict
604  );
605 
606  //- No copy construct
607  STDMD(const STDMD&) = delete;
608 
609  //- No copy assignment
610  void operator=(const STDMD&) = delete;
611 
612 
613  //- Destructor
614  virtual ~STDMD() = default;
615 
616 
617  // Member Functions
618 
619  //- Read STDMD settings
620  virtual bool read(const dictionary&);
621 
622  //- Execute STDMD
623  virtual bool execute();
624 
625  //- Write STDMD output
626  virtual bool write();
627 };
628 
629 
630 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
631 
632 } // End namespace functionObjects
633 } // End namespace Foam
634 
635 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
636 
637 #ifdef NoRepository
638  #include "STDMDTemplates.C"
639 #endif
640 
641 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
642 
643 #endif
644 
645 // ************************************************************************* //
runTime
engineTime & runTime
Definition: createEngineTime.H:13
writeFile.H
EigenMatrix.H
Foam::Enum< modeSorterType >
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
fvMeshFunctionObject.H
Foam::DynamicList< label >
QRMatrix.H
A
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
Foam::functionObjects::STDMD::read
virtual bool read(const dictionary &)
Read STDMD settings.
Definition: STDMD.C:1135
Foam::functionObjects::fvMeshFunctionObject
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Definition: fvMeshFunctionObject.H:64
Foam::DiagonalMatrix< scalar >
Foam::functionObjects::STDMD::TypeName
TypeName("STDMD")
Runtime type information.
Foam::RectangularMatrix< scalar >
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam::functionObjects::STDMD::~STDMD
virtual ~STDMD()=default
Destructor.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::complex
A complex number, similar to the C++ complex type.
Definition: complex.H:82
Foam::functionObjects::STDMD::operator=
void operator=(const STDMD &)=delete
No copy assignment.
Foam::OFstream
Output to file stream, using an OSstream.
Definition: OFstream.H:53
Foam::functionObjects::STDMD::STDMD
STDMD()=delete
No default construct.
Foam::SquareMatrix< scalar >
Foam::functionObject::name
const word & name() const
Return the name of this functionObject.
Definition: functionObject.C:131
Foam::functionObjects::STDMD
(Beta Release) STDMD (i.e. Streaming Total Dynamic Mode Decomposition) is a variant of a data-driven ...
Definition: STDMD.H:409
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: BitOps.H:63
Foam::UList< label >
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
STDMDTemplates.C
Foam::functionObjects::writeFile
Base class for writing single files from the function objects.
Definition: writeFile.H:119
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::functionObjects::STDMD::execute
virtual bool execute()
Execute STDMD.
Definition: STDMD.C:1207
Foam::functionObjects::STDMD::write
virtual bool write()
Write STDMD output.
Definition: STDMD.C:1251