import scipy.integrate as integrate
import numpy as np


def Debye(x):
    return (((x**4) * np.e**x) / (np.e**x - 1)**2)

integral_value = integrate.quad(Debye, 0, O_D/300)