FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

fftRenumber.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 Description
00025     Multi-dimensional renumbering used in the Numerical Recipes
00026    fft routine. This version is recursive, so works in n-d :
00027    determined by the length of array nn
00028 
00029 \*---------------------------------------------------------------------------*/
00030 
00031 #include <randomProcesses/fftRenumber.H>
00032 
00033 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00034 
00035 namespace Foam
00036 {
00037 
00038 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00039 
00040 // recursively evaluate the indexing necessary to do the folding of
00041 // the fft data. We recurse until we have the indexing ncessary for
00042 // the folding in all directions.
00043 
00044 void fftRenumberRecurse
00045 (
00046     List<complex>& data,
00047     List<complex>& renumData,
00048     const labelList& nn,
00049     label nnprod,
00050     label ii,
00051     label l1,
00052     label l2
00053 )
00054 {
00055     if (ii == nn.size())
00056     {
00057         // we've worked out the renumbering scheme. Now copy
00058         // the components across
00059 
00060         data[l1] = complex(renumData[l2].Re(),renumData[l2].Im());
00061     }
00062     else
00063     {
00064         // do another level of folding. First work out the
00065         // multiplicative value of the index
00066 
00067         nnprod /= nn[ii];
00068         label i_1(0);
00069 
00070         for (label i=0; i<nn[ii]; i++)
00071         {
00072             // now evaluate the indices (both from array 1 and to
00073             // array 2). These get multiplied by nnprod to (cumulatively)
00074             // find the real position in the list corresponding to
00075             // this set of indices.
00076 
00077             if (i<nn[ii]/2)
00078             {
00079                 i_1 = i + nn[ii]/2;
00080             }
00081             else
00082             {
00083                 i_1 = i - nn[ii]/2;
00084             }
00085 
00086 
00087             // go to the next level of recursion.
00088 
00089             fftRenumberRecurse
00090             (
00091                 data,
00092                 renumData,
00093                 nn,
00094                 nnprod,
00095                 ii+1,
00096                 l1+i*nnprod,
00097                 l2+i_1*nnprod
00098             );
00099         }
00100     }
00101 }
00102 
00103 
00104 // fftRenumber : fold the n-d data array to get the fft components in
00105 // the right places.
00106 
00107 #include <randomProcesses/fftRenumber.H>
00108 
00109 void fftRenumber
00110 (
00111     List<complex>& data,
00112     const labelList& nn
00113 )
00114 {
00115     List<complex> renumData(data);
00116 
00117     label nnprod(1);
00118     for (label i=0; i<nn.size(); i++)
00119     {
00120         nnprod *= nn[i];
00121     }
00122 
00123     label ii(0), l1(0), l2(0);
00124 
00125     fftRenumberRecurse
00126     (
00127         data,
00128         renumData,
00129         nn,
00130         nnprod,
00131         ii,
00132         l1,
00133         l2
00134     );
00135 }
00136 
00137 
00138 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00139 
00140 } // End namespace Foam
00141 
00142 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines