FreeFOAM The Cross-Platform CFD Toolkit
UOprocess.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 <OpenFOAM/error.H>
27 
28 #include "UOprocess.H"
29 #include <randomProcesses/Kmesh.H>
30 #include <OpenFOAM/dictionary.H>
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
38 
39 complexVector UOprocess::WeinerProcess()
40 {
41  return RootDeltaT*complexVector
42  (
43  complex(GaussGen.GaussNormal(), GaussGen.GaussNormal()),
44  complex(GaussGen.GaussNormal(), GaussGen.GaussNormal()),
45  complex(GaussGen.GaussNormal(), GaussGen.GaussNormal())
46  );
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
52 // from components
54 (
55  const Kmesh& kmesh,
56  const scalar deltaT,
57  const dictionary& UOdict
58 )
59 :
60  GaussGen(label(0)),
61  Mesh(kmesh),
62  DeltaT(deltaT),
63  RootDeltaT(sqrt(DeltaT)),
64  UOfield(Mesh.size()),
65 
66  Alpha(readScalar(UOdict.lookup("UOalpha"))),
67  Sigma(readScalar(UOdict.lookup("UOsigma"))),
68  Kupper(readScalar(UOdict.lookup("UOKupper"))),
69  Klower(readScalar(UOdict.lookup("UOKlower"))),
70  Scale((Kupper - Klower)*pow(scalar(Mesh.size()), 1.0/vector::dim))
71 {
72  const vectorField& K = Mesh;
73 
74  scalar sqrKupper = sqr(Kupper);
75  scalar sqrKlower = sqr(Klower) + SMALL;
76  scalar sqrK;
77 
78  forAll(UOfield, i)
79  {
80  if ((sqrK = magSqr(K[i])) < sqrKupper && sqrK > sqrKlower)
81  {
82  UOfield[i] = Scale*Sigma*WeinerProcess();
83  }
84  else
85  {
86  UOfield[i] = complexVector
87  (
88  complex(0, 0),
89  complex(0, 0),
90  complex(0, 0)
91  );
92  }
93  }
94 }
95 
96 
97 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
98 
100 {
101  const vectorField& K = Mesh;
102 
103  label count = 0;
104  scalar sqrKupper = sqr(Kupper);
105  scalar sqrKlower = sqr(Klower) + SMALL;
106  scalar sqrK;
107 
108  forAll(UOfield, i)
109  {
110  if ((sqrK = magSqr(K[i])) < sqrKupper && sqrK > sqrKlower)
111  {
112  count++;
113  UOfield[i] =
114  (1.0 - Alpha*DeltaT)*UOfield[i]
115  + Scale*Sigma*WeinerProcess();
116  }
117  }
118 
119  Info<< " Number of forced K = " << count << nl;
120 
121  return UOfield;
122 }
123 
124 
125 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
126 
127 } // End namespace Foam
128 
129 // ************************ vim: set sw=4 sts=4 et: ************************ //
130