Maple
Son bon usage en mathématiques

Philippe Dumas, Xavier Gourdon

Springer-Verlag, 1997

Code Maple du Chapitre 1, Livre premier.

# section 1, page 3.

1+1;

sqrt(2);

evalf(");

evalf("",20);

seq(exp(-k/10*x)*sin(x),k=1..4);

plot({"},x=0..6*Pi,title=`sinusoides amorties`);
#plot(["],x=0..6*Pi,title=`sinusoides amorties`);

plot([exp(t)*cos(t),exp(t)*sin(t),t=-Pi..2*Pi/3],
              title=`spirale logarithmique`,scaling=constrained);

f:=x^4+x^2+1:
solve(f,x);

solve(f=0,x);

fsolve(f,x,complex);

series(tan(x),x,20);

expr:=sin(omega*t)*exp(lambda*t);

diff(expr,t);
factor(");

diff(expr,t,t);    
diff(expr,t$2);

int(expr,t);

Int(expr,t);

value(");

?plot

# section 2, page 7.

restart;

x:=2:
int(tan(x),x);

x:='x':

for k to 10 do
  diff(1/x,x$k)
od;

plot(sin(x,x=0..Pi);

?printlevel

series(1/convert([seq(1-x^k,k=1..10)],`*`),x,10);

series(1/mul(1-x^k,k=1..10),x,10);

[seq(1-x^k,k=1..10)];
convert(",`*`);
series(1/convert([seq(1-x^k,k=1..10)],`*`),x,10);

z:=exp(I*Pi/8):
bound:=20:
s:=0:
x[0]:=Re(s):
y[0]:=Im(s):
Z:=1:
for n to bound do 
  Z:=evalf(z*Z);
  s:=s+Z/n;
  x[n]:=Re(s);
  y[n]:=Im(s);
od:
pic[brokenspiral]:=plot([seq([x[n],y[n]],n=0..bound)],color=red):
pic[limitpoint]:=plot([[Re(-ln(1-z)),Im(-ln(1-z))]],color=blue,
                                      style=point,symbol=circle):
plots[display]({pic[brokenspiral],pic[limitpoint]},
                                            scaling=constrained);


# section 3, page 10.

p:=x^4+x+1;

solve(p,x);

evalf(");

sys:={x^3-y^3+x*y,x^2+y^2-1};
inc:={x,y};

fsolve(sys,inc);

expr:=int(1/(1+x^3),x);

nops(expr);

op(0,expr),[op(expr)],op(2,expr),op([2,2,1],expr);

S:=series(sin(x)^x,x):

op(0,S);op(S);

equivalent:=table([sin(x)=x,1-cos(x)=x^2/2,exp(x)-1=x]);

equivalent[ln(1+x)]:=x;

op(0,eval(equivalent));
op(eval(equivalent));

# exercice 1.5, page 12.
int(1/(x^7-x),x);
# fin exercice 1.5.

L[i];

# exercice 1.7, page 13.
L:=[seq(factor(diff(exp(-x^2),x$k)),k=1..4)];
op(L);
L[2];
# fin exercice 1.7.


convert(rationalexpr,parfrac,name);

convert(seriesexpr,polynom);

# exercice 1.8, page 13.
ci[1]:=int(1/cos(x),x);
ci[2]:=convert(ci[1],tan);
ci[3]:=convert(ci[1],exp);
ci[4]:=convert(ci[2],exp);
ci[5]:=convert(ci[1],sincos);
ci[6]:=convert(ci[3],trig);

cI[1]:=Int(1/cos(x),x);
value(");
cI[2]:=convert(cI[1],tan);
value(");
cI[3]:=convert(cI[1],exp);
value(");
cI[4]:=convert(cI[2],exp);
value(");
cI[5]:=convert(cI[1],sincos);
value(");
# fin exercice 1.8.

evalf([seq([cos(k*Pi/5),sin(k*Pi/5)],k=0..10)],4);

n:=3:
A:=matrix(3,3,[[1,2,3],[2,3,4],[4,5,6]]);
B:=linalg[transpose](A);

X:=matrix(n,n):
AX:=evalm(A&*X):
XB:=evalm(X&*B):
S:=solve({seq(seq(AX[i,j]=XB[i,j],j=1..n),i=1..n)},
                               {seq(seq(X[i,j],j=1..n),i=1..n)});

Xsol:=subs(S,evalm(X));

# exercice 1.9, page 15.
Eseq:=entries(Xsol);

Eset:=map(op,{Eseq});

indets(Eset,name);
# fin exercice 1.9.

# exercice 1.10, page 16.
";
print(x);
";
x:=2;
";
x$0;
";
# fin exercice 1.10.

# exercice 1.11, page 16.
equivalent:=table([sin(x)=x,1-cos(x)=x^2/2,exp(x)-1=x]):
# fin exercice 1.11.

Euler:=x^2+x+41;

[seq(subs(x=i,Euler),i=1..39)];

T:=[cos(u)*(2+cos(v)),sin(u)*(2+cos(v)),sin(v)];
C:=subs(u=sqrt(2)*t,v=t,T);

pic[T]:=  plot3d(T,u=0..2*Pi,v=0..2*Pi):
pic[C]:=plots[spacecurve](C,t=0..100,numpoints=1000,thickness=3):
plots[display]({pic[T],pic[C]},scaling=constrained);

subs(x=y,y=x,[x,y]);

subs({x=y,y=x},[x,y]);

# exercice 1.13, page 18.
T:=[cos(u)*(2+rho*cos(v)),sin(u)*(2+rho*cos(v)),rho*sin(v)];
pic[T]:=plot3d(subs(rho=1,T),u=0..2*Pi,v=0..2*Pi,color=black):
pic[C]:=plots[spacecurve](subs(rho=1.1,u=sqrt(2)*t,v=t,T),
                                        t=0..100,numpoints=1000,
                                        thickness=3,color=black):
plots[display]({pic[T],pic[C]},scaling=constrained);
# fin exercice 1.13.

# section 4, page 19.

whattype(exp(-x^2/2));

type(expr,t);

seq([i,type(exp(-x^2/2),i)],i=[function,algebraic,algfun,ratpoly,
                               freeof(x),exp(ratpoly),anything]);

type(2+3/2*x+x^15,polynom(rational,x));

type(Pi*x+x^5,polynom(numeric));

type(Pi*x+x^5,polynom(realcons));

# exercice 1.14, page 20.
type(sqrt(2),numeric);
type(fsolve(x^3-x+1,x),numeric);
type(fsolve(x^3-3*x+1,x),numeric);
type([fsolve(x^3-3*x+1,x)],list(numeric));
type(x:=2,equation);
type(x<=2,relation);
S:=series(sin(x),x):
type(S,`+`);
P:=convert(S,polynom):
type(P,`+`);
type(toto,name);
type(toto,string);
# fin exercice 1.14.

# exercice 1.15, page 20.
trace(type):
convert(327,base,2);
sin(1+I);
evalc(sin(1+I));
# fin exercice 1.15.

# section 5, page 21.

do od;

f:=x^3-x+1;

u:=0.5;
for n to 4 do u:=subs(x=u,f) od;

n:='n';    
n:=evaln(n);

# exercice 1.16, page 22.
for theta from 0 by Pi/10 to Pi/2 do sin(theta) od;
for k from 0 to 5 do sin(k*Pi/10) od;
for n to 5 do
  T[n]:=subs(cos(theta)=x,expand(cos(n*theta)))
od;
for n to 5 do 
  P[n]:=collect((-1)^n/2^n/n!*diff((1-x^2)^n,x$n),x) 
od;
# fin exercice 1.16.

expr:=diff((z^2+1)/(z^2-1),z);

for i in expr do antider[i]:=int(i,z) od;

i;

n:=10:
sigma:=combstruct[draw](Permutation(n));

e:=[seq(i,i=1..n)];

tau:=sigma:
for k while tau <> e do 
  tau:=subs({seq(i=sigma[i],i=1..n)},tau)
 #tau:=[seq(sigma[i],i=tau)]
od;

k;

V:=[r*cos(omega*t),r*sin(omega*t),h];

M:=map(int,V,t);

R:=matrix(3,3):
R[1,1]:=cos(phi)*cos(psi)-cos(theta)*sin(phi)*sin(psi):
R[1,2]:=subs(psi=psi+Pi/2,R[1,1]):
R[2,1]:=subs(phi=phi-Pi/2,R[1,1]):
R[2,2]:=subs(phi=phi-Pi/2,R[1,2]):
R[1,3]:=sin(phi)*sin(theta):
R[2,3]:=-cos(phi)*sin(theta):
R[3,1]:=sin(theta)*sin(psi):
R[3,2]:=sin(theta)*cos(psi):
R[3,3]:=cos(theta):
R:=map(expand,R);

diffRphi:=diff(eval(R),phi);

diffRphi:=map(diff,R,phi);

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

map2(diff,peano,[seq([x$(2-k),y$k],k=0..2)]):
normal(");

[seq([x$(2-k),y$k],k=0..2)];

pf:=convert(1/z^3/(1-z)^3,parfrac,z);

select(has,pf,-1+z);

remove(has,pf,-1+z);

# exercice 1.19, page 26.
complexlist:=[seq(evalc(((1+I*3^(1/2))/(3/2))^k),k=0..8)];
# fin exercice 1.19.

add(x[i],i=0..5);

add(1./n^2,n=1..100);

add(k,k=1..n);

sum(k,k=1..n):
factor(");

# exercice 1.20, page 27.
iterate:=proc(path)
  local directionlist,n,direction,newpoint,T;
  directionlist:=[[1,0],[0,1],[-1,0],[0,-1]];
  n:=nops(path);
  for direction in directionlist do
    newpoint:=[path[n][1]+direction[1],path[n][2]+direction[2]];
    if member(newpoint,path) 
    then T[direction]:=NULL
    else T[direction]:=[op(path),newpoint]
    fi;
  od;
  seq(T[direction],direction=directionlist)
end: # iterate

PATH:=[[0,0],[1,0],[1,1],[1,2],[0,2],[0,1],[-1,1],[-1,0]]:
iterate(PATH);
# fin exercice 1.20.

# section 6, page 29.

maxi:=1000:
s:=0:
for n to maxi do
  s:= s+1./n;
  if (n mod 100)=0 then print(n,`--> `,evalf(s-ln(n))) fi;
od: 
n:='n':

# exercice 1.24, page 30.
asympt(1/(n^2+1),n);
# fin exercice 1.24.

# exercice 1.25, page 30.
maxi:=1000:
s:=0:
  for n to maxi do
    s:= s+1/n;
  if (n mod 100)=0 then print(n,`--> `,evalf(s-ln(n))) fi;
  od: n:='n':
# fin exercice 1.25.

# section 7, page 31.

expr_to_seq:=proc(expr::anything) 
  expr_to_seq_aux([],expr) 
end: # expr_to_seq
expr_to_seq_aux:=proc(path,expr)
  local i,n;
  if type(expr,{name,integer})            # termination condition
  then [path,expr]               
  else 
    n:=nops(expr);                        # recursive call
    [path,op(0,expr)],
      seq(expr_to_seq_aux([op(path),i],op(i,expr)),i=1..n) 
  fi 
end: # expr_to_seq_aux

expr:=factor(diff(exp(x^2)*(x^2+1),x));

expr_to_seq(expr);

trace(expr_to_seq_aux):
expr_to_seq(expr);
untrace(expr_to_seq_aux):

# exercice 1.26, page 33.
maple_to_lisp((1-4*x)^(-1/2));
# fin exercice 1.26.

interface(verboseproc=2);
print(student[intparts]);

# exercice 1.27, page 34.
trace(member);
student[intparts](Int(x*exp(x),x),x);
# fin exercice 1.27.

convert(series(tan(x),x,10),confrac);

interface(verboseproc=2);
readlib(`convert/confrac`);

interface(verboseproc=2);
print(evalm);

evalm(x,y);

convert(13,base,2);

nu:=proc(n::nonnegint)
  if n=0 then 0
  elif type(n,even) then nu(n/2)
  else nu((n-1)/2)+1
  fi
end: # nu
nu(13);

trace(nu):
nu(13):

countnu:=proc(n::nonnegint)
  global T;
  if assigned(T[n]) then T[n]:=T[n]+1 else T[n]:=1 fi;
  if n=0 then 0
  elif type(n,even) then countnu(n/2)
  else countnu((n-1)/2)+1
  fi
end: # countnu
for n from 0 to 10 do countnu(n) od:
seq([n,T[n]],n=0..10);

remembernu:=proc(n::nonnegint)
  option remember;
  if n=0 then 0
  elif type(n,even) then remembernu(n/2)
  else remembernu((n-1)/2)+1
  fi
end: # remembernu
for n from 0 to 10 do remembernu(n) od:
R:=op(4,eval(remembernu)):
seq([n,R[n]],n=0..10);

remembernu:=subsop(4=NULL,eval(remembernu)):

# exercice 1.30, page 39.
series(exp(exp(x)-1),x,11);
add(bell(n)/n!*x^n,n=0..10);
# fin exercice 1.30.

lucas:=proc(n::nonnegint)
  option remember;
  lucas(n-1)+lucas(n-2)
end: # lucas
lucas(0):=2:
lucas(1):=1:
op(4,eval(lucas));

f(0):=1:

# exercice 1.32, page 40.
for k from 0 to 4 do
  convert(series(add(fibonacci(i)^k*t^i,i=0..20),t,20),ratpoly)
od;
# fin exercice 1.32.

# section 8, page 41.

T:=table([seq(c=evalf(c,10),c=constants)]):
T;

print(T);

op(4,ln);

op(4,eval(ln));

p[1]:=proc()
  local U,V,W;
  U:=V;
  V:=array(1..2,1..2,[[1,2],[3,4]]);
  linalg[transpose](U)
end: # p[1]
p[1]();

p[2]:=proc()
  local U,V,W;
  U:=V;
  V:=W;
  W:=array(1..2,1..2,[[1,2],[3,4]]);
  linalg[transpose](U)
end: # p[2]
p[2]();

p[3]:=proc()
  local U,V,W;
  U:=V;
  V:=W;
  W:=array(1..2,1..2,[[1,2],[3,4]]);
  linalg[transpose](eval(U))
end: # p[3]
p[3]();

# exercice 1.33, page 43.
trace(linalg[transpose]);
# fin exercice 1.33.

A:=matrix(2,2,[[0,-1],[1,0]]);

B:=subs(x=A,x^2+1);

evalm(B);

z:=(1-I)*(5+I)^4;

arctan(Im(z),Re(z));

evala((5-sqrt(2))/(3+sqrt(2)));

evala((5-RootOf(x^2-2))/(3+RootOf(y^2-2)));

alias(theta=RootOf(_Z^3+_Z+1)):
trace(evala):
evala(1/add(theta^k,k=0..5));

# section 9, page 45.

F[0]:=Int(x^(2*m+1)/(x^2-a^2)^(n+1/2),x);

G:=subs(n=3,m=2,F[0]);

value(G);

F[1]:=student[intparts](F[0],x^(2*m)/(x^2-a^2)^n);

F[2]:=map(combine,F[1],power,symbolic);

F[3]:=map(combine,F[2],power,symbolic);

F[4]:=map(expand,F[3]);

F[5]:=map(combine,F[4],power,symbolic);

evalf(int(exp(-sec(phi)),phi=0..Pi/4));

evalf(Int(exp(-sec(phi)),phi=0..Pi/4));


2 &^ (10^6) mod 3^5;

# exercice 1.37, page 48.
binomialdecomposition(22,5);
# fin exercice 1.37.

# section 10, page 48.

# exercice 1.1, page 48.
dsolve(x*diff(y(x),x)+y(x)=sin(x),y);

dsolve(x*diff(y(x),x)+y(x)=sin(x),y(x));
# fin exercice 1.1.

# exercice 1.2, page 48.
plot([sin(x),x,x-x^3/6],x=-Pi/2..Pi);

plot([sin(x),x,x-x^3/6],x=-Pi/2..Pi,color=[red,blue,green]);

plot([sin(x),x,x-x^3/6],x=-Pi/2..Pi,color=black);

L:=[sin(x),x,x-x^3/6]:
C:=[red,blue,green]:
n:=nops(L):
for i to n do
  pic[i]:=plot(L[i],x=-Pi/2..Pi,color=C[i])
od:
plots[display]({seq(pic[i],i=1..n)});
# fin exercice 1.2.

# exercice 1.3, page 49.
series(sin(x),x,10);

series(sin(x)/x,x,10);

series(1/sin(x),x,10);

series(1/sin(x)-1/x,x,10);

s[1]:=series(1/sin(x),x);

s[2]:=series(sin(x),x);

series(s[1]*s[2],x);
# fin exercice 1.3.

# exercice 1.4, page 50.
expand(cos(5*phi));

subs(cos(phi)=x,");

orthopoly[T](5,x);

solve(",x);

fsolve("",x);
# fin exercice 1.4.

# exercice 1.6, page 51.
expr:=3/2; 
op(0,expr),op(expr);

expr:=I; 
op(0,expr),op(expr);

expr:=z=sin(z/2); 
op(0,expr),op(expr);
# fin exercice 1.6.

# exercice 1.7, page 52.
L:=[seq(factor(diff(exp(-x^2),x$k)),k=1..4)];

L[2];
# fin exercice 1.7.

# exercice 1.8, page 52.
ci[1]:=int(1/cos(x),x);

ci[5]:=convert(ci[1],sincos);

cI[1]:=Int(1/cos(x),x);

cI[2]:=convert(cI[1],tan);

value(");

cI[4]:=convert(cI[2],exp);

value(");
# fin exercice 1.8.

# exercice 1.9, page 53.
indices(Xsol);

Eset:=map(op,{Eseq});

indets((ln(x+sqrt(exp(x)+1))+x)/(x+sin(x*y)));

indets(convert(Xsol,set));
# fin exercice 1.9.

# exercice 1.11, page 53.
seq([i],i=eval(equivalent));
# fin exercice 1.11.

# exercice 1.13, page 54.
K:=[rho*u*cos(v),rho*u*sin(v),2*u];
pic[K]:=plot3d(subs(rho=1,K),u=-1..1,v=0..2*Pi,color=black):
pic[H]:=plots[spacecurve](subs(rho=1.1,u=t,v=10*t,K),t=-1..1,
                           numpoints=50,thickness=3,color=black):
plots[display]({pic[K],pic[H]},orientation=[70,85],
                                            scaling=constrained);
# fin exercice 1.13.

# exercice 1.14, page 54.
type(sqrt(2),numeric);

fsolve(x^3-x+1,x),type(fsolve(x^3-x+1,x),numeric);

fsolve(x^3-3*x+1,x),type(fsolve(x^3-3*x+1,x),numeric);

type([fsolve(x^3-3*x+1,x)],list(numeric));

type(x:=2,equation);

type(x<=2,relation);

S:=series(sin(x),x):
op(0,S);
op(S);

type(S,`+`);

P:=convert(S,polynom):
type(P,`+`);

series((1+x)^x,x);
lprint(");
series((1+sqrt(x))^x,x);
lprint(");
type(",`+`);
# fin exercice 1.14.

# exercice 1.15, page 56.
untrace(type);
# fin exercice 1.15.

# exercice 1.17, page 56.
Digits:=25:
olda:=1.: oldb:=sqrt(2.):
while abs(olda-oldb)>10^(-20) do
  newa:=(olda+oldb)/2;
  newb:=sqrt(olda*oldb);
  olda:=newa;
  oldb:=newb;
od:
M:=(olda+oldb)/2:
1/M;

evalf(2/Pi*Int(1/sqrt(1-t^4),t=0..1));
# fin exercice 1.17.

# exercice 1.18, page 57.
H:=0:
for n while H<=10 do H:=evalf(H+1/n) od:
n,H;
# fin exercice 1.18.

# exercice 1.19, page 57.
complexlist:=[seq(evalc(((1+I*3^(1/2))/(3/2))^k),k=0..8)];
pointlist:=map(proc(z) [Re(z),Im(z)] end,complexlist);
plot(pointlist,scaling=constrained);
# fin exercice 1.19.

# exercice 1.20, page 57.
T[0]:=1:
pathlist:=[[[0,0]]]:
nmax:=10:
for n to nmax do 
  pathlist:=map(iterate,pathlist);
  T[n]:=nops(pathlist)
od:
seq(T[n],n=0..nmax);

for n from 1 to 9 do evalf(T[n+1]/T[n]) od;
# fin exercice 1.20.

# exercice 1.21, page 58.
for p to 10 do 
    factor(sum(k^p,k=1..n))
od;
# fin exercice 1.21.

# exercice 1.22, page 58.
m:=5:
for i to m do
    p[i]:=1/mul(1-lambda[i]/lambda[j],j=1..i-1)
           /mul(1-lambda[i]/lambda[j],j=i+1..m)
od:
for k to m-1 do s[k]:=add(lambda[i]^k*p[i],i=1..m) od:

seq(normal(s[k]),k=1..m-1);

restart:
p:=1/product(1-lambda[i]/lambda[j],j=1..i-1)
           /product(1-lambda[i]/lambda[j],j=i+1..m):
s:=sum(lambda[i]^k*p[i],i=1..m):

for k to m-1 do indets(s[k]) od;

indets(s);
# fin exercice 1.22.

# exercice 1.23, page 59.
for n to 50 do
    convert(map(proc(pair) 1/pair[1]/pair[2] end,
                select(proc(pair)
                              evalb((pair[1]n) 
                         and (igcd(pair[1],pair[2])=1)) end,
                       [seq(seq([p,q],q=2..n),p=2..n)])),
            `+`)
od;

for n to 50 do
     convert(map(proc(pair) 1/pair[1]/pair[2] end,
                 select(proc(pair)
                          evalb((igcd(pair[1],pair[2])=1))
                        end,
                  [seq(seq([p,q],q=max(p+1,n+1-p)..n),p=2..n)])),
             `+`)
od;

for n to 50 do
  s[n]:=0;
  for p from 2 to n do 
    for q from max(p+1,n+1-p) to n do 
      if igcd(p,q)=1 then s[n]:=s[n]+1/p/q fi
    od
  od;
od:
seq([n,s[n]],n=1..50);
# fin exercice 1.23.

# exercice 1.24, page 60.
da:=asympt(Psi(n+1),n);

subs(n=1000,da);

op(0,da),op(da);

subs(n=1000,O=0,[op(da)]);

evalf(");
# fin exercice 1.24.

# exercice 1.26, page 61.
maple_to_lisp:=proc(expr::anything)
  local i;
  if type(expr,{name,integer}) then expr
  else [op(0,expr),seq(maple_to_lisp(op(i,expr)),i=1..nops(expr))]
  fi
end: # maple_to_lisp

n:=evalc((1-2*I)^3);
expr_to_seq(n);
maple_to_lisp(n);

absn:=abs(n);
expr_to_seq(absn);
maple_to_lisp(absn);

exprtrig:=combine(cos(x)^3,trig);
expr_to_seq(exprtrig);
maple_to_lisp(exprtrig);

A:=array(symmetric,1..2,1..2,[(1,1)=x,(1,2)=y,(2,2)=z]);
expr_to_seq(A);
expr_to_seq(eval(A));
maple_to_lisp(eval(A));

B:=array(1..2,1..2,[[a,b],[c,d]]);
expr_to_seq(eval(B));
maple_to_lisp(eval(B));
# fin exercice 1.26.

# exercice 1.27, page 61.
untrace(member);
# fin exercice 1.27.

# exercice 1.28, page 61.
permutation_order_1:=proc(sigma::list(posint))
  local n,e,tau,k,i;
  n:=nops(sigma);
  e:=[seq(i,i=1..n)];
  if {op(sigma)}={op(e)} then
    tau:=sigma;
    for k while tau <> e do tau:=[seq(sigma[i],i=tau)] od;
    k;
  else ERROR(`expected a permutation of [1,`,n,`],
                                             but received`,sigma)
  fi
end: # permutation_order_1

`type/permutation`:=proc(sigma)
  local i;
  type(sigma,list(posint)) and
    {op(sigma)}={seq(i,i=1..nops(sigma))}
end: # type/permutation
permutation_order_2:=proc(sigma::permutation)
  local n,e,tau,k,i;
  n:=nops(sigma);
  e:=[seq(i,i=1..n)];
  tau:=sigma;
    for k while tau <> e do tau:=[seq(sigma[i],i=tau)] od;
  k;
end: # permutation_order_2

sigma:=[1,2,7,4,5,6];
permutation_order_1(sigma);

permutation_order_2(sigma);
# fin exercice 1.28.

# exercice 1.29, page 62.
op(4,eval(evalf));
# fin exercice 1.29.

# exercice 1.30, page 63.
bell:=proc(n::nonnegint)
  local minus1,k;
  option remember;
  if n=0 then 1
  else 
    minus1:=n-1;
    add(binomial(minus1,k)*bell(k),k=0..minus1)
  fi
end: # bell
seq(bell(n),n=0..10);

bell:=proc(n::nonnegint)
  local minus1,k;
  option remember;
  minus1:=n-1;
  add(binomial(minus1,k)*bell(k),k=0..minus1)
end: # bell
bell(0):=1:
# fin exercice 1.30.

# exercice 1.31, page 63.
directnu:=proc(n::nonnegint)
  #local i;
  convert(convert(n,base,2),`+`)
  #add(i,i=convert(n,base,2));
end: # directnu
directnu(13);
# fin exercice 1.31.

# exercice 1.32, page 64.
fibonacci:=proc(n::nonnegint)
  option remember;
  fibonacci(n-1)+fibonacci(n-2)
end: # fibonacci
fibonacci(0):=1:
fibonacci(1):=1:

coeff_a:=proc(k::nonnegint,j::nonnegint)
  local m;
  option remember;
  if j=0 then 
    if k=0 then 1 else 0 fi
  else
    add(coeff_a(m,j-1)*fibonacci(k-m-2),m=0..k-2)
  fi
end: # coeff_a

Phi:=proc(k::nonnegint)
  option remember;
  normal((1+k*t*add((-1)^j*coeff_a(k,j)/j*
                     subs(t=(-1)^j*t,Phi(k-2*j)),j=1..iquo(k,2)))
                                      /(1-lucas(k)*t+(-1)^k*t^2))
end: # Phi
Phi(0):=1/(1-t):
Phi(1):=1/(1-t-t^2):

seq(Phi(k),k=0..3);
# fin exercice 1.32.

# exercice 1.33, page 64.
for i to m do for j to n do B[i,j]:=A[j,i] od od;
# fin exercice 1.33.

# exercice 1.34, page 65.
evalc((7+I)^5*(79+3*I)^2);
# fin exercice 1.34.

# exercice 1.35, page 65.
T[0]:=Int(cos(x)^n*sin(x)^m,x);
T[1]:=student[intparts](T[0],cos(x)^n*sin(x)^(m-1));

T[2]:=expand(T[1]);
T[3]:=combine(T[2],power,symbolic);
# fin exercice 1.35.

# exercice 1.36, page 65.
R:=Int(arctan(t)/t,t=0..x);

evalf(subs(x=1,R));

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

series(R,x,10);
evalf(Sum((-1)^n/(2*n+1)^2,n=0..infinity));
# fin exercice 1.36.

# exercice 1.37, page 66.
binomialdecomposition:=proc(N::nonnegint,l::posint)
  local k,counter;
  if N=0 then add(Binomial(k-1,k),k=1..l)
  elif l=1 then Binomial(N,1)
  else for counter while binomial(counter,l)<=N do od;
    Binomial(counter-1,l)
      +binomialdecomposition(N-binomial(counter-1,l),l-1)
  fi
end: # binomialdecomposition

`value/Binomial`:=proc(expr)
  subs(Binomial=binomial,expr)
end: # value/Binomial
# fin exercice 1.37.