QOJ.ac

QOJ

IDProblemSubmitterResultTimeMemoryLanguageFile sizeSubmit timeJudge time
#91166#6126. Sequence and SequenceGolovanov399AC ✓4588ms225116kbPython37.8kb2023-03-27 15:50:212023-03-27 15:50:23

Judging History

你现在查看的是最新测评结果

  • [2023-08-10 23:21:45]
  • System Update: QOJ starts to keep a history of the judgings of all the submissions.
  • [2023-03-27 15:50:23]
  • 评测
  • 测评结果:AC
  • 用时:4588ms
  • 内存:225116kb
  • [2023-03-27 15:50:21]
  • 提交

answer

#!/usr/bin/env python3

import math
from fractions import Fraction

def add(p, q):
	res = [0] * max(len(p), len(q))
	for i in range(len(p)):
		res[i] += p[i]
	for i in range(len(q)):
		res[i] += q[i]
	return res

def sub(p, q):
	res = [0] * max(len(p), len(q))
	for i in range(len(p)):
		res[i] += p[i]
	for i in range(len(q)):
		res[i] -= q[i]
	return res

def mult(p, q):
	res = [0] * ((len(p) - 1) + (len(q) - 1) + 1)
	for i in range(len(p)):
		for j in range(len(q)):
			res[i + j] += p[i] * q[j]
	return res

def subst(p, q):
	# P(Q(x))
	res = [0]
	cur = [1]
	for c in p:
		res = add(res, [x * c for x in cur])
		cur = mult(cur, q)
	return res

def apply_poly(p, x):
	res = 0
	for c in p[::-1]:
		res = res * x + c
	return res

def norm_to_bin(p):
	res = [0] * len(p)
	bins = [[1]]
	for i in range(len(p) - 1):
		q = mult(bins[-1], [i, 1])
		bins.append(q)
	for i in range(len(p) - 1, -1, -1):
		c = p[i]
		res[i] = c
		p = add(p, [-c * x for x in bins[i]])
	return res

def lift(p):
	p = [0] + p
	for i in range(1, len(p)):
		p[i] = Fraction(p[i]) / i
	return p

def bin_to_norm(p):
	res = [0]
	q = [1]
	for i in range(len(p)):
		res = add(res, [x * p[i] for x in q])
		q = mult(q, [i, 1])
	return res

def _poly_to_sum(p):
	return bin_to_norm(lift(norm_to_bin(p)))

# N = 3_000
# D = 37

N = 184_000
D = 17
pss = [_poly_to_sum([0] * i + [1]) for i in range(D)]

def poly_to_sum(p):
	res = [0]
	for i in range(len(p)):
		res = add(res, [x * p[i] for x in pss[i]])
	return res

def sumpoly(p, x):
	return apply_poly(poly_to_sum(p), x)

def varsumpoly(p, q):
	return subst(poly_to_sum(p), q)


def isqrt(n: int) -> int: # pylint: disable=too-many-branches
    # pylint: disable=line-too-long # Accommodate long link URLs on docstring.
    """
    Returns the largest root such that ``root**2 <= n (root + 1)**2 > n``.
    When using Python 3.8 or later, this function acts as a wrapper for the
    built-in :obj:`math.isqrt` function.
    For all other supported versions of Python, this function reverts to a
    pure Python algorithm that is adapted from an
    `implementation by Alexander Gosselin <https://gist.github.com/castle-bravo/e841684d6bad8e0598e31862a7afcfc7>`__,
    which is based on a `Stack Overflow answer by Tobin Fricke <http://stackoverflow.com/a/23279113/2738025>`__.
    >>> isqrt(4)
    2
    >>> isqrt(16)
    4
    >>> list(map(isqrt, range(16, 26)))
    [4, 4, 4, 4, 4, 4, 4, 4, 4, 5]
    >>> from random import randint
    >>> all([isqrt(r**2 + randint(0, r)) == r  for r in range(0, 1000)])
    True
    >>> r = randint(2**511, 2**512 - 1)
    >>> isqrt(r**2) == r
    True
    >>> isqrt(2**30000) == 2**15000
    True
    The type and sign of the input are checked.
    >>> isqrt(-2)
    Traceback (most recent call last):
      ...
    ValueError: input must be a non-negative integer
    >>> isqrt('abc')
    Traceback (most recent call last):
      ...
    TypeError: input must be an integer
    Test scenarios in which the :obj:`math.isqrt` function is not available.
    >>> if hasattr(math, 'isqrt'):
    ...     delattr(math, 'isqrt')
    >>> isqrt(16)
    4
    >>> isqrt(2**30000) == 2**15000
    True
    >>> isqrt(-2)
    Traceback (most recent call last):
      ...
    ValueError: input must be a non-negative integer
    >>> isqrt('abc')
    Traceback (most recent call last):
      ...
    TypeError: input must be an integer
    """
    # Try using built-in integer square root function (available in Python 3.8 or later).
    # To ensure this implementation is backwards-compatible with previous versions while
    # not introducing performance overheads for the most common case, only convert
    # exceptions if the initial call to the built-in function raises an exception.
    if hasattr(math, 'isqrt'):
        try:
            return math.isqrt(n)
        except ValueError as e:
            if str(e) == 'isqrt() argument must be nonnegative':
                raise ValueError('input must be a non-negative integer') from None
            # Continue to default implementation to ensure backwards-compatible behavior.
        except TypeError as e:
            if str(e).endswith('object cannot be interpreted as an integer'):
                raise TypeError('input must be an integer') from None
            # Continue to default implementation to ensure backwards-compatible behavior.
        except: # pylint: disable=bare-except # pragma: no cover
            pass # Continue to default implementation to ensure backwards-compatible behavior.

    try: # Attempt to use the :obj:`math.sqrt` function.
        root = int(math.sqrt(n))
        if root * root == n: # No error from floating point conversion.
            return root
    except ValueError as e:
        if str(e) == 'math domain error':
            raise ValueError('input must be a non-negative integer') from None
        # Continue to default implementation to ensure backwards-compatible behavior.
    except TypeError as e:
        if str(e).startswith('must be real number'):
            raise TypeError('input must be an integer') from None
        # Continue to default implementation to ensure backwards-compatible behavior.
    except OverflowError:
        pass # Use the integer-only bit-wise algorithm.

    if n is None or (not isinstance(n, int)):
        raise TypeError('input must be an integer') # pragma: no cover

    if n < 0:
        raise ValueError('input must be a non-negative integer') # pragma: no cover

    root = 0 # Running result.
    rmdr = 0 # Running remainder n - root**2.
    for s in reversed(range(0, n.bit_length(), 2)): # Shift n by s bits.
        bits = n >> s & 3 # The next two most significant bits of n.
        rmdr = rmdr << 2 | bits # Increase the remainder.
        cand = root << 2 | 1 # Shifted candidate root value to try.
        bit_next = int(rmdr >= cand) # The next bit in the remainder.
        root = root << 1 | bit_next # Add next bit to running result.
        rmdr -= cand * bit_next # Reduce the remainder if bit was added.
    return root

def S(n):
	return (n + 1) * (n + 2) // 2 - 1

def P(n):
	# max k: (k + 1) * (k + 2) / 2 - 1 >= n
	# max k: (k + 1) * (k + 2) >= 2 * n + 2
	# max k: k^2 + 3k >= 2 * n
	# max k: 4k^2 + 12k >= 8 * n
	# max k: 4k^2 + 12k + 9 >= 8 * n + 9
	k = (isqrt(8 * n + 9) - 3) // 2
	while S(k) < n:
		k += 1
	return k

Q = [1] * N
prefQ = [[1] * N for i in range(D)]

for i in range(2, N):
	Q[i] = Q[i - 1] + Q[P(i)]
	cur = Q[i]
	for j in range(D):
		prefQ[j][i] = prefQ[j][i - 1] + cur
		cur *= i

denom = 107520
for v in pss:
	for i in range(len(v)):
		if isinstance(v[i], Fraction) and denom % v[i].denominator == 0:
			v[i] = int(v[i] * denom)
		else:
			v[i] *= denom

def f(n, p):
	global find_Q
	# returns sum for k = 1 .. n of Q(k) * p(k)
	if n < N:
		return sum(prefQ[i][n] * p[i] for i in range(len(p)))
	
	mx = P(n)
	fr = S(mx - 1) + 1

	p = poly_to_sum(p)
	q = [-x for x in subst(p, [-1, 1])]
	q[0] += apply_poly(p, n)
	q = poly_to_sum(q)

	# ans += sum_{fr <= i <= n} (sumpoly(p, n) - sumpoly(p, i - 1))
	ans = find_Q(mx) * (apply_poly(q, n) - apply_poly(q, fr - 1))

	# ans += sum_{j < mx} sum_{S(j - 1) < i <= S(j)} (sumpoly(p, n) - sumpoly(p, i - 1))
	mult = 2 ** (len(q) - 1)
	ans *= mult
	for i in range(len(q)):
		q[i] *= 2 ** (len(q) - 1 - i)
	# q = sub(subst(q, [0, Fraction(3, 2), Fraction(1, 2)]), subst(q, [-1, Fraction(1, 2), Fraction(1, 2)]))
	q = sub(subst(q, [0, 3, 1]), subst(q, [-2, 1, 1]))
	ans += f(mx - 1, q)

	return ans // (denom**2 * mult)


def find_Q(n):
	if n < N:
		return Q[n]

	mx = P(n)
	return f(mx, [1, 1]) - find_Q(mx) * (S(mx) - n)

if __name__ == "__main__":
	t = int(input())
	for _ in range(t):
		n = int(input())
		# n = 10**40
		print(find_Q(n))

Details

Tip: Click on the bar to expand more detailed information

Test #1:

score: 100
Accepted
time: 1214ms
memory: 224600kb

input:

4
10
100
1000
987654321123456789

output:

30
2522
244274
235139898689017607381017686096176798

result:

ok 4 lines

Test #2:

score: 0
Accepted
time: 1940ms
memory: 225004kb

input:

10000
613939207402503646
408978283345333197
976677512231308716
629053280913850466
148712339
236220313279945487
590396556995977994
9226
215693877607285701
649702683896705
343173887453826567
847003949499596615
867133040287550291
159928123569892354
864534948175618654
209739383170746721
4295456752378791...

output:

91030728117067063595428375419346402
40469246710473908695676971059074376
229951682342450470520349294486964970
95558039501640054006660579352994194
5418340556433412
13536357243714243974772906966693552
84197207203086091733385317631619504
20934656
11291075621624104092841708040651034
104019777231815308683...

result:

ok 10000 lines

Test #3:

score: 0
Accepted
time: 1538ms
memory: 224828kb

input:

1000
6403632579734162001008235137760132245297
1307698664787972023762442022139627469666
668870338048562416595095770565441759482
5092270
8806864498496723812760785099973409980711
2178464116010775202899984038946879469187
204381824371
8638495456004772442511643693521120926431
45954412333082528168092594892...

output:

9822905445826021159291522774438593145331066315784567505849706373529921001845336
409552844852728078841401625660602682494169687360338453221088647649526339235330
107160056181509722327918304399871120167096186888511567354513496578559803370274
6354453295964
185817323129718525790559571287806776421589155460...

result:

ok 1000 lines

Test #4:

score: 0
Accepted
time: 1848ms
memory: 224804kb

input:

2000
2587627816607030340103003175959756184662
7728753705064569253253827797978613582938
6507847628603052973674714721
6989857636824717968431061258300652290539
4734281027640913533293237760425416005062
9411123735455625690768098631226366597446
8309753930995253536640660321476246470149
63565157427098067709...

output:

1603610451790269237852635641930301658193367441164307312552842461667027137613454
14309943493171496191506053530749878345271155973702143153225815632926701434642842
10414697803791950572309383355053080338229465674803757518
11704102206894264239190418551798088411458157423624785417336589284698838535371266
5...

result:

ok 2000 lines

Test #5:

score: 0
Accepted
time: 2940ms
memory: 225116kb

input:

5000
6701025283447923273597553918313900029215
1618190467906965189844764312914748628527
2135261797018456059451326428589353332126
8027429917075086154217136768450383650445
8263301530632040969183919589286944799002
376886964626613356031779878
1191561726595055348898524031899625958102
453561433135467095374...

output:

10756647971303093856509780939435968749671842310025383400207261624750784751725876
627115945498078452730193858129037650151913122482133071938783012929533045529940
1091921493223233296724616654299408127388059493196845968361252389122720100379070
154375707400033063668088306493199269136326791073281723548116...

result:

ok 5000 lines

Test #6:

score: 0
Accepted
time: 4588ms
memory: 225032kb

input:

10000
260865660295317278841
5232352637620496342310336202478387251106
7108789244285764135987032973912918380
12766535319519586095540974948550152061
5138306300154916887884637635208203681949
7603163140266839847784708214115317398585
149590294591155205497765668766786424787
63283557694500045989299147454323...

output:

16323111740957704392106241109891718054228
6557703685144914472554701877112177422100848067214049191882271294010013621817762
12143115079716078114619105501427985631361994195400195527663921137836417758
39139456824156526604158618001888125076786817219954316014947704612553450312
6324051018379978443719363340...

result:

ok 10000 lines