Numerical solution of a Gamma - integral equation using a higher order composite Newton-Cotes formulas

The paper aims at solving a complex equation with Gamma - integral. The solution is the infected size (p) at equilibrium. The approaches are both numerical and analytical methods. As a numerical method, the higher-order composite Newton-Cotes formula is developed and implemented. The results show that the infected size ( p ) increases along with the shape parameter ( k ). But the increase has two phases: an increasing rate phase and a decreasing rate phase; both phases can be explained by the instantaneous death rate characteristics of the Gamma distribution hazard function. As an analytical method, the Extreme Value Theory consolidates the numerical solutions of the infected size ( p ) when k ≥ 1 and provides a solution limit ( p = 1 − 1 2 R ) as k goes to +∞.