import scipy.constants as c
import numpy as np
a=c.physical_constants['Bohr radius'][0]

# Put all other functions here

def radial_wave_func(n,l,r):
	# Write here
    if l == 0 and n == 1:
        return np.round(2.0*(1/a)**(3/2.0)*np.exp(-r/a),5)
    if l == 0 and n == 2:
        return np.round((1/(np.sqrt(2)))*(1/a)**(3/2.0)*(2-(r/a))*np.exp(-r/(2.0*a)),5)
    if l == 1 and n == 2:
        return np.round((1/(2.0*np.sqrt(6)))*(r/a)*np.exp(-r/(2.0*a)),5)
    if l == 0 and n == 3:
        return np.round((2/(81.0*np.sqrt(3)))*(1/a)**(3/2.0)*(2.0*(r/a)**2-18.0*(r/a)+27)*np.exp(-r/(3.0*a)),5)
    if l == 1 and n == 3:
        return np.round((4/(81.0*np.sqrt(6)))*(1/a)**(3/2.0)*(6.0*(r/a)-(r/a)**2)*np.exp(-r/(3.0*a)),5)
    if l == 2 and n == 3:
        return np.round((4/(81.0*np.sqrt(30)))*(1/a)**(3/2.0)*(r/a)**2*np.exp(-r/(3.0*a)),5)
    if l == 0 and n == 4:
        return np.round((1/384.0)*(1/a)**(3/2.0)*np.exp(-r/(4.0*a))*(96.0-(72.0*r/a)+((12.0*r**2)/a**2)-((r**3)/(2.0*a**3))),5)
    if l == 1 and n == 4:
        return np.round((1/(256.0*np.sqrt(15)))*(1/a)**(3/2.0)*np.exp(-r/(4.0*a))*(r/a)(80.0-((20.0*r)/a)+((r**2)/(a**2))),5)
    if l == 2 and n == 4:
        return np.round((1/(768.0*np.sqrt(5)))*(1/a)**(3/2.0)*np.exp(-r/(4.0*a))*(r/a)**2*(12.0-(r/a)),5)
    if l == 3 and n == 4:
        return np.round((1/(768.0*np.sqrt(35)))*(1/a)**(3/2.0)*(r/a)**3*np.exp(-r/(4.0*a)),5)