Trying to solve the following integral, with $n,m \in \mathbb{Z^+}$, $\alpha>1$, $0 < \epsilon < 1$, and $\Gamma(.)$ and $\Gamma(.,.)$ the gamma and incomplete gamma functions, respectively: $$I(\alpha,m,n) = \frac{(\alpha (n-1))^n }{\Gamma (n)-\Gamma \left(n,\frac{(n-1) \alpha }{\epsilon +1}\right)}\int_1^{1/\epsilon}e^{\frac{(y-1) (\alpha -\alpha n)}{y}}(y-1)^{n-1} y^{m-n-1} dy$$ My approach has been as follows, leading to expressing the integral with gamma functions, which causes precision problems. Since $(y-1)^{n-1}= \sum _{i=0}^{n-1} (-1)^i \binom{n-1}{i} y^{n-1-i}$, and since $$\int_1^{1/\epsilon } y^{m-2-i} e^{\frac{(y-1) }{y}(\alpha -\alpha n)} \, dy = \frac{e^{\alpha (1- n)} \left(\frac{1}{\alpha -\alpha n}\right)^{i-m} (\Gamma (i-m+1,\alpha -n \alpha )-\Gamma (i-m+1,-(n-1) \alpha \epsilon ))}{\alpha (n-1)},$$ we get as a possible solution a summation showing ratios of differences of various gamma functions:
$I(\alpha, m,n)=e^{\alpha (1- n)} (\alpha (n-1))^{n-1} $
$$\sum _{i=0}^{n-1} (-1)^i \binom{n-1}{i}\frac{ \left(\frac{1}{\alpha -\alpha n}\right)^{i-m} (\Gamma (i-m+1,\alpha -n \alpha )-\Gamma (i-m+1,-(n-1) \alpha \epsilon ))}{\Gamma (n)-\Gamma \left(n,\frac{(n-1) \alpha }{\epsilon +1}\right)}$$
Now the hitch. I am computing (for control) the expression of $I(\alpha,m, n)$ integral using high precision numerical integration. With n=10, the numerical integration matches the gamma function. Beyond, it, the gamma functions (using Wolfram's Mathematica) produce different results (and from familiarity with the problem, I know that the gamma is wrong compared to the numerical). Yet I need $n$ around $10^4$. Further, other special functions such as the exponential integral give even worse results.
The results are in the attached picture. With $\epsilon=\frac{1}{100}$, $I(\frac{3}{2},1,n)$:
The questions are:
1) Is there a possible closed form for $I$not entailing special functions?
2) Is there a better way to solve the integral by avoiding the summation?
With gratitude.
I don't think closed form is possible, but here is (possibly) a more controlled series formula for the function: $$ I(\alpha, m, n) = \sum_{j=0}^{\infty} {}^{m+j-1}C_{j} \,\,(\alpha(n-1))^{-j} \left[\frac{\Gamma (j+n) - \Gamma (j+n, \alpha(n-1)(1-\epsilon))}{\Gamma(n) - \Gamma(n, (n-1)\alpha/(1+\epsilon))}\right] $$
I got this using the variable substitution $z=(y-1)/y$ and carrying out a Taylor expansion in powers of $z$.