Calculate the Euler-Mascheroni constant without the math module The below code was written to generate γ, for educational purposes.

Single threaded, no functional zeroes required, no binary splitting(which can all be used to compute competitively like y-cruncher, that version in works).
Uses the Arithmetic Geometric Mean to calculate large logarithms quickly.
Uses decimal module for precision management. 

I've computed 3000 digits in a few hours with it, and 200 in about a minute. I am happy calculating.

import decimal

D = decimal.Decimal

def agm(a, b):      #Arithmetic Geometric Mean
    a, b = D(a),D(b)
    for x in range(prec):
        a, b = (a + b) / 2, (a * b).sqrt()
    return a

def pi_agm():       #Pi via AGM and lemniscate
    a, b, t, p, pi, k = 1, D(2).sqrt()/2, 1/D(2), 2, 0, 0
    while 1:
        an    = (a + b) / 2
        b     = (a * b).sqrt()
        t    -= p * (a - an)**2
        a, p  = an, 2**(k+2)
        piold = pi
        pi    = (a + b) * (a + b) / (2*t)
        k    += 1
        if pi == piold:
            break
    return pi

def factorial(x):       #factorial fast loop
    x = int(x)
    factorial = D(1)
    for i in range(1, x+1):
        factorial *= i
    return factorial

def lntwo():        #Fast converging Ln 2
    logsum, logold, n = D(0), D(0), 0
    while 1:
        logold = logsum
        logsum += D(1/((D(961**n))*((2*n)+1)))
        n += 1
        if logsum == logold:
            logsum1 = (D(14)/D(31))*logsum
            break
    logsum, logold, n = D(0), D(0), 0
    while 1:
        logold = logsum
        logsum += D(1/((D(25921**n))*((2*n)+1)))
        n += 1
        if logsum == logold:
            logsum2 = (D(6)/D(161))*logsum
            break
    logsum, logold, n = D(0), D(0), 0
    while 1:
        logold = logsum
        logsum += D(1/((D(2401**n))*((2*n)+1)))
        n += 1
        if logsum == logold:
            logsum3 = (D(10)/D(49))*logsum
            break
    ln2 = logsum1 + logsum2 + logsum3
    return ln2

def lnagm(x):   #Natural log via AGM,
    try:
        if int(x) == 1:
            return 0
        if int(x) == 2:
            return lntwo()
    except:
        pass
    m = prec*2
    ln2 = lntwo()
    decimal.getcontext().prec = m
    pi = D(pi_agm())
    twoprec = D(2**(2-D(m)))/D(x)
    den = agm(1, twoprec)*2
    diff = m*ln2
    result = (D(pi/den) - D(diff))
    logr = D(str(result)[:m//2])
    decimal.getcontext().prec = prec
    return logr

def gamma():   #Compute Gamma from Digamma Expansion
    print('Computing Gamma!')
    k = D(prec/2)
    print('Calculating Logarithms...')
    lnk = lnagm(k)
    logsum = D(0)
    upper = int((12*k)+2)
    print('Summing...')
    for r in range(1, upper):
        logsum += D((D(-1)**D(r-1))*D(k**D(r+1)))/D(factorial(r-1)*D(r+1))
        if r%1000==0:
            print(str((D(r)/D(upper))*100)[:5], '% ; Sum 1 of 2')
    logsum1 = D(0)
    print('...')
    for r in range(1, upper):
        logsum1 += D((D(-1)**D(r-1))*(k**D(r+1)))/D(factorial(r-1)*D(D(r+1)**2))
        if r%1000==0:
            print(str((D(r)/D(upper))*100)[:5], '% ; Sum 2 of 2')
    twofac = D(2)**(-k)
    gammac = str(D(1)-(lnk*logsum)+logsum1+twofac)
    return D(gammac[:int(prec//6.66)])

#Calling Gamma
prec = int(input('Precision for Gamma: '))*8
decimal.getcontext().prec = prec
gam = gamma()
print(gam)
Computation of logsum and logsum1 in gamma() are suboptimal. You do costly operations of raising to power, and recompute factorial on each iteration (the latter invokes the quadratic time complexity BTW). Notice that in the \$\sum \dfrac{(-1)^{r-1} k^{r+1}}{(r+1)(r-1)!}\$ a consecutive term can be expressed via the previous one, as\$T_{r+1} = -k\dfrac{r+1}{(r+2)r} T_n\$. Instead of computing each term from scratch, use this recurrence, and enjoy a significant performance boost.

Converting the summation into a Horner schema would likely improve the accuracy, that is to achieve the desired number of digits you'd need less number of terms.
Tree loops in lntwo cry to become a function.