Stable $C^1$-conforming finite element methods for a class of nonlinear fourth-order evolution equations

Resumen Técnico: Métodos de Elementos Finitos Estables con Conformidad $C^1$ para una Clase de Ecuaciones de Evolución Cuarto Orden No Lineales

Formulación del Problema
El artículo aborda el análisis numérico de una clase general de ecuaciones diferenciales parciales (EDP) de evolución de cuarto orden no lineales definidas en un dominio acotado $\Omega \subset \mathbb{R}^d$ ($d \in \{1, 2, 3\}$). La ecuación gobernante viene dada por:
$$\frac{\partial u}{\partial t} - \beta_1 \Delta u + \beta_2 \Delta^2 u = \beta_3(1 - |u|^2)u - \beta_4 u \times \Delta u + \beta_5 \Delta(|u|^2 u) + \beta_6 (j \cdot \nabla)u$$
sujeta a condiciones de contorno de Neumann ($\frac{\partial u}{\partial n} = 0, \frac{\partial (\Delta u)}{\partial n} = 0$) y una condición inicial $u(0) = u_0$. Aquí, $u$ puede ser escalar ($m=1$) o vectorial ($m=3$). La ecuación abarca varios modelos físicos dependiendo de los coeficientes $\beta_i$ y la dimensión $d$, incluyendo:

• Micromagnetismo: La ecuación de Landau–Lifshitz–Baryakhtar (LLB) (con o sin torque de espín) y regularizaciones de la ecuación de Landau–Lifshitz–Bloch.
• Formación de Patrones: La ecuación de Swift–Hohenberg.
• Separación de Fases y Difusión: Ecuaciones de Cahn–Hilliard convectivas con términos fuente, modelos de difusión de biflujo anómalo y modelos de dispersión de población generalizados.

Los autores señalan que, si bien existen métodos numéricos para instancias específicas de estas ecuaciones (por ejemplo, Swift–Hohenberg escalar o Cahn–Hilliard estándar), faltaba en la literatura un esquema numérico unificado y general para el sistema completo (1.1) con todos los coeficientes potencialmente distintos de cero.

Metodología
Los autores proponen y analizan métodos de elementos finitos totalmente discretos utilizando elementos con conformidad $C^1$. Esta elección es necesaria debido a la derivada espacial de cuarto orden ($\Delta^2 u$), que requiere que el espacio de elementos finitos sea un subespacio de $H^2(\Omega)$.

• Discretización Espacial: El método emplea polinomios por trozos de grado $r$ (donde $r=3$ para $d=1,2$ y $r=9$ para $d=3$ para mantener la continuidad $C^1$). Los elementos específicos mencionados incluyen elementos Hermite $C^1$-P3 (1D), macroelementos de Hsieh–Clough–Tocher (HCT) (2D) y elementos de Zenísek (3D).
• Discretización Temporal: Se proponen dos esquemas de paso de tiempo linealizados:
  1. Un método de Euler semi-implícito (primer orden en el tiempo).
  2. Un método BDF2 semi-implícito (segundo orden en el tiempo).
     En ambos esquemas, los términos no lineales se linealizan evaluándolos en el paso de tiempo anterior, lo que resulta en un sistema de ecuaciones lineales en cada paso. Esto evita el costo computacional de resolver sistemas no lineales mediante el método de Newton.

Contribuciones Clave y Resultados Teóricos
El artículo proporciona un análisis de error riguroso y una prueba de estabilidad para los esquemas propuestos:

1. Estabilidad:

   • Se demuestra que los esquemas son incondicionalmente estables en la norma $H^2$ para las dimensiones $d=1$ y $d=2$. Esto significa que la estabilidad se mantiene sin restricciones en el tamaño del paso de tiempo $k$ respecto al tamaño de la malla $h$, y sin asumir triangulaciones cuasi-uniformes.
   • Para $d=3$, los esquemas son condicionalmente estables. La estabilidad está garantizada ya sea si los datos iniciales son suficientemente pequeños o si el paso de tiempo satisface una restricción específica relativa a la energía inicial. Los autores reconocen que los elementos con conformidad $C^1$ en 3D (que requieren polinomios de grado 9) son computacionalmente prohibitivos, sugiriendo que los métodos no conformes son más prácticos para $d=3$, aunque el análisis teórico se incluye por completitud.

2. Estimaciones de Error y Convergencia:

   • Los autores establecen tasas de convergencia óptimas en las normas $L^\infty(0, T; H^\beta)$ para $\beta = 0, 1, 2$.
   • Para el esquema semi-discreto, el error está acotado por $O(h^{r+1-\beta})$.
   • Para los esquemas totalmente discretos:
     • El método de Euler semi-implícito logra un error de $O(h^{r+1-\beta} + k)$.
     • El método BDF2 semi-implícito logra un error de $O(h^{r+1-\beta} + k^2)$.
   • El análisis se basa en el desarrollo de una proyección elíptica especializada adaptada al problema no lineal y en la estimación cuidadosa de los términos no lineales utilizando inserciones de Sobolev e desigualdades inversas.

3. Experimentos Numéricos:

   • Las tasas de convergencia teóricas fueron verificadas utilizando el paquete de código abierto FreeFEM con elementos HCT en 2D.
   • Los experimentos se realizaron sobre la ecuación de Landau–Lifshitz–Baryakhtar (valor vectorial) y la ecuación de Cahn–Hilliard con un término fuente (valor escalar).
   • Los resultados confirmaron los órdenes de convergencia predichos: 4º orden en $L^2$, 3º orden en $H^1$ y 2º orden en $H^2$ para la discretización espacial, consistente con el uso de elementos cúbicos ($r=3$).

Significancia y Reivindicaciones
El artículo afirma proporcionar el primer esquema numérico capaz de aproximar la solución al sistema general de cuarto orden no lineal (1.1) en toda su generalidad. Las principales ventajas destacadas son:

• Linealidad: Al linealizar los términos no lineales, el método evita la complejidad y la potencial inestabilidad de resolver sistemas no lineales en cada paso de tiempo, lo que resulta ventajoso para aplicaciones como las simulaciones micromagnéticas.
• Estabilidad: La estabilidad incondicional en $H^2$ para $d=1,2$ asegura un rendimiento robusto sin restricciones estrictas del paso de tiempo.
• Precisión: El uso de elementos con conformidad $C^1$ permite una mayor precisión en la captura de los gradientes de la solución en comparación con los métodos no conformes, a pesar del mayor costo de la continuidad $C^1$.

Los autores mantienen un tono modesto respecto al caso 3D, declarando explícitamente que, aunque el análisis teórico está completo, la aplicación práctica de elementos con conformidad $C^1$ en 3D es actualmente computacionalmente inviable debido al alto grado polinómico requerido. El trabajo sirve como un marco teórico fundacional para resolver ecuaciones de evolución de cuarto orden no lineales complejas que surgen en la física, la ciencia de materiales y la biología.