Maple
Son bon usage en mathématiques

Philippe Dumas, Xavier Gourdon

Springer-Verlag, 1997

Code Maple du Chapitre 2, Livre premier.



# section 1, page 67.

f:=(x^2-y^2)/(x-y)-(x+y);

normal(f);

r:=(1-2*sin(theta)^2)/(2*cos(theta)^2-1)-1;

normal(r);

combine(r,trig);

convert(r,tan);

normal(");

# exercice 2.1, page 68.
f:=expand((x-2*sin(3*theta))^3-(x+2*cos(2*theta))^3);
# fin exercice 2.1.

# exercice 2.2, page 69.
normal((x^2-1)/((x+2)*(x-2)));
normal((x-1)*(x+1)/(x^2-4));
# fin exercice 2.2.

# exercice 2.4, page 69.
A:=array(1..2,1..2,[[cos(theta),-sin(theta)],
                    [sin(theta),cos(theta)]]);
B:=subs(x=A,x^2-2*x*cos(theta)+1);
evalm(B);
# fin exercice 2.4.

# exercice 2.5, page 70.
alpha:=(1+3^(1/2)+(12)^(1/4))/2:
kappa[1]:=(alpha^3+1)/alpha/(alpha-1);
kappa[2]:=(9+6*3^(1/2))^(1/2);
# fin exercice 2.5.

f:=proc(x) 1/x end:

f1:=proc(x) diff(f(x),x) end:

f1(x);

f1(2);

f1:=proc(x) local w; subs(w=x,diff(f(w),w)) end:

g:=proc(x)
  if x<0 then x-1
  else 1
  fi
end: # g

diff(g(x),x);

diff('g'(x),x);

g:=proc(x)
  if type(x,numeric) then
    if x<0 then x-1
    else 1
    fi    
  else 'g'(x)
  fi
end: # g

diff(g(x),x);

f:=1/x;
f1:=diff(f,x);

subs(x=2,f1);

g:=piecewise(x<0,x-1,1);

g1:=diff(g,x);

plot(g,x=-1..1);

# section 2, page 73.


_Z:=1:
int(x,x);

p:=x^5-5*x+2:
lambda:=-0.03: 
f:=x+lambda*p:
Digits:=20:
w:=-1.7: 
for n to 20 do w:=subs(x=w,f) od:
w;

equ:=4*x*diff(y(x),x$2)+2*diff(y(x),x)-y(x);

sol:=dsolve(equ,y(x));

Y[1,0]:=subs(_C1=1,_C2=0,op(2,sol));

Y[0,1]:=subs(_C1=0,_C2=1,op(2,sol));

equ:=4*x*diff(y(x),x$2)+10*diff(y(x),x)-y(x);

restart;

''diff(sin(x)/x,x)'';

";

";

for n to 100 do od:
n:='n';

nextpartition:=proc(lambda)
  local s,r,c,d,alpha,i;
  s:=nops(lambda);
  for r from 0 to s-1 while lambda[s-r]=1 do od;
  if r=s then FAIL
  else
    c:=lambda[s-r];
    alpha:=iquo(c+r,c-1,d);
    if d=0 then [seq(lambda[i],i=1..s-r-1),seq(c-1,i=1..alpha)]
    else [seq(lambda[i],i=1..s-r-1),seq(c-1,i=1..alpha),d]
    fi
  fi
end: # nextpartition

N:=10:
lambda:=[N]:                                                   #1
lambdaseq:=[N]:
for k while lambda<>FAIL do 
  lambda:=nextpartition(lambda);
  lambdaseq:=lambdaseq,lambda; 
od:
lambda:=[N]:                                                   #2
lambdatable[1]:=[N]:
for k from 2 while lambda<>FAIL do 
  lambda:=nextpartition(lambda);
  lambdatable[k]:=lambda; 
od:
k:=k-2:
seq(lambdatable[i],i=1..k);

for i from 0 to 5 do
  x[i]:=cos(t[i]);y[i]:=sin(t[i]) 
od:
for i from 0 to 5 do
  equD[i]:=collect(linalg[det](
                        array(1..3,1..3,[[x[i],x[i+1 mod 6],X],
                                         [y[i],y[i+1 mod 6],Y],
                                         [1,1,1]])),
                       {X,Y},distributed,readlib(`combine/trig`))
od:
for i from 0 to 2 do
  M[i]:=combine(subs(solve({equD[i],equD[i+3 mod 6]},{X,Y}),
                                                     [X,Y]),trig)
od:
combine(linalg[det](array(1..3,1..3,
                   [seq([seq(M[i][j],i=0..2)],j=1..2),[1,1,1]])),
                                                           trig);

quantity[1]:=int(x^2/sqrt(x^2-x+1),x);

quantity[2]:=convert(quantity[1],ln);

quantity[3]:=subs(sqrt(3*x^2-3*x+3)=sqrt(3)*sqrt(x^2-x+1),
                                    sqrt(x^2-x+1)=y,quantity[2]);

quantity[4]:=collect(quantity[3],y);

quantity[5]:=map(factor,quantity[4]);

cp:=proc(pentagon::list(list(numeric)))          # main procedure
  if nops(pentagon) <> 5 or
     map(nops,{op(pentagon)})<>{2} then
     ERROR(`expected a list of 5 pairs, but received`,pentagon)
  fi;
  convexprop:=proc(pentagon)                     # auxiliary 
    cvxhull:=simplex[convexhull](pentagon);      # procedure
    if {op(cvxhull)}<>{op(pentagon)}             # convexprop
    then false
    else 
      for i to 5 do
        member(pentagon[i],cvxhull,'index');
        permutation[i]:=index;
      od;
      convert([seq(permutation[(i mod 5)+1]
                        =(permutation[i] mod 5)+1,i=1..5)],`and`)
      or       
      convert([seq(permutation[(i mod 5)+1]
                     =(permutation[i]-2  mod 5)+1,i=1..5)],`and`)
    fi
  end; # convexprop                      
  transfpentagon:=proc(pentagon)                # auxiliary 
    [seq([seq((pentagon[j][coord]+              # procedure
               pentagon[(j mod 5)+1][coord])/2, # transfpentagon
              coord=1..2)],j=1..5)]
  end; # transfpentagon               
  localpentagon:=pentagon;                      # main part
  for n from 0 to nmax 
  while not convexprop(localpentagon) do
    localpentagon:=transfpentagon(localpentagon);

  od;
  if n> nmax then FAIL else n fi
end: # cp 

# exercice 2.7, page 81.
transfpentagon:=proc(pentagon::list(list(numeric)))
  local j,coord;
  if nops(pentagon) <> 5 or
     map(nops,{op(pentagon)})<>{2} then
     ERROR(`expected a list of 5 pairs, but received`,pentagon)
  fi;
  [seq([seq((pentagon[j][coord]+pentagon[(j mod 5)+1][coord])/2,
                                            coord=1..2)],j=1..5)]
end: # transfpentagon
# fin exercice 2.7.

xrandom:=rand(0..100):
yrandom:=rand(0..100):
seq([xrandom(),yrandom()],i=1..5);

random:=rand(0..100);
seq([random(),random()],i=1..5);

# exercice 2.8, page 82.
_seed:=0;
r:=rand(0..10^6-1);
seq(r(),i=1..10);
# fin exercice 2.8.

# exercice 2.9, page 83.
solve(sin(x)=cos(x),x);
_EnvAllSolutions:=true:
solve(sin(x)=cos(x),x);
# fin exercice 2.9.

pattern:=proc(n::posint)
  local l,k;
  l:=length(n);
  Digits:=2*l;               # print(`inside, Digits = `,Digits);
  for k from 0 while 
     iquo(op(1,2.^k),10^l)<>n do od;
  k
end: # pattern

Digits:=20:
pattern(100);
print(`outside, Digits = `,Digits);

# exercice 2.11, page 84.
fib_exact:=proc(n::nonnegint)
    local oldold,old,new,i;
    oldold:=1;
    old:=1;
      for i from 2 to n do
        new:=oldold+old;     
        oldold:=old;         
        old:=new;            
      od;
    new;
  end: # fib_exact
  fib_exact(0):=1:
  fib_exact(1):=1:   

phi:=(1+sqrt(5))/2:
seq([fib_exact(n),evalf(phi^(n+1)/sqrt(5),5)],n=0..10);

fib_float(10^3);
fib_mod(10^3,10^50);
# fin exercice 2.11.

eval(plots);

with(plots):

with(linalg):

# section 3, page 87.

A:=1:B:=1:
alpha:=2/5:beta:=1/3:
phi:=0:psi:=Pi/2:

n:=lcm(denom(alpha),denom(beta)):
lissajous:=[A*cos(alpha*t+phi),B*cos(beta*t+psi),t=0..2*n*Pi];

plot(lissajous,scaling=constrained,color=black);

f:=piecewise(x<1/2,2*x,2*(1-x)):
interval:=0..1:
u0:=1/sqrt(2):
Number:=5:

window:=x=interval,y=interval:
bissect:=plot([t,t,t=interval],window):
graph:=plot([x,f,x=interval],window,thickness=3):
u[0]:=evalf(u0):
pt[0]:=[u[0],0]:
for k to Number do
  u[k]:=evalf(subs(x=u[k-1],f));
  pt[2*k-1]:=[u[k-1],u[k]]:
  pt[2*k]:=[u[k],u[k]]:
od:
brokenline:=plot([seq(pt[i],i=0..2*Number)]):
plots[display]({bissect,graph,brokenline},scaling=constrained);

astroid:=[8*a*cos(t)^3,8*a*sin(t)^3]:
diffastroid:=diff(astroid,t):
equ[1]:=linalg[det](array(1..2,1..2,
                               [[x-astroid[1],diffastroid[1]],
                                [y-astroid[2],diffastroid[2]]])):
equ[2]:=factor(equ[1]);

op(equ[2]);
equ[3]:=op(5,equ[2]):

equ[4]:=select(has,equ[2],{x,y}):
equ[5]:=map(combine,equ[4],trig);

sol:=solve({equ[5],x=0},{x,y}):
pt[Ox]:=expand(subs(sol,[x,y])):
sol:=solve({equ[5],y=0},{x,y}):
pt[Oy]:=expand(subs(sol,[x,y])):
combine(add((pt[Ox][i]-pt[Oy][i])^2,i=1..2),trig);

pointset:={[2,0],[0,1],[-1,0],[0,-1]}:
equ:=A*x^2+2*B*x*y+C*y^2+2*D*x+2*E*y+F:
inc:={A,B,C,D,E,F}:
sys:={seq(subs(x=i[1],y=i[2],equ),i=pointset)};

sol:=solve(sys,inc):
conic:=subs(sol,equ);

sys:={seq(subs(x=i[1],y=i[2],equ),i=pointset)};

sol:=solve(sys,inc):
conic:=subs(sol,equ);

# section 4, page 92.

x:=cos(t)^2+ln(sin(t)):
y:=sin(t)*cos(t):
diffx:=combine(diff(x,t),trig):
diffy:=combine(diff(y,t),trig):
solve({diffx,diffy},t);

eval(subs(t=Pi/4,[x,y,diffx,diffy]));

f:=ln(x+sqrt(x^2+1));

p[1]:=f+subs(x=-x,f);

combine(p[1],ln);

p[2]:=combine(p[1],ln,symbolic);

expand(p[2]);

p[3]:=map(expand,p[2]);

p[3];

q[1]:=(f+subs(x=-x,f))/2;

q[2]:=combine(q[1],ln,symbolic);

q[3]:=map(combine,q[2],radical,symbolic);

int(1/(u-cos(x)),x=0..Pi);

assume(u^2>1):
int(1/(u-cos(x)),x=0..Pi);

int(sin(n*theta)*sin(m*theta),theta=-Pi..Pi);

assume(n,integer):
assume(m,integer):
int(sin(n*theta)*sin(m*theta),theta=-Pi..Pi);

# section 6, page 96.

# exercice 2.1, page 96.
f:=expand((x-2*sin(3*theta))^3-(x+2*cos(2*theta))^3):
collect(f,x);

collect(f,x,readlib(`combine/trig`));
# fin exercice 2.1.

# exercice 2.2, page 97.
normal((x^2-1)/((x+2)*(x-2))),
normal((x-1)*(x+1)/(x^2-4)),
normal((x^2-1)/((x+2)*(x-2)),expanded),
normal((x-1)*(x+1)/(x^2-4),expanded);
# fin exercice 2.2.

# exercice 2.3, page 97.
h:=-1/(x-1):
hh:=h:
for k to 10 while hh<>x do hh:=normal(subs(x=h,hh)) od:
k;

hh:=h:
for k to 10 while hh<>x do hh:=subs(x=h,hh) od;
k;
# fin exercice 2.3.

# exercice 2.4, page 98.
map(combine,",trig);

A:=array(1..2,1..2,[[a,b],[c,d]]);
evalm(A^2-(a+d)*A+(a*d-b*c));
map(expand,");
# fin exercice 2.4.

# exercice 2.5, page 98.
alpha:=(1+3^(1/2)+(12)^(1/4))/2:
kappa[1]:=(alpha^3+1)/alpha/(alpha-1):
kappa[2]:=(9+6*3^(1/2))^(1/2):
simplify(kappa[1]-kappa[2],radical);

radnormal(kappa[1]-kappa[2]);

radnormal(kappa[1]),radnormal(kappa[2]);

radnormal(kappa[1],rationalized);
# fin exercice 2.5.

# exercice 2.6, page 99.
oddproperty:=proc(lambda)
  local i;
  for i to nops(lambda) do
    if type(lambda[i],even) then RETURN(false) fi
  od;
  true
end: # oddproperty
distinctproperty:=proc(lambda)
  local i;
  for i from 2 to nops(lambda) do
    if lambda[i-1]=lambda[i] then RETURN(false) fi
  od;
  true
end: # distinctproperty

oddproperty:=proc(lambda)
  type(lambda,list(odd))
end: # oddproperty
distinctproperty:=proc(lambda)
  evalb(nops(lambda)=nops({op(lambda)}))
end: # distinctproperty

niceproperty:=eval(oddproperty):
for N to Nmax do
  lambda:=[N];
  lambdatable[0]:=lambda;
  for i while lambda<>FAIL do
    lambda:=nextpartition(lambda);
    lambdatable[i]:=lambda
  od;
  nicecounter[N]:=nops(select(niceproperty,
                                [seq(lambdatable[j],j=0..i-2)]));
od:

seq(nicecounter[N],N=1..20);
# fin exercice 2.6.

# exercice 2.7, page 100.
nmax:=15:                    # maximal number of iterations in cp
rp:=rand(-1000..1000):
counter:=0:
to 100 do 
  P:=[seq([rp(),rp()],i=1..5)]; 
  if cp(P)=FAIL then 
    counter:=counter+1;
    example[counter]:=P;
    print(P)   
  fi 
od: 

invtransfpentagon:=proc(P::list(list(numeric)))
  local tP,x,y,i,S,newP;
  if nops(P) <> 5 or                                           #1
    (not convert([seq(evalb(nops(i)=2),i=P)],`and`)) then
     ERROR(`expected a list of 5 pairs, but received`,P) 
  fi;
  tP[0]:=[x,y];         
  for i to 5 do         
    tP[i]:=[2*P[i][1]-tP[i-1][1],2*P[i][2]-tP[i-1][2]] 
  od;
  S:=solve({tP[5][1]=x,tP[5][2]=y},{x,y});                     #2
  newP[0]:=subs(S,[x,y]);                 
  for i to 5 do 
    newP[i]:=subs(x=newP[0][1],y=newP[0][2],tP[i]) 
  od;
  [seq(newP[i],i=0..4)];
end: # invtransfpentagon

`type/pentagon`:=proc(P::list(list(numeric)))
  local i;
  evalb(nops(P)=5 and convert([seq(evalb(nops(i)=2),i=P)],`and`))
end: # type/pentagon

transfpentagon:=proc(P::pentagon)
#  ...
end: # transfpentagon
invtransfpentagon:=proc(P::pentagon)
#  ...
end: # invtransfpentagon
# fin exercice 2.7.

# exercice 2.10, page 101.
pattern:=proc(n::posint)
  local l,k,twotok,mantissa,m;
  l:=length(n);
  Digits:=2*l+5;
  twotok:=1.;
  mantissa:=op(1,twotok);
  m:=length(mantissa);
  for k from 0 while iquo(mantissa,10^(max(m-l,0)))<>n do 
    twotok:=2.*twotok;
    mantissa:=op(1,twotok);
    m:=length(mantissa);
  od:
  k
end: # pattern
plot([seq([n,pattern(n)],n=1..1000)]);
# fin exercice 2.10.

# exercice 2.11, page 102.
fib_float:=proc(n::nonnegint)
  global global_log10_phi;
  global_log10_phi:=evalf(log[10]((1+sqrt(5))/2),10):
  Digits:=round((n+1)*global_log10_phi*1.1+5);
  round(evalf(((1+sqrt(5))/2)^(n+1)/sqrt(5)))  
end: # fib_float

fib_mod:=proc(n::nonnegint,m::integer)
  local oldold,old,new,i;
  oldold:=1 mod m;
  old:=1 mod m;
    for i from 2 to n do
      new:=oldold+old mod m;
      oldold:=old;
      old:=new;
    od;
  new;
end: # fib_mod
fib_mod(0):=1 mod m:
fib_mod(1):=1 mod m:   
# fin exercice 2.11.

# exercice 2.12, page 103.
pointset:={[2,0],[0,1],[-1,0],[0,-1]}:
equ:=A*x^2+2*B*x*y+C*y^2+2*D*x+2*E*y+F:
inc:={A,B,C,D,E,F}:
sys:={seq(subs(x=i[1],y=i[2],equ),i=pointset)}:
sol:=solve(sys,inc):
conic:=subs(sol,equ);

hom2:=coeff(subs(x=t*x,y=t*y,conic),t,2);

SYS:=map(combine,{subs(x=cos(theta),y=sin(theta),hom2),
                  subs(x=-sin(theta),y=cos(theta),hom2)},trig);

SOL:=solve(SYS,inc union {theta});

for i to 3 do eqhyper[i]:=subs(SOL[i],conic) od:
seq(eqhyper[i],i=1..3);

pointset:={[2,1],[0,1],[-1,0],[0,-1]};
#...
conic:=subs(sol,equ);

hom2:=coeff(subs(x=t*x,y=t*y,conic),t,2);
#...
SOL:=[solve(SYS,inc union {theta})];
for i in SOL do 
  eqhyper[i]:=factor(subs(i,conic))
od:
EH:={seq(eqhyper[i],i=SOL)};

EH:=map2(select,has,EH,{x,y});
EH:=remove(type,EH,`*`);
# fin exercice 2.12.

# exercice 2.13, page 104.
J[1]:=Int(1/sqrt((1-t^2)*(1-m*t^2)),t=0..1);

J[2]:=student[changevar](t=sin(theta),J[1],theta);

J[3]:=subs(1-sin(theta)^2=cos(theta)^2,J[2]);

J[4]:=subs(cos(theta)=1,J[3]);

J[4]:=simplify(J[3],radical,symbolic):
# fin exercice 2.13.

# exercice 2.14, page 105.
E[1]:=evalc((1+I*tan(alpha))/(1-I*tan(alpha))):
E[2]:=convert(E[1],sincos):
E[3]:=normal(E[2]):
E[4]:=combine(E[3],trig);

C[1]:=evalc((1+I)^n/(1-I)^(n-2)):
C[2]:=combine(C[1],trig):
C[3]:=simplify(C[2],exp):
C[4]:=factor(C[3]):
C[5]:=simplify(C[4],power):
C[6]:=2*convert(normal(C[5]/2),exp):
C[7]:=expand(C[6]);
# fin exercice 2.14.

# exercice 2.15, page 105.
product(k*(k+2)/(k+1),k=1..n):
simplify(",GAMMA);

convert(",factorial);
expand(");
# fin exercice 2.15.

# exercice 2.16, page 106.
sum(cos((2*k+1)*Pi/(4*l+1)),k=0..l-1):
combine(",trig);

subs(l=10,");
# fin exercice 2.16.

# exercice 2.17, page 106.
K[1]:=Int(1/sqrt(a^2*cos(theta)^2+b^2*sin(theta)^2),
                                                  theta=0..Pi/2):
K[2]:=student[changevar](b*tan(theta)=t,K[1],t):
trace(simplify):
assume(b>0):
K[3]:=simplify(K[2],radical,symbolic);
assume(a>0):
assume(t,real):
K[4]:=simplify(K[3],power,symbolic);

K[3]:=simplify(K[2],radical,symbolic);
K[4]:=simplify(K[3],power,symbolic);
# fin exercice 2.17.