Sum of Prod of Mod of Linear_abc402G

最新推荐文章于 2025-12-19 15:31:04 发布

原创最新推荐文章于 2025-12-19 15:31:04 发布 · 953 阅读

7 ·

CC 4.0 BY-SA版权

文章标签：

#算法

G - Sum of Prod of Mod of Linear

这个问题需要前置知识：下取整和（floor sum），广义下取整和（generalized floor sum）。这属于类欧几里得算法，因为他们的时间复杂度和欧几里得算法一样在log级别。

下取整和的推导和结论，结论在最后一张，可以直接去看结论：

广义下取整和的推导和下取整和思路类似，但需要引入Sp(x)和sp(x)两个函数：

错误指正：以上写的程序的终止条件是在归约步骤，a%m=0，则程序结束，这里有误，应该是在交换步骤，当a=0时，还要计算一个常数（第一项的值），计算完常数后不再继续递归，将其作为结果返回，故判断a=0是在交换步骤进行。

其实这样存在大量的重复计算，程序的效率还比较低，以上我们已经有的优化效率的方法有：
1.当q=0时，直接返回S(p,N)作为答案

2.当n<=0时，返回0作为答案

针对大量重复计算，可以采用记忆化的方式，发现，针对同一组n,m,a,b，当递归层数相同的时候，当前的n,m,a,b是相同的，因此记忆化可以用递归层数代替前四个参数，同时这也是因为递归层数不会很大，用memo[dep][p][q]来记忆递归到第dep层，p，q值如何时，f(n,m,a,b,p,q)的结果。

递归层数可以通过实验试探出，当n,m,a,b上限是1e6时，大概是60层，当他们上限是1e9时，会变为85层左右。

接下来给出广义下取整和 (generalized floor sum,gfs) 的程序，可以在P5170 【模板】类欧几里得算法 - 洛谷验证程序是否正确。

#include<bits/stdc++.h>
using namespace std;
using ll=long long;

const ll maxn=1e6+5,maxpq=10,mod=998244353,maxdep=85;
ll fact[maxpq+5],inv[maxpq+5],invden[maxpq+5];
ll s[maxpq+5][maxpq+5]={
    {0,1},
    {0,-1,1},
    {0,1,-3,2},
    {0,0,1,-2,1},
    {0,-1,0,10,-15,6},
    {0,0,-1,0,5,-6,2},
    {0,1,0,-7,0,21,-21,6},
    {0,0,2,0,-7,0,14,-12,3},
    {0,-3,0,20,0,-42,0,60,-45,10},
    {0,0,-3,0,10,0,-14,0,15,-10,2},
    {0,5,0,-33,0,66,0,-66,0,55,-33,6}
};
ll den[maxpq+5]={1,2,6,4,30,12,42,24,90,20,66};
ll memo[maxdep][4][4];  //由于p+q+1最大是3，故将memo的p,q维度定义为4

inline void init_memo(){
    memset(memo,-1,sizeof(memo));
}

inline ll ipow(ll a,ll n){
    ll res=1;
    a=a%mod;
    while(n){
        if(n&1) res=res*a%mod;
        a=a*a%mod;
        n>>=1;
    }
    return res;
}


inline ll com(ll n,ll m){
    if(m<0 || m>n) return 0;
    if(m==n || m==0) return 1;
    ll res=fact[n]*inv[n-m]%mod*inv[m]%mod;
    return res;
}


inline ll get_S(ll p,ll x){
    ll res=0;
    for(ll i=0;i<=p+1;i++) res=(res+(s[p][i]*ipow(x,i)%mod))%mod;
    res=res*invden[p]%mod;
    return res;
}

inline ll getinv(ll x){
    return ipow(x,mod-2);
}

inline ll gfs(ll n,ll m,ll a,ll b,ll p,ll q,ll dep){
    if(memo[dep][p][q]!=-1) {
        return memo[dep][p][q];
    }
    if(n<=0) return memo[dep][p][q]=0;
    if(q==0) return memo[dep][p][q]=get_S(p,n);
    if(a>=m || b>=m || b<0) {  //归约步骤
        ll a1=a/m,b1=b/m;
        ll a2=a%m;
        ll b2=b%m;
        ll ans=0;
        for(ll j=0;j<=q;j++){
            for(ll t=0;t<=j;t++){
                ans=(ans+com(q,j)*com(j,t)%mod*ipow(a1,t)%mod*ipow(b1,j-t)%mod*gfs(n,m,a2,b2,p+t,q-j,dep+1)%mod)%mod;
            }
        }
        memo[dep][p][q]=ans;
        return ans;
    }else { //交换步骤
        ll K=(a*(n-1)+b)/m;
        ll ans=ipow(K,q)*get_S(p,n)%mod;
        if(a<=0) return memo[dep][p][q]=ans;
        ll delta=0;
        for(ll i=0;i<=q-1;i++){
            for(ll j=0;j<=p+1;j++){
                if(s[p][j]!=0){     //因很多s[x][y]=0，故进行此优化终止递归
                    ll t=gfs(K,a,m,m-b+a-1,i,j,dep+1);
                    delta=(delta+com(q,i)*t%mod*s[p][j]%mod)%mod;
                }
            }
        }
        delta=delta*invden[p]%mod;
        ans=(ans-delta)%mod;
        memo[dep][p][q]=ans;
        return ans;
    }
    return 0;       //guard
}

int main()
{
   // freopen("D:/in.txt","r",stdin);
   // freopen("D:/out.txt","w",stdout);
    ios::sync_with_stdio(0);cin.tie(0);
    //预计算一些值
    fact[0]=1;
    for(ll i=1;i<=maxpq;i++) fact[i]=fact[i-1]*i%mod;
    for(ll i=0;i<=maxpq;i++) {
        inv[i]=ipow(fact[i],mod-2); //阶乘和阶乘的逆元用于快速计算组合数
        invden[i]=ipow(den[i],mod-2);   //经测试，预处理小s公分母的逆元能很大程度上提高效率
    }

    ll T;
    cin>>T;
    while(T--)
    {
        ll n,m,a,b,p,q;
        cin>>n>>a>>b>>m;
        n++;
        init_memo();    //以下几组n,m,a,b相同，故他们共用memo
        ll ans2=gfs(n,m,a,b,0,2,0);
        ll ans3=gfs(n,m,a,b,1,1,0);
        ll ans1=gfs(n,m,a,b,0,1,0);

        cout<<(ans1+mod)%mod<<" ";
        cout<<(ans2+mod)%mod<<" ";
        cout<<(ans3+mod)%mod<<" ";
        cout<<"\n";
    }

    return 0;
}

标题题目解析

建议先去看官方题解，因为给出的代码实现是基于官方题解的，以下的题解思路与官方完全相同。

在实现过程中，需要注意：
1.对同一组(n,m,a,b)使用gfs的时候，共用memo以提高效率。

2.这里的b可能出现为负数的情况，需要略微修改归约步骤来应对b<0的情况，注意：C++负数取模的方法。

3.中间运算过程可能超过long long 范围，使用__int128进行计算，但最终结果在long long范围内，可以选择将__int128转成long long输出。

以下的实现是基于官方题解的，和上面所写实现方式略有不同，但思路一样。

#include<bits/stdc++.h>
using namespace std;
using ll=__int128;
using ull=unsigned long long ;

const ll maxn=1e6+5,maxpq=10,maxdep=60;
ll fact[maxpq+5];
ll s[maxpq+5][maxpq+5]={
    {0,1},
    {0,-1,1},
    {0,1,-3,2},
    {0,0,1,-2,1},
    {0,-1,0,10,-15,6},
    {0,0,-1,0,5,-6,2},
    {0,1,0,-7,0,21,-21,6},
    {0,0,2,0,-7,0,14,-12,3},
    {0,-3,0,20,0,-42,0,60,-45,10},
    {0,0,-3,0,10,0,-14,0,15,-10,2},
    {0,5,0,-33,0,66,0,-66,0,55,-33,6}
};
ll den[maxpq+5]={1,2,6,4,30,12,42,24,90,20,66};

ll memo[maxdep][4][4];

inline void init_memo(){
    memset(memo,-1,sizeof(memo));
}

inline ll ipow(ll a,ll n){
    ll res=1;
    while(n){
        if(n&1) res=res*a;
        a=a*a;
        n>>=1;
    }
    return res;
}

inline ll com(ll n,ll m){
    if(m<0 || m>n) return 0;
    if(m==n || m==0) return 1;
    ll res=fact[n]/(fact[n-m]*fact[m]);
    return res;
}


inline ll get_S(ll p,ll x){
    ll res=0;
    for(ll i=0;i<=p+1;i++) res+=s[p][i]*ipow(x,i);
    res/=den[p];
    return res;
}

inline ll gfs(ll n,ll m,ll a,ll b,ll p,ll q,ll dep){
    if(memo[dep][p][q]!=-1) {
        return memo[dep][p][q];
    }
    if(n<=0) return memo[dep][p][q]=0;
    if(q==0) return memo[dep][p][q]=get_S(p,n);
    if(a>=m || b>=m || b<0) {  //归约步骤
        ll a1=a/m,b1=b/m;
        if(b<0 && b%m!=0) b1--;
        ll a2=a%m;
        ll b2=(b%m+m)%m;
        ll ans=0;
        for(ll j=0;j<=q;j++){
            for(ll t=0;t<=j;t++){
                ans+=com(q,j)*com(j,t)*ipow(a1,t)*ipow(b1,j-t)*gfs(n,m,a2,b2,p+t,q-j,dep+1);
            }
        }
        memo[dep][p][q]=ans;
        return ans;
    }else { //交换步骤
        ll K=(a*(n-1)+b)/m;
        ll ans=ipow(K,q)*get_S(p,n);
        if(a<=0) return memo[dep][p][q]=ans;
        ll delta=0;
        for(ll i=0;i<=q-1;i++){
            for(ll j=0;j<=p+1;j++){
                if(s[p][j]!=0){
                    ll t=gfs(K,a,m,m-b+a-1,i,j,dep+1);
                    delta+=com(q,i)*t*s[p][j];
                }
            }
        }
        delta/=den[p];
        ans-=delta;
        memo[dep][p][q]=ans;
        return ans;
    }
    return 0;       //guard
}

istream& operator >> (istream &in,ll &x) {
    long long t;
    in>>t;
    x=t;
    return in;
}

ostream& operator << (ostream &out, const ll &x) {
        long long t=x;
        out<<t;
        return out;
}

int main()
{
    ios::sync_with_stdio(0);cin.tie(0);
    //预计算一些值
    fact[0]=1;
    for(ll i=1;i<=maxpq;i++) fact[i]=fact[i-1]*i;

    ll T;cin>>T;
    while(T--)
    {
        ll n,m,a,b1,b2;
        cin>>n>>m>>a>>b1>>b2;
        if(b1>b2) swap(b1,b2);
        if(a==0) {
            cout<<n*b1*b2<<"\n";
            continue;
        }
        ll res1=a*a*(n-1)*n/2*(2*n-1)/3+a*b2*n*(n-1)/2+a*b1*n*(n-1)/2+b1*b2*n;

        init_memo();
        ll res2=a*gfs(n,m,a,b2,1,1,0)+b1*gfs(n,m,a,b2,0,1,0);
        ll res3=gfs(n,m,a,b2,0,2,0);

        init_memo();
        res2+=a*gfs(n,m,a,b1,1,1,0)+b2*gfs(n,m,a,b1,0,1,0);
        res2*=m;
        ll delta;
        ll X=(a*(n-1)+b2)/m;
        ll dX_1=min((X*m-b1+a-1)/a,n);
        ll dX_2=min((X*m-b2+a-1)/a,n);

        init_memo();
        delta=-gfs(X,a,m,-b1+a-1,1,1,0);
        init_memo();
        delta+=gfs(X,a,m,-b2+a-1,1,1,0)-X*(dX_1-dX_2);
        res3+=delta;
        res3*=m*m;
        ll ans=res1-res2+res3;
        cout<<ans<<"\n";
    }

    return 0;
}