FreeFOAM The Cross-Platform CFD Toolkit
fvFieldDecomposer.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 \*---------------------------------------------------------------------------*/
25 
26 #include "fvFieldDecomposer.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 (
37  const unallocLabelList& addressingSlice,
38  const label addressingOffset
39 )
40 :
41  directAddressing_(addressingSlice)
42 {
43  forAll (directAddressing_, i)
44  {
45  // Subtract one to align addressing.
46  directAddressing_[i] -= addressingOffset + 1;
47  }
48 }
49 
50 
51 fvFieldDecomposer::processorVolPatchFieldDecomposer::
52 processorVolPatchFieldDecomposer
53 (
54  const fvMesh& mesh,
55  const unallocLabelList& addressingSlice
56 )
57 :
58  directAddressing_(addressingSlice.size())
59 {
60  const labelList& own = mesh.faceOwner();
61  const labelList& neighb = mesh.faceNeighbour();
62 
63  forAll (directAddressing_, i)
64  {
65  // Subtract one to align addressing.
66  label ai = mag(addressingSlice[i]) - 1;
67 
68  if (ai < neighb.size())
69  {
70  // This is a regular face. it has been an internal face
71  // of the original mesh and now it has become a face
72  // on the parallel boundary.
73  // Give face the value of the neighbour.
74 
75  if (addressingSlice[i] >= 0)
76  {
77  // I have the owner so use the neighbour value
78  directAddressing_[i] = neighb[ai];
79  }
80  else
81  {
82  directAddressing_[i] = own[ai];
83  }
84  }
85  else
86  {
87  // This is a face that used to be on a cyclic boundary
88  // but has now become a parallel patch face. I cannot
89  // do the interpolation properly (I would need to look
90  // up the different (face) list of data), so I will
91  // just grab the value from the owner cell
92 
93  directAddressing_[i] = own[ai];
94  }
95  }
96 }
97 
98 
99 fvFieldDecomposer::processorSurfacePatchFieldDecomposer::
100 processorSurfacePatchFieldDecomposer
101 (
102  const unallocLabelList& addressingSlice
103 )
104 :
105  addressing_(addressingSlice.size()),
106  weights_(addressingSlice.size())
107 {
108  forAll (addressing_, i)
109  {
110  addressing_[i].setSize(1);
111  weights_[i].setSize(1);
112 
113  addressing_[i][0] = mag(addressingSlice[i]) - 1;
114  weights_[i][0] = sign(addressingSlice[i]);
115  }
116 }
117 
118 
119 fvFieldDecomposer::fvFieldDecomposer
120 (
121  const fvMesh& completeMesh,
122  const fvMesh& procMesh,
123  const labelList& faceAddressing,
124  const labelList& cellAddressing,
125  const labelList& boundaryAddressing
126 )
127 :
128  completeMesh_(completeMesh),
129  procMesh_(procMesh),
130  faceAddressing_(faceAddressing),
131  cellAddressing_(cellAddressing),
132  boundaryAddressing_(boundaryAddressing),
133  patchFieldDecomposerPtrs_
134  (
135  procMesh_.boundary().size(),
136  static_cast<patchFieldDecomposer*>(NULL)
137  ),
138  processorVolPatchFieldDecomposerPtrs_
139  (
140  procMesh_.boundary().size(),
141  static_cast<processorVolPatchFieldDecomposer*>(NULL)
142  ),
143  processorSurfacePatchFieldDecomposerPtrs_
144  (
145  procMesh_.boundary().size(),
146  static_cast<processorSurfacePatchFieldDecomposer*>(NULL)
147  )
148 {
149  forAll (boundaryAddressing_, patchi)
150  {
151  if (boundaryAddressing_[patchi] >= 0)
152  {
153  patchFieldDecomposerPtrs_[patchi] = new patchFieldDecomposer
154  (
155  procMesh_.boundary()[patchi].patchSlice(faceAddressing_),
156  completeMesh_.boundaryMesh()
157  [
158  boundaryAddressing_[patchi]
159  ].start()
160  );
161  }
162  else
163  {
164  processorVolPatchFieldDecomposerPtrs_[patchi] =
165  new processorVolPatchFieldDecomposer
166  (
167  completeMesh_,
168  procMesh_.boundary()[patchi].patchSlice(faceAddressing_)
169  );
170 
171  processorSurfacePatchFieldDecomposerPtrs_[patchi] =
172  new processorSurfacePatchFieldDecomposer
173  (
174  static_cast<const unallocLabelList&>
175  (
176  procMesh_.boundary()[patchi].patchSlice
177  (
178  faceAddressing_
179  )
180  )
181  );
182  }
183  }
184 }
185 
186 
187 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
188 
189 fvFieldDecomposer::~fvFieldDecomposer()
190 {
191  forAll (patchFieldDecomposerPtrs_, patchi)
192  {
193  if (patchFieldDecomposerPtrs_[patchi])
194  {
195  delete patchFieldDecomposerPtrs_[patchi];
196  }
197  }
198 
199  forAll (processorVolPatchFieldDecomposerPtrs_, patchi)
200  {
201  if (processorVolPatchFieldDecomposerPtrs_[patchi])
202  {
203  delete processorVolPatchFieldDecomposerPtrs_[patchi];
204  }
205  }
206 
207  forAll (processorSurfacePatchFieldDecomposerPtrs_, patchi)
208  {
209  if (processorSurfacePatchFieldDecomposerPtrs_[patchi])
210  {
211  delete processorSurfacePatchFieldDecomposerPtrs_[patchi];
212  }
213  }
214 }
215 
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 } // End namespace Foam
220 
221 // ************************ vim: set sw=4 sts=4 et: ************************ //