STDMDTemplates.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 | 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 \*---------------------------------------------------------------------------*/
27 
28 #include "volFields.H"
29 #include "surfaceFields.H"
30 
31 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
32 
33 template<class GeoFieldType>
34 bool Foam::functionObjects::STDMD::getSnapshot()
35 {
36  if (!initialised_)
37  {
38  init();
39  }
40 
41  // Move previous-time snapshot into previous-time slot in z_
42  // Effectively moves the lower half of z_ to its upper half
43  std::rotate(z_.begin(), z_.begin() + nSnap_, z_.end());
44 
45  // Copy new current-time snapshot into current-time slot in z_
46  // Effectively copies the new field elements into the lower half of z_
47  const GeoFieldType& Field = lookupObject<GeoFieldType>(fieldName_);
48  const label nField = Field.size();
49 
50  for (direction dir = 0; dir < nComps_; ++dir)
51  {
52  z_.subColumn(0, nSnap_ + dir*nField, nField) = Field.component(dir);
53  }
54 
55  return true;
56 }
57 
58 
59 template<class Type>
60 bool Foam::functionObjects::STDMD::getSnapshotType()
61 {
62  typedef GeometricField<Type, fvPatchField, volMesh> VolFieldType;
63  typedef GeometricField<Type, fvsPatchField, surfaceMesh> SurfaceFieldType;
64 
65  if (foundObject<VolFieldType>(fieldName_))
66  {
67  return getSnapshot<VolFieldType>();
68  }
69  else if (foundObject<SurfaceFieldType>(fieldName_))
70  {
71  return getSnapshot<SurfaceFieldType>();
72  }
73 
74  return false;
75 }
76 
77 
78 template<class Type>
79 bool Foam::functionObjects::STDMD::getComps()
80 {
81  typedef GeometricField<Type, fvPatchField, volMesh> VolFieldType;
82  typedef GeometricField<Type, fvsPatchField, surfaceMesh> SurfaceFieldType;
83 
84  if (foundObject<VolFieldType>(fieldName_))
85  {
86  nComps_ = pTraits<typename VolFieldType::value_type>::nComponents;
87  return true;
88  }
89  else if (foundObject<SurfaceFieldType>(fieldName_))
90  {
91  nComps_ = pTraits<typename SurfaceFieldType::value_type>::nComponents;
92  return true;
93  }
94 
95  return false;
96 }
97 
98 
99 template<class Type>
100 void Foam::functionObjects::STDMD::filterIndexed
101 (
102  List<Type>& lst,
103  const UList<label>& indices
104 )
105 {
106  // Elems within [a, b]
107  List<Type> lstWithin(indices.size());
108 
109  // Copy if frequency of elem is within [a, b]
110  label j = 0;
111  for (const auto& i : indices)
112  {
113  lstWithin[j] = lst[i];
114  ++j;
115  }
116  lst.transfer(lstWithin);
117 }
118 
119 
120 template<class MatrixType>
121 void Foam::functionObjects::STDMD::filterIndexed
122 (
123  MatrixType& mat,
124  const UList<label>& indices
125 )
126 {
127  // Elems within [a, b]
128  MatrixType matWithin(labelPair(mat.m(), indices.size()));
129 
130  // Copy if frequency of elem is within [a, b]
131  label j = 0;
132  for (const auto& i : indices)
133  {
134  matWithin.subColumn(j) = mat.subColumn(i);
135  ++j;
136  }
137  mat.transfer(matWithin);
138 }
139 
140 
141 // ************************************************************************* //
volFields.H
surfaceFields.H
Foam::surfaceFields.
Foam::labelPair
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:54
Foam::Matrix::begin
iterator begin()
Return an iterator to begin traversing a Matrix.
Definition: MatrixI.H:513
Foam::Matrix::end
iterator end()
Return an iterator to end traversing a Matrix.
Definition: MatrixI.H:521
Foam::direction
uint8_t direction
Definition: direction.H:52
Foam::Matrix::subColumn
ConstMatrixBlock< mType > subColumn(const label colIndex, const label rowIndex=0, label len=-1) const
Return const column or column's subset of Matrix.
Definition: MatrixI.H:267