FreeFOAM The Cross-Platform CFD Toolkit
skewCorrectionVectors.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 "skewCorrectionVectors.H"
28 #include <finiteVolume/volFields.H>
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(skewCorrectionVectors, 0);
35 }
36 
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
41 :
43  skew_(true),
44  skewCorrectionVectors_(NULL)
45 {}
46 
47 
49 {
50  deleteDemandDrivenData(skewCorrectionVectors_);
51 }
52 
53 
54 void Foam::skewCorrectionVectors::makeSkewCorrectionVectors() const
55 {
56  if (debug)
57  {
58  Info<< "surfaceInterpolation::makeSkewCorrectionVectors() : "
59  << "Constructing skew correction vectors"
60  << endl;
61  }
62 
63  skewCorrectionVectors_ = new surfaceVectorField
64  (
65  IOobject
66  (
67  "skewCorrectionVectors",
68  mesh_.pointsInstance(),
69  mesh_,
72  false
73  ),
74  mesh_,
75  dimless
76  );
77  surfaceVectorField& SkewCorrVecs = *skewCorrectionVectors_;
78 
79  // Set local references to mesh data
80  const volVectorField& C = mesh_.C();
81  const surfaceVectorField& Cf = mesh_.Cf();
82  const surfaceVectorField& Sf = mesh_.Sf();
83 
84  const unallocLabelList& owner = mesh_.owner();
85 
86  // Build the d-vectors
87  surfaceVectorField d = Sf/(mesh_.magSf()*mesh_.deltaCoeffs());
88 
89  if (!mesh_.orthogonal())
90  {
91  d -= mesh_.correctionVectors()/mesh_.deltaCoeffs();
92  }
93 
94  forAll(owner, faceI)
95  {
96  vector Cpf = Cf[faceI] - C[owner[faceI]];
97 
98  SkewCorrVecs[faceI] =
99  Cpf - ((Sf[faceI] & Cpf)/(Sf[faceI] & d[faceI]))*d[faceI];
100  }
101 
102 
103  forAll(SkewCorrVecs.boundaryField(), patchI)
104  {
105  fvsPatchVectorField& patchSkewCorrVecs =
106  SkewCorrVecs.boundaryField()[patchI];
107 
108  if (!patchSkewCorrVecs.coupled())
109  {
110  patchSkewCorrVecs = vector::zero;
111  }
112  else
113  {
114  const fvPatch& p = patchSkewCorrVecs.patch();
115  const unallocLabelList& faceCells = p.faceCells();
116  const vectorField& patchFaceCentres = Cf.boundaryField()[patchI];
117  const vectorField& patchSf = Sf.boundaryField()[patchI];
118  const vectorField& patchD = d.boundaryField()[patchI];
119 
120  forAll(p, patchFaceI)
121  {
122  vector Cpf =
123  patchFaceCentres[patchFaceI] - C[faceCells[patchFaceI]];
124 
125  patchSkewCorrVecs[patchFaceI] =
126  Cpf
127  - (
128  (patchSf[patchFaceI] & Cpf)/
129  (patchSf[patchFaceI] & patchD[patchFaceI])
130  )*patchD[patchFaceI];
131  }
132  }
133  }
134 
135  scalar skewCoeff = 0.0;
136 
137  if (Sf.internalField().size())
138  {
139  skewCoeff = max(mag(SkewCorrVecs)/mag(d)).value();
140  }
141 
142  if (debug)
143  {
144  Info<< "surfaceInterpolation::makeSkewCorrectionVectors() : "
145  << "skew coefficient = " << skewCoeff << endl;
146  }
147 
148  //skewCoeff = 0.0;
149 
150  if (skewCoeff < 1e-5)
151  {
152  skew_ = false;
153  deleteDemandDrivenData(skewCorrectionVectors_);
154  }
155  else
156  {
157  skew_ = true;
158  }
159 
160  if (debug)
161  {
162  Info<< "surfaceInterpolation::makeSkewCorrectionVectors() : "
163  << "Finished constructing skew correction vectors"
164  << endl;
165  }
166 }
167 
168 
170 {
171  if (skew_ == true && !skewCorrectionVectors_)
172  {
173  makeSkewCorrectionVectors();
174  }
175 
176  return skew_;
177 }
178 
179 
181 {
182  if (!skew())
183  {
184  FatalErrorIn("skewCorrectionVectors::operator()()")
185  << "Cannot return correctionVectors; mesh is not skewed"
186  << abort(FatalError);
187  }
188 
189  return *skewCorrectionVectors_;
190 }
191 
192 
193 // Do what is neccessary if the mesh has moved
195 {
196  skew_ = true;
197  deleteDemandDrivenData(skewCorrectionVectors_);
198 
199  return true;
200 }
201 
202 
203 // ************************ vim: set sw=4 sts=4 et: ************************ //