Problem:
 Score : 500 points    Problem Statement Takahashi has come to a party as a special guest.
There are N ordinary guests at the party. The i -th ordinary guest has a power of A_i .  Takahashi has decided to perform M  handshakes to increase the happiness of the party (let the current happiness be 0 ).
A handshake will be performed as follows:   Takahashi chooses one (ordinary) guest x for his left hand and another guest y for his right hand ( x and y can be the same).  Then, he shakes the left hand of Guest x and the right hand of Guest y simultaneously to increase the happiness by A_x+A_y .   However, Takahashi should not perform the same handshake more than once. Formally, the following condition must hold:   Assume that, in the k -th handshake, Takahashi shakes the left hand of Guest x_k and the right hand of Guest y_k . Then, there is no pair p, q  (1 \leq p < q \leq M) such that (x_p,y_p)=(x_q,y_q) .   What is the maximum possible happiness after M handshakes?      Constraints  1 \leq N \leq 10^5  1 \leq M \leq N^2  1 \leq A_i \leq 10^5  All values in input are integers.         Input Input is given from Standard Input in the following format:  N  M  A_1  A_2  ...  A_N       Output Print the maximum possible happiness after M handshakes.        Sample Input 1 5 3
10 14 19 34 33      Sample Output 1 202  Let us say that Takahashi performs the following handshakes:   In the first handshake, Takahashi shakes the left hand of Guest 4 and the right hand of Guest 4 .  In the second handshake, Takahashi shakes the left hand of Guest 4 and the right hand of Guest 5 .  In the third handshake, Takahashi shakes the left hand of Guest 5 and the right hand of Guest 4 .   Then, we will have the happiness of (34+34)+(34+33)+(33+34)=202 .  We cannot achieve the happiness of 203 or greater, so the answer is 202 .       Sample Input 2 9 14
1 3 5 110 24 21 34 5 3      Sample Output 2 1837       Sample Input 3 9 73
67597 52981 5828 66249 75177 64141 40773 79105 16076      Sample Output 3 8128170

Solution (Python):
def solve(m, a) -> int:
    from bisect import bisect_left
    from itertools import accumulate

    *a, = sorted(a)
    n = len(a)

    acc = (0,) + tuple(accumulate(a))

    def is_ok(h):
        ret = 0
        for x in a:
            j = bisect_left(a, h - x)
            ret += (n - j)
            if ret >= m:
                return True
        return False

    def binary_search():
        """M通り以上の握手ができる幸福度の上限"""
        ok = 0
        ng = a[-1] * 2 + 1
        while abs(ng - ok) > 1:
            mid = (ok + ng) // 2
            if is_ok(mid):
                ok = mid
            else:
                ng = mid
        return ok

    sup = binary_search()  

    ret = 0
    rest = m
    for x in a:
        j = bisect_left(a, sup + 1 - x)  
        
        use = n - j
        rest -= use

        ret += use * x + (acc[n] - acc[n - use])

    ret += rest * sup

    return ret


def main():
    n, m = map(int, input().split())
    *a, = map(int, input().split())
    print(solve(m, a))


if __name__ == '__main__':
    main()

Solution (Python):
def solve(m, a) -> int:
    from bisect import bisect_left
    from itertools import accumulate

    *a, = sorted(a)
    n = len(a)

    acc = (0,) + tuple(accumulate(a))

    def is_ok(h):
        ret = 0
        for x in a:
            j = bisect_left(a, h - x)
            ret += (n - j)
            if ret >= m:
                return True
        return False

    def binary_search():
        """M通り以上の握手ができる幸福度の上限"""
        ok = 0
        ng = a[-1] * 2 + 1
        while abs(ng - ok) > 1:
            mid = (ok + ng) // 2
            if is_ok(mid):
                ok = mid
            else:
                ng = mid
        return ok

    sup = binary_search()  

    ret = 0
    rest = m
    for x in a:
        j = bisect_left(a, sup - x)
        
        use = n - j
        rest -= use

        ret += use * x + (acc[n] - acc[n - use])

    ret += rest * sup

    return ret


def main():
    n, m = map(int, input().split())
    *a, = map(int, input().split())
    print(solve(m, a))


if __name__ == '__main__':
    main()

Solution (Python):
def ints():
  return [int(x) for x in input().split()]
def ii():
  return int(input())

N, M = ints()
A = ints()
A.sort()
A.reverse()

def combinations_and_kofukudo(x):
  c = 0
  k = 0
  si = 0
  i = 0
  for j in reversed(range(N)):
    while i<N and A[i]+A[j]>=x:
      si += A[i]
      i += 1
    c += i
    k += si + A[j]*i
    
  return (c, k)

def bsearch(l, u):
  m = (l+u)//2
  c, k = combinations_and_kofukudo(m)
  if c<M:
    return bsearch(l, m)
  else:
    if l==m:
      print(k-(c-M)*l)
      exit()
    return bsearch(m, u)

bsearch(0, A[0]*2+1)

Solution (Python):
import bisect


n, m = map(int, input().split())
a = list(map(int, input().split()))
a = sorted(a)

ruiseki = [0]*(n+1)
for i in range(n):
    ruiseki[i+1] = ruiseki[i] + a[i]

def solve(min_power):
    ans = 0
    for pos in range(n):
        tmp = min_power - a[pos]
        next_pos = max(bisect.bisect_left(a, tmp), pos)
        times = n - pos - (next_pos - pos)
        if times == n - pos:
            ans += 2 * (times - 1) + 1 
        else:
            ans += 2 * times
    if ans >= m:
        return True
    else:
        return False

ok = 0
ng = 2*10**5 + 1
while abs(ok - ng) > 1:
    mid =  (ok + ng) // 2 
    if solve(mid):
        ok = mid
    else:
        ng = mid

cnt = 0
ans = 0

for pos in range(n):
    tmp = ok - a[pos]
    next_pos = max(bisect.bisect_left(a, tmp), pos)
    times = n - pos - (next_pos - pos)
    if times == n - pos:
        cnt += 2 * (times - 1) + 1 
        ans += 2 * (ruiseki[n] - ruiseki[n-times+1] + (times-1)*a[pos]) + 2*a[pos]
    else:
        cnt += 2 * times
        ans += 2 * (ruiseki[n] - ruiseki[n-times] + times*a[pos])

print(ans - (cnt-m)*ok)

Solution (Python):
def ints():
  return [int(x) for x in input().split()]
def ii():
  return int(input())

N, M = ints()
A = ints()
A.sort()
A.reverse()

def combinations(x):
  s = 0
  i = 0
  for j in reversed(range(N)):
    while i<N and A[i]+A[j]>=x:
      i += 1
    s += i
    
  return s

def koufukudo(x):
  s = 0
  si = 0
  i = 0
  for j in reversed(range(N)):
    while i<N and A[i]+A[j]>=x:
      si += A[i]
      i += 1
    s += si + A[j]*i
    
  return s


def bsearch(l, u):
  m = (l+u)//2
  c = combinations(m)
  if c<M:
    return bsearch(l, m)
  else:
    if l==m:
      return (l, c-M)
    return bsearch(m, u)

x, dm = bsearch(0, A[0]*2+1)
print(koufukudo(x)-dm*x)

Solution (Python):
import sys
import numpy as np

read = sys.stdin.read
N, M, *A = map(int, read().split())

A.sort()
A_as = np.array(A, np.int64)
A_rev = np.array(A[::-1], np.int64)
a = A_rev.cumsum()
a_max = a[-1]

left = 0
right = A[-1] * 2 + 1
while right - left > 1:
    mid = (right + left) // 2
    n = N - np.searchsorted(A_as, mid - A_as)
    n[n < 0] = 0
    m = n.sum()
    if m >= M:
        left = mid
    else:
        right = mid

threshold = left

n = N - np.searchsorted(A_as, threshold - A_as, side='right')
n[n < 0] = 0
answer = a[n[n > 0] - 1].sum() + (A_as * n).sum()
if n.sum() != M:
    answer -= threshold * (n.sum() - M)

print(answer)

Solution (Python):
from bisect import bisect_left
from bisect import bisect_right

n,m = map(int, input().split())
a = list(map(int, input().split()))
a = sorted(a)

a_count = [0]*100010
a_sum = [0]*100010

for i in range(n):
    a_count[a[i]] += 1
    a_sum[a[i]] += a[i]

for i in range(100005,-1,-1):
    a_count[i] += a_count[i+1]
    a_sum[i] += a_sum[i+1]

def chk(x):
    case = 0
    for i in range(n):
        if 0<= x-a[i]<= 10**5+5:
            case += a_count[x-a[i]]
        elif x-a[i] < 0:
            case += n
    if case >= m:
        ch = True
    else:
        ch = False
    return ch

l = 2*a[-1]+1    
r = 1 

while abs(r-l) > 1:
    kk = (r+l)//2
    if chk(kk):
        r = kk
    else:
        l = kk

ans = 0
sm = sum(a)

l = r

case = 0
for i in range(n):
    num = n-bisect_right(a,l-a[i]-1)
    case += num
    ans += num*a[i]*2

ans -= (case-m)*l

print(ans)

Solution (Python):
def ints():
  return [int(x) for x in input().split()]
def ii():
  return int(input())

N, M = ints()
A = ints()
A.sort()

def combinations(x):
  s = 0
  i = 0
  j = N-1
  while j>=0:
    while i<N and A[i]+A[j]<x:
      i += 1
    s += N-i

    j -= 1
    
  return s

def koufukudo(x):
  s = 0
  si = 0
  j = 0
  i = N-1
  while j<N:
    while i>=0 and A[i]+A[j]>=x:
      si += A[i]
      i -= 1
    s += si + A[j]*(N-1-i)

    j += 1
    
  return s


def bsearch(lower, upper):
  l = lower
  u = upper
  m = (l+u)//2
  c = combinations(m)
  if c<M:
    return bsearch(l, m)
  else:
    if l==m:
      return (l, c-M)
    return bsearch(m, u)

x, dm = bsearch(0, A[-1]*2+1)
print(koufukudo(x)-dm*x)

Solution (Python):
import sys
import bisect
import itertools


def solve():
    readline = sys.stdin.buffer.readline
    mod = 10 ** 9 + 7
    n, m = list(map(int, readline().split()))
    a = list(map(int, readline().split()))
    a.sort()
    cor_v = a[0] * 2
    inc_v = a[-1] * 2 + 1
    while - cor_v + inc_v > 1:
        bin_v = (cor_v + inc_v) // 2
        cost = 0
        
        for i in range(n-1, -1, -1):
            cost += n - bisect.bisect_left(a, bin_v - a[i])
            if cost >= m:
                break
        
        if cost >= m:
            cor_v = bin_v
        else:
            inc_v = bin_v
    acm = list(itertools.accumulate([0] + a))
    t = 0
    c = 0
    for i in range(n):
        p = bisect.bisect_left(a, cor_v - a[i])
        t += acm[-1] - acm[p] + a[i] * (n - p)
        c += n - p
    print(t - (c - m) * cor_v)


if __name__ == '__main__':
    solve()

Solution (Python):
import sys
import bisect
import itertools


def solve():
    readline = sys.stdin.buffer.readline
    mod = 10 ** 9 + 7
    n, m = list(map(int, readline().split()))
    a = list(map(int, readline().split()))
    a.sort()
    cor_v = a[0] * 2
    inc_v = a[-1] * 2 + 1
    while - cor_v + inc_v > 1:
        bin_v = (cor_v + inc_v) // 2
        cost = 0
        
        for i in range(n-1, -1, -1):
            cost += n - bisect.bisect_left(a, bin_v - a[i])
            if cost >= m:
                break
        
        if cost >= m:
            cor_v = bin_v
        else:
            inc_v = bin_v
    acm = list(itertools.accumulate([0] + a))
    t = 0
    c = 0
    for i in range(n):
        p = bisect.bisect_left(a, cor_v - a[i])
        t += acm[-1] - acm[p] + a[i] * (n - p)
        c += n - p
    print(t - (c - m) * cor_v)


if __name__ == '__main__':
    solve()

Solution (Python):
import sys
import bisect
import itertools


def solve():
    readline = sys.stdin.buffer.readline
    mod = 10 ** 9 + 7
    n, m = list(map(int, readline().split()))
    a = list(map(int, readline().split()))
    a.sort()
    cor_v = a[0] * 2
    inc_v = a[-1] * 2 + 1
    while - cor_v + inc_v > 1:
        bin_v = (cor_v + inc_v) // 2
        cost = 0
        
        for i in range(n-1, -1, -1):
            cost += n - bisect.bisect_left(a, bin_v - a[i])
        
        if cost >= m:
            cor_v = bin_v
        else:
            inc_v = bin_v
    acm = list(itertools.accumulate([0] + a))
    t = 0
    c = 0
    for i in range(n):
        p = bisect.bisect_left(a, cor_v - a[i])
        t += acm[-1] - acm[p] + a[i] * (n - p)
        c += n - p
    print(t - (c - m) * cor_v)


if __name__ == '__main__':
    solve()

Solution (Python):
from bisect import bisect_left
from bisect import bisect_right

n,m = map(int, input().split())
a = list(map(int, input().split()))
a = sorted(a)

a_count = [0]*110000
a_sum = [0]*110000

for i in range(n):
    a_count[a[i]] += 1
    a_sum[a[i]] += a[i]

for i in range(100005,0,-1):
    a_count[i] += a_count[i+1]
    a_sum[i] += a_sum[i+1]

def chk(x):
    case = 0
    for i in range(n):
        case += n-bisect_right(a,x-a[i]-1)
    if case >= m:
        ch = True
    else:
        ch = False
    return ch


l = 2*a[-1]+1    
r = 1 

while abs(r-l) > 1:
    kk = (r+l)//2
    if chk(kk):
        r = kk
    else:
        l = kk

ans = 0
sm = sum(a)

l = r

case = 0
for i in range(n):
    num = n-bisect_right(a,l-a[i]-1)
    case += num
    ans += num*a[i]*2

ans -= (case-m)*l

print(ans)

Solution (Python):
import heapq
import bisect
n,m=map(int,input().split())
a=list(map(int,input().split()))
a.sort()
s=[0]*(n+1)
dp=[0]*(10**5+1)
h=0
for i in range(10**5+1):
    if h==n:
        dp[i]=h
        continue
    if i<=a[h]:
        dp[i]=h
    else:
        while h<n and i>a[h]:
            h+=1
        dp[i]=h



for i in range(n):
    s[i+1]=s[i]+a[i]
def test(x):
    count=0
    for i in range(n):
        buf=x-a[i]
        
        if buf<0:
            idx=0
        elif buf>=10**5+1:
            idx=n
        else:
            idx=dp[buf]
        if i>=idx:
            count+=(i-idx)*2+1
    return count


def search():
    mn=-1
    mx=m*max(a)*2
    buf=0
    while mx-mn>1:
        buf=(mx+mn)//2
        count=test(buf)
        if count>=m:
            mn=buf
        else:
            mx=buf
    return mn
x=search()
y=test(x)

count=0
for i in range(n):
    idx=bisect.bisect_left(a,x-a[i])
    if i>=idx:
        count+=(2*i-2*idx)*a[i]+(s[i+1]-s[idx])*2
        
print(count-(y-m)*x)

Solution (Python):
from bisect import bisect_left
n, m = map(int, input().split())
a = list(map(int, input().split()))
a.sort()


cum = [a[-1]]
for i in range(2, n+1):
    cum.append(cum[-1] + a[n-i])


def judge(val):
    cnt = 0
    sum = 0
    for ai in a:
        if ai + a[-1] < val:
            continue
        d = val - ai
        i = n - bisect_left(a, d)
        sum += cum[i-1] + ai * i
        cnt += i
    return cnt, sum

l = 0
r = 2 * a[-1]


while True:
    if r - l == 1:
        cnt, sum = judge(r)
        if cnt > m:
            print(m * r)
        else:
            print(sum + (m-cnt) * l)
        break
    mid = (l + r) // 2
    cnt, sum = judge(mid)
    if cnt < m:
        r = mid
    else:
        l = mid

Solution (Python):

  
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines




import numpy as np

N,M = map(int,readline().split())
A = np.array(read().split(),np.int64)

A.sort()

def shake_cnt(x):
    
    
    X = np.searchsorted(A,x-A)
    return N * N - X.sum()

left = 0 
right = 10 ** 6 
while left + 1 < right:
    x = (left + right) // 2
    if shake_cnt(x) >= M:
        left = x
    else:
        right = x

left,right



X = np.searchsorted(A,right-A) 
Acum = np.zeros(N+1,np.int64) 
Acum[1:] = np.cumsum(A)

shake = N * N - X.sum()
happy = (Acum[-1] - Acum[X]).sum() + (A * (N - X)).sum()

happy += (M - shake) *left
print(happy)

Solution (Python):
k = 34
K = 1<<k
nu = lambda L: int("".join([bin(K+a)[-k:] for a in L[::-1]]), 2)
st = lambda n: bin(n)[2:] + "0"
li = lambda s: [int(a, 2) if len(a) else 0 for a in [s[-(i+1)*k-1:-i*k-1] for i in range(len(B)*2-1)]]

N, M = map(int, input().split())
A = [int(a) for a in input().split()]
B = [0] * 100001
for a in A:
    B[a] += 1
C = li(st(nu(B) ** 2))
ans = 0
for i in range(200001)[::-1]:
    a = min(M, C[i])
    M -= a
    ans += a * i
print(ans)

Solution (Python):
from bisect import bisect_left, bisect_right
n, m = map(int, input().split())
a = list(map(int, input().split()))
a.sort()



a_cnt = [0] * (2 * a[-1] + 5)
for ai in a:
    a_cnt[ai] += 1
for i in reversed(range(len(a_cnt) - 1)):
    a_cnt[i] += a_cnt[i + 1]


a_cs = [0] * (n + 1)
for i in range(n):
    a_cs[i + 1] = a_cs[i] + a[n - 1 - i]


def check(x):
    
    
    cnt = 0
    for ai in a:
        if x >= ai:
            idx = a_cnt[x - ai]
        else:
            idx = n
        cnt += idx
    return cnt >= m


lb = 1  
ub = 2 * a[-1] + 3  
while ub - lb > 1:
    mid = (ub + lb) // 2
    if check(mid):
        lb = mid
    else:
        ub = mid



ans = 0
cnt = 0
for ai in a:
    if lb >= ai:
        idx = a_cnt[lb - ai]
    else:
        idx = n
    ans += ai * idx + a_cs[idx]
    cnt += idx
ans -= lb * (cnt - m)
print(ans)

Solution (Python):
import sys

sys.setrecursionlimit(10 ** 6)
int1 = lambda x: int(x) - 1
p2D = lambda x: print(*x, sep="\n")
def II(): return int(sys.stdin.readline())
def MI(): return map(int, sys.stdin.readline().split())
def LI(): return list(map(int, sys.stdin.readline().split()))
def LLI(rows_number): return [LI() for _ in range(rows_number)]

def main():
    def ok(c):
        cnt = 0
        i = 0
        for a in aa[::-1]:
            while i < n and a + aa[i] >= c: i += 1
            cnt += i
            if cnt >= m: return True
        return False

    n, m = MI()
    aa = LI()
    aa.sort(reverse=True)
    
    l = aa[-1] * 2 - 1
    r = aa[0] * 2 + 1
    while l + 1 < r:
        c = (l + r) // 2
        if ok(c):
            l = c
        else:
            r = c
    min_once = l
    
    ans = 0
    i = 0
    s = 0
    cnt = 0
    for a in aa[::-1]:
        while i < n and a + aa[i] >= min_once:
            s += aa[i]
            i += 1
        ans += s + i * a
        cnt += i
    
    if cnt - m:
        ans -= min_once * (cnt - m)
    print(ans)

main()

Solution (Python):
from bisect import bisect_left, bisect_right
n, m = map(int, input().split())
a = list(map(int, input().split()))
a.sort()


a_cnt = [0] * (2 * a[-1] + 5)
for ai in a:
    a_cnt[ai] += 1
for i in reversed(range(len(a_cnt) - 1)):
    a_cnt[i] += a_cnt[i + 1]
a_cs = [0] * (n + 1)
for i in range(n):
    a_cs[i + 1] = a_cs[i] + a[n - 1 - i]


def check(x):
    
    
    cnt = 0
    for ai in a:
        if x >= ai:
            cnt += a_cnt[x - ai]
        else:
            cnt += n
    return cnt >= m


lb = 1  
ub = 2 * a[-1] + 3  
while ub - lb > 1:
    mid = (ub + lb) // 2
    if check(mid):
        lb = mid
    else:
        ub = mid



ans = 0
cnt = 0
for ai in a:
    if lb >= ai:
        idx = a_cnt[lb - ai]
    else:
        idx = n
    ans += ai * idx + a_cs[idx]
    cnt += idx
ans -= lb * (cnt - m)
print(ans)

Solution (Python):
from itertools import accumulate
N, M = map(int, input().split())
A = list(map(int, input().split()))

A.sort(reverse = True)

B = list(accumulate(A))

def pnum(X,N,A,B):
    l = 0
    r = N-1
    num = 0
    s = 0
    while l != N and r != -1:
        if A[l] + A[r] >= X:
            num += r+1
            s += (A[l] * (r+1)) + B[r] 
            l += 1
        else:
            r -= 1
    return num, s


l = 0
r = A[0] * 2 + 1

while r-l > 1:
    m = (r + l) // 2
    a, b = pnum(m, N, A, B)
    if a < M:
        r = m
    else:
        l = m
        num = a
        ans = b

if num > M:
    ans -= l * (num - M)
print(ans)

Solution (Python):


import numpy as np

N,M = map(int,input().split())
A = np.array(list(map(int,input().split())))
A.sort()

def shake_cnt(x):
    X = np.searchsorted(A,x-A)  
    return N*N - X.sum()
    


L = 0
R = 10**6   
while L+1 < R:
    mid = (L+R) // 2
    if shake_cnt(mid) >= M:
        L = mid
    else:
        R = mid


Acum = np.zeros(N+1)
Acum[1:] = np.cumsum(A) 
X = np.searchsorted(A,R-A)
shake = N*N - X.sum()
ans = (Acum[-1] - Acum[X]).sum() + (A*(N-X)).sum()

ans += (M-shake)*L
print(int(ans))

Solution (Python):
from bisect import bisect_left, bisect_right
n, m = map(int, input().split())
a = list(map(int, input().split()))
a.sort()


a_cs = [0] * (n + 1)
for i in range(n):
    a_cs[i + 1] = a_cs[i] + a[n - 1 - i]


def check(x):
    
    
    cnt = 0
    for ai in a:
        idx = bisect_left(a, x - ai)
        cnt += n - idx
    return cnt >= m


lb = 1  
ub = 2 * 10**5 + 5  
while ub - lb > 1:
    mid = (ub + lb) // 2
    if check(mid):
        lb = mid
    else:
        ub = mid



ans = 0
cnt = 0
for ai in a:
    idx = bisect_left(a, lb - ai)
    ans += ai * (n - idx) + a_cs[n - idx]
    cnt += n - idx
ans -= lb * (cnt - m)
print(ans)

Solution (Python):
import sys
readline = sys.stdin.buffer.readline
import numpy as np


N,M = map(int,readline().split())
A = np.array(readline().split(),np.int64)
A.sort()

def shake_cnt(x):
    
    X = np.searchsorted(A,x-A)  
    return N * N - X.sum()


left = 0        
right = 10 ** 6 
while left + 1 < right:
    x = (left + right) // 2
    if shake_cnt(x) >= M:
        left = x
    else:
        right = x


X = np.searchsorted(A,right-A) 
Acum = np.zeros(N+1,np.int64)  
Acum[1:] = np.cumsum(A)

shake = N * N - X.sum()     
happy = (Acum[-1] - Acum[X]).sum() + (A * (N - X)).sum()    

happy += (M - shake) * left     
print(happy)

Solution (Python):
n, m = map(int, input().split())
A = list(map(int, input().split()))
A.sort(reverse=True)
cumsumA = []

total = 0
for a in A:
    total += a
    cumsumA.append(total)
cumsumA.append(0)

l = 0
r = 10 ** 6
while l + 1 < r:
    num = (l+r)//2
    pair_idx = n-1
    count = 0
    for a in A:
        while a + A[pair_idx] < num and pair_idx >= 0:
            pair_idx -= 1

        count += pair_idx+1

    if count < m:
        r = num
        continue
    elif count >= m:
        l = num


ans = 0
count = 0
num = l
pair_idx = n-1
for a in A:
    while a + A[pair_idx] < num and pair_idx >= 0:
        pair_idx -= 1

    count += pair_idx+1
    ans += cumsumA[pair_idx] + a*(pair_idx+1)

diff = count - m
ans -= num * diff
print(ans)

Solution (Python):
n, m = map(int, input().split())
A = list(map(int, input().split()))
A.sort(reverse=True)
cumsumA = []

total = 0
for a in A:
    total += a
    cumsumA.append(total)
cumsumA.append(0)

l = 0
r = 10 ** 6
while l + 1 < r:
    num = (l+r)//2
    pair_idx = n-1
    count = 0
    for a in A:
        while a + A[pair_idx] < num and pair_idx >= 0:
            pair_idx -= 1

        count += pair_idx+1

    if count < m:
        r = num
        continue
    elif count >= m:
        l = num


ans = 0
count = 0
num = l
pair_idx = n-1
for a in A:
    while a + A[pair_idx] < num and pair_idx >= 0:
        pair_idx -= 1

    count += pair_idx+1
    ans += cumsumA[pair_idx] + a*(pair_idx+1)

diff = count - m
ans -= num * diff
print(ans)

Solution (Python):
import bisect

n, m = map(int, input().split())
a = list(map(int, input().split()))
a.sort()


left = a[0] * 2
right = a[-1] * 2 + 1
mid = (left + right) // 2
while right - left > 1:
    mid = (left + right) // 2
    count = 0
    for i in range(n):
        j = bisect.bisect_left(a, mid - a[i])
        count += n - j

    if count >= m:
        left = mid
    else:
        right = mid

acum = [0] * (n+1)
for i in range(n):
    acum[i+1] = acum[i] + a[i]

ans = 0
count = 0
for i in range(n):
    j = bisect.bisect_left(a, left - a[i])
    ans += acum[n] - acum[j] + a[i] * (n-j)
    count += n-j
ans -= left * (count - m)
print(ans)

Solution (Python):
from bisect import bisect_left as bisect
N,M = map(int,input().split())
A = list(map(int,input().split()))
A.sort()
def combination(X):
    res = 0
    for i in range(N):
        j = bisect(A,X-A[i])
        res += N-j
    return res
l = 0
r = 2*10**5+1
m = (l+r)//2
for i in range(20):
    c = combination(m)
    if c < M:
        r = m
    else:
        l = m
    m = (l+r)//2
    
if combination(m) >= M:
     m += 1

S = [0]*(N+1)
for i in range(N):
    S[i]=S[i-1]+A[i]

ans = 0
for i in range(N):
    j = bisect(A,m-A[i])
    ans += A[i]*(N-j)
    ans += S[N-1]-S[j-1]
ans += (M-combination(m))*(m-1)
print(ans)

Solution (Python):
from bisect import bisect_left, bisect_right
n, m = map(int, input().split())
a = list(map(int, input().split()))
a.sort()


a_cs = [0] * (n + 1)
for i in range(n):
    a_cs[i + 1] = a_cs[i] + a[n - 1 - i]


def check(x):
    
    
    cnt = 0
    for ai in a:
        idx = bisect_left(a, x - ai)
        cnt += n - idx
    return cnt >= m


lb = 1  
ub = 2 * 10**5 + 5  
while ub - lb > 1:
    mid = (ub + lb) // 2
    if check(mid):
        lb = mid
    else:
        ub = mid



ans = 0
cnt = 0
for ai in a:
    idx = bisect_left(a, lb - ai)
    ans += ai * (n - idx) + a_cs[n - idx]
    cnt += n - idx
ans -= lb * (cnt - m)
print(ans)

Solution (Python):
from numpy.fft import rfft, irfft
import numpy as np

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

a = np.bincount(A, minlength=2 ** 18)

fa = rfft(a)
a2 = irfft(fa ** 2)
a2 = np.round(a2).astype(np.int)

cum_a2 = np.cumsum(a2)

index = np.searchsorted(cum_a2, N ** 2 - M)

a2[index] -= cum_a2[index] - (N ** 2 - M)

ans = sum(A) * 2 * N - (a2[: index + 1] * np.arange(index + 1)).sum()
print(ans)

Solution (Python):
N, M, *A = map(int, open(0).read().split())


def func(X):
    cnt = 0
    j = N
    for i in range(1, N + 1):
        while A[i] + A[j] < X:
            j -= 1
        cnt += j
    return cnt

A = [10 ** 11] + sorted(A, reverse=True)

lb = 0
ub = 10 ** 11
while ub - lb > 1:
    mid = (lb + ub) // 2
    if func(mid)>= M:
        lb = mid
    else:
        ub = mid

acc = A[:]
acc[0] = 0
for i in range(N):
    acc[i + 1] += acc[i]

ans = 0
j = N
for i in range(1, N + 1):
    while A[i] + A[j] < lb + 1:
        j -= 1
    ans += A[i] * j + acc[j]
    M -= j
    
if M:
    ans += lb * M
    
print(ans)

Solution (Python):
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import numpy as np

N,M = map(int,readline().split())
A = np.array(read().split(),np.int64)

f = np.bincount(A)

fft_len = 1<<18
fft = np.fft.rfft; ifft = np.fft.irfft
Ff = fft(f,fft_len)
G = np.rint(ifft(Ff * Ff,fft_len)).astype(np.int64) 

Gcum = G.cumsum()
remove_cnt = N * N - M
i = np.searchsorted(Gcum,remove_cnt)


x = remove_cnt - Gcum[i-1]
remove_sum = (G[:i] * np.arange(i,dtype=np.int64)).sum() + i * x

answer = A.sum() * 2 * N - remove_sum
print(answer)

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 LS(): return list(map(list, input().split()))
def S(): return list(input().rstrip())
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, M = LI()
    a = LI()
    a.sort()
    l = 2 * (10 ** 5) + 1
    r = 2
    def f(m):
        tmp = 0
        for i in range(n):
            tmp += bisect_left(a, m - a[i])

        return n ** 2 - tmp >= M
    while l - r > 1:
        m = (l + r) // 2
        if f(m):
            r = m
        else:
            l = m
    ans = 0
    tmp = 0
    for i in range(n):
        b = n - bisect_left(a, r + 1 - a[i])
        tmp += b
        ans += a[i] * b * 2
    tmp = M - tmp
    ans += tmp * r
    print(ans)
    return



if __name__ == '__main__':
    solve()

Solution (Python):
from numpy.fft import rfft, irfft
import numpy as np

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

a = np.bincount(A, minlength=2 ** 18)

fa = rfft(a)
a2 = irfft(fa ** 2)
a2 = np.round(a2).astype(np.int)

cum_a2 = np.cumsum(a2)

index = np.searchsorted(cum_a2, N ** 2 - M)

a2[index] -= cum_a2[index] - (N ** 2 - M)

ans = sum(A) * 2 * N - (a2[: index + 1] * np.arange(index + 1)).sum()
print(ans)

Solution (Python):
from bisect import bisect, bisect_left
from itertools import accumulate


def check(x):
    count = 0
    for a in aaa:
        if x <= a:
            count += n
        else:
            count += n - bisect_left(aaa, x - a)
    
    return count >= m


n, m = map(int, input().split())
aaa = list(map(int, input().split()))
aaa.sort()
acc = [0] + list(accumulate(aaa))
sum_a = acc[-1]

l, r = 0, 200001
while l + 1 < r:
    mid = (l + r) // 2
    if check(mid):
        l = mid
    else:
        r = mid



ans = 0
count = 0
for a in aaa:
    if l <= a:
        ans += sum_a + a * n
        count += n
    else:
        omit = bisect(aaa, l - a)
        if omit < n:
            ans += sum_a - acc[omit] + a * (n - omit)
            count += n - omit
        
ans += l * (m - count)
print(ans)

Solution (Python):
from bisect import bisect_left, bisect_right

maxA = 10**5

N, M = map(int, input().split())
As = list(map(int, input().split()))

As.sort()

def isOK(v):
    num = 0
    for A in As:
        i = bisect_left(As, v-A)
        num += N-i
    return num >= M

ng, ok = 2*maxA+1, 0
while abs(ok-ng) > 1:
    mid = (ng+ok) // 2
    if isOK(mid):
        ok = mid
    else:
        ng = mid

accRevAs = [As[-1]]
for A in reversed(As[:-1]):
    accRevAs.append(accRevAs[-1]+A)
accRevAs.reverse()

v = ok
ans = 0
numGt = 0

for A in As:

    iR = bisect_right(As, v-A)
    if iR < N:
        ans += accRevAs[iR] + A*(N-iR)
    numGt += N-iR


ans += v * (M-numGt)
print(ans)

Solution (Python):



import sys
sys.setrecursionlimit(10**6)
readline = sys.stdin.readline


n,m = [int(i) for i in readline().split()]
a = [int(i) for i in readline().split()]



a.sort()



def check(mid):
    num = 0
    p = 0
    for i in range(n-1,-1,-1):
        c = 0
        while p<n and a[p]+a[i] < mid:
            p += 1
        num += n-p
    
    return num


ok = a[0]*2-1
ng = a[-1]*2+1

while ng - ok > 1:
    mid = (ng+ok)//2
    c = check(mid)
    if c<m:
        ng = mid
    else:
        ok = mid
    
    
    

    



ans = 0
p = n-1
num = 0

res = [0]*n
for i in range(n):
    c = 0
    while p>=0 and a[p]+a[i] >= ng:
        p -= 1
    res[i] = n-1-p
    num += n-1-p

for i in range(n):
    ans += res[i]*a[i]
ans *= 2

print(ans+ok*(m-num))

Solution (Python):


from collections import defaultdict
from bisect import bisect_left
def main():
    N, M = map( int, input().split())
    A = list( map( int, input().split()))
    A.sort()
    aaa = [0]*(10**5+1)
    d = defaultdict(int)
    appA = [0]*(N+1)
    for i in range(N,0,-1):
        appA[i-1] = A[i-1] + appA[i]

    for a in A:
        aaa[a] += 1
        d[a] += 1
    for i in range(10**5, 0, -1):
        aaa[i-1] += aaa[i]
    l = 0
    r = 2*10**5+1
    while r - l > 1:
        m = (l+r)//2
        check = True
        cnt = 0
        for k, v in d.items():
            if m - k > 10**5:
                continue
            cnt += v*aaa[max(m-k, 0)]

        if cnt >= M:
            l = m
        else:
            r = m
    ans = 0
    m = 0
    for k, v in d.items():
        if l - k > 10**5:
            continue
        if k + A[-1] < l:
            continue
        ans += k*v*(N - bisect_left(A, l-k))+appA[ bisect_left(A, l-k)]*v

        m += v*(N - bisect_left(A, l-k))
    if m > M:
        ans -= l*(m-M)
    print(ans)
if __name__ == '__main__':
    main()

Solution (Python):
import os
import sys

import numpy as np

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

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



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

A.sort()
A = np.array(A, dtype=int)
cs = A[::-1].cumsum()[::-1]


def test(k):
    
    counts = N - np.searchsorted(A, k - A)
    return counts.sum() >= M


def count(k):
    counts = N - np.searchsorted(A, k - A)
    cnt = counts.sum()

    pos = np.searchsorted(A, k - A)
    ret = ((A * counts)[pos < N] + cs[pos[pos < N]]).sum()
    ret -= abs(M - cnt) * k
    return ret


ok = 0
ng = max(A) * 2 + 1
while abs(ok - ng) > 1:
    mid = (ok + ng) // 2
    if test(mid):
        ok = mid
    else:
        ng = mid

ans = count(ok)
print(ans)