'''This file contains gaussian quadrature formula of order 2, 4 and 6
These functions are called by the script python13_1.py.
Input:
x... grid points of equidistant partition of interval.
fcn... the function to be integrated
Output:
approx... approximation to the integral by means of quadrature formula
'''
import numpy as np
def gauss_ord2(x, fcn):
n = len(x)-1 # number of subintervals of partition
h = x[1] - x[0] # length of subintervals (stepsize)
approx = 0
c1 = .5
for j in range(n):
Y1 = fcn(x[j] + h*c1)
approx += h*Y1
return approx
def gauss_ord4(x, fcn):
n = len(x)-1 # number of subintervals of partition
h = x[1] - x[0] # length of subintervals (stepsize)
approx = 0
c1 = .5 - np.sqrt(3)/6
c2 = .5 + np.sqrt(3)/6
for j in range(n):
Y1 = fcn(x[j] + h*c1)
Y2 = fcn(x[j] + h*c2)
approx += 0.5*h*(Y1 + Y2)
return approx
def gauss_ord6(x, fcn):
n = len(x)-1 # number of subintervals of partition
h = x[1] - x[0] # length of subintervals (stepsize)
approx = 0
c1 = .5 - np.sqrt(15)/10
c2 = .5
c3 = .5 + np.sqrt(15)/10
for j in range(n):
Y1 = fcn(x[j] + h*c1)
Y2 = fcn(x[j] + h*c2)
Y3 = fcn(x[j] + h*c3)
approx += h/18*(5*Y1 + 8*Y2 + 5*Y3)
return approx