import java.util.Arrays;public class FFT {private static final int N = 1024 * 1024;private static final double[] SIN = new double[N];private static final double[] COS = new double[N];static {if (Integer.highestOneBit(N) != N)throw new IllegalArgumentException("N is not a power of 2");for (int i = 0; i < N; i++) {double kth = i * -Math.PI / N;SIN[i] = Math.sin(kth);COS[i] = Math.cos(kth);}}private static void mult(double[] a, int k, double r, double i) {double ar = a[k], ai = a[k + 1];a[k] = ar * r - ai * i;a[k + 1] = ar * i + ai * r;}private static double[] fft(double[] x, int offset, int stride, int count) {if (count == 2)return new double[]{x[offset], x[offset + 1]};double[] q = fft(x, offset, stride * 2, count >> 1);double[] r = fft(x, offset + stride, stride * 2, count >> 1);double[] y = new double[count];count >>= 1;for (int i = 0, j = N * 2 / count, k = 0; k < count; k++, i += j) {mult(r, k, COS[i], SIN[i]);y[k] = q[k] + r[k];y[k + count] = q[k] - r[k];k++;y[k] = q[k] + r[k];y[k + count] = q[k] - r[k];}return y;}public static double[] fft(double[] ds) {if (Integer.highestOneBit(ds.length) != ds.length)throw new IllegalArgumentException("x.length is not a power of 2");return fft(ds, 0, 2, ds.length);}public static double[] ifft(double[] ds) {if (Integer.highestOneBit(ds.length) != ds.length)throw new IllegalArgumentException("x.length is not a power of 2");int n = ds.length;for (int i = 1; i < n; i += 2)ds[i] = -ds[i];ds = fft(ds);double r = 2.0 / n;for (int i = 0; i < n; i += 2) {ds[i] *= r;ds[i + 1] *= -r;}return ds;}public static double[] cconvolve(double[] da, double[] db) {int n = da.length;if (n != db.length)throw new IllegalArgumentException("x.length != y.length");double[] a = fft(da);double[] b = fft(db);for (int i = 0; i < n; i += 2)mult(a, i, b[i], b[i + 1]);return ifft(a);}public static double[] convolve(double[] da, double[] db) {return cconvolve(Arrays.copyOf(da, da.length * 2), Arrays.copyOf(db, db.length * 2));}public static double warmup() {double[] x = new double[N * 2];for (int i = 0; i < N * 2; i += 2)x[i] = (double)i / -N + 1;double[] y = fft(x);double[] z = ifft(y);double[] c = cconvolve(x, x);double[] d = convolve(x, x);double r = 0;for (int i = 0; i < N * 2; i += 2)r += y[i] + z[i] + c[i] + d[i];return r;}public static void main(String[] args) {for (int i = 0; i < 5; i++)warmup();for (int j = 0; j < 5; j++) {long begintime = System.currentTimeMillis();double r = warmup();System.out.println(System.currentTimeMillis() - begintime);System.out.println(r); // 输出结果,避免计算过程被编译器优化掉}}}
此处可能存在不合适展示的内容,页面不予展示。您可通过相关编辑功能自查并修改。
如您确认内容无涉及 不当用语 / 纯广告导流 / 暴力 / 低俗色情 / 侵权 / 盗版 / 虚假 / 无价值内容或违法国家有关法律法规的内容,可点击提交进行申诉,我们将尽快为您处理。