Problem:
 Score : 600 points    Problem Statement Given are a sequence of N integers A_1 , A_2 , \ldots , A_N and a positive integer S . For a pair of integers (L, R) such that 1\leq L \leq R \leq N , let us define f(L, R) as follows:   f(L, R) is the number of sequences of integers (x_1, x_2, \ldots , x_k) such that L \leq x_1 < x_2 < \cdots < x_k \leq R and A_{x_1}+A_{x_2}+\cdots +A_{x_k} = S .   Find the sum of f(L, R) over all pairs of integers (L, R) such that 1\leq L \leq R\leq N . Since this sum can be enormous, print it modulo 998244353 .      Constraints  All values in input are integers.  1 \leq N \leq 3000  1 \leq S \leq 3000  1 \leq A_i \leq 3000         Input Input is given from Standard Input in the following format:  N  S  A_1  A_2  ...  A_N       Output Print the sum of f(L, R) , modulo 998244353 .        Sample Input 1 3 4
2 2 4      Sample Output 1 5  The value of f(L, R) for each pair is as follows, for a total of 5 .   f(1,1) = 0  f(1,2) = 1 (for the sequence (1, 2) )  f(1,3) = 2 (for (1, 2) and (3) )  f(2,2) = 0  f(2,3) = 1 (for (3) )  f(3,3) = 1 (for (3) )        Sample Input 2 5 8
9 9 9 9 9      Sample Output 2 0       Sample Input 3 10 10
3 1 4 1 5 9 2 6 5 3      Sample Output 3 152

Solution (Python):
def main():
  N,S = map(int,input().split())
  A = list(map(int,input().split()))
  MOD = 998244353

  
  
  
  
  
  
  
  
  
  ans = 0
  q = [0 for _ in range(S+1)]
  for i in range(len(A)):
    q[0] += 1
    q[0] %= MOD
    q2 = q[:]
    for j in range(S):
      if j+A[i] <= S:
        q2[j+A[i]] += q[j]
        q2[j+A[i]] %= MOD
    q = q2[:]
    ans += q[S]
  
  print(ans%MOD)
  
main()

Solution (Python):
def main():
  N,S = map(int,input().split())
  A = list(map(int,input().split()))
  MOD = 998244353

  ans = 0
  q = [0 for _ in range(S+1)]
  for i in range(len(A)):
    q[0] += 1
    q[0] %= MOD
    q2 = q[:]
    for j in range(S):
      if j+A[i] <= S:
        q2[j+A[i]] += q[j]
        q2[j+A[i]] %= MOD
    q = q2[:]
    ans += q[S]
  
  print(ans%MOD)
  
main()

Solution (Python):
mod = 998244353


N, S, *A = map(int, open(0).read().split())

dp = [0] * (S + 1)
ans = 0
for i, a in enumerate(A, 1):
    if a <= S:
        T = dp[:]
        for j in range(a, S + 1):
            dp[j] += T[j - a]
        dp[a] += i

    ans += dp[-1]
    ans %= mod

print(ans)

Solution (Python):

import sys
MOD = 998244353  

def solve(N: int, S: int, A: "List[int]"):
    dp = [[0]*(S+1) for _ in range(N+1)]

    for i in range(N):
        for s in range(S+1):
            if s >= A[i]:
                dp[i+1][s] = dp[i][s]+dp[i][s-A[i]]+int(A[i]==s or s==0)
            else:
                dp[i+1][s] = dp[i][s]+int(A[i]==s or s==0)
            dp[i+1][s]%=MOD
    
    answer = 0
    for i in range(N):
        answer += dp[i+1][S]
    print(answer%MOD)
    return


def main():
    def iterate_tokens():
        for line in sys.stdin:
            for word in line.split():
                yield word
    tokens = iterate_tokens()
    N = int(next(tokens))  
    S = int(next(tokens))  
    A = [int(next(tokens)) for _ in range(N)]  
    solve(N, S, A)

if __name__ == '__main__':
    main()

Solution (Python):
import sys
import numpy as np

sr = lambda: sys.stdin.readline().rstrip()
ir = lambda: int(sr())
lr = lambda: list(map(int, sr().split()))

MOD = 998244353
N, S = lr()
A = lr()
dp = np.zeros(S+1, np.int64)
answer = 0
for a in A:
    dp[0] += 1 
    prev = dp.copy()
    dp[a:] += prev[:-a]
    dp %= MOD
    answer += dp[-1] 

print(answer % MOD)

Solution (Python):
import sys
import numpy as np

sr = lambda: sys.stdin.readline().rstrip()
ir = lambda: int(sr())
lr = lambda: list(map(int, sr().split()))

MOD = 998244353
N, S = lr()
A = lr()
dp = np.array([0] * (S+1), dtype=np.int64)
answer = 0
for a in A:
    dp[0] += 1 
    prev = dp.copy()
    dp[a:] += prev[:-a]
    answer += dp[-1] 
    dp %= MOD

print(answer % MOD)

Solution (Python):
import sys

sr = lambda: sys.stdin.readline().rstrip()
ir = lambda: int(sr())
lr = lambda: list(map(int, sr().split()))

MOD = 998244353
N, S = lr()
A = lr()
dp = [0] * (S+1)
answer = 0
for a in A:
    dp[0] += 1 
    prev = dp[:]
    for i in range(S-a+1):
        dp[i+a] += prev[i]
    answer += dp[-1] 

print(answer % MOD)

Solution (Python):
import sys

sr = lambda: sys.stdin.readline().rstrip()
ir = lambda: int(sr())
lr = lambda: list(map(int, sr().split()))

MOD = 998244353
N, S = lr()
A = lr()
dp = [0] * (S+1)
answer = 0
for a in A:
    dp[0] += 1 
    prev = dp.copy()
    for i in range(S-a+1):
        dp[i+a] += prev[i]
    answer += dp[-1] 

print(answer % MOD)

Solution (Python):
n,s = map(int,input().split())
a = list(map(int,input().split()))
c = [0]*(s+1)
ans = 0
mod = 998244353
for i in range(n):
    c_=c.copy()
    if a[i] <= s:
        for j in range(a[i],s+1):
            c[j] += c_[j-a[i]]
        c[a[i]] += i+1
    ans += c[s]
    ans %= mod

print(ans)

Solution (Python):
N,S =map(int,input().split())
A=list(map(int,input().split()))
mod=998244353
ans=0
q=[0]*(S+1)
for i in range(N):
    q2=q.copy()
    if A[i]<=S:
        for j in range(A[i],S+1):
            q[j]+=q2[j-A[i]]
        q[A[i]]+=i+1
    ans+=q[S]
    ans%=mod
print(ans)

Solution (Python):
def main():
	n, s, *a = map(int, open(0).read().split())
	ans = 0
	mod = 998244353
	dp = [0] * (s + 1)
	for i, x in enumerate(a):
		dp[0] += 1
		if x > s:continue
		ans += dp[s - x] * (n - i)
		ans %= mod
		dp = dp[:x] + [(dp[j] + dp[j - x]) % mod for j in range(x, s + 1)]
	print(ans)
	
if __name__=="__main__":
	main()

Solution (Python):
n,s,*a=map(int,open(0).read().split())
ans=0
mod=998244353
dp=[0]*(s+1)
for i in range(n):
	dp[0]+=1
	if a[i]>s:
		continue
	ans+=dp[s-a[i]]*(n-i)
	ans%=mod
	for j in range(s,a[i]-1,-1):
		dp[j]+=dp[j-a[i]]
		dp[j]%=mod
print(ans)

Solution (Python):
def main():
    n,s = map(int,input().split())
    A = list(map(int,input().split()))
    dp = [[0 for j in range(s+1)] for i in range(n+1)]

    mod = 998244353

    ans = 0

    dp[0][0] = 1

    for i in range(n):
        dp[i+1][0] = dp[i][0]+1
        for j in range(1,s+1):
            if j-A[i] >= 0:
                dp[i+1][j] += dp[i][j-A[i]]
            dp[i+1][j] = (dp[i][j]+dp[i+1][j])%mod

        ans = (ans+dp[i+1][-1])%mod
    
    print(ans%mod)
if __name__ == "__main__":
    main()

Solution (Python):
N, S = map(int, input().split())
A = list(map(int, input().split()))

mod = 998244353

"""
多項式係数を使った数え上げでリベンジしてみる。
maspyさんのを参考に実装する。

（注意）多項式係数全てを書き下すと、d = sum(A)次式になるが、
これは Ai*N = 3000*3000 次になりTLEする。。。
S次の項までで打ち切る事を忘れない事。

やってみてわかったけど、実装自体はdpでやったのとほぼ同じ。
多項式係数を使う方法は実質flexibleなdpということ？
"""

f = [0] * (1 + S) 
ans = 0

for i in range(N):
    g = [0] * (1 + S) 
    f[0] += 1
    for j in range(1 + S):
        if j - A[i] < 0:
            g[j] = f[j]
        else:
            g[j] = f[j] + f[j - A[i]]
    ans += g[S]
    ans %= mod
    f = g

print(ans)

Solution (Python):
import numpy as np
n,s = map(int,input().split())
A = list(map(int,input().split()))
mod = 998244353
ans = 0
dp = np.zeros(s+1,np.int64)
for i in A:
    dp[0] += 1
    dp[i:] += dp[:-i].copy()
    dp %= mod
    ans += dp[-1]
print(ans%mod)

Solution (Python):
N, S = map(int, input().split())
A = [0] + list(map(int, input().split()))
mod = 998244353

dp = [0] * (S + 1)
dp2 = [0] * (S + 1)

"""
dp[i][j] = (i個までの数を使ってj点とる組み合わせ)
としたはいいが、f(L, R)をf(1,L)とf(1, R)に分解できずに爆死したのが敗因。

各(L, R)に対してdpするとTLEになるから手詰まりになると思ったわけだけど、
上のdpテーブルを各iに対して使い回しつつ、新たなdpテーブル[1, 0, 0, ...]を
作成して「元のdpテーブルに加算しておけば」dpテーブル一個で1番目スタート、2番目スタート、
..., L番目スタートのdpテーブル全てを考えられる！
（解説放送前半を参考）
"""


dp[0] = 1
ans = 0


for i in range(1, N + 1):
    for j in range(S + 1):
        if j - A[i] < 0:
            dp2[j] = dp[j]
        else:
            dp2[j] = dp[j] + dp[j - A[i]]
    ans += dp2[S]
    ans %= mod

    dp = dp2 
    dp[0] += 1
    dp2 = [0] * (S + 1)

print(ans)

Solution (Python):
import sys
input = sys.stdin.readline

mod = 998244353

N, S = map(int, input().split())
A = list(map(int, input().split()))

dp = [[0]*(S+1) for _ in range(N+1)]
dp[0][0] = 1
for i, a in enumerate(A):
    dp[i+1][0] = i+2
    for s in range(1, S+1):
        if s-a >= 0:
            dp[i+1][s] = (dp[i][s] + dp[i][s-a]) % mod
        else:
            dp[i+1][s] = dp[i][s]

ans = 0
for r in range(N+1):
    ans = (ans + dp[r][S]) % mod

print(ans)

Solution (Python):
def main():
  import sys
  read = sys.stdin.buffer.read
  r = 998244353
  U, S, *A = map(int, read().split())
  import copy
  ans = 0
  f = [0]*3010
  for a in A:
      c = copy.copy(f)
      f = [0]*a+c[:-a]
      f = [x+y for x,y in zip(f,c)]
      f[0] += 1
      f[a] += 1
      ans += f[S]
  ans %= r
  print(ans)
if __name__ == '__main__':
    main()

Solution (Python):
import sys
input = sys.stdin.readline

N,S=map(int,input().split())
A=list(map(int,input().split()))
mod=998244353

DP=[0]*(S+1)






DP[0]=1
ANS=0
for a in A:
    for i in range(S-a,-1,-1):
        DP[i+a]+=DP[i]
    DP[0]+=1
    ANS+=DP[S]
    
print(ANS%mod)

Solution (Python):
def examA():
    N, M = LI()
    ans = N*(N-1)//2 + M*(M-1)//2
    print(ans)
    return

def examB():
    def judge(A):
        
        n = len(A)
        for i in range(n//2):
            if A[i]!=A[n-i-1]:
                return False
        return True
    S = SI()
    N = len(S)
    if not judge(S):
        print("No")
        return
    if judge(S[:(N-1)//2]) and judge(S[(N+1)//2:]):
        print("Yes")
    else:
        print("No")
    return

def examC():
    L = I()
    ans = (L/3)**3
    print(ans)
    return

def examD():
    N = I()
    A = LI()
    D = defaultdict(bool)
    C = Counter(A)
    M = 0
    for key,c in C.items():
        M += c*(c-1)//2
    
    
    ans = [0]*N
    for i in range(N):
        a = A[i]
        c = C[a]
        ans[i] = M -(c*(c-1)//2) + ((c-1)*(c-2)//2)
    for v in ans:
        print(v)
    return

def examE():
    H, W, K = LI()
    S = [SI()for _ in range(H)]
    SH = [[0] * (W + 1) for _ in range(H + 1)]
    for i in range(H):
        now = 0
        for j in range(W):
            now += int(S[i][j])
            SH[i + 1][j + 1] = SH[i][j + 1] + now

    def calc(a,b,c,d):
        return SH[c][d] - SH[a][d] - SH[c][b] + SH[a][b]

    ans = H+W-2
    mask = 2**(H-1)-1
    for i in range(mask+1):
        cut = [0]
        cur = 0
        for j in range(H-1):
            if i&(1<<j)==(1<<j):
                cut.append(j+1)
                cur += 1
        cut.append(H)
        L = [0]*(len(cut)-1)
        flag = True
        for w in range(W):
            now_L = [0]*(len(cut)-1)
            for s in range(len(cut) - 1):
                n = calc(cut[s], w, cut[s + 1], w + 1)
                now_L[s] = n
                if n+L[s]> K:
                    if n>K:
                        flag = False
                        cur = inf
                        break
                    cur += 1
                    L = [0] * (len(cut))
            if not flag:
                break
            for s in range(len(cut) - 1):
                L[s] += now_L[s]
        if not flag:
            continue
        ans = min(ans,cur)

    print(ans)
    return
"""
2 1 1
1
1
"""

def examF():
    N, S = LI()
    A = LI()
    dp = [[0]*(S+1) for _ in range(N+1)]
    for i in range(N):
        dp[i][0] = 1
    for i in range(N):
        a = A[i]
        for j in range(S+1):
            dp[i+1][j] += dp[i][j]
            dp[i+1][j] %= mod2
            if j+a<=S:
                dp[i+1][j+a] += dp[i][j]
                dp[i+1][j+a] %= mod2
    ans = 0
    for i in range(N+1):
        ans += dp[i][S]
        ans %= mod2
    print(ans)
    return

import sys,bisect,itertools,heapq,math,random
from copy import deepcopy
from heapq import heappop,heappush,heapify
from collections import Counter,defaultdict,deque

def I(): return int(sys.stdin.readline())
def LI(): return list(map(int,sys.stdin.readline().split()))
def LSI(): return list(map(str,sys.stdin.readline().split()))
def LS(): return sys.stdin.readline().split()
def SI(): return sys.stdin.readline().strip()
global mod,mod2,inf,alphabet,_ep
mod = 10**9 + 7
mod2 = 998244353
inf = 10**18
_ep = 10**(-12)
alphabet = [chr(ord('a') + i) for i in range(26)]

sys.setrecursionlimit(10**6)

if __name__ == '__main__':
    examF()

"""

"""

Solution (Python):
import numpy as np
from sys import stdin

'''
Reference: https://codeforces.com/blog/entry/75046

(0-based index)

For a sequence A[l]..A[k]..A[r] with sum S contributes l * (N - r) to answer.
    0, 1, 2, ..., l....k....r, ..., N - 2, N - 1
    x  x  x   x   x         v   v     v      v
There are (l + 1) `x` and (N - r) `v`.

Enumerate all possible sequence is impossible. 
Consider all sequeuces with sum S and **ends** at specific r:
    A[l1]...A[r] -> l1 * (N - r)
    A[l2]...A[r] -> l2 * (N - r)
    A[l3]...A[r] -> l3 * (N - r)
    .
    .
    .
which contributes (l1 + l2 + ... + lm) * (N - r) in total.

We can find *sum of left index for all valid sequence* for each r using DP.

    dp[i][s] = the sum of left index for all sequence that
                    1. ends at i
                    2. has sum s
    
    dp[i][:A[i]] = 0  since any sequence ends at i has sum >= A[i]

    dp[i][A[i]] = i + 1
    dp[i][s] = sum(dp[j][s - A[i]] for j in range(i))

    Take A = [2, 3, 4] for example, sequences that end at A[2] are:
    
            Base case:
                   [4] -> dp[2][4] = 3

            Transition from dp[0][s - 4]:
             [2,    4]

            Transition from dp[1][s - 4]:
                [3, 4]
             [2, 3, 4]
    
    The answer is sum(dp[r][S] * (N - r) for r in range(N))

So the naive implementation is:

    dp = np.zeros((N, max(S, A.max()) + 1), dtype=np.int64)
    for i, a in enumerate(A):
        dp[i, a] = i + 1
        for s in range(a + 1, S + 1):
            dp[i, s] = dp[:i, s - a].sum() % M  
    
    ans = 0
    for i in range(N):
        ans += (N - i) * dp[i, S]
        ans = ans % M

has time O(N * S * N) which results in TLE.

Using another array to store the cumsum (prefix sum) of dp solves the problem.
One can compress the dp matrix to vector to reduce memory usage.
'''

N, S = map(int, input().split())
A = np.int32(input().split())
M = 998244353
V = max(A.max(), S)

ans = 0
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp = np.zeros(V + 1, dtype=np.int64)
    dp[a] = i + 1
    dp[a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp) % M
    ans = (ans + (N - i) * dp[S]) % M
print(ans)

Solution (Python):
import numpy as np
from sys import stdin

'''
Reference: https://codeforces.com/blog/entry/75046

(0-based index)

For a sequence A[l]..A[k]..A[r] with sum S contributes l * (N - r) to answer.
    0, 1, 2, ..., l....k....r, ..., N - 2, N - 1
    x  x  x   x   x         v   v     v      v
There are (l + 1) `x` and (N - r) `v`.

Enumerate all possible sequence is impossible. 
Consider all sequeuces with sum S and **ends** at specific r:
    A[l1]...A[r] -> l1 * (N - r)
    A[l2]...A[r] -> l2 * (N - r)
    A[l3]...A[r] -> l3 * (N - r)
    .
    .
    .
which contributes (l1 + l2 + ... + lm) * (N - r) in total.

We can find *sum of left index for all valid sequence* for each r using DP.

    dp[i][s] = the sum of left index for all sequence that
                    1. ends at i
                    2. has sum s
    
    dp[i][:A[i]] = 0  since any sequence ends at i has sum >= A[i]

    dp[i][A[i]] = i + 1
    dp[i][s] = sum(dp[j][s - A[i]] for j in range(i))

    Take A = [2, 3, 4] for example, sequences that end at A[2] are:
    
                   [4] -> dp[2][4] = 3
             [2,    4] ↘
                [3, 4] -> dp[2][s] = dp[2][s - 4] for all s
             [2, 3, 4] ↗
    
    The answer is sum(dp[r][S] * (N - r) for r in range(N))

So the naive implementation is:

    dp = np.zeros((N, max(S, A.max()) + 1), dtype=np.int64)
    for i, a in enumerate(A):
        dp[i, a] = i + 1
        for s in range(a + 1, S + 1):
            dp[i, s] = dp[:i, s - a].sum() % M  
    
    ans = 0
    for i in range(N):
        ans += (N - i) * dp[i, S]
        ans = ans % M

has time O(N * S * N) which results in TLE.

Using another array to store the cumsum (prefix sum) of dp solves the problem.
One can compress the dp matrix to vector to reduce memory usage.
'''

N, S = map(int, input().split())
A = np.int32(input().split())
M = 998244353
V = max(A.max(), S)

ans = 0
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp = np.zeros(V + 1, dtype=np.int64)
    dp[a] = i + 1
    dp[a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp) % M
    ans = (ans + (N - i) * dp[S]) % M
print(ans)

Solution (Python):
import numpy as np
from sys import stdin

'''
Reference: https://codeforces.com/blog/entry/75046

(0-based index)

For a sequence A[l]..A[k]..A[r] with sum S contributes l * (N - r) to answer.
    0, 1, 2, ..., l....k....r, ..., N - 2, N - 1
    x  x  x   x   x         v   v     v      v
There are l `x` and (N - r) `v`.

Enumerate all possible sequence is impossible. 
Consider all sequeuces with sum S and **ends** at specific r:
    A[l1]...A[r] -> l1 * (N - r)
    A[l2]...A[r] -> l2 * (N - r)
    A[l3]...A[r] -> l3 * (N - r)
    .
    .
    .
which contributes (l1 + l2 + ... + lm) * (N - r) in total.

We can find *sum of left index for all valid sequence* for each r using DP.

    dp[i][s] = the sum of left index for all sequence that
                    1. ends at i
                    2. has sum s
    
    dp[i][:A[i]] = 0  since any sequence ends at i has sum >= A[i]

    dp[i][A[i]] = i + 1
    dp[i][s] = sum(dp[j][s - A[i]] for j in range(i))

    Take A = [2, 3, 4] for example, sequences that end at A[2] are:
    
                   [4] -> dp[2][4] = 3
             [2,    4] ↘
                [3, 4] -> dp[2][s] = dp[2][s - 4] for all s
             [2, 3, 4] ↗
    
    The answer is sum(dp[r][S] * (N - r) for r in range(N))

So the naive implementation is:

    dp = np.zeros((N, max(S, A.max()) + 1), dtype=np.int64)
    for i, a in enumerate(A):
        dp[i, a] = i + 1
        for s in range(a + 1, S + 1):
            dp[i, s] = dp[:i, s - a].sum() % M  
    
    ans = 0
    for i in range(N):
        ans += (N - i) * dp[i, S]
        ans = ans % M

has time O(N * S * N) which results in TLE.

Using another array to store the cumsum (prefix sum) of dp solves the problem:
'''

N, S = map(int, input().split())
A = np.int32(input().split())
M = 998244353
V = max(A.max(), S)

ans = 0
dp = np.zeros((N, V + 1), dtype=np.int64)
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp[i, a] = i + 1
    dp[i, a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp[i]) % M
    ans = (ans + (N - i) * dp[i, S]) % M
print(ans)

Solution (Python):
import numpy as np
from sys import stdin

'''
Reference: https://codeforces.com/blog/entry/75046

(0-based index)

For a sequence A[l]..A[k]..A[r] with sum S contributes l * (N - r) to answer.
    0, 1, 2, ..., l....k....r, ..., N - 2, N - 1
    x  x  x   x   x         v   v     v      v
There are l `x` and (N - r) `v`.

Enumerate all possible sequence is impossible. 
Consider all sequeuces with sum S and **ends** at specific r:
    A[l1]...A[r] -> l1 * (N - r)
    A[l2]...A[r] -> l2 * (N - r)
    A[l3]...A[r] -> l3 * (N - r)
    .
    .
    .
which contributes (l1 + l2 + ... + lm) * (N - r) in total.

We can find *sum of left index for all valid sequence* for each r using DP.

    dp[i][s] = the sum of left index for all sequence that
                    1. ends at i
                    2. has sum s
    
    dp[i][:A[i]] = 0  since any sequence ends at i has sum >= A[i]

    dp[i][A[i]] = i + 1
    dp[i][s] = sum(dp[j][s - A[i]] for j in range(i))

    Take A = [2, 3, 4] for example, sequences that end at A[2] are:
    
                   [4] -> dp[2][4] = 3
             [2,    4] ↘
                [3, 4] -> dp[2][s] = dp[2][s - 4] for all s
             [2, 3, 4] ↗
    
    The answer is dp[r][S] * (N - r) for r in range(N)

So the naive implementation is:

    dp = np.zeros((N, max(S, A.max()) + 1), dtype=np.int64)
    for i, a in enumerate(A):
        dp[i, a] = i + 1
        for s in range(a + 1, S + 1):
            dp[i, s] = dp[:i, s - a].sum() % M  
    
    ans = 0
    for i in range(N):
        ans += (N - i) * dp[i, S]
        ans = ans % M

has time O(N * S * N) which results in TLE.

Using another array to store the cumsum (prefix sum) of dp solves the problem:
'''

N, S = map(int, input().split())
A = np.int32(input().split())
M = 998244353
V = max(A.max(), S)

ans = 0
dp = np.zeros((N, V + 1), dtype=np.int64)
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp[i, a] = i + 1
    dp[i, a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp[i]) % M
    ans = (ans + (N - i) * dp[i, S]) % M
print(ans)

Solution (Python):
import numpy as np
from sys import stdin

'''
Reference: https://codeforces.com/blog/entry/75046

(0-based index)

For a sequence A[l]..A[k]..A[r] with sum S contributes l * (N - r) to answer.
    0, 1, 2, ..., l....k....r, ..., N - 2, N - 1
    x  x  x   x   x         v   v     v      v
There are l `x` and (N - r) `v`.

Enumerate all possible sequence is impossible. 
Let's consider all sequeuces with sum S and **ends** at r:
    A[l1]...A[r] -> l1 * (N - r)
    A[l2]...A[r] -> l2 * (N - r)
    A[l3]...A[r] -> l3 * (N - r)
    .
    .
    .
which contributes (l1 + l2 + ... + lm) * (N - r) in total.

We can find *sum of left index for all valid sequence* for each r using DP.

    dp[i][s] = the sum of left index for all sequence that
                    1. ends at i
                    2. has sum s
    
    dp[i][:A[i]] = 0  since any sequence ends at i has sum >= A[i]

    dp[i][A[i]] = i + 1
    dp[i][s] = sum(dp[j][s - A[i]] for j in range(i))

    Take A = [2, 3, 4] for example, sequences that end at A[2] are:
    
             [      4] -> dp[2][4] = 3
             [2,    4] ↘
             [   3, 4] -> dp[2][s] = dp[2][s - 4] for all s
             [2, 3, 4] ↗
    
    The answer is dp[r][S] * (N - r) for r in range(N)

So the naive implementation is:

    dp = np.zeros((N, max(S, A.max()) + 1), dtype=np.int64)
    for i, a in enumerate(A):
        dp[i, a] = i + 1
        for s in range(a + 1, S + 1):
            dp[i, s] = dp[:i, s - a].sum() % M  
    
    ans = 0
    for i in range(N):
        ans += (N - i) * dp[i, S]
        ans = ans % M

has time O(N * S * N) which results in TLE.

Using another array to store the cumsum (prefix sum) of dp solves the problem:
'''

N, S = map(int, input().split())
A = np.int32(input().split())
M = 998244353
V = max(A.max(), S)

ans = 0
dp = np.zeros((N, V + 1), dtype=np.int64)
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp[i, a] = i + 1
    dp[i, a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp[i]) % M
    ans = (ans + (N - i) * dp[i, S]) % M
print(ans)

Solution (Python):
from heapq import heappush,heappop,heapify
from collections import deque,defaultdict,Counter
from itertools import permutations,combinations,groupby
import sys,bisect,string,math,time,functools,random
def Golf():*a,=map(int,open(0))
def I():return int(input())
def S_():return input()
def IS():return input().split()
def LS():return [i for i in input().split()]
def LI():return [int(i) for i in input().split()]
def LI_():return [int(i)-1 for i in input().split()]
def NI(n):return [int(input()) for i in range(n)]
def NI_(n):return [int(input())-1 for i in range(n)]
def StoLI():return [ord(i)-97 for i in input()]
def ItoS(n):return chr(n+97)
def LtoS(ls):return ''.join([chr(i+97) for i in ls])
def GI(V,E,Directed=False,index=0):
    org_inp=[];g=[[] for i in range(n)]
    for i in range(E):
        inp=LI();org_inp.append(inp)
        if index==0:inp[0]-=1;inp[1]-=1
        if len(inp)==2:
            a,b=inp;g[a].append(b)
            if not Directed:g[b].append(a)
        elif len(inp)==3:
            a,b,c=inp;aa=(inp[0],inp[2]);bb=(inp[1],inp[2]);g[a].append(bb)
            if not Directed:g[b].append(aa)
    return g,org_inp
def GGI(h,w,boundary=1,search=[],replacement_of_found='.',mp_def={'

    mp=[boundary]*(w+2);found={}
    for i in range(h):
        s=input()
        for char in search:
            if char in s:
                found[char]=((i+1)*(w+2)+s.index(char)+1)
                mp_def[char]=mp_def[replacement_of_found]
        mp+=[boundary]+[mp_def[j] for j in s]+[boundary]
    mp+=[boundary]*(w+2)
    return h+2,w+2,mp,found
def TI(n):return GI(n,n-1)
def bit_combination(k,n=2):return [[tb//(n**bt)%n for bt in range(k)] for tb in range(n**k)]
def show(*inp,end='\n'):
    if show_flg:print(*inp,end=end)
def show2d(g,h,w):
    for i in range(h):show(g[i*w:i*w+w])

YN=['YES','NO'];Yn=['Yes','No']
mo=10**9+7
inf=float('inf')
l_alp=string.ascii_lowercase

input=lambda: sys.stdin.readline().rstrip()

def ran_input():
    n=random.randint(4,16)
    rmin,rmax=1,10
    a=[random.randint(rmin,rmax) for _ in range(n)]
    return n,a

show_flg=False
show_flg=True


mo=998244353

ans=0
n,s=LI()
a=LI()

def dyn(n,A,a):
    dp=[[0]*(A+1) for i in range(n+1)]
    for i in range(n):
        x=a[i]
        dp[i+1]=dp[i][:]
        dp[i+1][x]+=1+i
        for j in range(x,A+1):
            dp[i+1][j]+=dp[i][j-x]
            dp[i+1][j]%=mo
    return dp

dp=dyn(n,3000,a)

for i in dp:
    ans+=i[s]
    ans%=mo

print(ans)

Solution (Python):
from heapq import heappush,heappop,heapify
from collections import deque,defaultdict,Counter
from itertools import permutations,combinations,groupby
import sys,bisect,string,math,time,functools,random
def Golf():*a,=map(int,open(0))
def I():return int(input())
def S_():return input()
def IS():return input().split()
def LS():return [i for i in input().split()]
def LI():return [int(i) for i in input().split()]
def LI_():return [int(i)-1 for i in input().split()]
def NI(n):return [int(input()) for i in range(n)]
def NI_(n):return [int(input())-1 for i in range(n)]
def StoLI():return [ord(i)-97 for i in input()]
def ItoS(n):return chr(n+97)
def LtoS(ls):return ''.join([chr(i+97) for i in ls])
def GI(V,E,Directed=False,index=0):
    org_inp=[];g=[[] for i in range(n)]
    for i in range(E):
        inp=LI();org_inp.append(inp)
        if index==0:inp[0]-=1;inp[1]-=1
        if len(inp)==2:
            a,b=inp;g[a].append(b)
            if not Directed:g[b].append(a)
        elif len(inp)==3:
            a,b,c=inp;aa=(inp[0],inp[2]);bb=(inp[1],inp[2]);g[a].append(bb)
            if not Directed:g[b].append(aa)
    return g,org_inp
def GGI(h,w,boundary=1,search=[],replacement_of_found='.',mp_def={'

    mp=[boundary]*(w+2);found={}
    for i in range(h):
        s=input()
        for char in search:
            if char in s:
                found[char]=((i+1)*(w+2)+s.index(char)+1)
                mp_def[char]=mp_def[replacement_of_found]
        mp+=[boundary]+[mp_def[j] for j in s]+[boundary]
    mp+=[boundary]*(w+2)
    return h+2,w+2,mp,found
def TI(n):return GI(n,n-1)
def bit_combination(k,n=2):return [[tb//(n**bt)%n for bt in range(k)] for tb in range(n**k)]
def show(*inp,end='\n'):
    if show_flg:print(*inp,end=end)
def show2d(g,h,w):
    for i in range(h):show(g[i*w:i*w+w])

YN=['YES','NO'];Yn=['Yes','No']
mo=10**9+7
inf=float('inf')
l_alp=string.ascii_lowercase

input=lambda: sys.stdin.readline().rstrip()

def ran_input():
    n=random.randint(4,16)
    rmin,rmax=1,10
    a=[random.randint(rmin,rmax) for _ in range(n)]
    return n,a

show_flg=False
show_flg=True


mo=998244353

ans=0

n,s=LI()
a=LI()

A=3000



def dyn(n,A,a):
    dp=[[0]*(A+1) for i in range(n+1)]
    
    for i in range(n):
        x=a[i]
        dp[i+1]=dp[i][:]
        dp[i+1][x]+=1+i
        for j in range(x,A+1):
            
            dp[i+1][j]+=dp[i][j-x]
            dp[i+1][j]%=mo
    return dp

dp=dyn(n,A,a)

ans=0
for i in dp:
    ans+=i[s]
    ans%=mo
    

print(ans)

Solution (Python):
import numpy as np

N, S = map(int, input().split())
A = np.int32(input().split())
M = 998244353
V = max(max(A), S)

ans = 0
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp = np.zeros(V + 1, dtype=np.int64)
    dp[a] = i + 1
    dp[a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp) % M
    ans = (ans + (N - i) * dp[S]) % M
print(ans)

Solution (Python):
import numpy as np
from sys import stdin

N, S = map(int, stdin.readline().split())
A = list(map(int, stdin.readline().split()))
M = 998244353
V = max(max(A), S)

ans = 0
cs = np.zeros(V + 1, dtype=np.int64)
for i, a in enumerate(A):
    dp = np.zeros(V + 1, dtype=np.int64)
    dp[a] = i + 1
    dp[a + 1 : S + 1] = cs[np.arange(a + 1, S + 1) - a]
    cs = (cs + dp) % M
    ans = (ans + (N - i) * dp[S]) % M
print(ans)

Solution (Python):
MOD = 998244353

N, S = map(int, input().split())
A = list(map(int, input().split()))

F = [0 for _ in range(S + 1)]
res = 0

for i in range(N):
    F[0] += 1
    for j in range(S + 1)[::-1]:
        if j - A[i] >= 0:
            F[j] += F[j - A[i]]
    res += F[S]
    res %= MOD

print(res)

Solution (Python):
from sys import stdin

N, S = map(int, stdin.readline().split())
A = list(map(int,stdin.readline().split()))
M = 998244353
V = max(max(A), S)

ans = 0
cs = [0 for _ in range(V + 1)]
for i, a in enumerate(A):
    dp = [0 for _ in range(V + 1)]
    dp[a] = i + 1
    for s in range(a + 1, S + 1):
        dp[s] = cs[s - a]
    cs = [(a + b) % M for a, b in zip(cs, dp)]
    ans = (ans + (N - i) * dp[S]) % M
print(ans)

Solution (Python):
N,S=map(int,input().split())
A=list(map(int,input().split()))
m=998244353
import numpy as np
dp=np.zeros((N+1,S+1),dtype=np.int64)
dp[:,0]=[1]*(N+1)
for i,a in enumerate(A):
    dp[i+1]+=dp[i]
    dp[i+1,a:]=(dp[i+1,a:]+dp[i,:-a])%m
print(np.sum(dp[:,-1])%m)

Solution (Python):
import os
import sys

if os.getenv("LOCAL"):
    sys.stdin = open("_in.txt", "r")

sys.setrecursionlimit(10 ** 9)
INF = float("inf")
IINF = 10 ** 18

MOD = 998244353


N, S = list(map(int, sys.stdin.buffer.readline().split()))
A = list(map(int, sys.stdin.buffer.readline().split()))


dp = [[0] * S for _ in range(N + 1)]
ans = 0
for i in range(N):
    for k in range(S):
        dp[i + 1][k] += dp[i][k]
        dp[i + 1][k] %= MOD
        if k == 0 and k + A[i] < S:
            
            
            dp[i + 1][A[i]] += i + 1
            dp[i + 1][A[i]] %= MOD
        elif k == 0 and k + A[i] == S:
            ans += (i + 1) * (N - i)
        elif k + A[i] < S:
            dp[i + 1][k + A[i]] += dp[i][k]
            dp[i + 1][k + A[i]] %= MOD
        elif k + A[i] == S:
            
            
            ans += dp[i][k] * (N - i)
        ans %= MOD
print(ans)