高速フーリエ変換(FFT)
高速フーリエ変換は新しい変換ではなく、フーリエ変換の高速アルゴリズムです.このアルゴリズムは通常の離散Fourier変換(DFT)の時間的複雑度O(n*n)をOに下げることができる.(n log n)は、フーリエ変換の速度を大幅に向上させた.高速フーリエ変換アルゴリズムの提案により、フーリエ変換は通信分野において極めて大きな運用と発展を遂げた.ACMでは、高速フーリエ変換は通常、大数の乗算に用いられる.2つの数がJAVAのBIにも耐えられないほど大きい場合には、高速フーリエアルゴリズムを用いる.このアルゴリズムは、2つの数を、いずれも多項式行列として書き,次いで離散フーリエ変換を行った.通信分野には、時間領域のボリューム、周波数領域の乗算という古典的な言葉がある.従って,2つの数の多項式係数行列を離散フーリエ変換した後,直接乗算することができる.乗算した結果をフーリエ逆変換(DTFT)する.
########私はこれまでどうしても高速フーリエアルゴリズムを理解できませんでした.当時は複雑すぎて難しいと思っていたので、見ても読めませんでした.しかし、私の専攻は電子情報工学なので、各専攻の授業でフーリエ変換と会うのは避けられない.先週、私たちのデジタル信号処理の先生は、いくつかの授業の時間を費やして、高速フーリエアルゴリズムを説明してくれました.私はまじめに聞いて、自分で理解したと思っていました.しかし、自分の頭の中の高速フーリエアルゴリズムで大数乗算の問題を解きたいと思ったとき、手がつけられないことに気づいた.どうやって運用するか全然分かりません.その时、私はこのアルゴリズムをマスターしなければならないと思っていました.ACMに現れるからだけでなく、私の専攻でも極めて重要な地位を占めているからです.もし私がフーリエ変換さえできなかったら、これから自分が電子情報専門だと言う顔はありません.そこで、私は2日かけて、資料を調べて、やっとこの強大なアルゴリズムを理解しました.私はいくつかの大神が書いたFFTに関する文章が確かにとても良くて、はっきりしていて分かりやすくて、私自身もコレクションして、みんなに推薦して、みんなで一緒に勉強します.
超分かりやすいフーリエ変換
高速フーリエ変換(FFT)
###JAVAコード
###C++コード
########私はこれまでどうしても高速フーリエアルゴリズムを理解できませんでした.当時は複雑すぎて難しいと思っていたので、見ても読めませんでした.しかし、私の専攻は電子情報工学なので、各専攻の授業でフーリエ変換と会うのは避けられない.先週、私たちのデジタル信号処理の先生は、いくつかの授業の時間を費やして、高速フーリエアルゴリズムを説明してくれました.私はまじめに聞いて、自分で理解したと思っていました.しかし、自分の頭の中の高速フーリエアルゴリズムで大数乗算の問題を解きたいと思ったとき、手がつけられないことに気づいた.どうやって運用するか全然分かりません.その时、私はこのアルゴリズムをマスターしなければならないと思っていました.ACMに現れるからだけでなく、私の専攻でも極めて重要な地位を占めているからです.もし私がフーリエ変換さえできなかったら、これから自分が電子情報専門だと言う顔はありません.そこで、私は2日かけて、資料を調べて、やっとこの強大なアルゴリズムを理解しました.私はいくつかの大神が書いたFFTに関する文章が確かにとても良くて、はっきりしていて分かりやすくて、私自身もコレクションして、みんなに推薦して、みんなで一緒に勉強します.
超分かりやすいフーリエ変換
高速フーリエ変換(FFT)
###JAVAコード
import java.util.Arrays;
import java.util.Scanner;
import java.util.Vector;
/**
* Created by jal on 2017/12/6 0006.
*/
public class DIT_FFT {
private static int maxn;
public static void main(String[] args) {
Scanner scanner = new Scanner(System.in);
String numstr1 = scanner.next();
String numstr2 = scanner.next();
int len =numstr1.length()+numstr2.length();
maxn = 0;
double temp = log2(len);
double floortemp = Math.floor(temp);
if(floortemp == temp){
maxn = len;
}else {
maxn = (int)Math.pow(2,floortemp+1);
}
Complex []arra = createArray(maxn);
Complex []arrb =createArray(maxn);
Complex []arrc =createArray(maxn);
Complex []arrA = createArray(maxn);
Complex []arrB = createArray(maxn);
Complex []arrC = createArray(maxn);
char []a = numstr1.toCharArray();
int j = 0;
for(int i = a.length - 1; i >= 0; i--){
arra[j++] = new Complex(a[i] - '0');
}
j = 0;
char []b = numstr2.toCharArray();
for(int i = b.length - 1; i >= 0; i--){
arrb[j++] = new Complex(b[i] - '0');
}
//System.out.println(Arrays.toString(arra));
//invert(arra);
arrA = fft(arra);
System.out.println(Arrays.toString(arrA));
//invert(arrb);
arrB = fft(arrb);
//System.out.println(Arrays.toString(arrB));
for(int i = 0; i < arrC.length; i++){
arrC[i] = arrA[i].times(arrB[i]);
}
//System.out.println(Arrays.toString(arrC));
arrc = ifft(arrC);
Vector vector = new Vector<>();
vector = toIntOfString(arrc);
String str = "";
char vectorCHar[] = new char[vector.size()];
for(int i = 0; i < vectorCHar.length; i++){
vectorCHar[i] = (char)(vector.get(i) + '0');
}
str =String.valueOf(vectorCHar);
str = new StringBuffer(str).reverse().toString();
//System.out.println("str:"+str);
str = trim(str);
System.out.println(str);
}
private static String trim(String value) {
int len = value.length();
int st = 0;
char[] val = value.toCharArray(); /* avoid getfield opcode */
while ((st < len) && (val[st] <= '0')) {
st++;
}
// while ((st < len) && (val[len - 1] <= ' ')) {
// len--;
// }
return value.substring(st, len);
}
private static Vector toIntOfString(Complex[] arrc) {
Vector result = new Vector<>();
int n = arrc.length;
int arrayTemp [] = new int[n+1];
for(int i = 0; i < arrc.length; i++){
arrayTemp[i] = (int)arrc[i].re();
}
int i;
arrayTemp[n] = 0;
for(i = 0; i < n; i++){
result.add(arrayTemp[i]%10);
arrayTemp[i+1] += arrayTemp[i] / 10;
}
int t = arrayTemp[n];
while(t > 0){
result.add(t % 10);
t/= 10;
}
return result;
}
private static Complex[] ifft(Complex[] arrC) {
int n = arrC.length;
Complex arrayResult[] = new Complex[n];
for(int i = 0; i < arrC.length; i++){
arrC[i] = arrC[i].conjugate();
}
arrayResult= fft(arrC);
for(int i = 0; i < arrayResult.length; i++){
arrayResult[i] = arrayResult[i].conjugate().divides(new Complex(n));
}
return arrayResult;
}
private static Complex[] fft(Complex[] arrA) {
int N = arrA.length;
if(N == 1){
return arrA;
}
Complex [] arrayEven = new Complex[N/2];
Complex [] arrayOdd = new Complex[N/2];
for(int i = 0; i < arrayEven.length; i++){
arrayEven[i] = arrA[2*i];
arrayOdd[i] = arrA[2*i+1];
}
arrayEven = fft(arrayEven);
arrayOdd = fft(arrayOdd);
Complex[]arrayResult = new Complex[N];
for(int i = 0; i < N/2; i++){
Complex W = new Complex(0,-2*Math.PI*i/N).exp();
arrayResult[i] = arrayEven[i].plus(arrayOdd[i].times(W));
arrayResult[i+N/2] = arrayEven[i].minus(arrayOdd[i].times(W));
}
return arrayResult;
}
private static void bit_reverse_swap(Complex [] a) {
int n = a.length;
for (int i = 1, j = n >> 1, k; i < n - 1; ++i) {
if (i < j) swap(a,i,j);
// tricky
for (k = n >> 1; j >= k; j -= k, k >>= 1) // inspect the highest "1"
;
j += k;
}
}
private static void invert(Complex[] arra) {
int n = (int)log2(maxn);
for(int i = 0; i < arra.length; i++){
String temp = Integer.toBinaryString(i);
temp = new StringBuffer(temp).reverse().toString();
int len = n - temp.length();
char [] zeros = new char[len];
Arrays.fill(zeros,'0');
temp+=String.valueOf(zeros);
int j = Integer.valueOf(temp,2);
System.out.println("i:" + i +"j:" + j);
if(j>i){
swap(arra,i,j);
}
}
}
private static void swap(Complex[] arra, int i, int j) {
Complex tmp = arra[i];
arra[i] = arra[j];
arra[j] = tmp;
}
private static Complex[] createArray(int maxn) {
Complex[] result =new Complex[maxn];
for(int i=0;i
* For additional documentation, see Section 9.9 of
* Algorithms, 4th Edition by Robert Sedgewick and Kevin Wayne.
*
* @author Robert Sedgewick
* @author Kevin Wayne
*/
class Complex {
private double re; // the real part
private double im; // the imaginary part
/**
* Initializes a complex number from the specified real and imaginary parts.
*
* @param real the real part
* @param imag the imaginary part
*/
public Complex(double real, double imag) {
re = real;
im = imag;
}
public Complex(double re) {
this.re = re;
this.im = 0;
}
/**
* Returns a string representation of this complex number.
*
* @return a string representation of this complex number,
* of the form 34 - 56i.
*/
public String toString() {
if (im == 0) return re + "";
if (re == 0) return im + "i";
if (im < 0) return re + " - " + (-im) + "i";
return re + " + " + im + "i";
}
/**
* Returns the absolute value of this complex number.
* This quantity is also known as the modulus or magnitude.
*
* @return the absolute value of this complex number
*/
public double abs() {
return Math.hypot(re, im);
}
/**
* Returns the phase of this complex number.
* This quantity is also known as the angle or argument.
*
* @return the phase of this complex number, a real number between -pi and pi
*/
public double phase() {
return Math.atan2(im, re);
}
/**
* Returns the sum of this complex number and the specified complex number.
*
* @param that the other complex number
* @return the complex number whose value is {@code (this + that)}
*/
public Complex plus(Complex that) {
double real = this.re + that.re;
double imag = this.im + that.im;
return new Complex(real, imag);
}
/**
* Returns the result of subtracting the specified complex number from
* this complex number.
*
* @param that the other complex number
* @return the complex number whose value is {@code (this - that)}
*/
public Complex minus(Complex that) {
double real = this.re - that.re;
double imag = this.im - that.im;
return new Complex(real, imag);
}
/**
* Returns the product of this complex number and the specified complex number.
*
* @param that the other complex number
* @return the complex number whose value is {@code (this * that)}
*/
public Complex times(Complex that) {
double real = this.re * that.re - this.im * that.im;
double imag = this.re * that.im + this.im * that.re;
return new Complex(real, imag);
}
/**
* Returns the product of this complex number and the specified scalar.
*
* @param alpha the scalar
* @return the complex number whose value is {@code (alpha * this)}
*/
public Complex scale(double alpha) {
return new Complex(alpha * re, alpha * im);
}
/**
* Returns the product of this complex number and the specified scalar.
*
* @param alpha the scalar
* @return the complex number whose value is {@code (alpha * this)}
* @deprecated Replaced by {@link #scale(double)}.
*/
@Deprecated
public Complex times(double alpha) {
return new Complex(alpha * re, alpha * im);
}
/**
* Returns the complex conjugate of this complex number.
*
* @return the complex conjugate of this complex number
*/
public Complex conjugate() {
return new Complex(re, -im);
}
/**
* Returns the reciprocal of this complex number.
*
* @return the complex number whose value is {@code (1 / this)}
*/
public Complex reciprocal() {
double scale = re*re + im*im;
return new Complex(re / scale, -im / scale);
}
/**
* Returns the real part of this complex number.
*
* @return the real part of this complex number
*/
public double re() {
return re;
}
/**
* Returns the imaginary part of this complex number.
*
* @return the imaginary part of this complex number
*/
public double im() {
return im;
}
/**
* Returns the result of dividing the specified complex number into
* this complex number.
*
* @param that the other complex number
* @return the complex number whose value is {@code (this / that)}
*/
public Complex divides(Complex that) {
return this.times(that.reciprocal());
}
/**
* Returns the complex exponential of this complex number.
*
* @return the complex exponential of this complex number
*/
public Complex exp() {
return new Complex(Math.exp(re) * Math.cos(im), Math.exp(re) * Math.sin(im));
}
/**
* Returns the complex sine of this complex number.
*
* @return the complex sine of this complex number
*/
public Complex sin() {
return new Complex(Math.sin(re) * Math.cosh(im), Math.cos(re) * Math.sinh(im));
}
/**
* Returns the complex cosine of this complex number.
*
* @return the complex cosine of this complex number
*/
public Complex cos() {
return new Complex(Math.cos(re) * Math.cosh(im), -Math.sin(re) * Math.sinh(im));
}
/**
* Returns the complex tangent of this complex number.
*
* @return the complex tangent of this complex number
*/
public Complex tan() {
return sin().divides(cos());
}
}
###C++コード
#include
#include
#include
#include
#include
#include
#include
using namespace std;
double log2(int n){
return log(n)/log(2);
}
const int MAXN = 1 << 20;
complex arra[MAXN],arrb[MAXN],arrC[MAXN];
complex* arrA;
complex* arrB;//,arrC[MAXN];
complex* arrc;
const double PI = 3.141592653;
complex* DIT_FFT(complex* arr,int len){
if(len == 1)return arr;
complex *arrayEven = new complex[len/2];
complex *arrayOdd = new complex[len/2];
for(int i = 0; i < len/2; i++){
arrayEven[i] = arr[2*i];
arrayOdd[i] = arr[2*i+1];
}
arrayEven = DIT_FFT(arrayEven,len/2);
arrayOdd = DIT_FFT(arrayOdd,len/2);
complex *arrayResult = new complex[len];
for(int i = 0; i < len/2; i++){
//TODO
//help me
//help me
// what is 'W' ?
complex W = exp(complex(0,-2*PI *i / len));//ecomplex(cos(-2*PI *i / len),sin(-2*PI *i / len));
cout<* IFFT(complex*arrC,int len){
int n = len;
complex* arrayResult = new complex[len];
for(int i = 0; i < len; i++){
arrC[i] = conj(arrC[i]);
}
arrayResult= DIT_FFT(arrC,len);
for(int i = 0; i < len; i++){
arrayResult[i] = conj(arrayResult[i]);
arrayResult[i]/=n;
}
return arrayResult;
}
int main(){
string str1,str2;
cin>>str1>>str2;
int len = str1.size() + str2.size();
int maxn = 0;
double floortemp =log2(len);
if(floortemp == floor(floortemp)){
maxn = len;
}
else{
maxn = pow(2,(int)(floor(floortemp) + 1) );
}
for(int i = str1.size() - 1; i >= 0; i--){
arra[i] = complex(str1[i]-'0',0);
}
for(int i = str1.size() - 1; i >= 0; i--){
arrb[i] = complex(str2[i]-'0',0);
}
arrA = DIT_FFT(arra,maxn);
arrB = DIT_FFT(arrb,maxn);
for(int i= 0; i < maxn; i++){
arrC[i] = arrA[i]*arrB[i];
}
arrc = IFFT(arrC,len);
vectorv;
for(int i = 0; i < maxn;i++){
v.push_back((int)arrc[i].real());
}
}