Cálculo fraccionario, caos y ABM–PECE

El cálculo fraccionario amplía la idea de derivar e integrar a órdenes no enteros, pero no consiste en sustituir mecánicamente un exponente. En una ecuación de Caputo, el cambio en $t$ depende de toda la historia previa mediante un núcleo de ley de potencia. Esa no localidad altera tanto la interpretación física como el costo numérico. Antes de simular conviene dominar tres elementos: la función gamma, las integrales de potencia y la diferencia entre Riemann–Liouville y Caputo.

La función gamma extiende el factorial:

$$\Gamma(z)=\int_0^\infty s^{z-1}e^{-s}\,ds, \qquad \Gamma(z+1)=z\Gamma(z).$$

La integral de Riemann–Liouville de orden $\alpha>0$ es

$$({}^{RL}I_{0+}^{\alpha}f)(t)= \frac1{\Gamma(\alpha)}\int_0^t(t-\tau)^{\alpha-1}f(\tau)\,d\tau.$$

Para $\beta>-1$,

$${}^{RL}I_{0+}^{\alpha}t^\beta= \frac{\Gamma(\beta+1)}{\Gamma(\beta+\alpha+1)}t^{\beta+\alpha}.$$

La fórmula se deriva insertando $\tau=t u$ en la integral:

$${}^{RL}I_{0+}^{\alpha}t^\beta =\frac{t^{\alpha+\beta}}{\Gamma(\alpha)} \int_0^1(1-u)^{\alpha-1}u^\beta\,du.$$

La integral restante es $B(\alpha,\beta+1)=\Gamma(\alpha)\Gamma(\beta+1)/\Gamma(\alpha+\beta+1)$. Así aparece el cociente gamma sin tratar el orden fraccionario como un factorial informal. Das desarrolla gamma, beta, operadores y ecuaciones en esa secuencia. Das, 2020 Shantanu Das (2020) Kindergarten of Fractional Calculus Cambridge Scholars Publishing Ubicación consultada: caps. 1–3

Cuando $0<q<1$, la derivada de Caputo se define por

$${}^{C}D_{0+}^{q}x(t)= \frac1{\Gamma(1-q)}\int_0^t\frac{x'(\tau)}{(t-\tau)^q}\,d\tau.$$

Caputo anula constantes y admite la condición inicial clásica $x(0)=x_0$. En cambio, la derivada de Riemann–Liouville de la constante vale $t^{-q}/\Gamma(1-q)$. Esta diferencia determina el tipo de datos iniciales y la interpretación física.

Para una función suficientemente regular y $0<q<1$, ambas derivadas se relacionan por

$${}^CD_{0+}^q x(t)= {}^{RL}D_{0+}^q\bigl[x(t)-x(0)\bigr].$$

Por eso coinciden sobre potencias que se anulan en el origen y difieren sobre la parte constante. En modelos dinámicos, Caputo resulta cómodo porque usa datos iniciales expresados mediante derivadas de orden entero. Esa conveniencia no vuelve intercambiables las definiciones: siempre se debe declarar operador, límite inferior y orden. Pham y colaboradores resumen las definiciones de Grünwald–Letnikov, Riemann–Liouville y Caputo antes de aplicarlas a sistemas biológicos fraccionarios. Pham, 2018 Viet-Thanh Pham, Sundarapandian Vaidyanathan, Christos Volos y Tomasz Kapitaniak, editores (2018) Nonlinear Dynamical Systems with Self-Excited and Hidden Attractors Springer Ubicación consultada: pp. 3–12, definiciones y formulación fraccionaria

La comparación de potencias se obtiene directamente de la identidad gamma. Para la constante, Riemann–Liouville conserva una singularidad $t^{-q}$ y Caputo da cero; para $t^2$, ambos operadores coinciden bajo las hipótesis usuales y producen $\Gamma(3)t^{2-q}/\Gamma(3-q)$. La escala logarítmica permite ver simultáneamente el comportamiento cerca del origen y a tiempos mayores.

Núcleo de memoria y dependencia histórica

En

$${}^CD_{0+}^qx(t)=\frac1{\Gamma(1-q)} \int_0^t(t-\tau)^{-q}x'(\tau)\,d\tau,$$

la edad de una contribución es $a=t-\tau$. El peso $a^{-q}$ disminuye al alejarse del presente, pero nunca se corta exactamente. Si $q$ se acerca a uno, el núcleo se concentra con mayor fuerza cerca del instante actual; para órdenes menores, el historial lejano conserva una participación relativa mayor. Herrmann interpreta esta no localidad mediante transformadas y memoria en modelos físicos. Herrmann, 2014 Richard Herrmann (2014) Fractional Calculus: An Introduction for Physicists 2.ª ed. · World Scientific Ubicación consultada: caps. 5–6 y 8

La singularidad en $\tau=t$ es integrable porque $0<q<1$. No conviene evaluarla como un valor puntual; los esquemas numéricos integran el núcleo por intervalos y producen pesos finitos. Tampoco debe decirse que «todo el pasado pesa igual»: cada tramo recibe un peso distinto y la normalización gamma cambia con $q$.

Si el orden fuera $m-1<q<m$, la definición de Caputo usaría la derivada entera $x^{(m)}$ y requeriría datos $x(0),x'(0),\ldots,x^{(m-1)}(0)$. Aquí nos concentramos en $0<q<1$ para que el problema inicial use un solo valor.

Relajación fraccionaria

Considera

$${}^{C}D_t^q x(t)=-\lambda x(t),\qquad x(0)=x_0,\quad 0<q\le1.$$

La transformada de Laplace de Caputo es

$$\mathcal L\{{}^{C}D_t^q x\}(s)=s^qX(s)-s^{q-1}x_0.$$

Entonces

$$X(s)=x_0\frac{s^{q-1}}{s^q+\lambda},$$

y la inversión produce

$$x(t)=x_0E_q(-\lambda t^q), \qquad E_q(z)=\sum_{k=0}^\infty\frac{z^k}{\Gamma(qk+1)}.$$

Para $q=1$, $E_1(z)=e^z$. Cuando $q<1$, la cola decae más lentamente y la historia se distribuye en escalas temporales de ley de potencia.

La inversión puede comprobarse término a término. Como

$$E_q(-\lambda t^q)= \sum_{k=0}^\infty\frac{(-\lambda)^kt^{qk}}{\Gamma(qk+1)},$$

la derivada de Caputo de cada potencia con $k\ge1$ es

$${}^CD_t^q\frac{t^{qk}}{\Gamma(qk+1)} =\frac{t^{q(k-1)}}{\Gamma(q(k-1)+1)}.$$

Al reindexar se obtiene $-\lambda E_q(-\lambda t^q)$, mientras la constante inicial queda anulada. Para $0<q<1$ y $\lambda>0$, la relajación no es una suma finita de exponenciales y su cola es mucho más larga. La solución exacta funciona a la vez como interpretación y como referencia para verificar el integrador.

Oliveira reúne ejercicios resueltos de potencias y transformadas que proporcionan controles algebraicos para estas identidades y sus límites de orden entero. Oliveira, 2019 Edmundo Capelas de Oliveira (2019) Solved Exercises in Fractional Calculus Springer Ubicación consultada: caps. 2–5

Formulación integral y esquema predictor–corrector

Aplicar la integral fraccionaria a ${}^CD_t^qy=f(t,y)$ produce la ecuación de Volterra

$$y(t)=y_0+\frac1{\Gamma(q)} \int_0^t(t-\tau)^{q-1}f(\tau,y(\tau))\,d\tau.$$

Esta forma explica por qué cada paso necesita el historial. El predictor aproxima $f$ por valores previos en cada subintervalo; el corrector añade la evaluación en el extremo nuevo y usa una cuadratura de orden superior.

Para ${}^{C}D_t^q y=f(t,y)$, $y(0)=y_0$, el predictor usa

$$y_{n+1}^{P}=y_0+\frac{h^q}{\Gamma(q+1)} \sum_{j=0}^{n}b_{j,n+1}f(t_j,y_j),$$

con $b_{j,n+1}=(n+1-j)^q-(n-j)^q$. El corrector evalúa $f$ en el estado predicho y combina la historia con pesos de Adams–Moulton, que son segundas diferencias de potencias de exponente $q+1$. La secuencia PECE significa predecir, evaluar, corregir y volver a evaluar. Cada paso consulta todos los anteriores; una aproximación de historia declarada controla el crecimiento de costo y memoria.

En la forma clásica para $0<q<1$, los pesos interiores del corrector pueden escribirse

$$a_{j,n+1}=(n-j+2)^{q+1}-2(n-j+1)^{q+1}+(n-j)^{q+1},$$

con un peso de borde tratado aparte. La actualización combina esos términos con $f(t_{n+1},y_{n+1}^P)$ y el factor $h^q/\Gamma(q+2)$. Los índices cambian entre convenciones de implementación; una prueba debe comparar la fórmula programada contra un problema de solución conocida, no sólo contra otra gráfica del mismo código.

Diethelm, Ford y Freed derivan y analizan un predictor–corrector de Adams para ecuaciones diferenciales fraccionarias. Es la referencia directa del algoritmo; Sun, He y Wang lo estudian dentro de sistemas caóticos; Das, Herrmann, Oliveira y Pham aportan operadores, memoria, cálculo y otros contextos dinámicos. Mantener esa distinción separa procedencia matemática, aplicación y material de práctica.

Para verificar el integrador se usa la relajación ${}^{C}D_t^q x=-x$, cuya solución $E_q(-t^q)$ puede evaluarse por serie. La figura superpone varias mallas ABM–PECE y calcula el error máximo contra esa referencia. La pendiente de la gráfica logarítmica aporta una estimación empírica del orden en el intervalo estudiado; una demostración de convergencia establece el resultado general.

La serie de Mittag–Leffler también necesita control numérico: se trunca cuando los términos dejan de cambiar el resultado bajo la precisión solicitada y se contrasta, si es posible, con una implementación independiente. Usar la misma rutina defectuosa para generar «referencia» y solución ocultaría el error.

Interpretación de la relajación y los pesos históricos

En el núcleo $(t-\tau)^{-q}$, los instantes recientes reciben más peso y el pasado lejano conserva una contribución. Compara curvas con la misma condición y escala. En relajación, una cola alta para $q<1$ indica decaimiento más lento que el exponencial. En la animación de pesos, comprueba que al avanzar un paso aparece una contribución nueva y se vuelven a ponderar todas las anteriores.

Con $N$ pasos y suma directa del historial, el trabajo total crece como $O(N^2)$ y el almacenamiento como $O(N)$. Reducir $h$ a la mitad no sólo duplica pasos: puede multiplicar aproximadamente por cuatro las operaciones de memoria. Existen técnicas de memoria corta, sumas rápidas y aproximaciones del núcleo, pero cambian el algoritmo y deben declarar tolerancia y ventana. En este cálculo la referencia es el historial completo.

Para una tabla de convergencia conserva $q$, intervalo, solución de referencia y norma. Calcula

$$e_h=\max_n|y_h(t_n)-y_{\mathrm{ref}}(t_n)|, \qquad p_{\mathrm{obs}}=\log_2(e_h/e_{h/2}).$$

No extrapoles el orden observado en una sola ecuación a todos los problemas: la regularidad de la solución cerca de $t=0$ influye especialmente en ecuaciones fraccionarias.

Equilibrio, estabilidad y trayectoria en caos fraccionario

Un sistema autónomo de Caputo puede escribirse componente a componente:

$${}^CD_t^{q_i}x_i=f_i(\mathbf x), \qquad 0<q_i\le1,\quad i=1,\ldots,d.$$

Es conmensurable si todos los órdenes son iguales, $q_1=\cdots=q_d=q$; es inconmensurable si al menos dos difieren. El punto de equilibrio satisface

$$\mathbf f(\mathbf x^*)=\mathbf0.$$

Esa ecuación algebraica coincide con la del modelo entero porque la derivada de Caputo de una constante es cero. Lo que cambia con los órdenes es la estabilidad y la evolución alrededor del equilibrio, no su localización.

Para el sistema lineal conmensurable

$${}^CD_t^q\boldsymbol\xi=A\boldsymbol\xi, \qquad 0<q<1,$$

el equilibrio es asintóticamente estable si cada autovalor satisface

$$|\arg\lambda_j|>\frac{q\pi}{2}.$$

Cuando $q=1$, la condición recupera el semiplano izquierdo. Para $q<1$, la región estable excluye un sector más estrecho alrededor del eje real positivo. Este criterio pertenece al sistema lineal con orden común. Aplicarlo sin modificación a órdenes inconmensurables, retardos o aproximaciones racionales del operador mezcla problemas distintos.

Sistema de Lorenz fraccionario

Considera

$$\begin{aligned} {}^CD_t^q x&=\sigma(y-x),\\ {}^CD_t^q y&=x(\rho-z)-y,\\ {}^CD_t^q z&=xy-\beta z. \end{aligned}$$

Los equilibrios son

$$O=(0,0,0)$$

y, si $\rho>1$,

$$C_\pm= \left( \pm\sqrt{\beta(\rho-1)}, \pm\sqrt{\beta(\rho-1)}, \rho-1 \right).$$

Son los mismos del sistema entero. En el origen,

$$J(O)= \begin{pmatrix} -\sigma&\sigma&0\\ \rho&-1&0\\ 0&0&-\beta \end{pmatrix},$$

con un autovalor $-\beta$ y dos raíces de

$$\lambda^2+(\sigma+1)\lambda+\sigma(1-\rho)=0.$$

Para $\sigma=10$, $\rho=28$ y $\beta=8/3$, el término constante de ese polinomio es negativo; una raíz es real positiva. Su argumento es cero, así que viola el criterio sectorial para todo $q>0$. El origen continúa inestable. Esta cuenta no demuestra caos: sólo clasifica localmente un equilibrio.

En $C_\pm$ se repite el proceso con el jacobiano evaluado en cada punto. Si algún autovalor cruza la frontera angular al variar $q$, cambia el diagnóstico lineal. Después todavía hay que estudiar trayectorias. Sun, He y Wang usan el principio de estabilidad fraccionaria como filtro y muestran que el «orden mínimo de caos» obtenido numéricamente depende del método, paso, parámetros y órdenes asignados Sun, 2022 Kehui Sun, Shaobo He y Huihai Wang (2022) Solution and Characteristic Analysis of Fractional-Order Chaotic Systems Springer Nature Singapore y Science Press Ubicación consultada: cap. 5, §5.3, pp. 67–74 Abrir fuente .

Protocolo de evidencia para barridos del orden

Una proyección con forma de mariposa no basta. Para cada valor de $q$ se conserva la siguiente cadena:

1. definición del operador, límite inferior, vector de órdenes y datos iniciales;
2. equilibrios, jacobianos y alcance exacto del criterio sectorial aplicado;
3. método numérico, paso, longitud de historial y tolerancias;
4. comparación con una malla más fina y, cuando sea posible, otro algoritmo;
5. transitorio descartado, prueba de acotación y horizonte retenido;
6. serie, sección o diagrama de bifurcación junto con un indicador cuantitativo;
7. espectro de Lyapunov calculado con un procedimiento compatible con memoria, incluyendo tiempo finito y prueba de convergencia;
8. repetición al cambiar paso, horizonte y condición inicial.

El resultado se informa como una clasificación bajo ese protocolo. Un exponente máximo positivo y estable al refinar aporta evidencia de sensibilidad; la coincidencia con bifurcaciones o una prueba 0–1 independiente refuerza la clasificación. Si el signo cambia con $h$ o con el horizonte, no se anuncia un umbral de caos. Los capítulos de dinámica y complejidad de Sun, He y Wang combinan espectros de Lyapunov, diagramas de bifurcación, prueba 0–1 y medidas de complejidad, y su comparación numérica documenta la dependencia respecto del algoritmo Sun, 2022 Kehui Sun, Shaobo He y Huihai Wang (2022) Solution and Characteristic Analysis of Fractional-Order Chaotic Systems Springer Nature Singapore y Science Press Ubicación consultada: cap. 6, §§6.1–6.2, pp. 77–114; cap. 7, §§7.1–7.2, pp. 117–140 Abrir fuente .

Afirmación	Comprobación numérica dentro del protocolo	Lo que todavía no prueba
equilibrio inestable	jacobiano y criterio sectorial aplicable	existencia de un atractor
órbita acotada	control de norma en horizonte y mallas refinadas	caos
sensibilidad numérica	Lyapunov finito positivo y convergente	resultado asintótico riguroso
ventana caótica	acuerdo entre al menos dos diagnósticos y refinamiento	universalidad fuera del protocolo

Comparación del orden fraccionario con un protocolo fijo

En el módulo de orden fraccionario, selecciona un sistema compatible y el método ABM predictor–corrector. Usa la configuración de Lorenz $(10,28,8/3)$, $(1,1,1)$, $h=0.005$, tiempo $80$ y el mismo historial inicial. Ejecuta primero $q=(1,1,1)$ y luego órdenes conmensurables $q=0.98$ y $q=0.95$.

Compara series, proyección $xz$ y costo de cálculo. Repite $q=0.98$ con $h=0.0025$. Atribuye diferencias al orden después de comprobar paso, duración y transitorio. Registra si el software conserva historia completa, la convención de Caputo y la forma de inicialización.

El caso $q=(1,1,1)$ funciona como límite entero y control de consistencia, pero no basta comparar su órbita punto a punto durante $80$ unidades: la sensibilidad puede separar ejecuciones. Primero verifica corto plazo contra RK4 con refinamiento; después compara acotación, proyecciones, estadísticas y costo. Para órdenes distintos por componente, declara la convención exacta del sistema y los datos iniciales de cada ecuación.

Una disminución de $q$ no implica automáticamente «más estabilidad» ni «menos caos». El orden modifica memoria, escalas y condiciones de estabilidad de manera dependiente del modelo. En los sistemas biológicos compilados por Pham y colaboradores, la formulación fraccionaria se estudia junto con equilibrios, bifurcaciones y control; las conclusiones se derivan para cada sistema, no de una regla universal sobre $q$. Pham, 2018 Viet-Thanh Pham, Sundarapandian Vaidyanathan, Christos Volos y Tomasz Kapitaniak, editores (2018) Nonlinear Dynamical Systems with Self-Excited and Hidden Attractors Springer Ubicación consultada: pp. 3–44, modelos fraccionarios y análisis dinámico

Verificación del método con la identidad gamma

Riemann–Liouville y Caputo fijan qué significa derivar; la relajación de Mittag–Leffler muestra la consecuencia dinámica de la memoria; ABM–PECE aproxima la forma integral mediante una suma histórica cuya procedencia, convergencia y costo deben quedar documentados.