Optimización de llamadas de cola en Mathematica?

Puedo resumir las conclusiones a las que me llevó mi experiencia personal, con una advertencia de que lo que sigue puede no ser la explicación completamente correcta. El anwer parece estar en las diferencias entre la pila de llamadas de Mathematica y las pilas de llamadas tradicionales, que se originan a partir de funciones definidas por patrones de Mathematica que son realmente reglas. Por lo tanto, no hay llamadas a funciones reales. Mathematica necesita una pila por una razón diferente: ya que la evaluación normal ocurre desde la parte inferior de una expresión árbol, debe mantener expresiones intermedias en caso de que partes cada vez más profundas de (sub)expresiones sean reemplazadas como resultado de la aplicación de reglas (algunas partes de una expresión crecen desde la parte inferior). Este es el caso, en particular, de las reglas que definen lo que llamaríamos funciones no recursivas de cola en otros lenguajes. Así que, una vez más, la pila en Mathematica es una pila de expresiones intermedias, no llamadas a funciones.

Esto significa que si, como resultado de la aplicación de la regla, un la (sub)expresión se puede reescribir en su totalidad, la rama de expresión no necesita mantenerse en la pila de expresiones. Esto es probablemente lo que se conoce como optimización de llamadas de cola en Mathematica - y esta es la razón por la que en tales casos tenemos iteración en lugar de recursión (este es un muy buen ejemplo de las diferencias entre las aplicaciones de reglas y las llamadas a funciones). Reglas como f[x_]:=f[x+1] son de este tipo. Sin embargo, si alguna subexpresión se reescribe, produciendo más estructura de expresión, entonces la expresión debe se almacenará en la pila. La regla f[x_ /; x < 5] := (n += 1; f[x + 1]) es de este tipo, que está un poco oculta hasta que recordemos que ()significa CompoundExpression[]. Esquemáticamente lo que sucede aquí es f[1] -> CompoundExpression[n+=1, f[2]] -> CompoundExpression[n+=1,CompoundExpression[n+=1,f[3]]]->etc. Aunque la llamada a f es la última cada vez, sucede antes de que se ejecute el CompoundExpression[] completo, por lo que debe mantenerse en la pila de expresiones. Tal vez se podría argumentar que este es un lugar donde se podría hacer la optimización, para hacer una excepción para CompoundExpression, pero esto probablemente no es fácil de implementar.

Ahora, para ilustrar el proceso de acumulación de pila que describí esquemáticamente anteriormente, limitemos el número de llamadas recursivas:

Clear[n, f, ff, fff];
n = 0;
f[x_ /; x < 5] := (n += 1; f[x + 1]);

ff[x_] := Null /; (n += 1; False);
ff[x_ /; x < 5] := ff[x + 1];

fff[x_ /; x < 5] := ce[n += 1, fff[x + 1]];

Rastreando la evaluación:

In[57]:= Trace[f[1],f]
Out[57]= {f[1],n+=1;f[1+1],{f[2],n+=1;f[2+1],{f[3],n+=1;f[3+1],{f[4],n+=1;f[4+1]}}}}

In[58]:= Trace[ff[1],ff]
Out[58]= {ff[1],ff[1+1],ff[2],ff[2+1],ff[3],ff[3+1],ff[4],ff[4+1],ff[5]}

In[59]:= Trace[fff[1],fff]
Out[59]= {fff[1],ce[n+=1,fff[1+1]],{fff[2],ce[n+=1,fff[2+1]],{fff[3],ce[n+=1,fff[3+1]],   
{fff[4],ce[n+=1,fff[4+1]]}}}}

Lo que puede ver de esto es que la pila de expresiones se acumula para f y fff (este último se usa solo para mostrar que este es un mecanismo general, con ce[] solo una cabeza arbitraria), pero no para ff, porque, para los propósitos de coincidencia de patrones, la definición para ff es una regla probada pero no coincidente, y la segunda definición reescribe ff[arg_] en su totalidad, y no genera sub-partes más profundas que necesitan más reescritura. Por lo tanto, la línea de fondo parece que debe analizar su función y ver si sus llamadas recursivas crecerán la expresión evaluada desde la parte inferior o no. En caso afirmativo, no es recursivo de cola en lo que respecta a Mathematica.

Mi respuesta no estaría completa sin mostrar cómo hacer la llamada de cola optimización manual. Como ejemplo, consideremos la implementación recursiva de Select. Trabajaremos con listas enlazadas de Mathematica para hacerlo razonablemente eficiente en lugar de un juguete. A continuación se muestra el código para la implementación no recursiva de cola:

Clear[toLinkedList, test, selrecBad, sel, selrec, selTR]
toLinkedList[x_List] := Fold[{#2, #1} &, {}, Reverse[x]];
selrecBad[fst_?test, rest_List] := {fst,If[rest === {}, {}, selrecBad @@ rest]};
selrecBad[fst_, rest_List] := If[rest === {}, {}, selrecBad @@ rest];
sel[x_List, testF_] := Block[{test = testF}, Flatten[selrecBad @@ toLinkedList[x]]]

La razón por la que uso Block y selrecBad es para facilitar el uso de Trace. Ahora, esto sopla la pila en mi máquina:

Block[{$RecursionLimit = Infinity}, sel[Range[300000], EvenQ]] // Short // Timing

Puede rastrear en pequeñas listas para ver por qué:

In[7]:= Trace[sel[Range[5],OddQ],selrecBad]

Out[7]= {{{selrecBad[1,{2,{3,{4,{5,{}}}}}],{1,If[{2,{3,{4,{5,{}}}}}==={},{},selrecBad@@{2,{3,{4, 
{5,{}}}}}]},{selrecBad[2,{3,{4,{5,{}}}}],If[{3,{4,{5,{}}}}==={},{},selrecBad@@{3,{4,{5, 
{}}}}],selrecBad[3,{4,{5,{}}}],{3,If[{4,{5,{}}}==={},{},selrecBad@@{4,{5,{}}}]},{selrecBad[4,
{5,{}}],If[{5,{}}==={},{},selrecBad@@{5,{}}],selrecBad[5,{}],{5,If[{}==={},{},selrecBad@@{}]}}}}}}

Lo que sucede es que el resultado obtiene acumulado más y más en la lista. La solución es no aumentar la profundidad de la expresión resultante, y una forma de lograrlo es hacer que selrecBad acepte un parámetro adicional, que es la lista (vinculada) de resultados acumulados:

selrec[{fst_?test, rest_List}, accum_List] := 
    If[rest === {}, {accum, fst}, selrec[rest, {accum, fst}]];
selrec[{fst_, rest_List}, accum_List] := 
    If[rest === {}, accum, selrec[rest, accum]]

Y modificar la función principal en consecuencia:

selTR[x_List, testF_] := Block[{test = testF}, Flatten[selrec[toLinkedList[x], {}]]]

Esto pasará nuestra prueba de poder muy bien:

In[14]:= Block[{$IterationLimit= Infinity},selTR[Range[300000],EvenQ]]//Short//Timing

Out[14]= {0.813,{2,4,6,8,10,12,14,16,18,20,
<<149981>>,299984,299986,299988,299990,299992,299994,299996,299998,300000}}

(tenga en cuenta que aquí tuvimos que modificar Iter IterationLimit, que es una buena señal). Y el uso de Trace revela la razón:

In[15]:= Trace[selTR[Range[5],OddQ],selrec]

Out[15]= {{{selrec[{1,{2,{3,{4,{5,{}}}}}},{}],If[{2,{3,{4,{5,{}}}}}==={},{{},1},selrec[{2,{3,{4, 
{5,{}}}}},{{},1}]],selrec[{2,{3,{4,{5,{}}}}},{{},1}],If[{3,{4,{5,{}}}}==={},{{},1},selrec[{3, 
{4,{5,{}}}},{{},1}]],selrec[{3,{4,{5,{}}}},{{},1}],If[{4,{5,{}}}==={},{{{},1},3},selrec[{4, 
{5,{}}},{{{},1},3}]],selrec[{4,{5,{}}},{{{},1},3}],If[{5,{}}==={},{{{},1},3},selrec[{5, 
{}},{{{},1},3}]],selrec[{5,{}},{{{},1},3}],If[{}==={},{{{{},1},3},5},selrec[{},{{{{},1},3},5}]]}}}

Es decir, esta versión no acumula la profundidad de la expresión intermedia, ya que los resultados se mantienen en una lista separada.