Problem:
 Score : 500 points    Problem Statement Let us denote by f(x, m) the remainder of the Euclidean division of x by m .  Let A be the sequence that is defined by the initial value A_1=X and the recurrence relation A_{n+1} = f(A_n^2, M) .
Find \displaystyle{\sum_{i=1}^N A_i} .      Constraints  1 \leq N \leq 10^{10}  0 \leq X < M \leq 10^5  All values in input are integers.         Input Input is given from Standard Input in the following format:  N  X  M       Output Print \displaystyle{\sum_{i=1}^N A_i} .        Sample Input 1 6 2 1001      Sample Output 1 1369  The sequence A begins 2,4,16,256,471,620,\ldots Therefore, the answer is 2+4+16+256+471+620=1369 .       Sample Input 2 1000 2 16      Sample Output 2 6  The sequence A begins 2,4,0,0,\ldots Therefore, the answer is 6 .       Sample Input 3 10000000000 10 99959      Sample Output 3 492443256176507

Solution (Python):
N, X, M = map(int, input().split())
flag = [False for _ in range(M)]
record = list()
record.append(X)
flag[X] = 1

An = X

for i in range(M + 1):
    An = pow(An, 2, M)
    if flag[An]:
        start = flag[An]
        cnt = i + 2 - start
        cost = record[-1] - record[start - 2] if start > 1 else record[-1]
        break

    else:
        record.append(An + record[-1])
        flag[An] = i + 2

if start >= N:
    print(record[N - 1])
else:
    print(((N - start) // cnt) * cost
                + record[(N - start) % cnt + start - 1])

Solution (Python):
N, X, M = map(int, input().split())
flag = [False for _ in range(M)]
record = list()
record.append(X)
flag[X] = 1

An = X

for i in range(M + 1):
    An = pow(An, 2, M)
    if flag[An]:
        start = flag[An]
        cnt = i + 2 - start
        cost = record[-1] - record[start - 2] if start > 1 else record[-1]
        break

    else:
        record.append(An + record[-1])
        flag[An] = i + 2

if start >= N:
    print(record[N - 1])
else:
    print(((N - start) // cnt) * cost
                + record[(N - start) % cnt + start - 1])

Solution (Python):
N,X,M = map(int, input().split())
A = [X]
S = set(A)
cnt = 0
i = 0
for i in range(1, N):
    A.append(A[i-1] ** 2 % M)
    if A[i] in S:
        break
    else:
        S.add(A[i])
roup_cycle = i - A.index(A[i])
before_roup_sums = sum(A[:A.index(A[i])])
roup_sums = sum(A[A.index(A[i]):i])

if roup_cycle != 0:
    roup_amount = (N - A.index(A[i])) // roup_cycle
    roup_mod = (N - A.index(A[i])) % roup_cycle
    print(before_roup_sums + roup_sums * roup_amount + sum(A[A.index(A[i]):(A.index(A[i]) + roup_mod)]))
else:
    print(sum(A))

Solution (Python):

import sys


def floyd_cycle_finding(f, x0):
    '''
    循環していない部分の長さをlam
    循環部分の長さをmuとして
    (lam, mu)
    を返す
    >>> floyd_cycle_finding(lambda x: x**2 % 3, 2)
    (1, 1)
    >>> floyd_cycle_finding(lambda x: x**2 % 1001, 2)
    (2, 4)
    >>> floyd_cycle_finding(lambda x: x**2 % 16, 2)
    (2, 1)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 3)
    (0, 4)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 0)
    (0, 1)
    '''
    tortoise = f(x0)
    hare = f(f(x0))
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(f(hare))
    lam = 0
    tortoise = x0
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(hare)
        lam += 1
    mu = 1
    hare = f(tortoise)
    while tortoise != hare:
        hare = f(hare)
        mu += 1
    return lam, mu

def solve(N: int, X: int, M: int):
    f = lambda x: x**2 % M

    def g(x):
        while True:
            yield x
            x = f(x)
            
    lam, mu = floyd_cycle_finding(f, X)
    ans = 0
    gg = g(X)
    k = max((N - lam) // mu, 0)
    for i in range(min(N, lam)):
        ans += next(gg)
    for i in range(mu):
        ans += k * next(gg)
    for i in range(max((N-lam)%mu, 0)):
        ans += next(gg)
    return ans


def main():
    def iterate_tokens():
        for line in sys.stdin:
            for word in line.split():
                yield word
    tokens = iterate_tokens()
    N = int(next(tokens))  
    X = int(next(tokens))  
    M = int(next(tokens))  
    print(solve(N, X, M))


def test():
    import doctest
    doctest.testmod()

if __name__ == '__main__':
    
    main()

Solution (Python):

import sys


def floyd_cycle_finding(f, x0):
    '''
    循環していない部分の長さをlam
    循環部分の長さをmuとして
    (lam, mu)
    を返す
    >>> floyd_cycle_finding(lambda x: x**2 % 3, 2)
    (1, 1)
    >>> floyd_cycle_finding(lambda x: x**2 % 1001, 2)
    (2, 4)
    >>> floyd_cycle_finding(lambda x: x**2 % 16, 2)
    (2, 1)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 3)
    (0, 4)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 0)
    (0, 1)
    '''
    
    
    
    
    
    
    
    
    
    
    
    
    
    tortoise = f(x0)
    hare = f(f(x0))
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(f(hare))
    lam = 0
    tortoise = x0
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(hare)
        lam += 1
    mu = 1
    hare = f(tortoise)
    while tortoise != hare:
        hare = f(hare)
        mu += 1
    return lam, mu

def solve(N: int, X: int, M: int):
    f = lambda x: x**2 % M

    def g(x):
        while True:
            yield x
            x = f(x)
            
    lam, mu = floyd_cycle_finding(f, X)
    ans = 0
    gg = g(X)
    k = max((N - lam) // mu, 0)
    for i in range(min(N, lam)):
        ans += next(gg)
    for i in range(mu):
        ans += k * next(gg)
    for i in range(max((N-lam)%mu, 0)):
        ans += next(gg)
    return ans


def main():
    def iterate_tokens():
        for line in sys.stdin:
            for word in line.split():
                yield word
    tokens = iterate_tokens()
    N = int(next(tokens))  
    X = int(next(tokens))  
    M = int(next(tokens))  
    print(solve(N, X, M))


def test():
    import doctest
    doctest.testmod()

if __name__ == '__main__':
    test()
    main()

Solution (Python):

import sys


def floyd_cycle_finding(f, x0):
    '''
    循環していない部分の長さをlam
    循環部分の長さをmuとして
    (lam, mu)
    を返す
    >>> floyd_cycle_finding(lambda x: x**2 % 3, 2)
    (1, 1)
    >>> floyd_cycle_finding(lambda x: x**2 % 1001, 2)
    (2, 4)
    >>> floyd_cycle_finding(lambda x: x**2 % 16, 2)
    (2, 1)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 3)
    (0, 4)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 0)
    (0, 1)
    '''
    
    
    
    
    
    
    
    
    
    
    
    
    
    tortoise = f(x0)
    hare = f(f(x0))
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(f(hare))
    lam = 0
    tortoise = x0
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(hare)
        lam += 1
    mu = 1
    hare = f(tortoise)
    while tortoise != hare:
        hare = f(hare)
        mu += 1
    return lam, mu

def solve(N: int, X: int, M: int):
    f = lambda x: x**2 % M

    def g(x):
        while True:
            yield x
            x = f(x)
            
    lam, mu = floyd_cycle_finding(f, X)
    ans = 0
    gg = g(X)
    k = max((N - lam) // mu, 0)
    for i in range(min(N, lam)):
        ans += next(gg)
    for i in range(mu):
        ans += k * next(gg)
    for i in range((N-lam)%mu):
        ans += next(gg)
    return ans


def main():
    def iterate_tokens():
        for line in sys.stdin:
            for word in line.split():
                yield word
    tokens = iterate_tokens()
    N = int(next(tokens))  
    X = int(next(tokens))  
    M = int(next(tokens))  
    print(solve(N, X, M))


def test():
    import doctest
    doctest.testmod()

if __name__ == '__main__':
    
    main()

Solution (Python):

import sys


def floyd_cycle_finding(f, x0):
    '''
    循環していない部分の長さをlam
    循環部分の長さをmuとして
    (lam, mu)
    を返す
    >>> floyd_cycle_finding(lambda x: x**2 % 3, 2)
    (1, 1)
    >>> floyd_cycle_finding(lambda x: x**2 % 1001, 2)
    (2, 4)
    >>> floyd_cycle_finding(lambda x: x**2 % 16, 2)
    (2, 1)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 3)
    (0, 4)
    >>> floyd_cycle_finding(lambda x: x*2 % 5, 0)
    (0, 1)
    '''
    
    
    
    
    
    
    
    
    
    
    
    
    
    tortoise = f(x0)
    hare = f(f(x0))
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(f(hare))
    lam = 0
    tortoise = x0
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(hare)
        lam += 1
    mu = 1
    hare = f(tortoise)
    while tortoise != hare:
        hare = f(hare)
        mu += 1
    return lam, mu

def solve(N: int, X: int, M: int):
    f = lambda x: x**2 % M

    def g(x):
        while True:
            yield x
            x = f(x)
            
    lam, mu = floyd_cycle_finding(f, X)
    ans = 0
    gg = g(X)
    k = max((N - lam) // mu, 0)
    for i in range(min(N, lam)):
        ans += next(gg)
    for i in range(mu):
        ans += k * next(gg)
    for i in range((N-lam)%mu):
        ans += next(gg)
    return ans


def main():
    def iterate_tokens():
        for line in sys.stdin:
            for word in line.split():
                yield word
    tokens = iterate_tokens()
    N = int(next(tokens))  
    X = int(next(tokens))  
    M = int(next(tokens))  
    print(solve(N, X, M))


def test():
    import doctest
    doctest.testmod()

if __name__ == '__main__':
    test()
    main()

Solution (Python):
import sys

input = sys.stdin.readline

N, X, M = map(int, input().split())
LIST = [X]
SET = {X}

a = X

while True:
    if a * a % M in SET:
        break
    else:
        a = a * a % M
        LIST.append(a)
        SET.add(a)



start = LIST.index(a * a % M)
LEN = len(LIST)

if N <= LEN:
    print(sum(LIST[:N]))
else:
    LoopTimes, Excess = divmod(N-LEN, LEN-start)
    LoopSum = sum(LIST[start:]) * LoopTimes
    ExSum = sum(LIST[start:start+Excess])
    print(sum(LIST[:]) + LoopSum + ExSum)

Solution (Python):
import sys

n,x,m = (int(i) for i in input().split())
l = [x]
ans = x
memo_array = [-1]*(m+1)

for i in range(1,n):
    ans = ans**2%m
    idx = memo_array[ans]
    if idx != -1:
        loop_num, mod = divmod((n-idx),(i-idx))
        print(sum(l[:idx]) + sum(l[idx:])*(loop_num)+sum(l[idx:idx+mod]))
        sys.exit()
    l.append(ans)
    memo_array[ans]=i
print(sum(l))

Solution (Python):
n, x, m = map(int, input().split())

tmp = [x]
while x**2%m not in tmp:
    x = x**2%m
    tmp.append(x)

cnt = len(tmp)
x = x**2%m

if n <= cnt:
    print(sum(tmp[:n]))
    exit()

for i in range(len(tmp)):
    if tmp[i] == x:
        st = sum(tmp[:i])
        rp = tmp[i:]
        div, mod = divmod(n-i, len(rp))
        en = sum(rp[:mod])
        print(st+sum(rp)*div+en)
        exit()

Solution (Python):
N, X, M = map(int, input().split())

G = {X: 0}
a = X
while True:
    na = pow(a, 2, M)
    if na in G:
        s = G[na]
        sa = na
        g = G[a]
        break
    else:
        G[na] = G[a] + 1
        a = na

if N <= s:
    ans1 = 0
    a = X
    for i in range(N):
        ans1 += a
        a = pow(a, 2, M)
    print(ans1)
else:
    ans1 = 0
    a = X
    for i in range(s):
        ans1 += a
        a = pow(a, 2, M)

    L = g - s + 1
    S = [0] * (L + 1)
    a = sa
    for i in range(L):
        S[i + 1] = S[i] + a
        a = pow(a, 2, M)

    rest = N - s
    x, y = divmod(rest, L)
    ans2 = S[-1] * x + S[y]

    print(ans1 + ans2)

Solution (Python):
import time
startTimeProblem=time.time()

import fileinput, sys, itertools, functools, copy, math
from bisect import *
from heapq import *
from collections import *


def lcm(a, b):  
    return (a*b)/gcd(a, b)

class InputHelper:
    def __init__(self):
        self.myinput = fileinput.input()

    def isLocal(self):
        return not fileinput.isstdin()

    def int(self):
        return int(self.myinput.readline().rstrip())

    def ints(self):
        return [int(_) for _ in self.myinput.readline().rstrip().split()]

    def str(self):
        return self.myinput.readline().rstrip()

    def strs(self):
        return [_ for _ in self.myinput.readline().rstrip().split()]

class OutputHelper:
    def int(self, a):
        print(a)    

    def ints(self, a):  
        print(" ".join([str(_) for _ in a]))
    
    def intsNL(self, a):
        for _ in a:
            print(_)
    
    def str(self, s):
        print(s)

    def strs(self, s):
        print(" ".join([_ for _ in s]))

    def strsNL(self, s):
        for st in s:
            print(st)

class ListNode:
    def __init__(self, val):
        self.val = val
        self.next = None
        self.prev = None

In = InputHelper()
Out = OutputHelper()



N, X, M = In.ints()

ans = 0
steps = 1
seen = dict()
stepsseen = dict()
stepsseen[0] = 0

while X != 0 and N>0 and X not in seen:    
    ans += X
    seen[X] = (ans, steps)
    stepsseen[steps] = ans
    X = (X**2)%M
    N-=1
    steps+=1


if X not in seen:
    Out.int(ans)
else:
    stepsDiff = steps-seen[X][1]
    ansDiff = ans-stepsseen[seen[X][1]-1]

    remloops = N//stepsDiff
    N = N%stepsDiff

    ans += remloops*ansDiff

    while X != 0 and N>0:    
        ans += X
        X = (X**2)%M
        N-=1
    
    Out.int(ans)




if len(sys.argv)>2 and sys.argv[2]=="TIMEIT":
    fin = (time.time()-startTimeProblem)*1000
    print("{:.2f}".format(fin) + "ms")

Solution (Python):
N, X, M = map(int, input().split())

existence = [False] * M
a = []
A = X
for i in range(N):
    if existence[A]:
        break
    existence[A] = True
    a.append(A)
    A = A * A % M

try:
    loop_start = a.index(A)
except:
    loop_start = len(a)

result = sum(a[:loop_start])
N -= loop_start
if N != 0:
    a = a[loop_start:]
    loops = N // len(a)
    remainder = N % len(a)
    result += sum(a) * loops + sum(a[:remainder])
print(result)

Solution (Python):
N, X, M = map(int, input().split())

existence = [False] * M
a = []
A = X
for i in range(N):
    if existence[A]:
        break
    existence[A] = True
    a.append(A)
    A = A * A % M

for i in range(len(a)):
    if a[i] == A:
        loop_start = i
        break
else:
    loop_start = len(a)

result = sum(a[:loop_start])
N -= loop_start
if N != 0:
    a = a[loop_start:]
    loops = N // len(a)
    remainder = N % len(a)
    result += sum(a) * loops + sum(a[:remainder])
print(result)

Solution (Python):
N, X, M = map(int, input().split())

existence = [False] * M
a = []
A = X
for i in range(N):
    if existence[A]:
        break
    existence[A] = True
    a.append(A)
    A = A * A % M

for i in range(len(a)):
    if a[i] == A:
        loop_start = i
        break
else:
    loop_start = len(a)

result = sum(a[:loop_start])
a = a[loop_start:]
N -= loop_start
if N != 0:
    loops = N // len(a)
    remainder = N % len(a)
    result += sum(a) * loops + sum(a[:remainder])
print(result)

Solution (Python):
n,x,m = map(int,input().split())

place = [0] * m
place[x] = 1
cur = [x]
ends = 0

for i in range(m):
    x = (x**2) % m
    if place[x] == 0:
        place[x] = i+2
        cur.append(x)
    else:
        ends = i+2
        break
        
roopsum = sum(cur[place[x]-1:ends-1])
rooplen = ends - place[x]
curlen = cur[place[x]-1:]

mm = min(n,ends-1)
ans = 0

for i in range(mm):
    ans += cur[i]
    
if n-mm == 0:
    print(ans)
else:
    roopamo = (n-mm) // rooplen
    roopama = (n-mm) % rooplen
    
    ans += roopamo * roopsum
    for i in range(roopama):
        ans += curlen[i]
        
    print(ans)

Solution (Python):
N,X,M=map(int,input().split())
modlist=[0]*M
modset=set()
ans=0
for i in range(M):
    if X in modset:
        startloop=modlist.index(X)
        sum_outofloop=sum(modlist[:startloop])
        inloop=ans-sum_outofloop
        numofloop=(N-startloop)//(i-startloop)
        ans=sum_outofloop+numofloop*inloop+sum(modlist[startloop:startloop+((N-startloop)%(i-startloop))])
        
        break
    else:
        modset.add(X)
        ans+=X
        modlist[i]=X
        X**=2
        X%=M
    if i==N-1:
        break
print(ans)

Solution (Python):
n, x, m = map(int, input().split())
arr = [x]
mod = x
memo = [-1 for _ in range(m+1)]

ans = 0
for i in range(1, n):
    mod = mod ** 2 % m
    idx = memo[mod]

    if idx != -1:
        a = i - idx
        q, _mod = divmod((n - idx), a)

        ans = sum(arr[0:idx]) + q * sum(arr[idx:idx+a]) + sum(arr[idx:idx+_mod])
        break
    arr.append(mod)
    memo[mod] = i

if ans == 0:
    ans = sum(arr)

print(ans)

Solution (Python):
def readinput():
    n,x,m=map(int,input().split())
    return n,x,m

def main(n,x,m):
    ans_list=[]
    hist=[0]*m
    val=x
    count=0
    while hist[val]==0 and count<n:
        ans_list.append(val)
        hist[val]+=1
        val=(val*val)%m
    if count==n:
        loop_start=0
        loop_len=n
    else:
        i=0
        while ans_list[i] != val:
            i+=1
        loop_start=i
        loop_len=len(ans_list)-i 
    
    
    
    ans=sum(ans_list[:loop_start])
    ans+=sum(ans_list[loop_start:])*((n-loop_start)//loop_len)
    ans+=sum(ans_list[loop_start:loop_start+(n-loop_start)%loop_len])
    
    return ans

if __name__=='__main__':
    n,x,m=readinput()
    ans=main(n,x,m)
    print(ans)

Solution (Python):
n,x,m=map(int,input().split())
l,s=[0]*m,[0]*m
t=p=0
while l[x]<1:
  t+=1
  l[x]=t
  s[x]=s[p]+x
  p=x
  x=p*p%m
k=l[x]
d,m=divmod(n-k,t+1-k)
print((s[p]+x-s[x])*d+s[l.index(k+m)])

Solution (Python):
import bisect,collections,copy,itertools,math,string
import sys
def I(): return int(sys.stdin.readline().rstrip())
def LI(): return list(map(int,sys.stdin.readline().rstrip().split()))
def S(): return sys.stdin.readline().rstrip()
def LS(): return list(sys.stdin.readline().rstrip().split())
def main():

    n, x, m = LI()
    start = -1
    goal = -1
    a = [0 for _ in range(m+2)]
    acc = [0 for _ in range(m+2)]
    st = {x}
    a[0] = -1
    a[1], acc[1] = x, x
    for i in range(1,m+1):
        a[i+1] = a[i]**2%m
        acc[i+1] = acc[i] + a[i+1]
        if a[i+1] in st:
            start = a.index(a[i+1])
            goal = i
            break
        else:
            st.add(a[i+1])

    ans = 0
    if n<=goal:
        ans = acc[n]
    else:
        ans = acc[goal] + (acc[goal]-acc[start-1])*((n-goal)//(goal-start+1)) + acc[(n-goal)%(goal-start+1)+start-1] - acc[start-1]

    print(ans)
main()

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

N,X,M = map(int,readline().split())
nex = [0] * M 
sums = [0] * M 

for i in range(M):
  nex[i] = i * i % M
  sums[i] = i

ans = 0
while N:
  if N & 1:
    ans += sums[X]
    X = nex[X]
  sums = [sums[i] + sums[nex[i]] for i in range(M)]
  nex = [nex[nex[i]] for i in range(M)]
  N >>= 1
  
print(ans)

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

N,X,M = map(int,readline().split())
nex = [0] * M
sums = [0] * M

for i in range(M):
  nex[i] = i * i % M
  sums[i] = i

ans = 0
while N:
  if N & 1:
    ans += sums[X]
    X = nex[X]
  sums = [sums[i] + sums[nex[i]] for i in range(M)]
  nex = [nex[nex[i]] for i in range(M)]
  N >>= 1
  
print(ans)

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

N,X,M = map(int,readline().split())
nex = [0] * M
sums = [0] * M

for i in range(M):
  nex[i] = i * i % M
  sums[i] = i

ans = 0
while N:
  if N & 1:
    ans += sums[X]
    X = nex[X]
  sums = [sums[i] + sums[nex[i]] for i in range(M)]
  nex = [nex[nex[i]] for i in range(M)]
  N >>= 1
  
print(ans)

Solution (Python):
N , X , M = map(int,input().split())
a = (X*X)%M
lis = [X]
while a not in lis:
    lis.append(a)
    a = a*a%M
i = lis.index(a)
lis1 = lis[:i]
lis2 = lis[i:]
ans=0
if N<=len(lis1):
    ans = sum(lis1[:N])
else:
    ans = sum(lis1)
    n = N - len(lis1)
    ans += (n//len(lis2))*sum(lis2)
    ans += sum(lis2[:n%len(lis2)])


print(ans)

Solution (Python):
def f(a):
    return a * a % m


n, x, m = map(int, input().split())

a = x
i = 1
li = [a]
num2index = {a: 1}
while f(a) not in num2index and i + 1 <= n:
    i += 1
    a = f(a)
    li.append(a)
    num2index[a] = i

if i == n:
    ans = sum(li)
else:
    before_loop = num2index[f(a)]
    loop, left = divmod(n - i, i - before_loop + 1)
    loop_sm = loop * sum(li[before_loop-1:])
    left_sm = sum(li[before_loop-1:before_loop-1+left])
    ans = sum(li) + loop_sm + left_sm

print(ans)

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

N,X,M = map(int,readline().split())
MAXD = 55






nex = [[-1] * M for i in range(MAXD + 1)]
sums = [[0] * M for i in range(MAXD + 1)]


for r in range(M):
  nex[0][r] = r * r % M
  sums[0][r] = r

for p in range(MAXD):
  for r in range(M):
    nex[p + 1][r] = nex[p][nex[p][r]]
    sums[p + 1][r] = sums[p][r] + sums[p][nex[p][r]]
    
    
    
    
    
res = 0
cur = X
for p in range(MAXD, -1, -1):
  if N & (1 << p):
    res += sums[p][cur]
    cur = nex[p][cur]
    
print(res)

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

def I(): return int(input())
def MI(): return map(int, input().split())
def LI(): return list(map(int, input().split()))

def main():
    mod=10**9+7
    N,X,M=MI()
    use=[0]*M
    loopF=0
    loopN=0
    loopS=0
    loopL=[]
    ans=0
    cnt=0
    
    now=X
    
    if N<=4*M:
        for i in range(N):
            ans+=now
            now=(now**2)%M
        print(ans)
        exit()
            
    
    for i in range(M*3):
        use[now]+=1
        if use[now]==2:
            loopF=1
        if use[now]==3:
            cnt=i
            break
        
        ans+=now
        if loopF:
            loopL.append(now)
            loopN+=1
            loopS+=now
        
        now=(now**2)%M
        
    N2=N-cnt
    loop=N2//loopN
    rem=N2%loopN
    
    ans+=loop*loopS
    ans+=sum(loopL[:rem])
    
    print(ans)
    

main()

Solution (Python):
N,X,M = map(int,input().split())
S = X
A = [X]
exists = [0]*M
exists[X] = 1

for i in range(2,N+1):
  A_i = (A[-1]*A[-1])%M
  A.append(A_i)
  S += A_i
  if exists[A_i] == 0:
    exists[A_i] = 1
  else:
    S -= A_i
    del A[-1]
    B = A[A.index(A_i):]
    S += ((N-len(A))//len(B))*sum(B) + sum(B[:((N-len(A))%len(B))])
    break

print(S)

Solution (Python):
import sys
def input(): return sys.stdin.readline().strip()
def mapint(): return map(int, input().split())
sys.setrecursionlimit(10**9)

N, X, mod = mapint()
doubling = [[i*i%mod for i in range(mod)]]
accu = [doubling[0][:]]

for i in range(40):
    tmp1 = [None]*mod
    tmp2 = [None]*mod
    for j in range(mod):
        tmp1[j] = doubling[-1][doubling[-1][j]]
        tmp2[j] = accu[-1][doubling[-1][j]] + accu[-1][j]
    doubling.append(tmp1)
    accu.append(tmp2)
now = X
ans = X
for i in range(40):
    if (N-1)>>i&1:
        ans += accu[i][now]
        now = doubling[i][now]
print(ans)

Solution (Python):
n, x, m = [int(x) for x in input().split()]
val = [0, x]
d = {x}
prev = x
ans = x
for i in range(2, n+1):
    cur = prev*prev%m
    if cur==0:
        break
    if cur in d:
        pos = val.index(cur)
        length = len(val)-pos
        left = n-i+1
        temp = 0
        for j in range(pos, len(val)):
            temp += val[j]
        ans += (left//length)*temp
        left -= (left//length)*length
        while left:
            ans += val[pos]
            pos += 1
            left-= 1
        break
    else:
        ans += cur
        prev = cur
        d.add(cur), val.append(cur)

print(ans)

Solution (Python):
n, x, m = [int(_) for _ in input().split()]
d = {}
d_rev = {}
d[x] = 0
d_rev[0] = x
c = 0
while True:
    c += 1
    x = (x ** 2) % m
    if x in d:
        roop_idx = d[x]
        break
    else:
        d[x] = c
        d_rev[c] = x

ans = 0
if roop_idx >= n:
    for i in range(n):
        ans += d_rev[i]
    print(ans)
else:
    for i in range(roop_idx):
        ans += d_rev[i]
    roop_num = len(d) - roop_idx
    roop_sum = 0
    for i in range(roop_idx, len(d)):
        roop_sum += d_rev[i]
    ans += ((n-roop_idx)//roop_num) * roop_sum
    roop_last = (n-roop_idx)%roop_num
    for i in range(roop_last):
        ans += d_rev[roop_idx + i]
    print(ans)

Solution (Python):
n,x,m = map(int, input().split())
b = x%m
add = 0


ans = 0
if n <= m:
    for i in range(n):
        ans += b    
        b = pow(b,2,m)
    print(ans)
    exit()


cnt = 0
while 2 ** cnt < m:
    cnt += 1


s = set()
tmp = 0
for i in range(m):
    if i < cnt:
        add += b
        tmp = b
    else:
        s.add(b)
    b = pow(b, 2, m)


rest = (n-cnt)%len(s)
div = (n-cnt)//len(s)

rep = []
for i in range(len(s)):
    tmp = pow(tmp, 2, m)
    rep.append(tmp)

ans = add + sum(rep)*div + sum(rep[:rest])
print(ans)

Solution (Python):
n, x, m = [int(_) for _ in input().split()]
d = {}
d_rev = {}
d[x] = 0
d_rev[0] = x
c = 0
while True:
    c += 1
    x = (x ** 2) % m
    if x in d:
        roop_idx = d[x]
        break
    else:
        d[x] = c
        d_rev[c] = x



ans = 0
if roop_idx >= n:
    for i in range(n):
        ans += d_rev[i]
    print(ans)
else:
    for i in range(roop_idx):
        ans += d_rev[i]
    roop_num = len(d) - roop_idx
    roop_sum = 0
    for i in range(roop_idx, len(d)):
        roop_sum += d_rev[i]
    ans += ((n-roop_idx)//roop_num) * roop_sum
    roop_last = (n-roop_idx)%roop_num
    for i in range(roop_last):
        ans += d_rev[roop_idx + i]
    print(ans)

Solution (Python):
import sys
def input(): return sys.stdin.readline().strip()
def mapint(): return map(int, input().split())
sys.setrecursionlimit(10**9)

N, X, mod = mapint()
first = set()
second = set()
in_loop = []
ans = 0
for i in range(min(10**5*2, N)):
    if X in first:
        if X in second:
            break
        in_loop.append(X)
        second.add(X)
    first.add(X)
    ans += X
    X = pow(X, 2, mod)
else:
    print(ans)
    exit()

loop_len = len(in_loop)
rest = N-len(first)-len(in_loop)
ans += (rest//loop_len)*sum(in_loop)
amari = rest%loop_len
for i in range(amari):
    ans += in_loop[i]
print(ans)

Solution (Python):
n,x,m = map(int,input().split())
mod = m
ans = 0
bit = [-1 for i in range(m)]
cycle = False
for i in range(n):
    if i == 0 :
        a = x
        bit[a] = i
        ans += a
    else:
        a = (a**2)% mod
        if bit[a] != -1:
            cy_st = bit[a]
            cy_fi = i -1
            cycle = True
            break
        else:
            bit[a] = i
            ans += a

if cycle:
    ans2 = 0
    b =  -1
    for j in range(cy_st):
        if j == 0 :
            b = x
            ans2 += b
        else:
            b = (b**2)% mod
            ans2 += b
    cy_num = ans - ans2
    cy_repe = (n-cy_st) // (cy_fi - cy_st + 1)
    ans3 = cy_num * cy_repe
    cy_amari = (n-cy_st) % (cy_fi - cy_st + 1)
    
    if b == -1:
        for j in range(cy_amari):
            if j == 0 :
                b = x
                ans3 += b
            else:
                b = (b**2)% mod
                ans3 += b
        
    else:
        for i in range(cy_amari):
            b = (b**2)% mod
            ans3 += b

    print(ans2+ans3)
else:
    print(ans)

Solution (Python):
n, x, m = map(int, input().split())

a = x
li = [a]
exists = dict()
i = 1
for i in range(2, n + 1):
    a = a * a % m
    if a in exists:
        break
    li.append(a)
    exists[a] = i

ans = sum(li)
if i != n:
    loop, left = divmod(n - i + 1, i - exists[a])
    ans += loop * sum(li[exists[a]-1:]) + sum(li[exists[a]-1:exists[a]-1+left])
print(ans)

Solution (Python):
n, x, m = map(int, input().split())
li = [x]
se = {x}
for i in range(n-1):
    x = x*x%m
    if x in se:
        idx = li.index(x)
        break
    li.append(x)
    se.add(x)
else:
    print(sum(li))
    exit(0)
ans = sum(li) + sum(li[idx:])*((n-len(li)) // (len(li)-idx)) + sum(li[idx:idx+(n-len(li)) % (len(li)-idx)])
print(ans)

Solution (Python):
n, x, m = map(int, input().split())
g = x * 1
ans = x
arr = []
res = set([])
loop = (0, 0)
k = 0
for i in range(m + 1):
	tmp = x ** 2 % m
	if tmp in res:
		loop = (i, tmp)
		break
	res.add(tmp)
	arr.append(tmp)
	ans += tmp
	x = tmp

for i, y in enumerate(arr):
	if y == loop[1]:
		k = i
		break

ini = g + sum(arr[:k])
mul = ans - ini
t, v = divmod(n - k - 1, loop[0] - k)
print(ini + mul * t + sum(arr[k:k + v]))

Solution (Python):
n, x, m = map(int, input().split())
li = [x]
seen = [False]*m
seen[x] = True
for i in range(n-1):
    x = x*x%m
    if seen[x]:
        idx = li.index(x)
        break
    seen[x] = True
    li.append(x)
else:
    print(sum(li))
    exit(0)
ans = sum(li) + sum(li[idx:])*((n-len(li)) // (len(li)-idx)) + sum(li[idx:idx+(n-len(li)) % (len(li)-idx)])
print(ans)

Solution (Python):
import math

N,X,M = map(int, input().split())
MAP = [-1]*(M+1)
TMAP = [0]*(M+1)
ans = 0

for i in range(N):
    if (X%M!=0):
        X = X%M
        ans += X
        if MAP[X]==-1:
            MAP[X] = i
            TMAP[X] = ans
        else:
            loopcnt = i-MAP[X]
            loopsum = ans-TMAP[X]
            ans += (((N-MAP[X])//loopcnt)-1)*loopsum
            jmax = N-(MAP[X]+((N-MAP[X])//loopcnt)*loopcnt)
            ans -= X
            for j in range(jmax):
                X = X%M
                ans += X
                X = X**2
            X=M
        X = X**2
    else:
        break

print(ans)

Solution (Python):
n, x, m = map(int,input().split())
a = x; s = a
k = {a}
memo = [a]
mode = 0
cnt = 1
for i in range(n-1):
	a = (a**2)%m
	cnt += 1
	memo.append(a)
	if a in k:
		mode = 1
		break
	s += a
	k.add(a)
if mode:
	targ = memo.index(a)
	nagasa = cnt-1-targ
	ruisekiwa = [0 for _ in range(nagasa)]
	ruisekiwa[0] = memo[targ]
	for i in range(1, nagasa):
		ruisekiwa[i] = ruisekiwa[i-1] + memo[i+targ]
	t = n - cnt
	s += ruisekiwa[-1]*(t//nagasa) + ruisekiwa[t%nagasa]
print(s)

Solution (Python):
n,x,m=map(int,input().split())
sm=[-1 for i in range(m)]+[0]
w=[-1 for i in range(m)]
d=[m]
a=x
t=0
s=0
while True:
  s+=a
  if w[a]!=-1:
    inter=t-w[a]
    fv=w[a]
    ls=d[(n-fv)%inter+fv]
    lls=d[w[a]]
    print(sm[lls]+(n-fv)//inter*(s-sm[a])+(sm[ls]-sm[lls]))
    break
  w[a]=t
  d.append(a)
  sm[a]=s
  t+=1
  a=(a*a)%m

Solution (Python):

from collections import defaultdict,deque
from heapq import heappush, heappop
from bisect import bisect_left, bisect_right
import sys, random, itertools, math
sys.setrecursionlimit(10**5)
input = sys.stdin.readline
sqrt = math.sqrt
def LI(): return list(map(int, input().split()))
def LF(): return list(map(float, input().split()))
def LI_(): return list(map(lambda x: int(x)-1, input().split()))
def II(): return int(input())
def IF(): return float(input())
def S(): return input().rstrip()
def LS(): return S().split()
def IR(n): return [II() for _ in range(n)]
def LIR(n): return [LI() for _ in range(n)]
def FR(n): return [IF() for _ in range(n)]
def LFR(n): return [LI() for _ in range(n)]
def LIR_(n): return [LI_() for _ in range(n)]
def SR(n): return [S() for _ in range(n)]
def LSR(n): return [LS() for _ in range(n)]
mod = 1000000007
inf = float("INF")


def solve():
    n, x, m = LI()
    if n <= 10 ** 5:
        ans = x
        for _ in range(n - 1):
            x = x ** 2 % m
            ans += x
        print(ans)
        return
    else:
        d = defaultdict(int)
        xd = x
        ans = 0
        for i in range(n + 1):
            if d[x]:
                break
            ans += x
            d[x] = (i + 1, ans)
            x = x ** 2 % m
        i += 1
        cycle = i - d[x][0]
        cycle_score = ans - d[x][1] + x
        n = n - d[x][0] + 1
        ans = d[x][1] - x
        diff = n % cycle
        tern = n // cycle
        ans += tern * cycle_score
        for j in range(diff):
            ans += x
            x = x ** 2 % m
        print(ans)
    return



if __name__ == '__main__':
    solve()