Select Git revision
test_quad_transformed.py
nguyed99 authored
Code owners
Assign users and groups as approvers for specific file changes. Learn more.
test_quad_transformed.py 1.11 KiB
import numpy as np
import quadrature
import polynomial_new as poly
""" Settings: """
# Parameters of the quadrature rule on reference interval:
x = np.array([-1.0, 0.0, 1.0])
w = np.array([1.0/3, 4.0/3, 1.0/3])
# Define linear transformation of [-1, 1] to [1, 5]:
a = 1.0
b = 5.0
T = lambda x: 2 * x + 3
dT = 2.0
# Number of tests:
ntest = 10
""" Function to integrate polynomial: """
def int_poly(P, a, b):
# Coefficients of the primitive function:
pc = P.coeff / np.arange(1, P.deg+2, dtype=float)
# Evaluate powers at both boundary points:
bx = np.array([b**i for i in range(1, P.deg+2)])
ax = np.array([a**i for i in range(1, P.deg+2)])
return (pc @ bx - pc @ ax)
""" Main Test: """
S = quadrature.LinearTransformed(x, w)
for ii in range(ntest):
# Generate random third order polynomial:
#c = np.array([1.0, 0.0, 0.0, 1.0])
c = np.random.randn(4)
print("Test %d: coeffs: "%ii, c)
P = poly.polynomial(c)
# Compute integral:
I = int_poly(P, a, b)
# Approximation with Simpson rule:
ISimp = S.integrate(P, T, dT)
print("Test %d: diff = "%ii, np.abs(I - ISimp))