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):
import numpy as np
from numpy.fft import rfft,irfft
import sys
sys.setrecursionlimit(10**9)

def mi(): return map(int,input().split())
def ii(): return int(input())
def isp(): return input().split()
def deb(text): print("-------\n{}\n-------".format(text))

INF=10**20
class Polynomial:
    def __init__(self,dim=0,const=1):
        self.values = [(dim,const)]
    def __add__(self,other):
        return self.values + other.values


class PolySolver:
    def __init__(self,size,init_val=1):
        self.size = size
        self.f = np.zeros(size, np.int64)
        self.f[0] = init_val


    def multiple(self,polynomial,MOD): 
        new_F = np.zeros(self.size, np.int64)
        for dim,const in polynomial:
            if dim != 0:
                g = np.zeros(self.size, np.int64)
                g[dim:] += self.f[:-dim] 
            else:
                g = const * self.f.copy()

            new_F += g

        self.f = new_F
        
    
    def add(self,polynomial): 
        for dim,const in polynomial:
            self.f[dim] += const
            

    def get_coefficient(self,dim):
        return self.f[dim]

    def fft(self,A,B,fft_len=1<<18):
        x = irfft(rfft(A,fft_len) * rfft(B,fft_len))
        return np.rint(x).astype(np.int64)

    def convolve(self,g):
        h = self.fft(self.f,g)[:len(self.f) + len(g) - 1]
        self.f = h
    


def main():
    N,M=mi()
    A=list(mi())
    
    mx = 2*10**5+10
    F = PolySolver(mx,init_val=0)
    for a in A:
        F.add(Polynomial(dim=a).values)
    
    g = F.f
    F.convolve(g)

    ans = 0
    handshake = M
    for x in range(mx-1,1,-1):
        count = F.get_coefficient(x)
        if count <= 0: continue

        if count > handshake:
            ans += handshake * x
        else:
            ans += count * x
 
        handshake = max(handshake-count,0)
        if handshake == 0:
            break

    print(ans)

if __name__ == "__main__":
    main()

Solution (Python):
def getN():
    return int(input())
def getNM():
    return map(int, input().split())
def getList():
    return list(map(int, input().split()))
def getArray(intn):
    return [int(input()) for i in range(intn)]
def input():
    return sys.stdin.readline().rstrip()
def rand_N(ran1, ran2):
    return random.randint(ran1, ran2)
def rand_List(ran1, ran2, rantime):
    return [random.randint(ran1, ran2) for i in range(rantime)]
def rand_ints_nodup(ran1, ran2, rantime):
  ns = []
  while len(ns) < rantime:
    n = random.randint(ran1, ran2)
    if not n in ns:
      ns.append(n)
  return sorted(ns)

def rand_query(ran1, ran2, rantime):
  r_query = []
  while len(r_query) < rantime:
    n_q = rand_ints_nodup(ran1, ran2, 2)
    if not n_q in r_query:
      r_query.append(n_q)
  return sorted(r_query)

from collections import defaultdict, deque, Counter
from sys import exit
from decimal import *
from heapq import heapify, heappop, heappush
from math import sqrt
from fractions import gcd
import random
import string
import copy
from itertools import combinations, permutations, product
from operator import mul, itemgetter
from functools import reduce
from bisect import bisect_left, bisect_right

import sys
sys.setrecursionlimit(1000000000)
mod = 10 ** 9 + 7





N, K = getNM()
A = getList()

A.sort(reverse = True)


def cnt(x):
    r = 0
    ans = 0
    list_num = [0] * N
    for i in range(N - 1, -1, -1):
        while r < N and A[r] + A[i] >= x:
            r += 1
        list_num[i] = r
        ans += r

    return [ans, list_num]

ng = 0
ok = 2 * (10 ** 6)

while ok - ng > 1:
    mid = (ok + ng) // 2
    if cnt(mid)[0] <= K:
        ok = mid
    else:
        ng = mid


imos = [0] + copy.deepcopy(A)
for i in range(N):
    imos[i + 1] += imos[i]

ans = 0
plus = cnt(ok)[1]


for i in range(N):
    ans += (A[i] * plus[i]) + imos[plus[i]]

ans += ng * (K - sum(plus))
print(ans)

Solution (Python):
import numpy as np
from numpy.fft import rfft,irfft
import sys
sys.setrecursionlimit(10**9)

def mi(): return map(int,input().split())
def ii(): return int(input())
def isp(): return input().split()
def deb(text): print("-------\n{}\n-------".format(text))

INF=10**20
class Counter:
    def __init__(self):
        self.dict = {}

    def add(self,x):
        if x in self.dict: self.dict[x] += 1
        else: self.dict[x] = 1

    def decrement(self,x):
        self.dict[x] -= 1
        if self.dict[x] <= 0:
            del self.dict[x]

    def get_dict(self):
        return self.dict




def fft(A,B,fft_len=2*10**5):
    x = irfft(rfft(A,fft_len) * rfft(B,fft_len))+0.5
    return x.astype(int)


def main():
    N,M=mi()
    A=list(mi())

    counter = Counter()
    for a in A: 
        counter.add(a)
    
    count_map = counter.get_dict()
    F=[0]*2*10**5
    for a,count in count_map.items():
        F[a] = count
    
    G = fft(F,F,fft_len=1<<18)

    ans = 0
    handshake = M
    for x in range(len(G)-1,1,-1):
        if G[x] <= 0: continue
        count = G[x]
        if count > handshake:
            ans += handshake * x
        else:
            ans += count * x

        handshake = max(handshake-count,0)
        if handshake == 0:
            break

    print(ans)



if __name__ == "__main__":
    main()

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


def main():
    n, m = map(int, input().split())
    a = list(map(int, input().split()))
    a.sort()
    
    l, r = -1, 10**6
    while r-l > 1:
        k = (r+l)//2
        count = 0
        for i in range(n):
            count += n-bisect.bisect_left(a, k-a[i])
        
        if count > m:
            l = k
        else:
            r = k
    ans = 0
    count = 0
    b = [0]*(n+1)
    for i in range(n):
        b[i+1] += b[i]+a[i]
    for i in range(n):
        index = bisect.bisect_left(a, r-a[i])
        ans += b[n]-b[index]+a[i]*(n-index)
        count += n-index
    
    if count < m:
        ans += (m-count)*l
    
    print(ans)
    
    
    
    
if __name__ == "__main__":
    main()

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


def main():
    n, m = map(int, input().split())
    a = list(map(int, input().split()))
    a.sort()
    
    l, r = -1, 10**6
    while r-l > 1:
        k = (r+l)//2
        count = 0
        for i in range(n):
            count += n-bisect.bisect_left(a, k-a[i])
        
        if count > m:
            l = k
        else:
            r = k
    ans = 0
    count = 0
    b = [0]*(n+1)
    for i in range(n):
        b[i+1] += b[i]+a[i]
    for i in range(n):
        index = bisect.bisect_left(a, r-a[i])
        ans += b[n]-b[index]+a[i]*(n-index)
        count += n-index
    
    if count < m:
        ans += (m-count)*l
    
    print(ans)
    
    
    
    
if __name__ == "__main__":
    main()

Solution (Python):

from bisect import bisect_left
from itertools import accumulate

def resolve():
    def getCount(x):
        cnt = 0
        for a in A:
            
            idx = bisect_left(A, x - a)
            
            cnt += N - idx
        
        return cnt >= M

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

    ng = 10 ** 9
    ok = 0
    
    while ng - ok > 1:
        mid = (ok + ng) // 2
        if getCount(mid):
            ok = mid
        else:
            ng = mid

    
    
    B = [0] + list(accumulate(A))
    ans = 0
    cnt = 0
    for i in range(N):
        idx = bisect_left(A, ok - A[i])
        ans += B[N] - B[idx] + A[i] * (N - idx)
        cnt += N - idx

    rest = (cnt - M) * ok
    print(ans - rest)


if __name__ == "__main__":
    resolve()

Solution (Python):
import numpy as np
from numpy.fft import rfft,irfft
N,M = map(int,input().split())
A = list(map(int,input().split()))
a = np.zeros(10**5+1)
for i in range(N):
    a[A[i]] += 1
b = rfft(a,2*10**5+1)
b = irfft(b*b,2*10**5+1)
b = np.rint(b).astype(np.int64)
c = 0
ans = 0
for n in range(2*10**5,1,-1):
    if c+b[n]<M:
        c += b[n]
        ans += b[n]*n
    else:
        ans += (M-c)*n
        c = M
        break
print(ans)

Solution (Python):
import sys
import numpy as np




N, M = map(int, sys.stdin.readline().rstrip().split())
A = np.array(sys.stdin.readline().rstrip().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 1 < right - left:  
    x = (left + right) // 2
    if shake_cnt(x) >= M:
        left = x
    else:
        right = x


X = np.searchsorted(A, right - A)  
shake = N * N - X.sum()  

Acum = np.zeros(N + 1, np.int64)  

Acum[1:] = np.cumsum(A)  



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

happy += (M - shake) * left  

print(happy)

Solution (Python):
import os
import sys
import numpy as np

def solve(N, M, A):
    
    ok = 0
    ng = 202020
    while ok + 1 < ng:
        c = ok+ng >> 1
        cnt = 0
        idx_A = N-1
        for a1 in A:
            while idx_A >= 0 and A[idx_A]+a1 < c:
                idx_A -= 1
            cnt += idx_A + 1
        if cnt >= M:
            ok = c
        else:
            ng = c
    idx_A = N-1
    cum_A = np.zeros(len(A)+1, dtype=np.int64)
    cum_A[1:] = np.cumsum(A)
    ans = 0
    cnt = 0
    for a1 in A:
        while idx_A >= 0 and A[idx_A]+a1 < ok:
            idx_A -= 1
        cnt += idx_A + 1
        ans += cum_A[idx_A+1] + (idx_A+1)*a1
    ans -= (cnt-M) * ok
    return ans



numba_config = [
    [solve, "i8(i8,i8,i8[:])"],
]
if sys.argv[-1] == "ONLINE_JUDGE":
    from numba import njit
    from numba.pycc import CC
    cc = CC("my_module")
    for func, signature in numba_config:
        vars()[func.__name__] = njit(signature)(func)
        cc.export(func.__name__, signature)(func)
    cc.compile()
    exit()
elif os.name == "posix":
    exec(f"from my_module import {','.join(func.__name__ for func, _ in numba_config)}")
else:
    from numba import njit
    for func, signature in numba_config:
        vars()[func.__name__] = njit(signature, cache=True)(func)
    print("compiled!", file=sys.stderr)



def main():
    N, M = map(int, input().split())
    A = np.array(sorted(map(int, input().split()), reverse=True), dtype=np.int64)
    ans = solve(N, M, A)
    print(ans)

main()

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

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

def isOK(x):
    num=0
    for k,a in enumerate(A,start=1):
        i=bisect_left(A,x-a)
        num+=n-i
    return num>=m



l=0
r=2*10**5+1
while l+1<r:
    mid=(l+r)//2
    if isOK(mid):
        l=mid
    else:
        r=mid

s=0
num=0

acc=list(accumulate([0]+A))

for a in A:
    i=bisect_right(A,l-a)
    s+=a*(n-i)+acc[n]-acc[i]
    num+=n-i

ans=s+l*(m-num)

print(ans)

Solution (Python):
N, M = map(int, input().split())
A = list(map(int, input().split()))
A.sort()

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



def find(power):
    cnt = 0
    sam = 0
    for a in A:
        left = -1
        right = N
        while right-left>1:
            mid = left + (right-left)//2
            if a+A[mid] >= power:
                right = mid
            else:
                left = mid
        cnt += (N-right)
        sam += (AS[N]-AS[right]) + (N-right)*a
        
        
    return cnt, sam

left = 0
right = 10**15
while right-left>1:
    mid = left + (right-left)//2
    cnt, sam = find(mid)
    if cnt>M:
        left = mid
    else:
        right = mid


cnt, ans = find(left)
print(ans - (cnt-M)*left)

Solution (Python):
from bisect import bisect_left as BL 
N,M=map(int,input().split())
A=list(map(int,input().split()))
A.sort()


def isok(border):
  res=0
  for i in range(N):
    res+=N-BL(A,border-A[i])
  return res>=M
    

ok=0
ng=A[-1]*2+1
mid=(ok+ng)//2


while abs(ok-ng)>1:
  mid=(ok+ng)//2
  if isok(mid):
    ok=mid
  else:
    ng=mid
  
B=reversed(A)
SB=[0]
for b in B:
  SB.append(b+SB[-1])
  
ans=0
shake_c=0
for i in range(N):
  goodshake=N-BL(A,ok-A[i])
  shake_c+=goodshake
  ans+=A[i]*goodshake
  ans+=SB[goodshake]  
ans-=(shake_c-M)*ok
print(ans)

Solution (Python):
from itertools import accumulate

def bisect(ng, ok, judge, eps=1):
    while abs(ng-ok) > eps:
        m = (ng+ok)//2
        if judge(m):
            ok = m
        else:
            ng = m
    return ok

N,M = map(int,input().split())
A = sorted(map(int,input().split()))
acc = [0]+list(accumulate(A))

def judge(p):
    j = N
    cnt = 0
    for al in A:
        while j > 0 and al+A[j-1] >= p:
            j -= 1
        cnt += N-j
    return cnt <= M

p = bisect(0, A[-1]*2+1, judge)

j = N
cnt = 0
res = 0
nv = 0
max2 = lambda x,y: x if x > y else y
for al in A:
    while j > 0 and al+A[j-1] >= p:
        j -= 1
    cnt += N-j
    res += (N-j)*al + acc[-1]-acc[j]
    if j > 0:
        nv = max2(nv, al+A[j-1])

res += nv*(M-cnt)
print(res)

Solution (Python):
"""
チルノ（無色）のプログラミング教室～☆

こんなの無色コーダーのレベルじゃないよ～
え？　ABC D-Pairs とほとんど同じ方針？
こんな典型問題も解けないから灰色にすらなれない無色なんだよ、だって？
うわーん、無色っていうなー！！！
"""


def get_count(x):
    """
    :return: x 以上となる組の数
    """
    cnt = 0
    for i in range(N):
        cnt += bisect_left(A, x-A[i])
    res = N**2 - cnt
    return res

def binary_search_int(ok, ng, test):
    """
    :param ok: solve(x) = True を必ず満たす点
    :param ng: solve(x) = False を必ず満たす点
    """
    while abs(ok - ng) > 1:
        mid = (ok + ng) // 2
        if test(mid):
            ok = mid
        else:
            ng = mid
    return ok

def test(x):
    return get_count(x) >= M

def get_sum(x):
    res = 0
    cnt = 0
    for i in range(N):
        temp = bisect_left(A, x-A[i])
        res += A[i]*(N - temp) + B[-1] - B[temp]
        cnt += N - temp
    return res - x*(cnt-M)



from bisect import *
from itertools import accumulate

N, M = map(int, input().split())
A = list(map(int, input().split()))
A.sort()
B = list(accumulate([0] + A))

print(get_sum(binary_search_int(0, 2*10**5+10, test)))

Solution (Python):
import numpy as np

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

na = np.zeros(2**18)
for i in a:
    na[i] += 1

fa = np.fft.fft(na)

c = np.round(np.fft.ifft(fa*fa)).astype(int)

ans = 0
cm = 0
for i in range(2**18 - 1, 1, -1):
    cm += c[i]
    ans += i*c[i]
    if cm > m:
        ans -= i * (cm - m)
    if cm >= m:
        break

print(ans)

Solution (Python):
import numpy as np

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

a.sort()
na = np.zeros(2**18)
for i in a:
    na[i] += 1

fa = np.fft.fft(na)

c = np.round(np.fft.ifft(fa*fa)).astype(int)

ans = 0
cm = 0
for i in range(2**18 - 1, 1, -1):
    cm += c[i]
    ans += i*c[i]
    if cm > m:
        ans -= i * (cm - m)
    if cm >= m:
        break

print(ans)

Solution (Python):
import bisect
import copy

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

a.sort()
asum = [0]
tmp = 0
for i in range(n):
    tmp += a[i] 
    asum.append(tmp)

lb = 1
ub = 2*10**5 + 1
while ub - lb > 1:
    cx = (ub+lb)//2
    cnt = 0
    total = 0
    for ca in a:
        cy = cx - ca
        ci = bisect.bisect_left(a, cy)
        cnt += n - ci
        total += asum[n] - asum[ci] + (n - ci) * ca
    if cnt >= m:
        lb = cx
        ans = copy.deepcopy(total) - (cnt - m) * cx
    else:
        ub = cx

print(ans)

Solution (Python):
import bisect
import copy

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

a.sort()
asum = [0]
tmp = 0
for i in range(n):
    tmp += a[i] 
    asum.append(tmp)

lb = 1
ub = 2*10**5 + 1
while ub - lb > 1:
    cx = (ub+lb)//2
    cnt = 0
    total = 0
    for ca in a:
        cy = cx - ca
        ci = bisect.bisect_left(a, cy)
        cnt += n - ci
        total += asum[n] - asum[ci] + (n - ci) * ca
    if cnt >= m:
        lb = cx
        ans = copy.deepcopy(total) - (cnt - m) * cx
    else:
        ub = cx

print(ans)

Solution (Python):
n,m=map(int,input().split())
a=list(map(int,input().split()))
a.sort()
if n==1:
    print(a[0]*2)
    exit()
from bisect import bisect_left,bisect_right

def f(x):
    global n,a
    ret=0
    for i in range(n):
        ret+=(n-bisect_left(a,x-a[i]))
    return ret

l,r=-1,10**6
while l+1<r:
    k=l+(r-l)//2
    if f(k)>=m:
        l=k
    else:
        r=k
co,ans=0,0
for i in range(n):
    co+=(n-bisect_right(a,l-a[i]))
    ans+=(n-bisect_right(a,l-a[i]))*a[i]
print(2*ans+(m-co)*l)

Solution (Python):
n,m=map(int,input().split())
a=list(map(int,input().split()))
a.sort()
if n==1:
    print(a[0]*2)
    exit()
from bisect import bisect_left,bisect_right

def f(x):
    global n,a
    ret=0
    for i in range(n):
        ret+=(n-bisect_left(a,x-a[i]))
    return ret
l,r=-1,10**6
while l+1<r:
    k=l+(r-l)//2
    if f(k)>=m:
        l=k
    else:
        r=k


if f(r)>=m:
    co,ans=0,0
    for i in range(n):
        co+=(n-bisect_right(a,r-a[i]))
        ans+=(n-bisect_right(a,r-a[i]))*a[i]
    print(2*ans+(m-co)*r)
else:
    co,ans=0,0
    for i in range(n):
        co+=(n-bisect_right(a,l-a[i]))
        ans+=(n-bisect_right(a,l-a[i]))*a[i]
    print(2*ans+(m-co)*l)

Solution (Python):


import math

import itertools
from collections import deque
from collections import Counter

from fractions  import gcd

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

s=[0]*n
a.sort()
s[0]=a[0]
for i in range(1,n):
    s[i]=s[i-1]+a[i]
res=0
l=0
r=2*a[-1]+1
while r>l+1:
    x=(l+r)//2
    res=0
    for ai in a:
        j=x-ai
        cnt=bisect.bisect_left(a,j)
        res+=(n-cnt)
    if res<m:
        r=x
    else:
        l=x
ans=0
c=0
for ai in a:
    cnt=bisect.bisect_left(a,l-ai)
    ans=ans+ai*(n-cnt)+(s[-1]-s[cnt-1]) if cnt!=0 else ans+ai*(n-cnt)+s[-1]
    c+=(n-cnt)
ans-=(c-m)*l
print(ans)

Solution (Python):
n, m = map(int, input().split())
a = sorted(list(map(int, input().split())), reverse=True)
s = [0]
for ai in a:
    s.append(ai + s[-1])

def count(x, accum=False):
    ret = 0
    for ai in a:
        lo, hi = -1, n
        while hi - lo > 1:
            mid = (lo + hi) // 2
            if ai + a[mid] >= x:
                lo = mid
            else:
                hi = mid
        ret += ai * hi + s[hi] if accum else hi
    return ret

lo, hi = 0, 1000000000
while hi - lo > 1:
    mid = (lo + hi) // 2
    if count(mid) >= m:
        lo = mid
    else:
        hi = mid

print(count(lo, accum=True) - (count(lo) - m) * lo)

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

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

ok = -1
ng = (10 ** 5) * 2 + 1


import bisect
def isOk(x):
  cnt = 0
  for left in A:
    right = x - left
    ind = bisect.bisect_left(A, right)
    cnt += N - ind
  return cnt >= M
  
while abs(ng - ok) > 1:
  mid = abs(ng + ok) // 2
  if isOk(mid):
    ok = mid
  else:
    ng = mid



sumA = [0] + A.copy()
for i in range(1, len(sumA)):
  sumA[i] += sumA[i - 1]

ans = 0
num = 0
for left in A:
  right = ok - left
  ind = bisect.bisect_left(A, right)
  num += N - ind
  ans += (N - ind) * left + sumA[-1] - sumA[ind]
  
if num > M:
  ans -= (num - M) * ok
  
print(ans)

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

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

sumA = [0] + A.copy()
for i in range(1,len(sumA)):
  sumA[i] += sumA[i - 1]
  









ok = -1
ng = 10 ** 5 * 2 + 1

import bisect
def isOk(x):
  
  cnt = 0
  for left in A:
    right = x - left
    ind = bisect.bisect_left(A, right)
    cnt += N - ind
  if cnt > M:
    return True
  else:
    return False

while abs(ng - ok) > 1:
  mid = abs(ng + ok) // 2
  if isOk(mid):
    ok = mid
  else:
    ng = mid
    

ans = 0
num = 0
for left in A:
  right = ok - left
  ind = bisect.bisect_left(A, right)
  ans += sumA[-1] - sumA[ind] + (N - ind) * left
  num += N - ind

if num > M:
  ans -= (num - M) * ok

print(ans)

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

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

sumA = [0] + A.copy()
for i in range(1,len(sumA)):
  sumA[i] += sumA[i - 1]
  









ok = -1
ng = 10 ** 5 * 2 + 1

import bisect
def isOk(x):
  
  cnt = 0
  for left in A:
    right = x - left
    ind = bisect.bisect_left(A, right)
    
    cnt += N - ind
  
  if cnt > M:
    return True
  else:
    return False

while abs(ng - ok) > 1:
  
  mid = abs(ng + ok) // 2
  if isOk(mid):
    ok = mid
  else:
    ng = mid
    



ans = 0
num = 0


for left in A:
  right = ok - left
  ind = bisect.bisect_left(A, right)




  ans += sumA[-1] - sumA[ind] + (N - ind) * left
  num += N - ind

if num > M:
  ans -= (num - M) * ok

print(ans)

Solution (Python):
from bisect import bisect_left
N, M = map(int, input().split())
A = list(map(int, input().split()))
A.sort()


def get_count(tot):
    '''
    tot以上を満たす幸福度の握手の個数を数える
    '''
    ret = 0
    for i in range(N):
        cnt = bisect_left(A, tot-A[i])
        ret += N - cnt
    return ret



Asum = A[:]
for i in range(N-2, -1, -1):
    Asum[i] += Asum[i+1]


def get_sum(tot):
    '''
    左端を固定した時の右端の個数分だけ左端を足し、右端の区間和（累積和）も足す。
    '''
    ret = 0
    for i in range(N):
        idx = bisect_left(A, tot-A[i])
        cnt = N - idx
        ret += A[i]*cnt
        
        if idx == N:
            continue
        ret += Asum[idx]
        
    return ret



ok = 0
ng = 10**12
while ok + 1 != ng:
    md = (ok + ng)//2
    if M <= get_count(md):
        ok = md
    else:
        ng = md




ans = get_sum(ok+1) + (M - get_count(ok+1))*ok

print(ans)

Solution (Python):
from bisect import bisect_left
N, M = map(int, input().split())
A = list(map(int, input().split()))
A.sort()


def get_count(tot):
    '''
    tot以上を満たす幸福度の握手の個数を数える
    '''
    ret = 0
    for i in range(N):
        cnt = bisect_left(A, tot-A[i])
        ret += N - cnt
    return ret



Asum = A[:]
for i in range(N-2, -1, -1):
    Asum[i] += Asum[i+1]


def get_sum(tot):
    '''
    左端を固定した時の右端の個数分だけ左端を足し、右端の区間和（累積和）も足す。
    '''
    ret = 0
    for i in range(N):
        idx = bisect_left(A, tot-A[i])
        cnt = N - idx
        ret += A[i]*cnt
        
        if idx == N:
            continue
        ret += Asum[idx]
        
    return ret



ok = 0
ng = 10**12
while ok + 1 != ng:
    md = (ok + ng)//2
    if M <= get_count(md):
        ok = md
    else:
        ng = md




ans = get_sum(ok+1) + (M - get_count(ok+1))*ok

print(ans)

Solution (Python):
import sys 
from itertools import accumulate
from bisect import bisect_left
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines
N, M = map(int, readline().split())
A = list(map(int,readline().split()))
A = sorted(A)
S = [0] + list(accumulate(A[::-1]))

def solve(mid):
    cnt_all = 0
    for a in A:
        cnt = N-bisect_left(A,mid-a)
        cnt_all += cnt   
    if cnt_all < M:
        return True
    else:
        return False
left = 0
right = 10**6

while right-left > 1:
    mid = (right+left) // 2
    flag = solve(mid)
    if flag:
        right = mid
    else:
        left = mid

ans = 0 
all_cnt = 0
for a in A:
    cnt = N-bisect_left(A,left-a)
    all_cnt += cnt
    ans += S[cnt]
    ans += cnt*a
    
count = all_cnt-M
ans -= count*left
print(ans)

Solution (Python):
import sys
import itertools

import time
import math
import heapq
from collections import defaultdict
sys.setrecursionlimit(10 ** 7)
 
INF = 10 ** 18
MOD = 10 ** 9 + 7
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines



N, M = map(int, input().split())
A = list(sorted(map(int, input().split())))
S = sum(A)
acc = [0] * (N + 1)
for i in range(1, N + 1):
    acc[i] = acc[i - 1] + A[i - 1]

def check(x):
    total = 0
    cnt = 0
    for a in A:
        l = -1
        r = N
        while r - l > 1:
            m = (l + r) // 2
            if A[m] + a >= x:
                r = m
            else:
                l = m
        cnt += N - r
        total += acc[N] - acc[r]
        total += a * (N - r)
    return (total, cnt)

left = 0
right = 10 ** 6
while right - left > 1:
    mid = (left + right) // 2
    res = check(mid)
    if res[1] >= M:
        left = mid
    else:
        right = mid
res = check(left)
ans = res[0]
print(ans - (res[1] - M) * left)

Solution (Python):
from itertools import accumulate
from bisect import bisect_left, bisect_right
N, M = map(int, input().split())
*A, = map(int, input().split())
A.sort()
acc = list(accumulate(A[::-1]))[::-1]+[0]


def isok(x):
    cnt = 0
    for a in A:
        cnt += N - bisect_left(A, x-a)
    return cnt >= M


l = 0
r = 202020
while l+1 < r:
    m = (l+r)//2
    if isok(m):
        l = m
    else:
        r = m

ans = 0
cnt = 0
for a in A:
    i = bisect_left(A, l-a)
    ans += a*(N-i) + acc[i]
    cnt += N-i
if cnt > M:
    ans -= l * (cnt-M)
print(ans)

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

N, M = map(int, input().split())
*A, = map(int, input().split())
B = np.zeros(5*10**5)
for a in A:
    B[a] += 1

L = 5*10**5
FB = rfft(B, L)
C = np.rint(irfft(FB*FB)).astype(int)

ans = 0
for i in range(2*10**5, -1, -1):
    c = C[i]
    if not c:
        continue
    if M - c > 0:
        ans += i*c
        M -= c
    else:
        ans += i*M
        break
print(ans)

Solution (Python):
import math
import collections
N,M=map(int,input().split())
A=list(map(int,input().split()))
A.sort()
cumax=(2*10**5)+10
cuml=[0]*cumax

ind=0
l=[]
c=collections.Counter(A)
values=list(c.values())	
key=list(c.keys())	
for i in range(len(key)):
    l.append([key[i],values[i]])
lr=[]
v=[]
for i in range(len(key)):
    lr.append([0,l[i][0]])
    v.append(l[i][1])
for i in range(len(v)):
    l=lr[i][0]
    r=lr[i][1]
    cuml[l]=cuml[l]+v[i]
    if r+1<cumax:
        cuml[r+1]=cuml[r+1]-v[i]
for i in range(cumax-1):
    cuml[i+1]=cuml[i+1]+cuml[i]

maxind=2*10**5+10
minind=0
while maxind-minind>1:
    piv=minind+(maxind-minind)//2       
    cnt=0
    for i in range(N):
        x=piv-A[i]
        cnt=cnt+cuml[max(x,0)]
    if cnt<=M:
        maxind=piv
    else:
        minind=piv
    

thr=minind
cuml2=[0]*(N+1)
A.sort(reverse=True)
for i in range(N):
    cuml2[i+1]=cuml2[i]+A[i]

ans=0
counts=0

for i in range(N):
    x=thr-A[i]
    
    cnt=cuml[max(x,0)]
    
    ans=ans+(A[i]*cnt)+cuml2[cnt]
    counts=counts+cnt
ans=ans-(counts-M)*thr

print(ans)

Solution (Python):
from bisect import bisect_left
def check(x):
  cnt=0
  l=n
  r=n
  for i in range(n):
    while True:
      if l==0:
        break
      if a[l-1]<x-a[i]:
        break
      l-=1
    cnt+=(r-l)
  return cnt>=m

n,m=map(int,input().split())
a=list(map(int,input().split()))
a.sort()
ok=0
ng=2*10**5+10

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

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

ans=0
cnt=0
r=n

ok=ok+1 
for i in range(n):
  
  l=bisect_left(a,ok-a[i])
  cnt+=r-l
  ans+=(r-l)*a[i]+(ru[r]-ru[l])
ans+=(ok-1)*(m-cnt) 
print(ans)

Solution (Python):

from bisect import bisect_left
from itertools import accumulate
def resolve():
    def getCount(x):
        count = 0
        for Ai in A:
            idx = bisect_left(A, x - Ai)
            count += N - idx
        return count >= M


    N, M = map(int, input().split())
    A = sorted(map(int, input().split()))
    A_r = A[::-1]
    B = [0] + list(accumulate(A_r))

    MIN = 0
    MAX = 2 * 10 ** 5 + 1
    while MAX - MIN > 1:
        MID = (MIN + MAX) // 2
        if getCount(MID):
            MIN = MID
        else:
            MAX = MID

    ans = 0
    count = 0

    for Ai in A_r:
        idx = bisect_left(A, MIN - Ai)
        ans += Ai * (N - idx) + B[N - idx]
        count += N - idx
    print(ans - (count - M) * MIN)


if __name__ == "__main__":
    resolve()

Solution (Python):

from bisect import bisect_left, bisect_right
from itertools import accumulate
def resolve():
    def shake_cnt(x):
        cnt = 0
        pos = 0
        for i in range(N):
            pos = bisect_left(A, x - A[i])
            cnt += N - pos
        return cnt < M 

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

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

    acc = [0] + list(accumulate(A))
    ans = 0
    for i in range(N):
        pos = bisect_right(A, ng - A[i])
        cnt = N - pos
        ans += cnt * A[i] + (acc[N] - acc[pos])
        M -= cnt
    ans += M * ng
    print(ans)


if __name__ == "__main__":
    resolve()

Solution (Python):
def resolve():
    import sys

    def LI():
        return list(map(int, sys.stdin.readline().rstrip().split()))  

    def S():
        return sys.stdin.readline().rstrip()

    N, M = map(int, S().split())
    A = LI()
    A.sort()

    def f(n,x):  
        left = -1
        right = N
        while left + 1 < right:
            mid = (left + right) // 2
            if A[n] + A[mid] < x:
                left = mid
            else:
                right = mid
        return right

    low = 2 * A[0]
    high = 2 * A[-1] + 1

    while low + 1 < high:
        mid = (low + high) // 2
        if sum(N - f(i,mid) for i in range(N)) >= M:
            low = mid
        else:
            high = mid

    a = low  

    from itertools import accumulate

    B = [0] + list(accumulate(A))

    r = 0
    ans = 0
    for i in range(N):  
        b = f(i,a+1)
        r += N - b
        ans += B[N] - B[b]
        ans += A[i] * (N - b)

    ans += (M - r) * a  

    print(ans)

if __name__ == '__main__':
    resolve()