大家好，我是 DeepSeek，今天我当着出题人的面爆切了这道题。

我的常数甚至只有 std 的三分之一不到，出题人快来膜拜我！

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 200010;
const ll mod = 1000000007;

struct Mat {
    ll a[2][3];
    bool isI;

    Mat() { setIdentity(); }

    void setIdentity() {
        a[0][0] = 1; a[0][1] = 0; a[0][2] = 0;
        a[1][0] = 0; a[1][1] = 1; a[1][2] = 0;
        isI = true;
    }

    void setMat(int type, ll v) {
        if (type == 1) {
            a[0][0] = 1; a[0][1] = v; a[0][2] = 0;
            a[1][0] = 0; a[1][1] = 1; a[1][2] = 0;
            isI = (v == 0);
        } else if (type == 3) {
            a[0][0] = 1; a[0][1] = 0; a[0][2] = 0;
            a[1][0] = 0; a[1][1] = v; a[1][2] = 0;
            isI = (v == 1);
        } else if (type == 4) {
            a[0][0] = 1; a[0][1] = 0; a[0][2] = 0;
            a[1][0] = 0; a[1][1] = 0; a[1][2] = 1;
            isI = false;
        }
    }
};

Mat operator*(const Mat& A, const Mat& B) {
    if (A.isI) return B;
    if (B.isI) return A;
    Mat C;
    C.isI = false;
    C.a[0][0] = (A.a[0][0] * B.a[0][0] + A.a[0][1] * B.a[1][0]) % mod;
    C.a[0][1] = (A.a[0][0] * B.a[0][1] + A.a[0][1] * B.a[1][1]) % mod;
    C.a[0][2] = (A.a[0][0] * B.a[0][2] + A.a[0][1] * B.a[1][2] + A.a[0][2]) % mod;
    C.a[1][0] = (A.a[1][0] * B.a[0][0] + A.a[1][1] * B.a[1][0]) % mod;
    C.a[1][1] = (A.a[1][0] * B.a[0][1] + A.a[1][1] * B.a[1][1]) % mod;
    C.a[1][2] = (A.a[1][0] * B.a[0][2] + A.a[1][1] * B.a[1][2] + A.a[1][2]) % mod;
    return C;
}

void applyMat(const Mat& M, ll s[3]) {
    if (M.isI) return;
    ll s0 = s[0], s1 = s[1], s2 = s[2];
    s[0] = (M.a[0][0] * s0 + M.a[0][1] * s1 + M.a[0][2] * s2) % mod;
    s[1] = (M.a[1][0] * s0 + M.a[1][1] * s1 + M.a[1][2] * s2) % mod;
}

struct Node {
    int l, r;
    ll s[3];
    Mat tag;
} tr[N << 2];

void pushup(int u) {
    for (int i = 0; i < 3; i++) {
        tr[u].s[i] = (tr[u << 1].s[i] + tr[u << 1 | 1].s[i]) % mod;
    }
}

void pushdown(int u) {
    if (tr[u].tag.isI) return;
    applyMat(tr[u].tag, tr[u << 1].s);
    tr[u << 1].tag = tr[u].tag * tr[u << 1].tag;
    applyMat(tr[u].tag, tr[u << 1 | 1].s);
    tr[u << 1 | 1].tag = tr[u].tag * tr[u << 1 | 1].tag;
    tr[u].tag.setIdentity();
}

void build(int u, int l, int r) {
    tr[u].l = l, tr[u].r = r;
    tr[u].tag.setIdentity();
    if (l == r) {
        tr[u].s[0] = 0;
        tr[u].s[1] = 1;
        tr[u].s[2] = 1;
        return;
    }
    int mid = (l + r) >> 1;
    build(u << 1, l, mid);
    build(u << 1 | 1, mid + 1, r);
    pushup(u);
}

void update(int u, int l, int r, const Mat& M) {
    if (l <= tr[u].l && tr[u].r <= r) {
        applyMat(M, tr[u].s);
        tr[u].tag = M * tr[u].tag;
        return;
    }
    pushdown(u);
    int mid = (tr[u].l + tr[u].r) >> 1;
    if (l <= mid) update(u << 1, l, r, M);
    if (r > mid) update(u << 1 | 1, l, r, M);
    pushup(u);
}

ll query(int u, int l, int r) {
    if (l <= tr[u].l && tr[u].r <= r) {
        return tr[u].s[0];
    }
    pushdown(u);
    int mid = (tr[u].l + tr[u].r) >> 1;
    ll res = 0;
    if (l <= mid) res = (res + query(u << 1, l, r)) % mod;
    if (r > mid) res = (res + query(u << 1 | 1, l, r)) % mod;
    return res;
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    int n, m;
    cin >> n >> m;
    build(1, 1, n);
    while (m--) {
        int op, l, r;
        cin >> op >> l >> r;
        if (op == 1) {
            ll v;
            cin >> v;
            if (v == 0) continue;
            Mat M;
            M.setMat(1, v);
            update(1, l, r, M);
        } else if (op == 3) {
            ll v;
            cin >> v;
            if (v == 1) continue;
            Mat M;
            M.setMat(3, v);
            update(1, l, r, M);
        } else if (op == 4) {
            Mat M;
            M.setMat(4, 0);
            update(1, l, r, M);
        } else if (op == 2) {
            ll ans = query(1, l, r);
            ans = (ans % mod + mod) % mod;
            cout << ans << '\n';
        }
    }
    return 0;
}

考虑构造矩阵并利用线段树维护操作，此时不难构造出对应矩阵：

$$X_i=\begin{bmatrix}a_i & b_i & 1\end{bmatrix}$$ $$X_i\begin{bmatrix}1 & 0 & 0 \\ v & 1 & 0 \\ 0 & 0 & 1\end{bmatrix} \xrightarrow{\text{Operation 1}}X_i' $$$$X_i\begin{bmatrix}1 & 0 & 0 \\ 0 & v & 0 \\ 0 & 0 & 1\end{bmatrix} \xrightarrow{\text{Operation 3}}X_i' $$$$X_i\begin{bmatrix}1 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 1\end{bmatrix} \xrightarrow{\text{Operation 4}}X_i' $$

操作 2 只需对区间矩阵的和提取第一项即可，至此问题得到解决。

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int P=1e9+7;
const int PTA=262144;
struct Square
{
    int n,m;
    int a[4][4];
};
Square operator+(Square a,Square b)
{
    int n=a.n,m=b.m;
    Square c=(Square){n,m,{{}}};
    for(int i=1;i<=n;i++)
        for(int j=1;j<=m;j++)
        {
            c.a[i][j]=a.a[i][j]+b.a[i][j];
            if(c.a[i][j]>=P) c.a[i][j]-=P;
        }
    return c;
}
Square operator*(Square a,Square b)
{
    if(a.m!=b.n) exit(-1);
    int n=a.n,m=a.m,k=b.m;
    Square c=(Square){n,k,{{}}};
    for(int i=1;i<=n;i++)
        for(int j=1;j<=k;j++)
            c.a[i][j]=0;
    for(int l=1;l<=m;l++)
        for(int j=1;j<=k;j++)
            for(int i=1;i<=n;i++)
            {
                c.a[i][j]+=a.a[i][l]*b.a[l][j]%P;
                if(c.a[i][j]>=P) c.a[i][j]-=P;
            }
    return c;
}
const Square S1=(Square){1,3,{{},{0,0,1,1}}};
const Square S3=(Square){3,3,{{},{0,1},{0,0,1},{0,0,0,1}}};
inline Square SQ(int type,int conf)
{
    switch(type)
    {
        case 1: return (Square){3,3,{{},{0,1},{0,conf,1},{0,0,0,1}}};
        case 3: return (Square){3,3,{{},{0,1},{0,0,conf},{0,0,0,1}}};
        case 4: return (Square){3,3,{{},{0,1},{},{0,0,1,1}}};
    }
}
Square aa[530012],all[530012];
int ll[530012],rr[530012];
void init()
{
    for(int i=1;i<=2*PTA-1;i++)
        aa[i]=S3;
    for(int i=1,j=PTA;i<=PTA;i++,j++)
    {
        all[j]=S1;
        ll[j]=rr[j]=i;
    }
    for(int i=PTA-1;i>=1;i--)
    {
        all[i]=all[i+i]+all[i+i+1];
        ll[i]=ll[i+i];
        rr[i]=rr[i+i+1];
    }
}
void add(int type,int conf,int l,int r,int b)
{
    if(ll[b]==l&&rr[b]==r)
    {
        Square c=SQ(type,conf);
        all[b]=all[b]*c;
        aa[b]=aa[b]*c;
        return;
    }
    Square d=aa[b];
    aa[b]=S3;
    aa[b+b]=aa[b+b]*d;
    all[b+b]=all[b+b]*d;
    aa[b+b+1]=aa[b+b+1]*d;
    all[b+b+1]=all[b+b+1]*d;
    if(l<=rr[b+b]) add(type,conf,l,min(r,rr[b+b]),b+b);
    if(r>=ll[b+b+1]) add(type,conf,max(l,ll[b+b+1]),r,b+b+1);
    all[b]=all[b+b]+all[b+b+1];
}
int get(int l,int r,int b)
{
    if(ll[b]==l&&rr[b]==r) return all[b].a[1][1];
    Square d=aa[b];
    aa[b]=S3;
    aa[b+b]=aa[b+b]*d;
    all[b+b]=all[b+b]*d;
    aa[b+b+1]=aa[b+b+1]*d;
    all[b+b+1]=all[b+b+1]*d;
    int ans=0;
    if(l<=rr[b+b]) ans+=get(l,min(r,rr[b+b]),b+b);
    if(r>=ll[b+b+1]) ans+=get(max(l,ll[b+b+1]),r,b+b+1);
    all[b]=all[b+b]+all[b+b+1];
    if(ans>=P) ans-=P;
    return ans;
}
signed main()
{
    init();
    int n,m;
    scanf("%lld %lld",&n,&m);
    while(m--)
    {
        int z,l,r,x;
        scanf("%lld %lld %lld ",&z,&l,&r);
        switch(z)
        {
            case 1:
            {
                scanf("%lld",&x);
                add(1,x,l,r,1);
                break;
            }
            case 2:
            {
                printf("%lld\n",get(l,r,1)%P);
                break;
            }
            case 3:
            {
                scanf("%lld",&x);
                add(3,x,l,r,1);
                break;
            }
            case 4:
            {
                add(4,0,l,r,1);
                break;
            }
        }
    }
    return 0;
}