How wolfram alpha calculates something like this?
Brute force, or has some kind of formula?
An integral representation derived by Euler gives H_n = Integrate[(1-x^n)/(1-x),{x,0,1)]
The limit n → infinity (H_n - log n) = gamma, the Euler-Mascheroni constant 0.5772…