I'm trying to use Gaussian quadrature to approximate the integral of a function. (More info here: ). The first function is on the interval [-1,1]. The second function is generalized to [a,b] by change of variable. The problem is that I keep getting the error "'numpy.ndarray' object is not callable". I assume (please correct me if I'm wrong) this means I've tried to call the arrays w and x as functions, but I'm not sure how to fix this.

This is the code

from __future__ import division
from pylab import *
from scipy.special.orthogonal import p_roots
def gauss1(f,n): [x,w] = p_roots(n+1) f = (1-x**2)**0.5 for i in range(n+1): G = sum(w[i]*f(x[i])) return G
def gauss(f,a,b,n): [x,w] = p_roots(n+1) f = (1-x**2)**0.5 for i in range(n+1): G = 0.5*(b-a)*sum(w[i]*f(0.5*(b-a)*x[i]+ 0.5*(b+a))) return G

These are the respective error messages

gauss1(f,4)
Traceback (most recent call last): File "<ipython-input-82-43c8ecf7334a>", line 1, in <module> gauss1(f,4) File "C:/Users/Me/Desktop/hw8.py", line 16, in gauss1 G = sum(w[i]*f(x[i]))
TypeError: 'numpy.ndarray' object is not callable
gauss(f,0,1,4)
Traceback (most recent call last): File "<ipython-input-83-5603d51e9206>", line 1, in <module> gauss(f,0,1,4) File "C:/Users/Me/Desktop/hw8.py", line 23, in gauss G = 0.5*(b-a)*sum(w[i]*f(0.5*(b-a)*x[i]+ 0.5*(b+a)))
TypeError: 'numpy.ndarray' object is not callable

4 Answers

As Will says you're getting confused between arrays and functions.

You need to define the function you want to integrate separately and pass it into gauss.

E.g.

def my_f(x): return 2*x**2 - 3*x +15
gauss(m_f,2,1,-1)

You also don't need to loop as numpy arrays are vectorized objects.

def gauss1(f,n): [x,w] = p_roots(n+1) G=sum(w*f(x)) return G
def gauss(f,n,a,b): [x,w] = p_roots(n+1) G=0.5*(b-a)*sum(w*f(0.5*(b-a)*x+0.5*(b+a))) return G

quadpy, a small project of mine, might help:

import numpy
import quadpy
def f(x): return numpy.exp(x)
scheme = quadpy.line_segment.gauss_legendre(10)
val = scheme.integrate(f, [0.0, 1.0])
print(val)

1.7182818284590464

There are many other quadrature schemes for 1D.

In gauss1(f,n), you are treating f as if it's a function when it is an array, since you are reassigning it;

def gauss1(f,n): [x,w] = p_roots(n+1) f = (1-x**2)**0.5 # This line is your problem. for i in range(n+1): G = sum(w[i]*f(x[i])) return G

You are doing something similar in the second function.

2

Example: solving using gaussian integral with n = 2 for integral 2+sinX with b = pi/2 and a = 0

import numpy as np
E = np.array([-0.774597, 0.000000, 0.774597])
A = np.array([0.555556, 0.888889, 0.555556])
def gauss(f, a, b, E, A): x = np.zeros(3) for i in range(3): x[i] = (b+a)/2 + (b-a)/2 *E[i] return (b-a)/2 * (A[0]*f(x[0]) + A[1]*f(x[1]) + A[2]*f(x[2]))
f = lambda x: 2 + np.sin(x)
a = 0.0; b =
areaGau = gauss(f, a, b, E, A)
print("Gaussian integral: ", areaGau)

Your Answer

Sign up or log in

Sign up using Google Sign up using Facebook Sign up using Email and Password Submit

Post as a guest

Post Your Answer Discard

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy