plot_trigonometric_interpolant.py

# function [done]  =  plot_trigonometric_interpolant (n, m, f)
# Accompanying program for the text
#
#    Classical and Modern Numerical Analysis:
#    Theory, Methods and Practice
# by Azmy S. Ackleh, Edward J. Allen,
#    R. Baker Kearfott, and Padmanabhan Seshaiyer
#
#    (Taylor and Francis / CRC Press, 2009)
#
# ok  =  plot_trigonometric_interpolant (n, m, f) plots a trigonometric
# interpolant to the function f over the interval [0, 2\pi], with
# 2n+1 interpolating points and m points in the plot.  The trigonometric
# polynomial is of the form
# a_0 + \sum_{k = 1}^{2n+1} (a_j \cos(jx) + b_j \sin(jx).
# For example
# ok  =  plot_trigonometric_interpolant(6, 1001, 'trig_interp_example')
# will plot the trigonometric interpolant to the function in
# trig_interp_example.m (supplied with these routines) [-\pi,\pi],
# using 13 interpolation points and 1001 points for the plot.

import numpy as np
import matplotlib.pyplot as plt

def plot_trigonometric_interpolant( n, m, f ):
        xi = np.linspace(0, 2*np.pi, 2*n+2)
        yi = f(xi)
        a0 = 0
        a0 += yi.sum()/(2*n+1)
        a = np.empty(n,dtype=complex)
        b = np.empty(n,dtype=complex)
        for j in range(n):
           s1 = (yi*np.exp(-1j*(j+1)*xi)).sum()/(2*n+1)
           s2 = (yi*np.exp( 1j*(j+1)*xi)).sum()/(2*n+1)
           a[j] = s1+s2
           b[j] = 1j*(s1-s2)

        x = np.linspace(0,2*np.pi,m)
        y = f(x)
        ya  =  trigonometric_interpolant_value (x, a0, a, b)
        plt.plot(x,y,'-',x,ya.real,'--')
        plt.title('Dashed curve is the trigonometric interpolant')
        return True

if __name__ == '__main__':

        from trigonometric_interpolant_value import trigonometric_interpolant_value

        def example(x): return np.sqrt(np.pi**2-(x-np.pi)**2) # from trig_interp_example.m

        plot_trigonometric_interpolant( 6, 1001, example )
        plt.show()