Potion Brewing Class (Codeforces 1654D)
Tags: Math, number theory

Alice's potion making professor gave the following assignment to his students: brew a potion using $$$n$$$ ingredients, such that the proportion of ingredient $$$i$$$ in the final potion is $$$r_i > 0$$$ (and $$$r_1 + r_2 + \cdots + r_n = 1$$$).He forgot the recipe, and now all he remembers is a set of $$$n-1$$$ facts of the form, "ingredients $$$i$$$ and $$$j$$$ should have a ratio of $$$x$$$ to $$$y$$$" (i.e., if $$$a_i$$$ and $$$a_j$$$ are the amounts of ingredient $$$i$$$ and $$$j$$$ in the potion respectively, then it must hold $$$a_i/a_j = x/y$$$), where $$$x$$$ and $$$y$$$ are positive integers. However, it is guaranteed that the set of facts he remembers is sufficient to uniquely determine the original values $$$r_i$$$.He decided that he will allow the students to pass the class as long as they submit a potion which satisfies all of the $$$n-1$$$ requirements (there may be many such satisfactory potions), and contains a positive integer amount of each ingredient.Find the minimum total amount of ingredients needed to make a potion which passes the class. As the result can be very large, you should print the answer modulo $$$998\,244\,353$$$.

Input: The first line contains a single integer $$$t$$$ ($$$1 \le t \le 10^4$$$) — the number of test cases.The first line of each test case contains a single integer $$$n$$$ ($$$2 \le n \le 2 \cdot 10^5$$$).Each of the next $$$n-1$$$ lines contains four integers $$$i, j, x, y$$$ ($$$1 \le i, j \le n$$$, $$$i\not=j$$$, $$$1\le x, y \le n$$$) — ingredients $$$i$$$ and $$$j$$$ should have a ratio of $$$x$$$ to $$$y$$$. It is guaranteed that the set of facts is sufficient to uniquely determine the original values $$$r_i$$$.It is also guaranteed that the sum of $$$n$$$ for all test cases does not exceed $$$2 \cdot 10^5$$$.

Output: For each test case, print the minimum total amount of ingredients needed to make a potion which passes the class, modulo $$$998\,244\,353$$$.

Note: In the first test case, the minimum total amount of ingredients is $$$69$$$. In fact, the amounts of ingredients $$$1, 2, 3, 4$$$ of a valid potion are $$$16, 12, 9, 32$$$, respectively. The potion is valid because   Ingredients $$$3$$$ and $$$2$$$ have a ratio of $$$9 : 12 = 3 : 4$$$;  Ingredients $$$1$$$ and $$$2$$$ have a ratio of $$$16 : 12 = 4 : 3$$$;  Ingredients $$$1$$$ and $$$4$$$ have a ratio of $$$16 : 32 = 2 : 4$$$. In the second test case, the amounts of ingredients $$$1, 2, 3, 4, 5, 6, 7, 8$$$ in the potion that minimizes the total amount of ingredients are $$$60, 60, 24, 48, 32, 60, 45, 30$$$.

Sample Inputs:
3
4
3 2 3 4
1 2 4 3
1 4 2 4
8
5 4 2 3
6 4 5 4
1 3 5 2
6 8 2 1
3 5 3 4
3 2 2 5
6 7 4 3
17
8 7 4 16
9 17 4 5
5 14 13 12
11 1 17 14
6 13 8 9
2 11 3 11
4 17 7 2
17 16 8 6
15 5 1 14
16 7 1 10
12 17 13 10
11 16 7 2
10 11 6 4
13 17 14 6
3 11 15 8
15 6 12 8

Sample Outputs:
69
359
573672453

Solution (Python):
import sys, os, io
input = io.BytesIO(os.read(0, os.fstat(0).st_size)).readline

def the_sieve_of_eratosthenes(n):
    s = [1] * (n + 1)
    x = []
    for i in range(2, n + 1):
        if s[i]:
            x.append(i)
            for j in range(i, n + 1, i):
                s[j] = 0
    return x

def fast_prime_factorization(n):
    d = [(i + 1) % 2 * 2 for i in range(n + 1)]
    (d[0], d[1]) = (0, 1)
    for i in range(3, n + 1):
        if d[i]:
            continue
        for j in range(i, n + 1, 2 * i):
            if not d[j]:
                d[j] = i
    return d

def pow(x, y, mod):
    if x == -1:
        return 1 - y % 2 * 2
    (ans, u) = (1, x)
    while y:
        if y & 1:
            ans = ans * u % mod
        y >>= 1
        u = u * u % mod
    return ans
t = int(input())
ans = []
mod = 998244353
inf = pow(10, 9, mod) + 1
for _ in range(t):
    n = int(input())
    G = [[] for _ in range(n + 1)]
    for _ in range(n - 1):
        (i, j, x, y) = map(int, input().split())
        G[i].append((j, x, y))
        G[j].append((i, y, x))
    d = the_sieve_of_eratosthenes(n)
    p = fast_prime_factorization(n)
    (cnt, mi) = ([inf] * (n + 1), [inf] * (n + 1))
    (st, e) = ([1], [0])
    (bx, by) = ([0] * (n + 1), [0] * (n + 1))
    for _ in range(3 * n):
        (i, j) = (st[-1], -1)
        if i == 1 and len(G[1]) == e[-1]:
            break
        if len(G[i]) == e[-1]:
            (j, x, y) = (st[-2], bx[i], by[i])
            st.pop()
            e.pop()
            e[-1] += 1
        else:
            (j, x, y) = G[i][e[-1]]
            if len(st) > 1 and st[-2] == j:
                (j, bx[i], by[i]) = (-1, x, y)
                e[-1] += 1
            else:
                st.append(j)
                e.append(0)
        if j ^ -1:
            u = x
            while p[u] ^ u:
                cnt[p[u]] -= 1
                u //= p[u]
            if u ^ 1:
                cnt[u] -= 1
            u = y
            while p[u] ^ u:
                cnt[p[u]] += 1
                u //= p[u]
            if u ^ 1:
                cnt[u] += 1
            u = x
            while p[u] ^ u:
                v = p[u]
                mi[v] = min(mi[v], cnt[v])
                u //= v
            if u ^ 1:
                mi[u] = min(mi[u], cnt[u])
    u = 1
    for i in d:
        u = u * pow(i, cnt[i] - mi[i], mod) % mod
    q = [1]
    z = [0] * (n + 1)
    z[1] = u
    visit = [0] * (n + 1)
    visit[1] = 1
    for k in range(n):
        i = q[k]
        for (j, x, y) in G[i]:
            if not visit[j]:
                visit[j] = 1
                q.append(j)
                z[j] = z[i] * y % mod * pow(x, mod - 2, mod) % mod
    ans0 = sum(z) % mod
    ans.append(ans0)
sys.stdout.write('\n'.join(map(str, ans)))

Solution (Python):
import os, sys
from io import BytesIO, IOBase

def main():
    mod = 998244353
    lim = 200001
    (sieve, primes) = ([[] for _ in range(lim)], [])
    for i in range(2, lim):
        if not len(sieve[i]):
            primes.append(i)
            for j in range(i, lim, i):
                xx = j
                while not xx % i:
                    sieve[j].append(i)
                    xx //= i
    for _ in range(int(input())):
        n = int(input())
        path = [[] for _ in range(n)]
        col = {}
        for _ in range(n - 1):
            (i, j, x, y) = map(int, input().split())
            path[i - 1].append((j - 1, x, y))
            path[j - 1].append((i - 1, y, x))
            col[i - 1, j - 1] = (x, y)
            col[j - 1, i - 1] = (y, x)
        curr = {i: 0 for i in primes if i <= n}
        mini = [0] * (n + 1)
        (st, poi, visi) = ([0], [0] * n, [1] + [0] * (n - 1))
        while len(st):
            (x, y) = (st[-1], path[st[-1]])
            if poi[x] != len(y) and visi[y[poi[x]][0]]:
                poi[x] += 1
            if poi[x] == len(y):
                st.pop()
                if len(st):
                    (co, co1) = col[x, st[-1]]
                    for xx in sieve[co]:
                        curr[xx] -= 1
                    for xx in sieve[co1]:
                        curr[xx] += 1
                    for xx in sieve[co] + sieve[co1]:
                        mini[xx] = min(mini[xx], curr[xx])
            else:
                (z, co, co1) = y[poi[x]]
                visi[z] = 1
                st.append(z)
                poi[x] += 1
                for xx in sieve[co]:
                    curr[xx] -= 1
                for xx in sieve[co1]:
                    curr[xx] += 1
                for xx in sieve[co] + sieve[co1]:
                    mini[xx] = min(mini[xx], curr[xx])
        ini = 1
        for i in primes:
            if i >= len(mini):
                break
            ini = ini * pow(i, -mini[i], mod) % mod
        (st, poi, visi) = ([0], [0] * n, [ini] + [0] * (n - 1))
        while len(st):
            (x, y) = (st[-1], path[st[-1]])
            if poi[x] != len(y) and visi[y[poi[x]][0]]:
                poi[x] += 1
            if poi[x] == len(y):
                st.pop()
            else:
                (z, co, co1) = y[poi[x]]
                visi[z] = visi[x] * pow(co, mod - 2, mod) * co1 % mod
                st.append(z)
                poi[x] += 1
        ans = 0
        for i in visi:
            ans = (ans + i) % mod
        print(ans)
BUFSIZE = 8192

class FastIO(IOBase):
    newlines = 0

    def __init__(self, file):
        self._fd = file.fileno()
        self.buffer = BytesIO()
        self.writable = 'x' in file.mode or 'r' not in file.mode
        self.write = self.buffer.write if self.writable else None

    def read(self):
        while True:
            b = os.read(self._fd, max(os.fstat(self._fd).st_size, BUFSIZE))
            if not b:
                break
            ptr = self.buffer.tell()
            (self.buffer.seek(0, 2), self.buffer.write(b), self.buffer.seek(ptr))
        self.newlines = 0
        return self.buffer.read()

    def readline(self):
        while self.newlines == 0:
            b = os.read(self._fd, max(os.fstat(self._fd).st_size, BUFSIZE))
            self.newlines = b.count(b'\n') + (not b)
            ptr = self.buffer.tell()
            (self.buffer.seek(0, 2), self.buffer.write(b), self.buffer.seek(ptr))
        self.newlines -= 1
        return self.buffer.readline()

    def flush(self):
        if self.writable:
            os.write(self._fd, self.buffer.getvalue())
            (self.buffer.truncate(0), self.buffer.seek(0))

class IOWrapper(IOBase):

    def __init__(self, file):
        self.buffer = FastIO(file)
        self.flush = self.buffer.flush
        self.writable = self.buffer.writable
        self.write = lambda s: self.buffer.write(s.encode('ascii'))
        self.read = lambda : self.buffer.read().decode('ascii')
        self.readline = lambda : self.buffer.readline().decode('ascii')
(sys.stdin, sys.stdout) = (IOWrapper(sys.stdin), IOWrapper(sys.stdout))
input = lambda : sys.stdin.readline().rstrip('\r\n')
if __name__ == '__main__':
    main()

Solution (Python):
from collections import deque, Counter

def sieve(n):
    isprime = [1] * (n + 1)
    for i in range(2, n + 1):
        if isprime[i] == 1:
            isprime[i] = i
            for j in range(i * i, n + 1, i):
                isprime[j] = i
    return isprime
prime = sieve(2 * 10 ** 5 + 1)
MOD = 998244353

def solve():
    n = int(input())
    g = [[] for i in range(n)]
    for _ in range(n - 1):
        (i, j, x, y) = map(int, input().split())
        i -= 1
        j -= 1
        g[i].append((j, x, y))
        g[j].append((i, y, x))
    queue = deque()
    queue.append(0)
    cnt = Counter()
    cnt1 = Counter()

    def f1(m):
        while m > 1:
            p = prime[m]
            m //= p
            cnt1[p] += 1
            cnt[p] = max(cnt[p], cnt1[p])

    def f2(m):
        while m > 1:
            p = prime[m]
            m //= p
            cnt1[p] -= 1
    v1 = [0] * n
    v2 = [0] * n
    r = [(1, 1)] * n
    st = [0]
    while st:
        i = st[-1]
        if v1[i] == 0:
            v1[i] = 1
            (x, y) = r[i]
            f2(y)
            f1(x)
            for (j, x, y) in g[i]:
                if v1[j]:
                    continue
                r[j] = (x, y)
                st.append(j)
        else:
            st.pop()
            v2[i] = 1
            (x, y) = r[i]
            f2(x)
            f1(y)
    w = 1
    for (p, k) in cnt.items():
        w = w * pow(p, k, MOD) % MOD
    ans = 0
    v = [0] * n
    v[0] = 1
    queue = deque()
    queue.append((0, w))
    while queue:
        (i, wi) = queue.popleft()
        ans = (ans + wi) % MOD
        for (j, x, y) in g[i]:
            if v[j] == 0:
                v[j] = 1
                queue.append((j, wi * y * pow(x, MOD - 2, MOD) % MOD))
    return ans
import sys
input = lambda : sys.stdin.readline().rstrip()
t = int(input())
for i in range(t):
    print(solve())