# Calculus
#
#  Derivatives To compute the derivative of the expression
#
#                            3
#                           x  - 2 x + 9
#  use

diff( x^3 - 2*x + 9, x);

#  To compute the second derivative, write the code

diff( x^3 - 2*x + 9, x, x );

#  or alternatively,

diff( x^3 - 2*x + 9, x$2 ); # This works because Maple translates the expression x$2 into the
#  sequence x, x.  By analogy, x\$3 would give the third derivative.  Thus
#  one can easily compute derivatives of any order.  Now suppose that we
#  are given a function g defined by

g := x -> x^3 - 2*x + 9;

#  It seems natural to use

diff( g, x );

#  to get the derivative of g. However, Maple expects an expression in x
#  and so interprets g as a constant, giving the wrong result.  The
#  command

diff( g(x), x );

# which uses the expression g(x), works correctly.  The subtlety here is
# an important one: diff operates on expressions, not on functions - g is
# a function while g(x) is an expression.  To define the derivative of a
# function, use Maple's D operator:

dg := D(g);

#  The result is a function.  You can work with it just as you worked with
#  g. Thus you can compute function values and make plots:

dg(1);
plot( { g(x), dg(x) }, x = -5..5,
title = "g(x) = x^3 - 2x + 9 and its derivative" );

# Partial Derivatives

#  To compute partial derivatives:

q := sin(x*y);     # expression with two variables
diff( q, x );      # partial with respect to x
diff( q, x, y );   # compute d/dy of dq/dx

#  As in the one variable case, there is an operator for computing
#  derivatives of functions (as opposed to expressions):

k := (x,y) -> cos(x) + sin(y);
D[1](k);        # partial with respect to x
D[2](k);        # partial with repsect to y
D[1,1](k);      # second partial with respect to x
D[1,2](k);      # partial with respect to y, then x
D[1](D[2](k));  # same as above;

# Integrals

#  To compute integrals, use int.  The indefinite integral
#  (antiderivative)
#                              /
#                             |   3
#                             |  x  dx
#                             |
#                            /
#
#  is given by int( x^3, x). The following examples illustrate maple
#  engine integration by parts, substitution, and partial fractions:

int( 1/x, x );
int( x*sin(x), x );
int( sin(3*x + 1), x );
int( x/ (x^2 - 5*x + 4), x );

#  Nonetheless, Maple can't do everything:

int( sin( sqrt(1-x^3)), x );

#  The last response echoed back the indefinite integral, a signal that
#  the Maple engine failed to find an antiderivative for
#
#                           /    /     3\\
#                        sin\sqrt\1 - x //
#
#  In fact, it can be proved that no such antiderivative exists. There are
#  cases when the engine fails but there is an antiderivative. And cases
#  when the engine gets the wrong answer. The lesson: don't trust maple;
#  use maple to assist solving a problem, do not use maple as an
#  authority. To compute definite integrals like
#
#                              /1
#                             |    3
#                             |   x  dx
#                             |
#                            /0
#   write the maple code

int( x^3, x = 0..1 ).

#   Note that the only difference is that we give an interval of
#   integration.

# Numerical Integration.

#                       /
#                      |     /    /     3\\
#                      |  sin\sqrt\1 - x // dx
#                      |
#                     /
#
#  which Maple did not evaluate symbolically.  Maple finds a numerical
#  value for the definite integral by this syntax:

int( sin( sqrt(1 - x^3) ), x = 0..1 );
evalf(%);      # Force numerical evaluation.
# Works only if there are no undefined symbols.
evalf( Int( sin( sqrt(1 - x^3) ), x = 0..1 ) );
# The inert integral Int(...) skips symbolic evaluation.

#  Another approach to numerical integration is to use the student
#  package.

with( student ):            # load the package

j := sin( sqrt(1 - x^3) );  # define the integrand

trapezoid( j, x = 0..1 );   # apply trapezoid rule

evalf(%);                   # put in decimal form

#  By default the trapezoid command approximates the area under the graph
#  of
#
#                           /    /     3\\
#                        sin\sqrt\1 - x //
#
#  with four trapezoidal panels. For greater accuracy use more panels,
#  i.e., a finer subdivision of the interval of integration:

panels:=10: evalf( trapezoid( j, x = 0..1, panels) );

#  Better yet, use a more sophisticated numerical method like Simpson's
#  rule:

evalf( simpson( j, x = 0..1 ));

#  Only an even number of subdivisions is allowed, as in simpson( j, x =
#  0..1, 10 ).

#  The student package is well worth exploring.  Among other things it has
#  tools for displaying figures which explain the meaning of integration:

leftbox( j, x = 0..1, 10 );

#  The area of the figure displayed by leftbox is computed by leftsum:

evalf( leftsum( j, x = 0..1, 10 ));

#  One can also experiment with rightbox and middlebox and their companion
#  functions rightsum and middlesum.

# Finding Exact Areas Using Limits
#
# Example: Determine the exact area of the region bounded by the curve y =
# 4 - x^2, the x-axis, and the vertical lines x= -2 and x=2.

y:=4-x^2; # Define curve expression
unassign('n'); # Restore symbol
ans1:=rightsum(y,x=-2..2,n); # n=panel count
ans2:= value(ans1); # Simplify sums
area:=limit(ans2,n=infinity);

# Multiple Integrals
#
#  Let R be the rectangular region defined by x between 0 and 1, and y
#  between 0 and 1. To integrate
#
#                        /   /
#                       |   |   2    2
#                       |   |  x  + y  dx dy
#                       |   |
#                      /   /
#
#  over R in Maple, we compute the repeated integral.  Here is one way to
#  do this:

int( x^2 + y^2, x = 0..1 );
int( %, y = 0..1);

#  The first command integrates with respect to the x variable, producing
#  an expression in y . The second command integrates the result with
#  respect to y to give a number.

# Other Calculus Tools
#
#  Limits:

g := x -> (x^3 - 2*x + 9)/(2*x^3 + x - 3);
limit( g(x), x = infinity );
limit( sin(x)/x, x = infinity );
limit( sin(x)/x, x = 0 );

#  Taylor expansions and sums:

taylor( exp(x), x = 0, 4 ); # expansion around x = 0
sum( i^2, i = 1..100 );    # be sure i is unassigned
sum( x^n, n = 5..10 );     # needs n to be unassigned
sum( 1/d^5, d=1..infinity );  evalf(%);