Sage MAxQ

emptwob=[]
def solve_pell(N,numT=100):
    cf=continued_fraction(sqrt( N ))
    for i in range (numT ):
        denom = cf.denominator (i )
        numer = cf.numerator (i )
        if numer^2-N*denom^2==1:
            return numer , denom
    return None , None
for v in range(2,10):
    for j in range(4,5):
        for h in range(1,j-1):
            if h==(j-1)/2:
                continue
            d=squarefree_part((h*(j^(v)-1)/(j-1))^2-1)
            g=solve_pell(d)[0]
            t=solve_pell(d)[1]
            kappaa=log(solve_pell(d)[0]+(solve_pell(d)[1])*sqrt(d),j)
            if g+t*sqrt(d)>j:

                if 6*10^(15)< denominator(continued_fraction(kappaa).convergent(49)):
                    emptwob.append(denominator(continued_fraction(kappaa).convergent(49)))
                else:
                    print(v,j,h)
            else:
                continue
print(max(emptwob))