Category: Integral Calculus

logabs, logarc: How to make integrate() return what you expect

logarc

In the back of my calculus book there is a table of famous integrals.  Here’s integral number 21 in that table:

20171004_113109

From Maxima integrate(), I get

integral21Maxima

 

What’s going on?

Both forms give a workable antiderivative for the original integrand:

integral21Maxima2

Furthermore, we believe that both forms are correct because of this helpful identity for hyperbolic sine:

 {\rm asinh}(z)=\ln(z+\sqrt{1+z^2}).  

Turns out (thanks to a Barton Willis for pointing me in the right direction) there’s a variable logarc that we can set to make Maxima return the logarithmic form instead of hyperbolic sine:

integral21Maxima3

I haven’t yet encountered cases where this would be a bad idea in general, but I’ll update this if I do.

logabs

In the first week of my differential equations course, we study methods of direct integration and separation of variables.  I like to emphasize that the absolute values can lend an extra degree of generality to solutions with antiderivatives of the form

\int \frac{1}{u}\;du = \ln |u|+ C.   

As an example, for the initial value problem

 y' = \frac{x}{1-x^2}    ,  y(0)=1    ,

it is convenient for treating all possible initial conditions (x \ne \pm 1) in one step to use the antiderivative

y = -\frac{1}{2} \ln | 1-x^2 | + C.   

However, Maxima omits the absolute values.

logabs1.PNG

For this case, we could consider only the needed interval -1<x<1, but still…

Turns out we can set the Maxima variable logabs to make integrate() include absolute values in cases like this:

logabs2.PNG

But then later in the course, I saw that logabs also impacts the Ordinary Differential Equation solver ode2().  I encountered an example for which Maxima, in particular solve() applied to expressions involving absolute value,  didn’t do what I wanted with logabs:true

For the logistic equation

\frac{dP}{dt} = kP\left ( 1-\frac{P}{P_c} \right )    ,  P(0)=P_0   

we expect that by separating variables we can obtain the solution

P(t) = \frac{P_0 P_c}{(P_c-P_0)e^{-kt}+ P_0}.   

Here’s what happens when we use ode2() with and without logabs:true:

logistice2

logistic3

 

 

Advertisement

Path Integrals

For my multivariable calculus class, I wanted an easy-to-call suite of symbolic integrators for path integrals of the form

\int_C f(x,y) dr,

\int_C {\bf F}(x,y)\cdot d{\bf r} = \int_C \langle P(x,y), Q(x,y) \rangle\cdot \langle dx, dy \rangle, or

\int_C {\bf F}(x,y,z)\cdot d{\bf r} = \int_C \langle P(x,y,z), Q(x,y,z), R(x,y,z) \rangle\cdot \langle dx, dy, dz \rangle.

My overarching design idea was that the input arguments needed to look the way they do when I teach the course:

  • a scalar field f(x,y):R^2 \rightarrow R or a vector field  {\bf F}(x,y): R^2 \rightarrow R^2 or {\bf F}(x,y,z): R^3 \rightarrow R^3
  • a curve C defined by a vector-valued function {\bf r}(t): R \rightarrow R^n, a\le t \le b where n=2,3 as appropriate.

It took me a while to work out how to evaluate the integrand along the path within my function.  Things that worked fine on the command line failed when embedded into a batch file to which I passed functions as arguments.  I ended up using subst, one variable at a time.  I’d like to be able to do this in a single command which can detect whether we’re in 2 or 3 dimensions so that I don’t need separate commands.

For now, here’s what I came up with along with some illustrative examples taken from Paul’s online math notes, that show how to call these new commands I, I2 and I3.

/* path integral of a scalar integrand f(x,y) on path r(t) in R^2, t from a to b */
 I(f,r,t,a,b):=block(
 [f1,f2,dr,Iout],
 f1:subst(x=r[1],f),
 f2:subst(y=r[2],f1),
 dr:sqrt(diff(r,t).diff(r,t)),
 Iout: integrate(f2*dr,t,a,b),
 Iout
 );
/* path integral of a vector integrand F(x,y) on path r(t) in R^2, t from a to b */
 I2(H,r,t,a,b):=block(
 [H1,H2,I],
 H1:subst(x=r[1],H),
 H2:subst(y=r[2],H1),
 I: integrate( H2.diff(r,t),t,a,b),
 I
 );
/* path integral of a vector integrand F(x,y,z) on path r(t) in R^3, t from a to b */
 I3(H,r,t,a,b):=block(
 [H1,H2,H3,I],
 H1:subst(x=r[1],H),
 H2:subst(y=r[2],H1),
 H3:subst(z=r[3],H2),
 I: integrate( H3.diff(r,t),t,a,b),
 I
 );

PathIntegrals

Here’s an update:  a related maxima function for evaluating a complex integral

\int_\Gamma f(z) dz

where f: C \rightarrow C and the curve \Gamma is given by r: R \rightarrow C.

/* path integral of a complex integrand f(z): C --> C, on path z(t): R --> C, t from a to b */
IC(f,r,t,a,b):=block(
[f1,dz,Iout],
f1:subst(z=r,f),
dz:diff(r,t),
Iout: integrate(f1*dz,t,a,b),
Iout
);

ComplexIntegral

Two Maxima Functions for Riemann Sums

Two early attempts at programming Maxima functions.  I’d love to hear your comments about how to make these work better.  Many thanks to the online community from whom I learned how to get this far!

The script is given at the bottom of this post.  Copy and paste into a file (in your working directory) called Riemann.mac, and then load into Maxima using

load(“Riemann.mac”);

Riemann

/* Two Maxima functions for introductory integral calculus
(c) 2016, themaximalist.org  */

/* RiemannSum
fn is the function to be integrated on the interval [a,b]
n rectangles
opt specifies:
0: left endpoints
1: midpoints
2: right endpoints
*/

RiemannSum(fn,a,b,n,opt):=
block([xx,s],
xx(i,n):= a+(i+opt/2)*(b-a)/n,
s: sum(ev(fn,x=xx(i,n)),i,0,n-1)*(b-a)/n,
float(s)

)$

/* RiemannRectangles to Draw rectangles
fn is the function to be integrated:
on the interval [a,b]
n rectangles
opt specifies:
0: left endpoints
1: midpoints
2: right endpoints
*/

load(draw)%

RiemannRectangles(fn,a,b,n,opt):=
block([xi,xx,rects,wxd],
xi(i,n):= a+i*(b-a)/n,
xx(i,n):= a+(i+opt/2)*(b-a)/n,
rects(n):=makelist(rectangle([xi(k,n),0], [xi(k+1,n),ev(fn,x=xx(k,n))]),k,0,n-1),
wxd(n):=apply(wxdraw2d,
append([xrange=[a-(b-a)/4, b+(b-a)/4],color=blue],
rects(n),
[transparent=true, explicit(fn,x,a, b)]
)
),
wxd(n)
)$