Maximum-area rectangle inscribed in a normal distribution

Here we find the maximum-area rectangle inscribed within the normal distribution (the Gaussian distribution around 0, normalized to have an area of 1).

import numpy as np
import math
from matplotlib import pyplot as plt

Some simple functions. The first is the PDF (probability density function) of the normal distribution, which returns the height of the normal distribution at any given x-value. The second returns the area of a rectangle inscribed from -x to x within the normal distribution. (It is obvious that the optimal rectangle will be symmetric, so the only goal we have is to find the optimal xvalue for which the rectangle area is maximized.)

def normal(x):
    return ((2*math.pi)**-0.5) * np.exp(-(x**2)/2)

def rect_area(x):
    return 2 * x * normal(x)

We now take iterative snapshots of the function for the area of the rectangle, starting with the range between 0 and 10:

arange = np.arange(0, 10, 0.05)
areas = rect_area(arange)
plt.title('$x$-offset vs. rect area ($x \\in {}..{}$)'.format(round(min(arange),5), round(max(arange),5)))
plt.plot(arange, areas)

And now the area between 0.5 and 1.5:

arange = np.arange(0.5, 1.5, 0.01)
areas = rect_area(arange)
plt.title('$x$-offset vs. rect area ($x \\in {}..{}$)'.format(round(min(arange),5), round(max(arange),5)))
plt.plot(arange, areas)

It looks like the optimal xis 1.0. We continue to verify:

arange = np.arange(0.9, 1.1, 0.001)
areas = rect_area(arange)
plt.title('$x$-offset vs. rect area ($x \\in {}..{}$)'.format(round(min(arange),5), round(max(arange),5)))
plt.plot(arange, areas)

arange = np.arange(0.999, 1.001, 0.00001)
areas = rect_area(arange)
plt.title('$x$-offset vs. rect area ($x \\in {}..{}$)'.format(round(min(arange),5), round(max(arange),5)))
plt.plot(arange, areas)
plt.plot([1,1],[0.4839409,0.4839415])
plt.xticks(arange[::25], rotation=20)
plt.gca().get_yaxis().get_major_formatter().set_useOffset(False)

Clearly, the optimal x-offset is 1, corresponding to the rectangle spanning between -1 and 1 with a height of normal(1) \approx 0.24197072451914337

arange = np.arange(-3, 3, 0.05)
plt.title('maximal rectangle inscribed in normal distribution')
plt.plot(arange, normal(arange))
plt.plot([-1, -1, 1, 1], [0, normal(-1), normal(1), 0])
plt.plot([min(arange), max(arange)], [0, 0], color='black')