import numpy as np

a = 1
n_vals = [n for n in range(1, 7)]

N = 1000
dx = 1 / N
x_vals = np.linspace(-a / 2, a / 2, N)

psi_evens = []
psi_odds = []
for n in n_vals:
    if n % 2 == 0:
        psi_evens.append(np.sqrt(2 / a) * np.sin(n * np.pi * x_vals / a))
    else:
        psi_odds.append(np.sqrt(2 / a) * np.cos(n * np.pi * x_vals / a))

for i in range(min(len(psi_evens), len(psi_odds))):
    print(f'psi_even_{i + 1} * psi_odd_{i + 1} = {np.dot(psi_evens[i], psi_odds[i]) * dx}')