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: ************************ //