Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
randomProcesses
fft
fftRenumber.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
Description
25
Multi-dimensional renumbering used in the Numerical Recipes
26
fft routine. This version is recursive, so works in n-d :
27
determined by the length of array nn
28
29
\*---------------------------------------------------------------------------*/
30
31
#include <
randomProcesses/fftRenumber.H
>
32
33
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34
35
namespace
Foam
36
{
37
38
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
39
40
// recursively evaluate the indexing necessary to do the folding of
41
// the fft data. We recurse until we have the indexing ncessary for
42
// the folding in all directions.
43
44
void
fftRenumberRecurse
45
(
46
List<complex>
& data,
47
List<complex>
& renumData,
48
const
labelList
& nn,
49
label nnprod,
50
label ii,
51
label l1,
52
label l2
53
)
54
{
55
if
(ii == nn.
size
())
56
{
57
// we've worked out the renumbering scheme. Now copy
58
// the components across
59
60
data[l1] =
complex
(renumData[l2].
Re
(),renumData[l2].
Im
());
61
}
62
else
63
{
64
// do another level of folding. First work out the
65
// multiplicative value of the index
66
67
nnprod /= nn[ii];
68
label i_1(0);
69
70
for
(label i=0; i<nn[ii]; i++)
71
{
72
// now evaluate the indices (both from array 1 and to
73
// array 2). These get multiplied by nnprod to (cumulatively)
74
// find the real position in the list corresponding to
75
// this set of indices.
76
77
if
(i<nn[ii]/2)
78
{
79
i_1 = i + nn[ii]/2;
80
}
81
else
82
{
83
i_1 = i - nn[ii]/2;
84
}
85
86
87
// go to the next level of recursion.
88
89
fftRenumberRecurse
90
(
91
data,
92
renumData,
93
nn,
94
nnprod,
95
ii+1,
96
l1+i*nnprod,
97
l2+i_1*nnprod
98
);
99
}
100
}
101
}
102
103
104
// fftRenumber : fold the n-d data array to get the fft components in
105
// the right places.
106
107
#include <
randomProcesses/fftRenumber.H
>
108
109
void
fftRenumber
110
(
111
List<complex>
& data,
112
const
labelList
& nn
113
)
114
{
115
List<complex>
renumData(data);
116
117
label nnprod(1);
118
for
(label i=0; i<nn.
size
(); i++)
119
{
120
nnprod *= nn[i];
121
}
122
123
label ii(0), l1(0), l2(0);
124
125
fftRenumberRecurse
126
(
127
data,
128
renumData,
129
nn,
130
nnprod,
131
ii,
132
l1,
133
l2
134
);
135
}
136
137
138
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
139
140
}
// End namespace Foam
141
142
// ************************ vim: set sw=4 sts=4 et: ************************ //