The important functions in my program are "matrixHN" and "matrixSN". I use the numerical integration method "NewtonCotesRule", because it gives the best results. Is there a better integration method?