lec23.mws

L-23 MCS 320 Wednesday 19 October 2005

> restart;

1. The assume facility

Maple knows the types of the objects.  With the assume we can add properties to the types.

> int(exp(a*t),t=0..infinity) assuming a > 0;

infinity

Depending on the sign of a, we get different values for the integral.

The assuming from above was a "local" assumption, it is only valid for the calculation of the integral.

> assume(a>0);

> a; # verify that there is an assumption on a

a

With the assume command we add a property to the variable a, which will be taken into account in all the calculations with a.

> about(a);

Originally a, renamed a~:
 is assumed to be: RealRange(Open(0),infinity)


With "additionally" we add assumptions:

> additionally(a < 4);

> about(a);

Originally a, renamed a~:
 is assumed to be: RealRange(Open(0),Open(4))


We see that the bounds of the interval (0,4) are not included.

> additionally(a <= 2);

> about(a);

Originally a, renamed a~:
 is assumed to be: RealRange(Open(0),2)


Now a lies in the half-open, half-closed interval (0,2] or ]0,2].

We can use these properties when we execute the bisection algorithm.

> plot(cos(x)-x^2,x=0..Pi);

[Plot]

We see from the plot that this function has a root between 0 and Pi.

> f := x -> cos(x)-x^2:

We will use r to denote a root of f, we know it lies between 0 and Pi:

> assume(r>0,r<Pi);

> about(r);

Originally r, renamed r~:
 is assumed to be: RealRange(Open(0),Open(Pi))


Exercise: add the attribute "is a root of cos(x) - x^2 = 0" to r.

The bisection algorithm will narrow the uncertainty interval for the root by taking the midpoint of the interval and then checking the sign of the function at the midpoint agains the signs of the end points.

> midpoint := (0+Pi)/2;

midpoint := 1/2*Pi

> is(f(midpoint)*f(0)<0);

true

> is(f(midpoint)*f(Pi)<0);

false

The use of the "is" command calculates with properties.  Without having to evaluate the product of the function values, Maple decides for us on the sign.

Here we keep the left bound on the uncertainty interval around the root of f, and the midpoint replaces the right bound of the interval:

> additionally(r<midpoint);

> about(r);

Originally r, renamed r~:
 is assumed to be: RealRange(Open(0),Open(1/2*Pi))


There is an analogy with the algebraic numbers.  We declared algebraic numbers as RootOf, as they were roots of polynomials.  Here we add properties to a variable and calculate with these properties, using the fact that we have an approximation for a root.

2. Simplification

We need assumptions to obtain simplifications of expressions.  We have seen this with regard to complex numbers.

> eln := expand(ln(x*y)) assuming x>0,y>0; # "local" assumption

eln := ln(x)+ln(y)

> combine(eln); # does not work because x and y can be anything

ln(x)+ln(y)

> assume(x>0,y>0); # make a "global" assumption

> combine(eln);

ln(x*y)

With the expand and combine we have verified that Maple knows the identity ln(x*y) = ln(x)+ln(y),

which is only valid for x and y both positive.

With the expand and combine we can transform general expressions.  What is the analogue to the "combine", or the opposite to "expand" for polynomials?  Remember the "factor".

In this lecture we deal with more general expressions than polynomials:

> x := 'x': y := 'y': # remove the assumptions

> expression := exp(x)*exp(y) + (sin(x))^2 + (cos(x))^2;

expression := exp(x)*exp(y)+sin(x)^2+cos(x)^2

> simplify(expression);

exp(x+y)+1

We can decide to simplify only the exponential:

> simplify(expression,exp);

exp(x+y)+sin(x)^2+cos(x)^2

> simplify(expression,trig); # or only simplify trigonometric

exp(x)*exp(y)+1

Some other arguments for the simplify are the "symbolic" and "delay":

> sqrtsquare := sqrt(x^2);

sqrtsquare := (x^2)^(1/2)

> simplify(sqrtsquare);

csgn(x)*x

The output of the simplify can be understood from knowing that the square root of a complex number consists of two values.  Depending on a complex sign function, we have either the positive or the negative branch of the double-valued square root function.  But that is rather cryptic, and we may prefer to simplify otherwise:

> simplify(sqrtsquare,symbolic); # sqrt is inverse of square

x

> simplify(sqrtsquare,delay); # do not simplify unless assumption are known

(x^2)^(1/2)

Maple has a whole set of rules to convert trigonometic expressions:

> trigsubs(cos(x+y));

[cos(x+y), cos(-x-y), 2*cos(1/2*x+1/2*y)^2-1, 1-2*sin(1/2*x+1/2*y)^2, cos(1/2*x+1/2*y)^2-sin(1/2*x+1/2*y)^2, 1/sec(x+y), 1/sec(-x-y), (1-tan(1/2*x+1/2*y)^2)/(1+tan(1/2*x+1/2*y)^2), 1/2*exp(I*(x+y))+1/...[cos(x+y), cos(-x-y), 2*cos(1/2*x+1/2*y)^2-1, 1-2*sin(1/2*x+1/2*y)^2, cos(1/2*x+1/2*y)^2-sin(1/2*x+1/2*y)^2, 1/sec(x+y), 1/sec(-x-y), (1-tan(1/2*x+1/2*y)^2)/(1+tan(1/2*x+1/2*y)^2), 1/2*exp(I*(x+y))+1/...

The output of trigsubs contains several alternatives to replace the cosine of a sum of two angles.

3. Side Relations

> simplify(x^3+y^2,{x^2+y^2 = 1});

y^2-x*y^2+x

This simplification is part of a general term rewriting mechanism.  We will illustrate this with Lagrange minimization.

> f := x^2 + 2*y*z*x - z^2: # target function

> g := x^2 + y^2 + z^2 - 1: # constraint

The problem is to find the extremal values (minima and maxima) of the function f, for points (x,y,z) satisfying the equation g(x,y,z) = 0.

> sys := {g,diff(g,x)-lambda*diff(f,x),diff(g,y)-lambda*diff(f,y),diff(g,z)-lambda*diff(f,z)};

sys := {x^2+y^2+z^2-1, 2*y-2*lambda*z*x, 2*x-lambda*(2*x+2*y*z), 2*z-lambda*(2*x*y-2*z)}

> gb := grobner[gbasis](sys,[x,y,z,lambda],plex);

gb := [6*x^2-2*lambda^3-2*lambda^2-lambda-1, 3*x*y-2*lambda^2*z-2*lambda*z, 3*z*x-2*lambda*y, -x+lambda*x+lambda*y*z, 3*y^2+2*lambda^2-2, -3*y+2*lambda^2*y, 2*lambda^3-2*lambda^2+6*z^2+lambda-1, -3*la...gb := [6*x^2-2*lambda^3-2*lambda^2-lambda-1, 3*x*y-2*lambda^2*z-2*lambda*z, 3*z*x-2*lambda*y, -x+lambda*x+lambda*y*z, 3*y^2+2*lambda^2-2, -3*y+2*lambda^2*y, 2*lambda^3-2*lambda^2+6*z^2+lambda-1, -3*la...gb := [6*x^2-2*lambda^3-2*lambda^2-lambda-1, 3*x*y-2*lambda^2*z-2*lambda*z, 3*z*x-2*lambda*y, -x+lambda*x+lambda*y*z, 3*y^2+2*lambda^2-2, -3*y+2*lambda^2*y, 2*lambda^3-2*lambda^2+6*z^2+lambda-1, -3*la...

The "gb" (the Groebner basis) is a simplification of the orginal problem.  We see that the last equation is an equation only in lambda.  This is always the case with a pure lexicographical ordering of the variables: we always obtain a triangular system:

> solve(gb[nops(gb)],lambda);

-1, 1, 1/2*6^(1/2), -1/2*6^(1/2)

We see that there are four candidate extremal values.  With backsubstition of these four values in the other equations, we obtain values for the corresponding coordinates of the extrema.