Using Trapz for Integral Evaluation

Hi,

I have an expression I want to integrate with respect to THETA between 0 and 2pi.

The expression to be integrated is: sin^2THETA / r^3 Where: r = sqrt(2*cosTHETA)

So this is my attempt to evaluate this integral.

    theta = 0:(2*pi/70):2*pi;
    r = sqrt(2*cos(theta));
    f = ((sin(theta)).^2) ./ ((r).^3);
    I = trapz(theta,f);

Is this correct??

Thanks for any advice/suggestion.

Answer 1

[Roger Stafford]

This integral is divergent, since theta crosses over 1/2*pi and 3/2*pi where the cos(theta) becomes zero. Having the cube of r in the denominator ensures that the integral will diverge to infinity at these places. This is therefore a decidedly inappropriate integration, not only for 'trapz', but for the other numerical matlab quadrature routines such as 'quad' and 'integral'. Moreover, over the range of theta from 1/2*pi to 3/2*pi the integrand becomes imaginary where cos(theta) is negative. Are you sure this is the integral you wish to solve? Perhaps if you tell us how you derived it, perhaps we could find a flaw in the reasoning.

[Walter Roberson]

The cases you list are isomorphic to the planar graph interpretation. Integration in polar coordinates is not isomorphic to those.

If S's question is whether the code was syntactically valid, that would be one matter. But that is not what S asked: S asked if the code was correct for finding that particular integral. Of course it could be that S just chose a very unfortunate set of names for a planar integral, but if it is intended to be an integral in polar coordinates, then I see no reason yet why trapz() would be appropriate.

[John D'Errico]

Sorry Walter, but this is a case where you are simply wrong. As long as the integral is formulated to integrate what was desired, trapz will be correct, to the limits of trapezoidal rule accuracy. You are still confused by your geometric interpretation of trapz, viewing it in Cartesian coordinates. But the fact is, had S called x the variable of integration, then you would have had no issue.

[Roger Stafford]

Well said, John. Thanks.

Answer 2

[Walter Roberson]

It does not look correct to me. You are attempting to do a polar integration with a technique that assumes rectangular coordinates.

[Roger Stafford]

@Walter: I see no indication of rectangular coordinates here. In theory it could be a perfectly valid integration of the given function of theta, except for the inconvenient fact that it is highly divergent, not to mention, complex-valued.

[Walter Roberson]

trapazoid method assumes that trapazoids can be built with right sides and only the top being angled. But when you are working in polar coordinates, the equi-theta lines are not right sides relative to the base -- indeed the base would be 0.

[Roger Stafford]

No, I very much disagree, Walter; 'trapz' makes no assumption about cartesian versus polar coordinates or any other coordinate system. It is simply a numerical approximation to an integral with respect to one variable, given in discrete form, for an integrand function, also given in discrete form, for the other variable. Trapezoids are merely a geometric visualization of this approximation method and have nothing to do with the actual significance of the variables involved.

[Walter Roberson]

http://en.wikipedia.org/wiki/Trapezoidal_rule

The trapezoidal rule works by approximating the region under the graph of the function f(x) as a trapezoid and calculating its area.

However it is obvious that the expression to be integrated is not in cartesian coordinates with variable names that just "happen" to be "r" and "theta": it is obvious that an angular integration is wanted, that the expression is a polar coordinate expression. And for such an expression, a change in the variable being integrated over, theta, does not correspond to a trapazoid. The approximation (b-a) * (f(a)+f(b))/2 does not apply.

The approximation formula, in polar terms, is like multiplying theta times the average r. But the area formula for sectors is theta/2 * r^2 and so in a/b terms would look like (b-a) * (f(a)+f(b))^2/2. So coordinate system is important.

[John D'Errico]

Walter - I think you are confusing the issue of converting an integral in cartesian form to an integral in polar form, with an integral already in polar form. Thus, you cannot simply append dr*dtheta to a function over a cartesian domain and perform an integral in polar coordinates.

However, IF the integral is properly formulated, with an appropriate differential element, then an integral over theta effectively DOES have parallel sides. Each trapezoidal side is indeed parallel in theta. Again, what matters is if the differential element is chosen properly, thus the transformation from cartesian to polar coordinates was done properly.

Once an integral is properly formulated as an integral, you could as easily call that variable x instead of theta, with limits of [0,2pi].

[Roger Stafford]

No, Walter, you are just flat out wrong on this issue! There is no other way to put it. The trapezoidal rule is so named because, if the two variables involved, namely the integrand and the variable of integration, are depicted in a graph, then the approximation that is made is that of a trapezoid approximating a section of an area beneath a curve. However, it does not follow that the actual interpretation of these variables must necessarily adhere to that graphical interpretation. The variable of integration might be, for example, time, and the other a varying current through a wire, with the integral therefore representing the accumulated charge on some condenser. These have nothing whatever to do with geometrical trapezoids. However, it is convenient for us mere mortals to think of that charge as the area under a curve whose height is the varying current level.

Just because "S" chose to name his variables r and theta does not mean that his integral is to be interpreted as having any relationship to polar coordinates, or even if they are polar coordinates, it does not mean that his integral referred to an area in the polar coordinate sense. He might have asked if 'trapz' was correct to use with the integral of sin(theta)^2*r^3, with r = 2*cos(theta), taken with respect to theta, or he might have asked for the integral of sin(x)^2*y^3 where y = 2*cos(x) taken with respect to x, and either way would have been a perfectly valid use of 'trapz'. 'Trapz' does not care about the names or the significances of the variables it uses in its computations. It just follows certain arithmetic rules to obtain a numeric answer which is designed to be an approximation to a desired definite integral value.