# How to minimize a numerical integration in python?

I need some help to minimize a numerical integration. It's about a classical problem in physics (hydrogen atom). It can be solved analytically but I need to solve it numerically in Python.

We have an electron's wave function $$\psi(r,a)=e^{-r/a}$$ and a Hamiltonian (in Ry units) $$H=-\nabla^2-\frac{2}{r}$$. The expected value of Energy is given by

$$\langle E(a)\rangle =\frac{\iiint \psi^{*}(r,a) H \psi(r,a) r^2 \sin(\theta) \mathrm dr \mathrm d\theta \mathrm d\phi}{\iiint |\psi(r,a)|^2 r^2 \sin(\theta) \mathrm dr \mathrm d\theta \mathrm d\phi}$$.

The integral it's over all space ($$0\leq r \leq \infty$$ ; $$0\leq \theta \leq \pi$$; $$0\leq \phi \leq 2\pi$$ )

Doing $$\langle E \rangle$$ numerically, how to find (also numerically) the value of $$a$$ which minimizes the value of the integral?

Edit 1- I know the value of '$$a$$' that minimizes the value of $$$$.

$$a=1$$ and $$=-1$$. But somehow I need to 'discovery' the value of $$'a'$$ that minimizes $$$$ numerically.

The code:

import numpy as np
import sympy as sp

# deffining the symbols
r, theta, phi, a = sp.symbols('r theta phi a')

# Deffining the function
f = sp.exp(-r/a)

# laplacian
laplacian = sp.diff(r**2 * sp.diff(f, r), r) / r**2 + \
sp.diff(sp.sin(theta) * sp.diff(f, theta), theta) /
(r**2 * sp.sin(theta)) + \
sp.diff(f, phi, phi) / (r**2 * sp.sin(theta)**2)

# changing the functions into a lambda functions
laplacian_numeric = sp.lambdify((r, theta, phi, a), laplacian,
'numpy')

f_numeric = sp.lambdify((r, theta, phi, a), f, 'numpy')

# deffining the kernel of integrations
def integrand(r, theta, phi, a):
jac = r**2 * sp.sin(theta)
return f_numeric(r, theta, phi, a) * laplacian_numeric(r,
theta, phi, a)*jac

def integrand2(r,theta,phi,a):
jac=r**2 * sp.sin(theta)
return 2*jac*f_numeric(r, theta, phi, a)**2 /r

def integrand3(r,theta,phi,a):
jac=r**2 * sp.sin(theta)
return jac*f_numeric(r, theta, phi, a)**2

# integrations limits
r_limits = (0, np.inf)
theta_limits = (0, np.pi)
phi_limits = (0, 2*np.pi)

# value of the parameter a
a_value = 1

# numerical integration
h1, _ = tplquad(integrand, *phi_limits, *theta_limits,*r_limits, args=(a_value,))

• If you want to minimize over only a single parameter (namely, $a$), why not just plot the quantity as a function of $a$? That will already give you intuition for what value of $a$ is approximately a minimum. Jan 22 at 18:34