¿Cómo se resuelve la siguiente ecuación diferencial?
$$2 x''' + xx'' =0$$
Sujeto a condiciones:
$$ x(0)=0$$ $$ x'(0)=0$$ $$ x'(\infty)=1$$
¿Existe algún método numérico para resolverlo o algún método general?
¿Cómo se resuelve la siguiente ecuación diferencial?
$$2 x''' + xx'' =0$$
Sujeto a condiciones:
$$ x(0)=0$$ $$ x'(0)=0$$ $$ x'(\infty)=1$$
¿Existe algún método numérico para resolverlo o algún método general?
Numéricamente no parece haber ningún problema, $T=20$ como aproximación de $T=\infty$ parece lo suficientemente grande, la solución puede calcularse a través del solucionador de valores límite, aquí python
's scipy.integrate.solve_bvp
:
T = 20
def x_ode(t,x): return [x[1], x[2], -0.5*x[0]*x[2]]
def x_bc(x0, xT): return [x0[0], x0[1], xT[1]-1]
s = np.linspace(0,1,11);
t_init = T*s
x_init = [ T*s , 1+0*s, 0+0*s ]
res = solve_bvp(x_ode, x_bc, t_init, x_init, tol=1e-5, max_nodes=80000)
print res.message
con el resultado
The algorithm converged to the desired accuracy.
El gráfico de la función y la derivada
se produce entonces a través de
if res.success:
plt.figure(figsize=(9,6))
plt.subplot(211); plt.plot(res.x, res.y[0], '-'); plt.grid();
plt.subplot(212); plt.plot(res.x, res.y[1], '-'); plt.grid();
plt.show()
Para transformar el problema a un intervalo finito las condiciones de contorno sugieren utilizar $s=x'$ , $s\in[0,1]$ como nuevo parámetro independiente, utilizando así la función inversa $u$ a $x'$ , $t=u(x')$ , $x=v_0(x')$ , $v_1(s)=x'(u(s))=s$ , $v_2(s)=x''(u(s))$ lo que lleva a las derivadas \begin{align} 1=\frac{dt}{dt}&=u'(x'(t))x''(t)\\ \frac{dt}{ds}&=u'(s)=\frac1{v_2(s)}\\ \frac{dv_0(s)}{ds}&=x'(u(s))u'(s)=\frac{s}{v_2(s)}\\ \frac{dv_2(s)}{ds}&=x'''(u(s))u'(s)=-\frac{v_0(s)}{2} \end{align} Para conseguir $v_0(s)\to\infty$ para $s\to 1$ necesitamos $v_2(1)=0$ junto con $u(0)=0$ y $v_0(0)=0$ . Para desingularizar $1/v_2$ el enfoque utilizando $v_2/(\epsilon^2+v_2^2)$ funciona mejor.
eps = 1e-8
def uv_ode(s,y):
u, v0, v2 = y;
v2inv = v2/(eps**2+v2**2);
return [v2inv, s*v2inv, -0.5*v0]
def uv_bc(u0, u1):
return [ u0[0], u0[1], u1[2]]
s = np.linspace(0,1,11);
y_init = [ 16*s**2 , 14*s**2, (1-s)**2/4 ]
res = solve_bvp(uv_ode, uv_bc, s, y_init, tol=1e-5, max_nodes=80000)
print res.message
if True or res.success:
plt.figure(figsize=(9,9))
plt.subplot(311); plt.plot(res.x, res.y[0], '-o',ms=1); plt.ylabel('$u=t$'); plt.grid();
plt.subplot(312); plt.plot(res.x, res.y[1], '-o',ms=1); plt.ylabel('$v_0=x$'); plt.grid();
plt.subplot(313); plt.plot(res.x, res.y[2], '-o',ms=1); plt.ylabel('$v_2=x\'\'$'); plt.grid(); plt.xlabel('$s=x\'$');
plt.show()
Demasiado largo para un comentario.
En su respuesta, JJacquelin terminó con $$t(x)=\int_0^x\frac{d\xi}{\sqrt{\text{erf}(\xi/2)}}$$ que no tiene solución explícita.
Sin embargo, podemos hacer una aproximación bastante buena construyendo el aproximante de Padé en $\xi=0$ $$\frac{1}{\sqrt{\text{erf}(\xi/2)}}\simeq \sqrt[4]\pi \, \frac{1+\frac{363 }{3680}\xi ^2+\frac{491 }{264960}\xi ^4-\frac{23789 }{445132800}\xi^6 }{\sqrt \xi \left(1+\frac{629 }{11040}\xi ^2 \right) }$$ que pueden integrarse para llegar a la desagradable (pero practicable) $$t(x)=-\frac{23789 \sqrt[4]{\pi } x^{9/2}}{114125760}+\frac{1085389 \sqrt[4]{\pi } x^{5/2}}{55389740}+\frac{3036269046 \sqrt[4]{\pi } \sqrt{x}}{1742007323}-\frac{111936400\ 2^{3/4} \sqrt[4]{\frac{345 \pi }{629}} \log \left(\sqrt{434010} x-4\ 345^{3/4} \sqrt[4]{1258} \sqrt{x}+2760\right)}{1742007323}+\frac{111936400\ 2^{3/4} \sqrt[4]{\frac{345 \pi }{629}} \log \left(\sqrt{434010} x+4\ 345^{3/4} \sqrt[4]{1258} \sqrt{x}+2760\right)}{1742007323}-\frac{223872800\ 2^{3/4} \sqrt[4]{\frac{345 \pi }{629}} \tan ^{-1}\left(1-\frac{\sqrt[4]{\frac{629}{345}} \sqrt{x}}{2^{3/4}}\right)}{1742007323}+\frac{223872800\ 2^{3/4} \sqrt[4]{\frac{345 \pi }{629}} \tan ^{-1}\left(\frac{\sqrt[4]{\frac{629}{345}} \sqrt{x}}{2^{3/4}}+1\right)}{1742007323}$$ que no es tan malo como se muestra a continuación $$\left( \begin{array}{ccc} x & \text{approximation} & \text{exact} \\ 0.00 & 0.00000 & 0.00000 \\ 0.25 & 1.33203 & 1.33203 \\ 0.50 & 1.88671 & 1.88671 \\ 0.75 & 2.31671 & 2.31671 \\ 1.00 & 2.68470 & 2.68470 \\ 1.25 & 3.01528 & 3.01528 \\ 1.50 & 3.32122 & 3.32122 \\ 1.75 & 3.61020 & 3.61020 \\ 2.00 & 3.88724 & 3.88725 \\ 2.25 & 4.15580 & 4.15584 \\ 2.50 & 4.41833 & 4.41843 \\ 2.75 & 4.67656 & 4.67681 \\ 3.00 & 4.93172 & 4.93227 \\ 3.25 & 5.18460 & 5.18576 \\ 3.50 & 5.43566 & 5.43794 \\ 3.75 & 5.68507 & 5.68926 \\ 4.00 & 5.93273 & 5.94004 \\ 4.25 & 6.17831 & 6.19048 \\ 4.50 & 6.42126 & 6.44074 \\ 4.75 & 6.66082 & 6.69087 \\ 5.00 & 6.89601 & 6.94094 \end{array} \right)$$
El primer paso en un "método general" para resolver este problema es reconocer que la ecuación diferencial ordinaria (EDO) de tercer orden dada es autónoma.
Usaremos subíndices para las derivadas en lo que sigue, porque vamos a usar derivadas con respecto a la variable independiente $t$ pero también con respecto a la variable dependiente $x$ .
A partir de $2 x_{ttt} + x x_{tt} = 0$ escribimos $x_t = w(x)$ con alguna función desconocida $w$ para obtener un problema de valores límite (BVP) con una EDO de segundo orden para $w$ : \begin{equation} w w_{xx} + w_x \left( w_x + \frac{x}{2} \right) = 0, \quad w(0) = 0, \quad w(\infty) = 1 \end{equation} (nótese que las derivadas de $w$ se toman con respecto a $x$ !).
Una vez que la función $w$ podemos resolver el siguiente problema de valor inicial (PIV) con una EDO de primer orden separable para $x$ : \begin{equation} x_t = w(x), \quad x(0) = 0. \end{equation} Así, hemos dividido el problema de valores iniciales (IBVP) con una EDO autónoma de tercer orden en un BVP con una EDO de segundo orden y un IVP con una EDO separable de primer orden.
Ir a la función inversa $t(x)$ como sugiere JJacquelin podemos escribir \begin{equation} t(x) = \int \limits_0^x \frac{1}{w(\xi)} \, \mathrm{d}\xi, \end{equation} con la función $w$ desde arriba. Por lo tanto, queda por resolver el BVP para $w$ y luego calcular la integral.
Ambas tareas pueden ser difíciles, y creo que hay un error en la solución proporcionada por JJacquelin (en el cálculo de la tercera derivada de la función inversa). Normalmente escribiría un comentario al respecto, pero aún no se me permite.
Sigue un script de MATHEMATICA que muestra una solución numérica. Aquí $t = 200 \approx\infty$
tmax = 200; solx = NDSolve[{x1'[t] == x2[t], x2'[t] == x3[t], 2 x3'[t] + x1[t] x3[t] == 0, x1[0] == 0, x2[0] == 0, x2[200] == 1}, {x1, x2, x3}, {t, 0, tmax}][[1]] Plot[Evaluate[{x1[t], x2[t], x3[t]} /. solx], {t, 0, 20}, PlotRange -> {0, 2}, PlotStyle -> {{Black, Thick}, {Blue, Thick}, {Red, Thick}}, PlotLegends -> {Subscript[x, 1][t], Subscript[x, 2][t], Subscript[x, 3][t]}]
$$ x_1(t) \to \mbox{black}\\ x_2(t) \to \mbox{blue}\\ x_3(t) \to \mbox{red} $$
I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.