What is the position of the maximal value of this bell-shaped function?

Question

Consider the following function : $$\tag{1} f(v) = v^{\frac{d}{2}} \int_v^{\infty} u^{\alpha \,-\, \smash{\frac{d}{2}} \,-\, 1} \; e^{-\, \alpha \, u} \; du, $$ where $d \le 6$ and $\alpha > 0$ are two positive constants (parameters). Notice the lower limit of the integral : $v$ is a variable. The plot of this function (for $0 \le v < \infty$) shows an almost bell-shaped curve, so it has a single maximal value.

Now, I would like to find the position $v = v_0(\alpha, d)$ of the maximal value of this function. I need an analytical expression for $v_0(\alpha, d)$, probably an approximation. $$\tag{2} v_0(\alpha, d) \approx \; ? $$ From the graph of $f(v)$, I know that $v_0 \propto d$ (maybe with some exponent).

Take note that the derivative of function (1), set to 0, give this relation : $$\tag{3} f(v_0) \equiv f_{\text{max}} = \frac{2}{d} \; v_0^{\alpha} \; e^{-\, \alpha \, v_0}, $$ where $v_0 \equiv v_0(\alpha, d)$ is the position of the max value of (1).

Someone knows a method to find the function (2) ?

---

EDIT : Function (1) describes the "deformation" of a black body luminosity caused by the expansion of space in a cosmology model. From Wien's and Planck's laws, $\alpha$ should be around 3 (depending on the presence of gaz and dust). $d$ describes the kind of fluid contained in the cosmological model. We have $d = 3$ for a dust filled universe, and $d = 4$ for a radiation universe. $d = 0$ for an empty universe with a cosmological constant. Since $\frac{d}{2} + 1$ gives 2.5 for dust and 3 for radiation, I suspect that $\alpha$ should be close to $\frac{d}{2} + 1$ (while it is an independant parameter). In this special case, it is easy to explicitely evaluate the integral in (1) and we get the special case $$\tag{4} v_0(d) = \frac{d}{d + 2}, \quad \text{if $\alpha = \frac{d}{2} + 1$}. $$

Answer 1

This is a partial solution only and needs a correction for more precision in the result.

I wrote a Mathematica code to show the curve of function (1), with sliders to interactively change the two parameters $\alpha$ and $d$. When $\alpha = \frac{d}{2} + 1$, the exact solution to the problem is $$\tag{1} v_0(d) = \frac{d}{d + 2}, \quad \text{if $\alpha = \frac{d}{2} + 1$}. $$ When $\alpha \ne \frac{d}{2} + 1$, the Mathematica graphics shows that we still have $v_0 \approx \frac{d}{d \,+\, 2}$, for all values in the ranges $0 \le \alpha \le 10$ and $0 \le d \le 6$.

When $\alpha$ is very different than $\frac{d}{2} + 1$, formula (1) misses the true value by a small amount only. So it needs a small correction like $$\tag{2} v_0(\alpha, d) \approx \frac{d}{d + 2} \; g(\alpha, d), $$ where $g(\alpha, d) \approx 1 + \delta(\alpha, d)$, or maybe something like $$\tag{3} v_0(\alpha, d) \approx \frac{d}{d + 2} + h(\alpha, d), $$ where $h(\alpha, d)$ is a small correction.

Any idea about how to find the correction $\delta(\alpha, d)$ or $h(\alpha, d)$?

Here's the Mathematica code for you to see how good is the partial solution (1) above :

Clear["Global`*"]

function[v_, a_, d_] := v^(d/2) NIntegrate[u^(a - d/2 - 1) Exp[-a u], {u, v, Infinity}]

Xpic[a_, d_] := d/(d + 2)
Ypic[a_, d_] := function[d/(d + 2), a, d]

Pic[a_, d_] := Graphics@{Black, PointSize -> 0.01, Point[{Xpic[a, d], Ypic[a, d]}]}

Manipulate[
    Show[{Plot[{function[v, a, d]}, {v, 0, 4}, PlotStyle -> Directive[Red, Thick]], Pic[a, d]},
    Frame -> True,
    AspectRatio -> 1,
    PlotRange -> Automatic,
    GridLines ->  Automatic,
    GridLinesStyle -> Directive[Dashed, GrayLevel[0.75]],
    ImageSize -> 500
    ],
    {{a, 3, Style["\[Alpha]", 10]}, 0.5, 10, 0.01, ImageSize -> Large},
    {{d, 3, Style["d", 10]}, 0, 6, 0.01, ImageSize -> Large},
    ControlPlacement -> Bottom,
    FrameMargins -> None
]

Answer 2

Taking logarithmic derivatives, the maximum is attained when $$ \frac{d}{2v} - \frac{\mathrm{e}^{-v \alpha} \alpha (v \alpha)^{\alpha - (d/2 + 1)}}{\Gamma(\alpha - d/2, v\alpha)} = 0 $$

The form of the first summand suggests finding the (Laurent) series expansion of the second term around $v = \infty$. (I.e., find a Taylor series in powers of $1/v$.) This series is $$ \alpha - \frac{\alpha - (d/2-1)}{v} + \frac{2\alpha -(d+2)}{2 \alpha v^2} + ... $$

Solving $\dfrac{d}{2v} - \left( \alpha - \frac{\alpha - (d/2-1)}{v} \right) = 0$ for $v$, we get $\frac{\alpha - 1}{\alpha}$. If we substitute $\frac{d}{2}+1$ for $\alpha$ in this, we get $\frac{d}{d+2}$, as expected from your prior work.

Solving $\dfrac{d}{2v} - \left( \alpha - \frac{\alpha - (d/2-1)}{v} + \frac{2\alpha -(d+2)}{2 \alpha v^2} \right) = 0$ for $v$, we get two roots: $$ v_\pm = \frac{\alpha - 1 \pm \sqrt{2d + (\alpha - 5)(\alpha-1)}}{2\alpha} \text{.} $$ If we substitute $\frac{d}{2}+1$ for $\alpha$ in these, we get $0$ for $v_-$ (corresponding to the zero at $v=0$ for certain parameter values) and $\frac{d}{d+2}$ for $v_+$, as expected. This suggests $v_+$ is an improved estimate for your $v_0(\alpha,d)$.

(Note: I haven't graphed, bounded, or sanity checked the above, and it looks like I won't have time for several days. I wouldn't be very surprised to find I have a sign error (or worse) in transcription. I'm considering making this Community Wiki so others can improve it while I'm otherwise engaged. If someone finds a meaningful improvement, correction, or unrecoverable defect, say so in comments and I'll flip the switch (or withdraw the answer).)