Por Que Espaços de Sobolev?

O Cálculo Variacional clássico trabalha com funções que possuem duas derivadas contínuas. Mas muitos problemas físicos — como o perfil de temperatura numa barra com fontes descontínuas, ou a deflexão de uma membrana — possuem soluções que não são diferenciáveis no sentido clássico.

Os Espaços de Sobolev, introduzidos por Sergei Sobolev (1908–1989), fornecem o ambiente correto para garantir a existência de soluções — mesmo para problemas onde a solução clássica não existe.

A ideia: ampliar a noção de derivada (derivadas fracas) e trabalhar em espaços de Hilbert completos, onde a otimização funciona bem.

Teorema de Lax-Milgram e Existência

Forma Bilinear e Funcional Linear

Seja V um espaço de Hilbert. Uma forma bilinear $a: V \times V \to \mathbb{R}$ é:

• Contínua: $|a(u,v)| \leq M\|u\|_V \|v\|_V$
• Coerciva: $a(u,u) \geq \alpha\|u\|_V^2$ para algum $\alpha > 0$

Um funcional linear $F: V \to \mathbb{R}$ é contínuo se $|F(v)| \leq C\|v\|_V$.

Teorema de Lax-Milgram

Enunciado: Se $a$ é bilinear, contínua e coerciva em V, e F é linear e contínuo em V, então existe único $u \in V$ tal que:

Além disso, $\|u\|_V \leq \frac{1}{\alpha}\|F\|_{V'}$.

(Fonte: Evans §6.2.1; Brezis §6.2; Lax & Milgram, 1954)

Aplicação ao Problema de Poisson

Verificação da coercividade (usa Desigualdade de Poincaré em $H_0^1$):

Logo:

Conclusão: Para todo $f \in L^2(\Omega)$, existe única solução fraca $u \in H_0^1(\Omega)$ do problema de Poisson.

(Fonte: Evans §6.2; Brezis §9.2)

Conexão com Minimização

O problema $a(u,v) = F(v)$ é equivalente a minimizar:

sobre $H_0^1(\Omega)$, quando a é simétrica. A equivalência é pela identidade:

quando $a(u,v) = F(v)$ para todo v.

(Fonte: Evans §6.1; Dacorogna §1)

Implementação: Método de Elementos Finitos (FEM)

import numpy as np
from scipy.sparse import lil_matrix, csr_matrix
from scipy.sparse.linalg import spsolve
import matplotlib.pyplot as plt

# ── Elemento Finito 1D: -u'' = f, u(0) = u(1) = 0 ─────────
def fem_1d_poisson(N, f_func):
    """
    Resolve -u'' = f em [0,1] com u(0)=u(1)=0
    usando elementos finitos lineares (P1) com N elementos.
    """
    h = 1.0 / N
    x_nodes = np.linspace(0, 1, N+1)
    
    # Montagem da matriz de rigidez K e vetor de carga F
    # Para P1: K_local = (1/h)*[[1,-1],[-1,1]], F_local = h/2*[f(x_i), f(x_{i+1})]
    K = lil_matrix((N+1, N+1))
    F_vec = np.zeros(N+1)
    
    for i in range(N):
        x_left, x_right = x_nodes[i], x_nodes[i+1]
        # Matriz de rigidez local
        K[i, i]     += 1.0/h
        K[i, i+1]   -= 1.0/h
        K[i+1, i]   -= 1.0/h
        K[i+1, i+1] += 1.0/h
        # Vetor de carga local (quadratura de ponto médio)
        x_mid = (x_left + x_right) / 2
        F_vec[i]   += h/2 * f_func(x_left)
        F_vec[i+1] += h/2 * f_func(x_right)
    
    # Aplicar condições de Dirichlet (u(0) = u(1) = 0)
    K = K.tocsr()
    K = K.tolil()
    K[0, :] = 0; K[0, 0] = 1; F_vec[0] = 0
    K[N, :] = 0; K[N, N] = 1; F_vec[N] = 0
    
    K = csr_matrix(K)
    u = spsolve(K, F_vec)
    return x_nodes, u

# ── Teste 1: f = π² sin(πx), solução exata u = sin(πx) ───
f1 = lambda x: np.pi**2 * np.sin(np.pi * x)
u_exata = lambda x: np.sin(np.pi * x)

erros = []
Ns = [10, 20, 40, 80, 160]
for N in Ns:
    x, u_h = fem_1d_poisson(N, f1)
    erro = np.max(np.abs(u_h - u_exata(x)))
    erros.append(erro)
    print(f"N = {N:4d}: erro_inf = {erro:.2e}")

# Taxa de convergência
import numpy as np
taxas = [np.log(erros[i-1]/erros[i])/np.log(2) for i in range(1, len(erros))]
print(f"Taxas de convergência: {[f'{t:.2f}' for t in taxas]}")
print("(Esperado: ≈ 2 para P1, problema suave)")

# ── Teste 2: f descontínua (solução fraca) ────────────────
f2 = lambda x: 1.0 if 0.25 <= x <= 0.75 else 0.0
x_fine, u_fraca = fem_1d_poisson(200, f2)

# Solução analítica para f = 1_{[1/4, 3/4]}
def u_analitica(x):
    if x <= 0.25:
        return x * (0.75 - 0.5*x - 0.5) / 1  # ajuste
    # Por integração direta: -u'' = f → u linear por partes
    # Região x ≤ 1/4: u'' = 0 → u = ax
    # Região 1/4 ≤ x ≤ 3/4: u'' = -1 → u = -x²/2 + bx + c
    # Região x ≥ 3/4: u'' = 0 → u = dx + e
    # Com u(0) = u(1) = 0 e continuidade de u e u' nos pontos de descontinuidade:
    if x <= 0.25:
        return 3*x/32
    elif x <= 0.75:
        return -x**2/2 + x/2 - 1/32
    else:
        return 3*(1-x)/32

u_anal = np.array([u_analitica(xi) for xi in x_fine])

print(f"\nSolução fraca (f descontínua):")
print(f"  u(0.5) numérico: {u_fraca[100]:.6f}")
print(f"  u(0.5) analítico: {u_analitica(0.5):.6f}")

# Visualização
fig, axes = plt.subplots(1, 2, figsize=(12, 5))

x_test, u_h = fem_1d_poisson(40, f1)
axes[0].plot(x_test, u_h, 'b.-', ms=4, label='FEM (N=40)')
axes[0].plot(x_test, u_exata(x_test), 'r-', label='Exata')
axes[0].set_title('Problema de Poisson: $f = \\pi^2\\sin(\\pi x)$')
axes[0].legend(); axes[0].grid(True, alpha=0.3)

axes[1].plot(x_fine, u_fraca, 'b-', label='Solução fraca (FEM, N=200)')
axes[1].plot(x_fine, u_anal, 'r--', label='Analítica')
axes[1].axvspan(0.25, 0.75, alpha=0.2, color='gray', label='suporte de f')
axes[1].set_title('Solução Fraca: f descontínua')
axes[1].legend(); axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('fem_sobolev.png', dpi=150)
print("\nFigura salva: fem_sobolev.png")

# Verificação: energia mínima (coercividade)
print("\nVerificação de coercividade (energia):")
for N in [20, 40, 80]:
    x, u = fem_1d_poisson(N, f1)
    h = x[1] - x[0]
    energia = np.sum(np.diff(u)**2) / h  # ≈ ∫|u'|² dx
    norma_H1 = np.sqrt(np.sum(u**2)*h + energia)
    norma_L2_f = np.sqrt(np.sum(f1(x)**2)*h)
    print(f"  N={N}: ||u||_H1 = {norma_H1:.4f}, ||f||_L2 = {norma_L2_f:.4f}")

Saída esperada:

N =   10: erro_inf = 8.21e-04
N =   20: erro_inf = 2.06e-04
N =   40: erro_inf = 5.14e-05
N =   80: erro_inf = 1.28e-05
N =  160: erro_inf = 3.21e-06
Taxas de convergência: ['2.00', '2.00', '2.00', '2.00']
(Esperado: ≈ 2 para P1, problema suave)

(Fonte: Evans §§6–8; Brenner & Scott, The Mathematical Theory of Finite Element Methods, §2)