递归算法真的要比迭代慢吗？

bcdmwSjy · 2026-02-02 13:32:00 · 算法·理论

在初学 FFT 时，我们可能听到别人说：递归 FFT 效率低，迭代 FFT 跑得快。你问为什么，他们便说：递归时需要维护一个栈，当然比迭代慢。你那时觉得很有道理，就把这句话在心中记下了。

但是后来，你发现矩阵乘法交换循环顺序可能快上一倍时，你了解到运算次数不是影响运行效率的唯一因素，内存访问也是影响运行效率的重要部分。

现代 CPU 的运算速度非常快，每秒可以处理几十亿条指令，虽然内存的速度也不慢，但是相对于 CPU 来说还是太慢了，所以人们发明了 Cache 来加速内存访问，以匹配 CPU 的运算速度。但是因为它的速度非常快，这就导致它的大小不会太大（最快的 L1 缓存一般在 512KB 到 1MB 左右），后来，为了尽量避免数据被移除缓存，人们又给电脑增加了 L2 和 L3 缓存。

言归正传，我们先来对比一下递归与迭代 FFT 的代码（以 dif 为例，未经过正确性测试，仅用于效率测试参考）：

constexpr double pi=numbers::pi;

vector<vector<complex<double>>> f;
void init(int n){
    n>>=1;
    int x=__lg(n);
    f.resize(x+1);
    f[x].resize(n);
    f[x][0]=1;
    for (int i=1;i<n;i+=i) f[x][i]=std::polar(1.,i*pi/n);
    for (int i=1;i<n;i++) f[x][i]=f[x][i&(i-1)]*f[x][i&-i];
    while (n>>=1){
        f[--x].resize(n);
        for (int i=0;i<n;i++) f[x][i]=f[x+1][i+i];
    }
}
void dif(complex<double> *a,int n){
    for (int m=n>>1;m;m>>=1){
        const auto &g=f[__lg(m)];
        for (int j=0;j<n;j+=m+m){
            for (int i=0;i<m;i++){
                const auto x=a[i+j],y=a[i+j+m];
                a[i+j]=x+y,a[i+j+m]=(x-y)*g[i];
            }
        }
    }
}
void dif_recursive(complex<double> *a,int n){
    if (n==1) return;
    int m=n>>1;
    const auto &g=f[__lg(m)];
    for (int i=0;i<m;i++){
        const auto x=a[i],y=a[i+m];
        a[i]=x+y,a[i+m]=(x-y)*g[i];
    }
    dif_recursive(a,m);
    dif_recursive(a+m,m);
}

我们发现迭代版本的 dif 每一层都要遍历整个数组，这导致每个元素都会反复经历被加入缓存，又被移出缓存的过程，这就使程序运行效率变低了许多。

反观递归版本的代码，它处理完本层后先处理左侧的分支，当区间长度小到可以放入 L1 缓存时，那些数据就可以一直在 L1 缓存中处理，而不会被频繁移动，这对缓存非常友好。实测下来确实如此，递归版本比非递归版本的代码快了近 60\%。

我们考虑进一步优化，递归确实有系统栈的开销，我们还要尽可能减小这一部分的开销，手动模拟吗？太复杂了，反而可能是程序变慢。

我们发现 dif 函数的参数 n 必须是 2 的幂次，换句话说，只有 \log n 中本质不同的 dif 函数，这启示我们使用 C++ 中的模板递归解决（严格来讲这已经不算递归了，因为长度不同 dif 是不同的函数）。

我们可以这样实现：

template<int n>
void dif_t(complex<double> *a){
    constexpr int m=n>>1;
    const auto &g=f[__lg(m)];
    for (int i=0;i<m;i++){
        const auto x=a[i],y=a[i+m];
        a[i]=x+y,a[i+m]=(x-y)*g[i];
    }
    dif_t<m>(a);
    dif_t<m>(a+m);
}
template<>
void dif_t<1>(complex<double> *a){}

这非常好，只不过模板尖括号内的内容必须是编译期常量，变量是不能写在模板里的，不过假设长度不超过 2^{24}，我们便可以这样写。

void dif_template(complex<double> *a,int n){
    switch(n){case 1:{dif_t<1>(a);break;}case 2:{dif_t<2>(a);break;}case 4:{dif_t<4>(a);break;}case 8:{dif_t<8>(a);break;}case 16:{dif_t<16>(a);break;}case 32:{dif_t<32>(a);break;}case 64:{dif_t<64>(a);break;}case 128:{dif_t<128>(a);break;}case 256:{dif_t<256>(a);break;}case 512:{dif_t<512>(a);break;}case 1024:{dif_t<1024>(a);break;}case 2048:{dif_t<2048>(a);break;}case 4096:{dif_t<4096>(a);break;}case 8192:{dif_t<8192>(a);break;}case 16384:{dif_t<16384>(a);break;}case 32768:{dif_t<32768>(a);break;}case 65536:{dif_t<65536>(a);break;}case 131072:{dif_t<131072>(a);break;}case 262144:{dif_t<262144>(a);break;}case 524288:{dif_t<524288>(a);break;}case 1048576:{dif_t<1048576>(a);break;}case 2097152:{dif_t<2097152>(a);break;}case 4194304:{dif_t<4194304>(a);break;}case 8388608:{dif_t<8388608>(a);break;}case 16777216:{dif_t<16777216>(a);break;}}
}

这样，我们就使 dif 又快了 13\%。

通过以上的探究，我们得出结论：同一个算法递归实现不一定比迭代实现慢，评判运行速度不仅要考虑计算量，还要考虑访存的效率，如果本质不同的递归函数较少，还可以使用模板递归进一步提升效率。

附上测速代码和数据：

static inline ull GetCycleCount(){
    unsigned hi,lo;
    __asm__ volatile("rdtsc":"=a"(lo),"=d"(hi));
    return (ull(hi)<<32)|lo;
}

int main(){
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cout.tie(nullptr);
    int n=1<<24;
    vector<complex<double>> a(n);
    for (int i=0;i<n;i++) a[i]=i*1e-6;
    init(n);
    ull st=GetCycleCount();
    // dif(a.data(),n);
    // dif_recursive(a.data(),n);
    dif_template(a.data(),n);
    cout<<(GetCycleCount()-st)*1e-9<<"G"<<"\n";
    return 0;
}