Problem:
 Score : 700 points    Problem Statement Takahashi has a tower which is divided into N layers.
Initially, all the layers are uncolored. Takahashi is going to paint some of the layers in red, green or blue to make a beautiful tower.
He defines the beauty of the tower as follows:   The beauty of the tower is the sum of the scores of the N layers, where the score of a layer is A if the layer is painted red, A+B if the layer is painted green, B if the layer is painted blue, and 0 if the layer is uncolored.   Here, A and B are positive integer constants given beforehand. Also note that a layer may not be painted in two or more colors.  Takahashi is planning to paint the tower so that the beauty of the tower becomes exactly K .
How many such ways are there to paint the tower? Find the count modulo 998244353 .
Two ways to paint the tower are considered different when there exists a layer that is painted in different colors, or a layer that is painted in some color in one of the ways and not in the other.      Constraints  1 ≤ N ≤ 3×10^5  1 ≤ A,B ≤ 3×10^5  0 ≤ K ≤ 18×10^{10}  All values in the input are integers.         Input Input is given from Standard Input in the following format:  N  A  B  K       Output Print the number of the ways to paint tiles, modulo 998244353 .        Sample Input 1 4 1 2 5      Sample Output 1 40  In this case, a red layer worth 1 points, a green layer worth 3 points and the blue layer worth 2 points. The beauty of the tower is 5 when we have one of the following sets of painted layers:   1 green, 1 blue  1 red, 2 blues  2 reds, 1 green  3 reds, 1 blue   The total number of the ways to produce them is 40 .       Sample Input 2 2 5 6 0      Sample Output 2 1  The beauty of the tower is 0 only when all the layers are uncolored. Thus, the answer is 1 .       Sample Input 3 90081 33447 90629 6391049189      Sample Output 3 577742975

Solution (Python):
n,a,b,k=map(int,input().split())
mod=998244353
ans=0
p=[1]*(n+1)
for i in range(n):
    p[i+1]=((i+1)*p[i])%mod
for ai in range(n+1):
    if k-a*ai>=0 and (k-a*ai)%b==0:
        bi=(k-a*ai)//b
        if bi<=n:
            ans+=p[n]*pow(p[ai],mod-2,mod)*pow(p[n-ai],mod-2,mod)*p[n]*pow(p[bi],mod-2,mod)*pow(p[n-bi],mod-2,mod)
            ans%=mod
print(ans)

Solution (Python):
N,A,B,K = map(int, input().split())

MOD = 998244353 
MAX = 10**6

fac = [0]*MAX
finv = [0]*MAX
inv = [0]*MAX

def init():
    fac[0] = fac[1] = 1
    finv[0] = finv[1] = 1
    inv[1] = 1
    for i in range(2,MAX):
        fac[i] = fac[i-1] * i % MOD
        inv[i] = MOD - inv[MOD%i] * (MOD // i) % MOD
        finv[i] = finv[i-1] * inv[i] % MOD

def C(n, k):
    if n<k:
        return 0
    if n<0 or k<0:
        return 0
    return fac[n]%MOD * (finv[k] * finv[n-k] % MOD) % MOD

res = 0
init()
for a in range(0,N+1):
    
    if (K - A*a)%B == 0 and (K -A*a)//B >= 0:
        b = (K -A*a)//B
        res += C(N,a)*C(N,b)
        res %= MOD
print(res)

Solution (Python):
def cmb1(n, r, mod):
    if ( r<0 or r>n ):
        return 0
    r = min(r, n-r)
    return g1[n] * g2[r] * g2[n-r] % mod

mod = 998244353 
N = 10**6
g1 = [1, 1] 
g2 = [1, 1] 
inverse = [0, 1] 

for i in range( 2, N + 1 ):
    g1.append( ( g1[-1] * i ) % mod )
    inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
    g2.append( (g2[-1] * inverse[-1]) % mod )


n, a, b, k = map(int, input().split())
mod = 998244353

ans = 0
for i in range(n+1):
    if k - a*i >= 0 and (k - a*i)%b == 0:
        j = (k - a*i)//b
        ans += cmb1(n, i, mod)*cmb1(n, j, mod)
        ans %= mod

print(ans)

Solution (Python):
def cmb1(n, r, mod):
    if ( r<0 or r>n ):
        return 0
    r = min(r, n-r)
    return g1[n] * g2[r] * g2[n-r] % mod

mod = 998244353 
N = 10**6
g1 = [1, 1] 
g2 = [1, 1] 
inverse = [0, 1] 

for i in range( 2, N + 1 ):
    g1.append( ( g1[-1] * i ) % mod )
    inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
    g2.append( (g2[-1] * inverse[-1]) % mod )


n, a, b, k = map(int, input().split())
mod = 998244353

ans = 0
for i in range(n+1):
    if k - a*i >= 0 and (k - a*i)%b == 0:
        j = (k - a*i)//b
        ans += cmb1(n, i, mod)*cmb1(n, j, mod)
        ans %= mod

print(ans)

Solution (Python):
mod = 998244353
rng = 300100
fctr = [1]
finv = [1]
for i in range(1,rng):
  fctr.append(fctr[-1]*i%mod)
for i in range(1,rng):
  finv.append(pow(fctr[i],mod-2,mod))
def cmb(n,k):
  if n<0 or k<0:
    return 0
  else:
    return fctr[n]*finv[n-k]*finv[k]%mod
n,a,b,k = map(int,input().split())
pair = []
for i in range(n+1):
  if (k-a*i)%b == 0 and (k-a*i)//b <= n:
    pair.append((i,(k-a*i)//b))
ans = 0
for x,y in pair:
  ans = (ans+cmb(n,x)*cmb(n,y))%mod
print(ans)

Solution (Python):
mod = 998244353
N, A, B, K = map(int, input().split())
def make_fact(n):
    fact = [1]*(n+1)
    for i in range(1, n+1):
        fact[i] = fact[i-1]*i%mod
    return fact
fact = make_fact(N) 
def make_fact_inv(n):
    fact_inv = [1]*(n+1)
    fact_inv[n] = pow(fact[n], mod-2, mod)
    for i in range(n, 0, -1):
        fact_inv[i-1] = fact_inv[i]*i%mod
    return fact_inv
fact_inv = make_fact_inv(N)
def comb(n, k):
    return fact[n]*fact_inv[k]*fact_inv[n-k]%mod
ans = 0
for k in range(N+1):
    r = K-k*A
    if r<0 or r%B != 0 or r//B>N:
        continue
    l = r//B
    ans+=comb(N, k)*comb(N, l)
    ans%=mod
print(ans)

Solution (Python):
mod = 998244353

n,a,b,k = map(int,input().split())

fact = [1]
fact_inv = [1]
for i in range(n+4):
    new_fact = fact[-1]*(i+1)%mod
    fact.append(new_fact)
    fact_inv.append(pow(new_fact,mod-2,mod))


def mod_comb_k(n,k,mod):
    return fact[n] * fact_inv[k] % mod * fact_inv[n-k] %mod


res = 0

for i in range(n+1):
    if k - a*i < 0:
        continue
    
    if (k-a*i)%b != 0:
        continue

    j = (k-a*i)//b

    if j > n:
        continue

    res = (res + mod_comb_k(n,j,mod)*mod_comb_k(n,i,mod)%mod)%mod

print(res%mod)

Solution (Python):
import sys
input = sys.stdin.readline
N, A, B, K = map(int, input().split())
mod = 998244353
class Factorial:
  def __init__(self, n, mod):
    self.f = [1]
    for i in range(1, n + 1):
      self.f.append(self.f[-1] * i % mod)
    self.i = [pow(self.f[-1], mod - 2, mod)]
    for i in range(1, n + 1)[: : -1]:
      self.i.append(self.i[-1] * i % mod)
    self.i.reverse()
  def factorial(self, i):
    return self.f[i]
  def ifactorial(self, i):
    return self.i[i]
  def combi(self, n, k):
    return self.f[n] * self.i[n - k] % mod * self.i[k] % mod

f = Factorial(N, mod)
res = 0
for r in range(N + 1):
  if K < r * A: break
  if (K - r * A) % B: continue
  b = (K - r * A) // B
  if b > N: continue
  res += f.combi(N, r) * f.combi(N, b)
  res %= mod
print(res)

Solution (Python):
n, a, b, k = map(int,input().split())

mod = 998244353
l = 3*10**5 + 3

M = [1]  
m = 1
for i in range(1, l):
    m = (m * i) % mod
    M.append(m)

def pow(x, y, mod):  
    ans = 1
    while y > 0:
        if y % 2 == 1:
            ans = (ans * x) % mod
        x = (x**2) % mod
        y //= 2
    return ans

def inv(x, mod):  
    return pow(x, mod-2, mod)



MI = [0] * (l-1) +[inv(M[l-1], mod)]  
for i in range(l-2, -1, -1):
    MI[i] = MI[i+1] * (i+1) % mod

def C(x, y):  
    if y < 0 or y > x:
        return 0
    ans = M[x]
    ans = (ans * MI[y]) % mod
    return (ans * MI[x-y]) % mod

ans = 0
for i in range(n+1):
    c0n = i
    c1 = k-i*a
    if c1 % b == 0 and 0 <= c1 // b <= n:
        c1n = c1 // b
        
        ans = (ans + C(n, c0n) * C(n, c1n)) % mod
print(ans%mod)

Solution (Python):
N, A, B, K = map(int, input().split())
maxAB = max(A, B)
minAB = min(A, B)
ans = 0


def cmb(n, r):
    if (r < 0 or r > n):
        return 0
    r = min(r, n-r)
    return fac[n] * finv[r] * finv[n-r] % mod


mod = 998244353  
N_cmb = 0
N_cmb += N  
fac = [1, 1]  
finv = [1, 1]  
inverse = [0, 1]  

for ZZZ in range(2, N_cmb + 1):
    fac.append((fac[-1] * ZZZ) % mod)
    inverse.append(mod-inverse[mod % ZZZ]*(mod//ZZZ) % mod)
    finv.append((finv[-1] * inverse[-1]) % mod)

for i in range(min(N+1, K//maxAB+1)):
    j = (K - maxAB * i) // minAB
    if (K - maxAB * i) % minAB == 0 and 0 <= j <= N:
        
        ans += cmb(N, i) * cmb(N, j) % mod
        ans %= mod
print(ans)

Solution (Python):
n,a,b,k=map(int,input().split())
m=998244353
s=[1]
for i in range(n):
    s.append((s[-1]*(n-i)*pow(i+1,m-2,m))%m)
p=0
for i in range(min(n+1,k//a+1)):
    if (k-i*a)%b==0:
        if (k-i*a)//b<=n:
            p=(p+s[i]*s[(k-i*a)//b])%m
print(p)

Solution (Python):



def modmal(a,b,mod): 

    return a * b % mod




def moddiv(a,b,mod): 

    return (a * pow(b,mod-2)) % mod




def inverse(a,mod): 
    return pow(a,mod-2)







def modfac(n, MOD):
 
    f = 1
    factorials = [1]
    for m in range(1, n + 1):
        f *= m
        f %= MOD
        factorials.append(f)
    inv = pow(f, MOD - 2, MOD)
    invs = [1] * (n + 1)
    invs[n] = inv
    for m in range(n, 1, -1):
        inv *= m
        inv %= MOD
        invs[m - 1] = inv
    return factorials, invs


def modnCr(n,r,mod,fac,inv): 

    return fac[n] * inv[n-r] * inv[r] % mod

    

N,A,B,K = map(int,input().split())

mod = 998244353
fac,inv = modfac(N,mod)
ans = 0

for i in range(K // A + 1):

    if i * A <= K and (K-i*A) % B == 0:
        j = (K-i*A)//B

        
        if i > N or j > N:
            continue
        
        ans += modnCr(N,i,mod,fac,inv) * modnCr(N,j,mod,fac,inv)
        ans %= mod

print (ans)

Solution (Python):
n,a,b,k = map(int,input().split())
mod = 998244353


def factorial(N, MOD, r=True):
    fact = [1] * (N + 1)
    rfact = [1] * (N + 1)
    r = 1
    for i in range(1, N + 1):
        fact[i] = r = r * i % MOD
    rfact[N] = r = pow(fact[N], MOD - 2, MOD)
    for i in range(N, 0, -1):
        rfact[i - 1] = r = r * i % MOD
    if r:
        return fact, rfact
    else:
        return fact, rfact

fact,rfact = factorial(n,mod)


def comb(n, k):  
    return fact[n] * rfact[k] * rfact[n - k] % mod

ans = 0
for x in range(n+1):
    if k-a*x >= 0 and (k-a*x)%b == 0:
        y = (k-a*x)//b
        if y <= n:
            ans += comb(n,x) * comb(n,y)
            ans %= mod

print(ans)

Solution (Python):


import sys

def input(): return sys.stdin.readline().strip()
def list2d(a, b, c): return [[c] * b for i in range(a)]
def list3d(a, b, c, d): return [[[d] * c for j in range(b)] for i in range(a)]
def list4d(a, b, c, d, e): return [[[[e] * d for j in range(c)] for j in range(b)] for i in range(a)]
def ceil(x, y=1): return int(-(-x // y))
def INT(): return int(input())
def MAP(): return map(int, input().split())
def LIST(N=None): return list(MAP()) if N is None else [INT() for i in range(N)]
def Yes(): print('Yes')
def No(): print('No')
def YES(): print('YES')
def NO(): print('NO')
sys.setrecursionlimit(10 ** 9)
INF = 10 ** 18
MOD = 998244353

class ModTools:
    """ 階乗たくさん使う時用のテーブル準備 """

    def __init__(self, MAX, MOD):
        """ MAX：階乗に使う数値の最大以上まで作る """
        
        MAX += 1
        self.MAX = MAX
        self.MOD = MOD
        
        factorial = [1] * MAX
        factorial[0] = factorial[1] = 1
        for i in range(2, MAX):
            factorial[i] = factorial[i-1] * i % MOD
        
        inverse = [1] * MAX
        
        inverse[MAX-1] = pow(factorial[MAX-1], MOD-2, MOD)
        for i in range(MAX-2, 0, -1):
            
            inverse[i] = inverse[i+1] * (i+1) % MOD
        self.fact = factorial
        self.inv = inverse
    
    def nCr(self, n, r):
        """ 組み合わせの数 (必要な階乗と逆元のテーブルを事前に作っておく) """

        if n < r: return 0
        
        r = min(r, n-r)
        
        numerator = self.fact[n]
        
        denominator = self.inv[r] * self.inv[n-r] % self.MOD
        return numerator * denominator % self.MOD

N, A, B, K = MAP()

ans = 0
mt = ModTools(N, MOD)
for x in range(N+1):
    
    y = (K -  A * x) / B
    if y < 0 or not y.is_integer():
        continue
    y = int(y)
    ans += mt.nCr(N, x) * mt.nCr(N, y)
    ans %= MOD
print(ans)

Solution (Python):
def examA():
    N = SI()
    ans = 0
    for i in N:
        ans += int(i)
    if ans==1:
        ans = 10
    print(ans)
    return

class combination():
    
    def __init__(self, n, mod):
        self.n = n
        self.fac = [1] * (n + 1)
        self.inv = [1] * (n + 1)
        for j in range(1, n + 1):
            self.fac[j] = self.fac[j - 1] * j % mod

        self.inv[n] = pow(self.fac[n], mod - 2, mod)
        for j in range(n - 1, -1, -1):
            self.inv[j] = self.inv[j + 1] * (j + 1) % mod

    def comb(self, n, r, mod):
        if r > n or n < 0 or r < 0:
            return 0
        return self.fac[n] * self.inv[n - r] * self.inv[r] % mod
def examB():
    N, A, B, K = LI()
    C = combination(N, mod2)
    Da = []
    for i in range(N+1):
        if (K-A*i)%B!=0 or (K-A*i)//B>N:
            continue
        Da.append((i, (K-A*i)//B))

    ans = 0
    for a,b in Da:
        cur = C.comb(N,a,mod2)*C.comb(N,b,mod2)
        ans +=cur
        ans %=mod2
    print(ans)
    return

import sys,copy,bisect,itertools,heapq,math
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
mod = 10**9 + 7
mod2 = 998244353
inf = 10**18
alphabet = [chr(ord('a') + i) for i in range(26)]

if __name__ == '__main__':
    examB()

Solution (Python):
def main():
    class Factorial:
        def __init__(self, n, mod):
            self.f = [1]
            self.mod = mod
            for j in range(1, n + 1):
                self.f.append(self.f[-1] * j % mod)
            self.i = [pow(self.f[-1], mod - 2, mod)]
            for j in range(n, 0, -1):
                self.i.append(self.i[-1] * j % mod)
            self.i.reverse()
        def factorial(self, j):
            return self.f[j]
        def ifactorial(self, j):
            return self.i[j]
        def comb(self, n, k):
            return self.f[n] * self.i[n - k] % self.mod * self.i[k] % self.mod if n >= k else 0

    N, A, B, K = map(int, input().split())
    MOD = 998244353

    c = Factorial(N + 1, MOD).comb

    ans = 0
    for a in range(N + 1):
        b = (K - (A * a)) // B
        if A * a + B * b == K and 0 <= b <= N:
            ans += (c(N, a) * c(N, b)) % MOD
            ans %= MOD

    print (ans)

if __name__ == '__main__':
    main()

Solution (Python):
class Factorial:
    def __init__(self, n, mod):
        self.f = [1]
        self.mod = mod
        for j in range(1, n + 1):
            self.f.append(self.f[-1] * j % mod)
        self.i = [pow(self.f[-1], mod - 2, mod)]
        for j in range(n, 0, -1):
            self.i.append(self.i[-1] * j % mod)
        self.i.reverse()
    def factorial(self, j):
        return self.f[j]
    def ifactorial(self, j):
        return self.i[j]
    def comb(self, n, k):
        return self.f[n] * self.i[n - k] % self.mod * self.i[k] % self.mod if n >= k else 0


N, A, B, K = map(int, input().split())
MOD = 998244353

c = Factorial(N + 1, MOD).comb

def calc(a, b):
    return (c(N, a) * c(N, b)) % MOD

ans = 0
for a in range(N + 1):
    b = (K - (A * a)) // B
    if A * a + B * b == K and 0 <= b <= N:
        ans += calc(a, b)
        ans %= MOD

print (ans)

Solution (Python):
import numpy as np
MOD = 998244353

class Factorial:
    def __init__(self, n, mod):
        self.f = [1]
        self.mod = mod
        for j in range(1, n + 1):
            self.f.append(self.f[-1] * j % mod)
        self.i = [pow(self.f[-1], mod - 2, mod)]
        for j in range(n, 0, -1):
            self.i.append(self.i[-1] * j % mod)
        self.i.reverse()
    def factorial(self, j):
        return self.f[j]
    def ifactorial(self, j):
        return self.i[j]
    def comb(self, n, k):
        return self.f[n] * self.i[n - k] % self.mod * self.i[k] % self.mod if n >= k else 0

def gcd(a, b):
    if a > b:
        a, b = b, a
    if b % a == 0:
        return a
    return gcd(b % a, a)

N, A, B, K = map(int, input().split())

if A < B:
    A, B = B, A


n = min((K + A - 1) // A, N) + 1
lcm = A * B // gcd(A, B)
step = lcm // A

c = Factorial(N + 1, MOD).comb 

def calc(a, b):
    if a > b:
        a, b = b, a 
    return (c(N, a) * c(N, b)) % MOD


def count(a, b):
    if a <= N and b <= N:
        return calc(a, b)
    return 0

def maike_satrt(step):
    for i in range(step):
        if (K - A * i) % B == 0:
            return i

start = maike_satrt(step)

ans = 0
for a in range(start, n, step): 
    if A * a > K:
        break
    
    
    b = (K - A * a) // B
    
    ans += count(a, b)
    ans %= MOD

print (ans)

Solution (Python):
import sys
mod=998244353
MAX=3*10**5+100
g1=[1,1]
g2=[1,1]
for i in range(2,MAX+1):
    num_1=g1[-1]*i%mod
    g1.append(num_1)
    g2.append(pow(num_1,mod-2,mod))
def cmb(n,r):
    return g1[n]*g2[r]*g2[n-r]%mod

N,A,B,K=map(int,input().split())
data=[]
for x in range(N+1):
    if (K-A*x)%B!=0:
        continue
    y=(K-A*x)//B
    if 0<=x<=N and 0<=y<=N:
        data.append([x,y])
ans=0
for x,y in data:
    ans=(ans+cmb(N,x)*cmb(N,y)%mod)%mod
print(ans)

Solution (Python):

MAX = 310000
MOD = 998244353
factrial = [0]*MAX
inverse = [0]*MAX
factrial_inverse = [0]*MAX



def COMinit():
    global factrial, inverse, factrial_inverse
    factrial[0] = 1
    factrial[1] = 1
    inverse[1] = 1
    factrial_inverse[0] = 1
    factrial_inverse[1] = 1
    for i in range(2, MAX):
        factrial[i] = factrial[i-1] * i % MOD
        inverse[i] = MOD - inverse[MOD % i] * (MOD//i) % MOD
        factrial_inverse[i] = factrial_inverse[i-1] * inverse[i] % MOD



def COM(n, k):
    global factrial, inverse, factrial_inverse
    if n < k:
        return 0
    if n < 0 or k < 0:
        return 0
    return factrial[n] * (factrial_inverse[k] * factrial_inverse[n-k] % MOD) % MOD



COMinit()
N, A, B, K = map(int, input().split())
ans = 0
for a in range(N+1):
    
    if (K - A*a) % B == 0:
        b = (K - A*a) // B
        ans += COM(N, a) * COM(N, b) % MOD

print(ans % MOD)

Solution (Python):
def main():
    MOD = 998244353

    def cmb(n, r):
        if n - r < r:
            return cmb(n, n - r)
        if r < 0:
            return 0
        if r == 0:
            return 1
        if r == 1:
            return n
        return f[n] * invs[r] * invs[n - r] % MOD

    N, A, B, K = map(int, input().split())

    f = [-1] * (N + 1)
    f[0] = 1
    t = 1
    for i in range(1, N + 1):
        t *= i
        t %= MOD
        f[i] = t

    inv = pow(t, MOD - 2, MOD)
    invs = [-1] * (N + 1)
    invs[N] = inv
    for i in range(N - 1, 0, -1):
        inv *= (i + 1)
        inv %= MOD
        invs[i] = inv

    ret = 0
    for i in range(K // A + 1):
        j, r = divmod(K - A * i, B)
        if r == 0:
            ret = (ret + cmb(N, i) * cmb(N, j)) % MOD

    print(ret)


if __name__ == '__main__':
    main()

Solution (Python):
mod = 998244353

class ModComb:
    def __init__(self, MAX):
        fac = [1, 1]
        finv = [1, 1]
        inv = [0, 1]
        for i in range(2, MAX):
            fac.append(fac[i - 1] * i % mod)
            inv.append(mod - inv[mod % i] * (mod // i) % mod)
            finv.append(finv[i - 1] * inv[i] % mod)
        self.fac, self.finv, self.mod = fac, finv, mod

    def nCk(self, n, k):
        if n < k or n < 0 or k < 0:
            return 0
        fac, finv, mod = self.fac, self.finv, self.mod
        return fac[n] * (finv[k] * finv[n - k] % mod) % mod

N, A, B, K = map(int, input().split())

mc = ModComb(3 * 10 ** 5 + 10)

ans = 0
for a in range(N + 1):
    b, r = divmod(K - A * a, B)
    if r != 0:
        continue
    ans += mc.nCk(N, a) * mc.nCk(N, b) % mod
    ans %= mod

print(ans)

Solution (Python):
import sys
sys.setrecursionlimit(2147483647)
INF=float("inf")
MOD=998244353
input=lambda:sys.stdin.readline().rstrip()

class modfact(object):
    def __init__(self,n):
        fact=[1]*(n+1); invfact=[1]*(n+1)
        for i in range(1,n+1): fact[i]=i*fact[i-1]%MOD
        invfact[n]=pow(fact[n],MOD-2,MOD)
        for i in range(n-1,-1,-1): invfact[i]=invfact[i+1]*(i+1)%MOD
        self.__fact=fact; self.__invfact=invfact

    def inv(self,n):
        assert(n>0)
        return self.__fact[n-1]*self.__invfact[n]%MOD

    def fact(self,n):
        return self.__fact[n]

    def invfact(self,n):
        return self.__invfact[n]

    def comb(self,n,k):
        if(k<0 or n<k): return 0
        return self.__fact[n]*self.__invfact[k]*self.__invfact[n-k]%MOD

    def perm(self,n,k):
        if(k<0 or n<k): return 0
        return self.__fact[n]*self.__invfact[n-k]%MOD

from fractions import gcd
def resolve():
    n,a,b,k=map(int,input().split())

    
    g=gcd(a,b)
    if(k%g):
        print(0)
        return
    k//=g; a//=g; b//=g

    
    for x in range(b):
        if((k-a*x)%b==0):
            y=(k-a*x)//b
            break
    
    
    ans=0
    mf=modfact(n)
    while(y>=0):
        ans+=mf.comb(n,x)*mf.comb(n,y)
        x,y=x+b,y-a

    print(ans%MOD)
resolve()

Solution (Python):
N, A, B, K = map(int, input().split())
mod = 998244353
m = N

fac = [1] * (m + 1)
inv = [1] * (m + 1)
for i in range(1, m + 1):
	fac[i] = fac[i - 1] * i % mod
inv[-1] = pow(fac[-1], mod - 2, mod)
for i in range(m - 1, -1, -1):
	inv[i] = inv[i + 1] * (i + 1) % mod

def cmb(n, r):
	assert n >= r >= 0
	return fac[n] * inv[n - r] * inv[r] % mod

ans = 0
for r in range(N + 1):
	if (K - A * r) % B == 0:
		b = (K - A * r) // B
		if 0 <= b <= N:
			ans += cmb(N, r) * cmb(N, (K - A * r) // B)
			ans %= mod
print(ans)

Solution (Python):
MOD = 998244353
n, a, b, k = map(int, input().split())
fac = [1] * (n + 1)
for i in range(n):
    fac[i + 1] = fac[i] * (i + 1) % MOD


def comb(n, k):
    return fac[n] * pow(fac[n - k], MOD - 2, MOD) * pow(fac[k], MOD - 2, MOD)


ans = 0
for x in range(n + 1):
    if k - x * a < 0 or (k - x * a) % b:
        continue
    y = (k - x * a) // b
    if n < y:
        continue
    ans += (comb(n, x) * comb(n, y)) % MOD
print(ans % MOD)

Solution (Python):
def main():
    n,a,b,k=map(int,input().split())
    fact=[1]
    mod=998244353
    for i in range(1,n+1):fact.append(fact[-1]*i%mod)
    fn=fact[n]
    ncr=[pow(fact[n-i]*fact[i],mod-2,mod)*fn%mod for i in range(n+1)]
    ans=0
    ans_list=[]
    for i in range(n+1):
        j=(k-a*i)//b
        if i*a+j*b==k:
            if 0<=j<=n:
                ans_list.append([i,j])
    for p,q in ans_list:
        ans=(ans+ncr[p]*ncr[q])%mod
    print(ans)
main()

Solution (Python):
def main():
    n,a,b,k=map(int,input().split())
    fact=[1]
    mod=998244353
    for i in range(1,n+1):fact.append(fact[-1]*i%mod)
    def inv_n(n):return pow(n,mod-2,mod)
    def nCr(n,r):return inv_n(fact[n-r]*fact[r]%mod)*fact[n]%mod
    ans=0
    ans_list=[]
    for i in range(n+1):
        j=(k-a*i)//b
        if i*a+j*b==k:
            if 0<=j<=n:
                ans_list.append([i,j])
    for p,q in ans_list:
        ans=(ans+nCr(n,p)*nCr(n,q))%mod
    print(ans)
main()

Solution (Python):
mod = 998244353
N, A, B, K = [int(_) for _ in input().split()]


def extgcd(a, b):
    if b:
        x, y = extgcd(b, a % b)
        x, y = y, x - (a // b) * y
    else:
        x, y = 1, 0
    return x, y


E = extgcd(A, B)
g = A * E[0] + B * E[1]
l = A * B // g
F = [E[0] * (K // g), E[1] * (K // g)]
if F[0] > F[1]:
    F = F[::-1]
    A, B = B, A
H = [F[0] + B * (F[1] // A), F[1] - A * (F[1] // A)]
L = []
for i in range(-100000, max(H[0] // (B // g), H[1] // (A // g)) + 100000):
    J = [H[0] - (B // g) * i, H[1] + (A // g) * i]
    if 0 <= J[0] <= N and 0 <= J[1] <= N:
        L += [J]
f = [1] * (N + 1)
for i in range(2, N + 1):
    f[i] = f[i - 1] * i % mod
fi = [pow(f[-1], mod - 2, mod)]
for i in range(N, 0, -1):
    fi += [fi[-1] * i % mod]
fi = fi[::-1]


def comb(x, y):
    return (f[x] * fi[y] * fi[x - y]) % mod


ans = 0
for pq in L:
    p, q = pq
    ans += comb(N, p) * comb(N, q)
    ans %= mod
print(ans)

Solution (Python):
N, A, B, K = map(int, input().split())
mod = 998244353

fac = [1]
inv = [1]
for i in range(1, N + 1) :
    fac.append(fac[-1] * i % mod)
    inv.append(pow(fac[-1], mod - 2, mod))

ret = 0
for na in range(N + 1) :
    if (K - na * A) % B != 0 :
        continue
        
    nb = (K - na * A) // B
    if 0 <= na <= N and 0 <= nb <= N :
        ret += fac[N] * inv[N-na] * inv[na] * fac[N] * inv[N-nb] * inv[nb]
        ret %= mod
        
print(ret)

Solution (Python):

from collections import defaultdict
con = 998244353
INF = float("inf")

def getlist():
	return list(map(int, input().split()))

N, A, B, K = getlist()

fac = [0] * (N + 1)
inv = [0] * (N + 1)
finv = [0] * (N + 1)
fac[0], fac[1] = 1, 1
inv[1] = 1
finv[0], finv[1] = 1, 1


for i in range(2, N + 1):
	fac[i] = fac[i - 1] * i % con
	inv[i] = inv[con % i] * (con - (con // i)) % con
	finv[i] = finv[i - 1] * inv[i] % con


def main():
	ans = 0
	for i in range(N + 1):
		if (K - i * A) % B == 0:
			j = int((K - i * A) // B)
			if 0 <= j and j <= N:
				ans += (fac[N] * finv[i] * finv[N - i] * fac[N] * finv[j] * finv[N - j]) % con
				ans %= con

	
	print(ans)

if __name__ == '__main__':
	main()

Solution (Python):

mod = 998244353 

def cmb(n, r, mod=mod):
    if ( r<0 or r>n ):
        return 0
    r = min(r, n-r)
    return g1[n] * g2[r] * g2[n-r] % mod

NN = 10**6 
g1 = [1, 1] 
g2 = [1, 1] 
inverse = [0, 1] 

for i in range( 2, NN + 1 ):
    g1.append( ( g1[-1] * i ) % mod )
    inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
    g2.append( (g2[-1] * inverse[-1]) % mod )

N, A, B, K = map(int, input().split())

ans = 0
for n in range(N+1):
    m = (K-n*A)//B
    if (K-n*A)%B == 0 and 0 <= m <= N:
        ans += cmb(N, m)*cmb(N, n) % mod
        ans %= mod
print(ans)

Solution (Python):
from heapq import heappush, heappop
from collections import deque,defaultdict,Counter
import itertools
from itertools import permutations,combinations
import sys
import bisect
import string



def I():
    return int(input())
def MI():
    return map(int,input().split())
def LI():
    return [int(i) for i in input().split()]
def LI_():
    return [int(i)-1 for i in input().split()]
def StoI():
    return [ord(i)-97 for i in input()]
def show(*inp,end='\n'):
    if show_flg:
        print(*inp,end=end)
YN=['Yes','No']
mo=10**9+7
inf=float('inf')


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

import sys
sys.setrecursionlimit(10**7)
mo=998244353

class PERM_COMB_MOD():
    
    def __init__(self, max_n=510000, mod=10**9+7):
        self.fac = [0]*max_n
        self.finv = [0]*max_n
        self.inv = [0]*max_n
        self.fac[0] = self.fac[1] = 1
        self.finv[0] = self.finv[1] = 1
        self.inv[1] = 1
        self.max = max_n
        self.mod = mod
        self._maesyori()
 
    def _maesyori(self):
        for i in range(2,self.max):
            self.fac[i] = self.fac[i-1] * i % self.mod
            self.inv[i] = self.mod - self.inv[self.mod % i] * (self.mod // i) % self.mod
            self.finv[i] = self.finv[i-1] * self.inv[i] % self.mod
 
    def comb(self, n, k):
        if n < k : return 0
        if n < 0 or k < 0:return 0
        return self.fac[n] * (self.finv[k] * self.finv[n-k] % self.mod) % self.mod

n,a,b,k=LI()
com=PERM_COMB_MOD(mod=mo)

show_flg=False
show_flg=True

if a<b:
    a,b=b,a
ans=0
for i in range(min(n+1,k//a+1)):
    j=(k-a*i)//b
    if (k-a*i)%b!=0 or j<0 or j>n:
        continue
    try:
        ans+=com.comb(n,i)*com.comb(n,j)
    except:
        0
    ans%=mo

print(ans%mo)

Solution (Python):
from fractions import gcd

class PERM_COMB_MOD():
    
    def __init__(self, max_n=510000, mod=10**9+7):
        self.fac = [0]*max_n
        self.finv = [0]*max_n
        self.inv = [0]*max_n
        self.fac[0] = self.fac[1] = 1
        self.finv[0] = self.finv[1] = 1
        self.inv[1] = 1
        self.max = max_n
        self.mod = mod
        self._maesyori()

    def _maesyori(self):
        for i in range(2,self.max):
            self.fac[i] = self.fac[i-1] * i % self.mod
            self.inv[i] = self.mod - self.inv[self.mod % i] * (self.mod // i) % self.mod
            self.finv[i] = self.finv[i-1] * self.inv[i] % self.mod

    def comb(self, n, k):
        if n < k : return 0
        if n < 0 or k < 0:return 0
        return self.fac[n] * (self.finv[k] * self.finv[n-k] % self.mod) % self.mod

mod = 998244353
n,a,b,k = map(int, input().split())
if k % gcd(a,b) != 0:
    print(0)
    exit()
hoge = PERM_COMB_MOD(max_n=n+10, mod=mod)
ans = 0
for i in range(n+1):
    tmp = k-a*i
    if tmp % b:
        continue
    j = tmp//b
    if j < 0 or n < j:
        continue
    ans += hoge.comb(n, i)*hoge.comb(n, j)
    ans %= mod
print(ans)

Solution (Python):
from fractions import gcd

class PERM_COMB_MOD():
    
    def __init__(self, max_n=510000, mod=10**9+7):
        self.fac = [0]*max_n
        self.finv = [0]*max_n
        self.inv = [0]*max_n
        self.fac[0] = self.fac[1] = 1
        self.finv[0] = self.finv[1] = 1
        self.inv[1] = 1
        self.max = max_n
        self.mod = mod
        self._maesyori()

    def _maesyori(self):
        for i in range(2,self.max):
            self.fac[i] = self.fac[i-1] * i % self.mod
            self.inv[i] = self.mod - self.inv[self.mod % i] * (self.mod // i) % self.mod
            self.finv[i] = self.finv[i-1] * self.inv[i] % self.mod

    def comb(self, n, k):
        if n < k : return 0
        if n < 0 or k < 0:return 0
        return self.fac[n] * (self.finv[k] * self.finv[n-k] % self.mod) % self.mod

mod = 998244353
hoge = PERM_COMB_MOD(mod=mod)
n,a,b,k = map(int, input().split())
if k % gcd(a,b) != 0:
    print(0)
    exit()
ans = 0
for i in range(n+1):
    tmp = k-a*i
    if tmp % b:
        continue
    j = tmp//b
    if j < 0 or n < j:
        continue
    ans += hoge.comb(n, i)*hoge.comb(n, j)
    ans %= mod
print(ans)

Solution (Python):

import sys
def input():
	return sys.stdin.readline().strip()

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

mod = 998244353

fact_list = [1]
for i in range(1, N + 1):
	fact_list.append((fact_list[-1]*i) % mod)

def pow(x, y, mod=1000000007):
	pow_list = [x]
	i = 1
	while 2**i <= y:
		a = pow_list[-1]**2
		if mod != 0:
			a = a % mod
		pow_list.append(a)
		i += 1
	ans = 1
	for bit in range(len(pow_list)):
		if (2**bit) & y != 0:
			ans *= pow_list[bit]
			if mod != 0:
				ans = ans % mod
	return ans

def fact(n, mod=1000000007):
	return fact_list[n]

inv_t = [0] + [1]
for i in range(2,N + 1):
	inv_t += [inv_t[mod % i] * (mod - int(mod / i)) % mod]

def mod_rev(x, mod):
	"""
	関数powが必要
	"""
	return inv_t[x]

inv_fact = [1] + [1]

for i in range(2, N + 1):
	inv_fact.append((inv_fact[-1] * inv_t[i]) % mod)


def comb(a, b, mod=1000000007):
	"""
	関数mod_rev, fact, powが必要
	"""
	if a < b or b < 0:
		return 0
	if a == b or b == 0:
		return 1
	a_fact = fact(a, mod)
	b_fact = fact(b, mod)
	a_b_fact = fact(a - b, mod)
	if mod != 0:
		return (a_fact * inv_fact[b] * inv_fact[a - b])%mod
	else:
		return a_fact//(b_fact*a_b_fact)

ans = 0
times = 0
for r in range(N + 1):
	rest = K - r * A
	if rest % B == 0 and rest // B >= 0 and rest // B <= N:
		b = rest // B
		ans += comb(N, r, mod)*comb(N, b, mod)
		ans %= mod
		times += 1

print(ans)