Philippe.Dumas@inria.fr


http://algo.inria.fr/dumas/Maple/

Retour en page principale Table des matières Index

Maple V.3 worksheet Maple V.4 worksheet, bis, ter, proc Maple V.5 worksheet, bis, proc

Structure de DAG. Le but de cet exercice est de mettre en valeur la structure de DAG (directed acyclic graph, autrement dit graphe orienté sans cycle) et l'emploi que l'on peut en faire pour augmenter l'efficacité des calculs. Cette notion de DAG a été présentée à la page 69 et illustrée à la page 70.

Considérons l'expression suivante et ses dérivées successives.

> f:=exp(sin(x));

> for k to 3 do f.k:=diff(f,x$k) od;

La dérivée première f1 peut être vue comme un arbre, mais elle est en fait stockée comme un DAG. Les deux objets sont illustrés ci-dessous.

Dans le DAG, le x qui figure deux fois dans l'expression est partagée, ce qui représente un gain en mémoire. On pourrait rétorquer que le gain est bien faible. Mais considérons la derivée seconde.

Cette fois-ci le partage est plus important, puisque exp(sin(x)) apparaît deux fois, sin(x) apparaît trois fois et x apparaît quatre fois. Pour la dérivée troisième le gain en mémoire devient encore plus important.

Insistons sur le fait qu'il s'agit bien d'un gain. En effet chaque objet Maple n'existe qu'en un exemplaire. Nous nous en convainquons en regardant les adresses mémoire de quelques expressions et sous-expressions que nous venons de rencontrer.

> addressof(f);

> op([2,2],f2);

> addressof(op([2,2],f2));

> addressof(cos(x));

> op([3,1,1],f3);

> addressof(op([3,1,1],f3));

Emploi des DAGs à travers l'option remember. Comment utiliser au mieux cette structure de DAG ? Une première solution consiste à employer l'option remember dans les procédures. Reprenons la fonction f3 et calculons en un développement limité au voisinage de 0.

> f:=exp(sin(x)):

> f3:=diff(f,x$3);

La procédure series possède une table de remember. Comme on le voit dans le calcul suivant, les deux procédures series/cos et series/sin ne sont appelées qu'une fois, alors que cos(x) et sin(x) apparaissent chacun plusieurs fois dans l'expression de f3. Ceci est dû à la présence de cette table de remember ; avant de démarrer le calcul, le système regarde si le résultat n'est pas déjà stocké dans la table.

> trace(`series/cos`,`series/sin`):

> series(f3,x);

{--> enter series/cos, args = x, x

<-- exit series/cos (now at top level) = series(1-1/2*x^2+1/24*x^4+O(x^6),x,6)}

{--> enter series/sin, args = x, x

<-- exit series/sin (now in series/exp) = series(1*x-1/6*x^3+1/120*x^5+O(x^6),x,6)}

> untrace(`series/cos`,`series/sin`):

On a bien ainsi exploité la structure de DAG, puisqu'on n'a pas parcouru tout l'arbre associé à f3 mais seulement le DAG. Cependant cette méthode a l'inconvénient de charger la mémoire car la table de remember se remplit de plus en plus au cours de la session.

> interface(verboseproc=2):

> op(4,eval(series));

Emploi efficace de la structure de DAG. Pour ne pas tomber dans cet écueil, l'exercice proposait une autre approche. Cependant l'énoncé est si obscur et confus que le lecteur qui l'aurait résolu mériterait toute notre considération pour ses capacités divinatoires. On y demandait d'abord d'employer optimize.

En Maple V.3 ou V.4, on commence par se renseigner sur la syntaxe d'optimize puis on passe à l'acte.

> ?optimize

> ?optimize/makeproc

> f:=2*x^2*sin(x*ln(1+x))+sin(x^3+x*sin(x*ln(1+x)));

> optimf:=readlib(optimize)(f);

> `optimize/makeproc`([optimf],parameters=[x]);

En Maple V.5, la syntaxe a changé, mais on arrive au même résultat. Plus précisément, il semble qu'il y ait deux syntaxes, mais il n'y a en fait qu'une seule procédure cachée derrière ces deux syntaxes.

> f:=2*x^2*sin(x*ln(1+x))+sin(x^3+x*sin(x*ln(1+x)));

> oldoptimf:=readlib(optimize)(f);

> newoptimf:=codegen[optimize](f);

> Optimf:=[oldoptimf]:

> F:=codegen[makeproc](Optimf,[x]);

> F(t);

La procédure optimize a fourni une chaîne de calculs, telle que chaque sous-expression n'est calculée qu'une fois. On est exactement dans l'optique des DAGs.

Pour utiliser cette chaîne de calcul, on peut procéder comme suit. On construit une table tosubs que l'on va faire évoluer. Les indices de la table sont les variables t1, ..., t12, qui servent à coder la chaîne de calcul. Les valeurs associées sont les résultats de l'application de series aux expressions correspondantes. Ces valeurs vont être utilisées pour opérer les substitutions qui permettent d'arriver au développement demandé. On démarre avec tosubs[t1]:=x^2. Ensuite on rencontre t3 qui vaut ln(1+x) ; on va donc affecter à tosubs[t3] la valeur renvoyée par series(ln(1+x),x). Cette valeur va être utilisée dans la chaîne de calcul pour être substituée à t3. Au bout du compte, on aura les développements des expressions associées à chacune des variables t1, ..., t12.

> for i to nops(Optimf) do

    tosubs[op(1,Optimf[i])]:=

      series(subs([seq(j=tosubs[j],

                       j=seq(op(1,Optimf[k]),k=1..i-1))],

                  op(2,Optimf[i])),x)

  od;

> series(F(x),x);

On obtient bien le résultat voulu. Cette façon de procéder est tout à fait satisfaisante pour une procédure qui ne possède pas de table de remember. En effet chaque sous-expression utile n'a été calculée qu'une fois et on n'a pas stocké de résultat indésirable. Hélas, dans le cas qui nous occupe nous employons la procédure series qui possède une table de remember.

Reprenons donc cette boucle en évitant au mieux les appels à series. Si une expression s'écrit f(s), nous appelons series/f sur s. D'habitude de telles procédures sont appelées par series. Ici il faut un readlib pour en disposer. Dans le cas d'une somme ou d'un produit nous sommes obligés d'appeler series lui-même.

> restart:

> f:=2*x^2*sin(x*ln(1+x))+sin(x^3+x*sin(x*ln(1+x))):

> Optimf:=[readlib(optimize)(f)];

> printlevel:=2:

> for i to nops(Optimf) do

    expr:=op(2,Optimf[i]);

    if type(expr,function) then

      tosubs[op(1,Optimf[i])]:=

        readlib(`series/`.(op(0,expr)))(

          op(subs([seq(j=tosubs[j],

                       j=seq(op(1,Optimf[k]),k=1..i-1))],

                                        expr)),x)

    elif nops(expr)>1 then

      tosubs[op(1,Optimf[i])]:=
        series(subs([seq(j=tosubs[j],

                         j=seq(op(1,Optimf[k]),k=1..i-1))],

                                        expr),x)

    else 

      subs([seq(j=tosubs[j],j=seq(op(1,Optimf[k]),k=1..i-1))],expr)

    fi

  od;

Emploi efficace de l'option remember. Passons maintenant à une autre approche du problème. Nous allons d'abord l'illustrer avec la procédure subs. Celle-ci n'a pas de table de remember ; assez curieusement, cette procédure de base n'exploite pas la structure de DAG. La procédure qui suit, au contraire, exploite cette structure. Cependant elle n'a pas de table de remember. Autrement dit, on a les avantages sans avoir les inconvénients. Ce petit miracle est réalisé en incluant dans la procédure dagsubs une procédure avec table de remember. Comme la procédure incluse est locale, elle est évacuée quand on ressort de dagsubs. On notera l'emploi de procname, indispensable en Maple V.4 parce qu'il n'y a jusqu'à la version V.4 qu'un niveau de variables locales. En Maple V.5 on pourrait écrire doit, le nom de la procédure au lieu de procname.

> dagsubs:=proc(u,expr)

    local doit;

    doit:=proc(expr)

      option remember;

      if nops(expr)>1 then map(procname,expr)

      elif type(expr,function) then op(0,expr)(procname(op(expr)))

      else expr

      fi

    end;

    doit(op(1,u)):=op(2,u);

    doit(expr)

  end:

Cette élégante solution apparaît dans [GiHaLeMaSa97].

Hélas, la mémoire ocupée par le code et la table de remember n'est pas récupérée ! On peut s'en rendre compte comme suit. On fait sortir de la procédure dagsubs, les adresses mémoire de la procédure incluse doit et de sa table de remember.

> dagsubs:=proc(u,expr)

    local doit;

    doit:=proc(expr)

      option remember;

      if nops(expr)>1 then map(procname,expr)

      elif type(expr,function) then op(0,expr)(procname(op(expr)))

      else expr

      fi

    end;

    doit(op(1,u)):=op(2,u);

    [doit(expr),addressof(eval(doit)),addressof(op(4,eval(doit)))]

  end:

> e:=diff(exp(sin(x)),x$3);

> `bloody hell`:=dagsubs(x=Pi,e);

La procédure incluse n'est pas accessible par son nom et l'accès à la table de remember par ce nom n'a pas de sens.

> eval(doit);

> op(4,eval(doit));

Error, improper op or subscript selector

Pourtant la procédure et sa table de remember sont toujours en mémoire.

> pointto(`bloody hell`[2]);

> pointto(`bloody hell`[3]);

Ce comportement est passablement étonnant : la mémoire est gérée comme une poubelle.

Il ne nous reste qu'une solution qui est de vider explicitement la table de remember. On notera d'ailleurs qu'il n'y a du coup plus de raison d'inclure doit dans dagsubs ; on pourrait aussi bien la mettre à l'extérieur de dagsubs.

> dagsubs:=proc(u,expr)

    local doit,result;

    doit:=proc(expr)

      option remember;

      if nops(expr)>1 then map(doit,expr)

      elif type(expr,function) then op(0,expr)(doit(op(expr)))

      else expr

      fi

    end;

    doit(op(1,u)):=op(2,u);

    result:=doit(expr);

    subsop(4=NULL,eval(doit));

    [result,addressof(eval(doit)),addressof(op(4,eval(doit)))]

end:

> `bloody hell`:=dagsubs(x=Pi,e);

> pointto(`bloody hell`[2]);

> pointto(`bloody hell`[3]);

Reprenons l'idée sous-jacente pour écrire une procédure dagseries.

> dagseries:=proc(f,x)

  local doit,result;

  if nargs>2 then Order:=args[3] fi;

  doit:=proc(f,x)

    option remember;

    if type(f,function) then

      readlib(`series/`.(op(0,f)))

        (op(map(doit,[op(f)],x)),x)

    elif nops(f)>1 then

      series(map(doit,f,x),x)

    else f

    fi

  end;

  result:=doit(f,x);

  subsop(4=NULL,eval(doit));

  result

end:

> dagseries(diff(exp(sin(x)),x$3),x,10);

L'énoncé de l'exercice commençait en évoquant des évaluateurs. Les différentes méthodes présentées ci-dessus peuvent être appliquées quand on emploie des évaluateurs dont l'action sur une expression ne dépend que de l'étiquette de la racine de l'arbre qui la représente.

Retour en page principale