Diferencias finitas para la ecuación de ondas - Historial de revisiones

Carlos Castro: /* Programa para OCTAVE o MATLAB */

2013-06-07T14:30:26Z

‎Programa para OCTAVE o MATLAB

Deprecated: The each() function is deprecated. This message will be suppressed on further calls in /home/mat/public_html/w/includes/diff/DairikiDiff.php on line 434

Carlos Castro: /* Programa para OCTAVE o MATLAB */

2013-06-07T14:21:40Z

‎Programa para OCTAVE o MATLAB

Herraiz: /* Ejemplo */

2013-05-28T13:25:48Z

‎Ejemplo

Carlos Castro en 10:07 28 may 2013

2013-05-28T10:07:01Z

Carlos Castro en 09:12 28 may 2013

2013-05-28T09:12:29Z

Carlos Castro en 09:07 28 may 2013

2013-05-28T09:07:14Z

Carlos Castro en 09:03 28 may 2013

2013-05-28T09:03:49Z

Carlos Castro en 08:56 28 may 2013

2013-05-28T08:56:48Z

Carlos Castro en 08:55 28 may 2013

2013-05-28T08:55:31Z

Carlos Castro en 08:52 28 may 2013

2013-05-28T08:52:41Z

@@ Línea 63: / Línea 63: @@
 {{matlab|codigo=
-%%%% Diferencias finitas para resolver la ecuacion de ondasclear allL1=0; L2=1;       % extremos del intervalo%% Discretizacion espacialN=10;             % numero de subintervalos dx=(L2-L1)/N;     % paso en espaciox=L1+dx:dx:L2-dx;U=(x.^2.*(1-x))'; % posicion inicialV=x';             % velocidad inicialTF=3;             % tiempo final%% discretización en tiempodt=dx/10;         % paso en tiempot=0:dt:TF;%% matriz que aproxima -u_xxK=(1/(dx^2))*(2*diag(ones(1,N-1)) -diag(ones(1,N-2),1)-diag(ones(1,N-2),-1)      );k(1,1)=1/dx^2;
+%%%% Diferencias finitas para resolver la ecuacion de ondas
-g=4./(1+x.^2)';  % segundo miembrosol(1,:)=[U(1),U',0]; % matriz que contiene la solucion para dibujarlaf=zeros(N-1,1);       % vector con la aportacion de la frontera for k=1:length(t)-1    f(1)=-(2*t(k)*exp(-t(k)^2)/dx);  % condicion neumann en x=L1    f(N-1)=(t(k)/((1+t(k)^2)*dx^2)); % condicion dirichlet en x=L2    U=U+dt*V;                         % método explícito    V=V+dt*(-K*U+g*exp(-t(k))+f);    sol(k+1,:)=[U(1)-2*t(k)*exp(-t(k))*dx , U', t(k)/(1+t(k)^2)];end%Dibujo de la solucionx1=L1:dx:L2[X,T]=meshgrid(x1,t)surf(X,T,sol)
+clear all
 }}

@@ Línea 63: / Línea 63: @@
 {{matlab|codigo=
-L1=0; L2=1;   % extremos del intervalo
+%%%% Diferencias finitas para resolver la ecuacion de ondasclear allL1=0; L2=1;       % extremos del intervalo%% Discretizacion espacialN=10;             % numero de subintervalos dx=(L2-L1)/N;     % paso en espaciox=L1+dx:dx:L2-dx;U=(x.^2.*(1-x))'; % posicion inicialV=x';             % velocidad inicialTF=3;             % tiempo final%% discretización en tiempodt=dx/10;         % paso en tiempot=0:dt:TF;%% matriz que aproxima -u_xxK=(1/(dx^2))*(2*diag(ones(1,N-1)) -diag(ones(1,N-2),1)-diag(ones(1,N-2),-1)      );k(1,1)=1/dx^2;
-% Discretizacion espacial
+g=4./(1+x.^2)';  % segundo miembrosol(1,:)=[U(1),U',0]; % matriz que contiene la solucion para dibujarlaf=zeros(N-1,1);       % vector con la aportacion de la frontera for k=1:length(t)-1    f(1)=-(2*t(k)*exp(-t(k)^2)/dx);  % condicion neumann en x=L1    f(N-1)=(t(k)/((1+t(k)^2)*dx^2)); % condicion dirichlet en x=L2    U=U+dt*V;                         % método explícito    V=V+dt*(-K*U+g*exp(-t(k))+f);    sol(k+1,:)=[U(1)-2*t(k)*exp(-t(k))*dx , U', t(k)/(1+t(k)^2)];end%Dibujo de la solucionx1=L1:dx:L2[X,T]=meshgrid(x1,t)surf(X,T,sol)
-N=10;         % número de subintervalos
-dx=(L2-L1)/N; % paso de discretizacion
-x=L1+dx:dx:L2-dx;  % malla espacial
-U=x.^2.*(1-x);     % posición inicial
-V=x;               % velocidad inicial
-% Discretización en tiempo
-TF=3;              % tiempo final
-dt=dx/10;          % paso de tiempo
-t=0:dt:TF;         % malla temporal
-% Matriz discretización de -u_xx
-K=(1/(dx^2))*(2*diag(ones(1,N-1)) -diag(ones(1,N-2),1)-diag(ones(1,N-2),-1));
-K(1,1)=1/dx^2;
-g=4./(1+x.^2); % solucion en x=0
-% Guardamos la solucion en la matriz sol
-sol(1,:)=[U(1),U,0];
-f=zeros(1,N-1);
-% Bucle temporal
 for k=1:length(t)-1
-    f(1)=-(2*t(k)*exp(-t(k)^2)/dx);
-    f(N-1)=(t(k)/((1+t(k)^2)*dx^2));
-    Z=U'+dt*V';
-    W=V'+dt*(-K*U'+g'*exp(-t(k))+f')+(4./(1+X.^2)).*exp(-T));
-    U=Z';
-    V=W';
-    sol(k+1,:)=[U(1)-2*t(k)*exp(-t(k))*dx , U, t(k)/(1+t(k)^2)];
 end
-% Dibujo de la solucion
-x1=L1:dx:L2;
-[X,T]=meshgrid(x1,t);
 surf(X,T,sol)
 }}

@@ Línea 99: / Línea 99: @@
 == Ejemplo ==
-<gallery>
+[[Archivo:Diferenciasfiinitaseuler.jpg|600px|thumb|centro|Aproximación numérica de la ecuación de ondas]]
-Archivo:Diferenciasfiinitaseuler.jpg|Aproximación numérica de la cuación de ondas
-</gallery>
 [[Categoría:Grado en Ingeniería Civil y Territorial]]
 [[Categoría:Ecuaciones Diferenciales]]

← Revisión anterior		Revisión del 10:07 28 may 2013
Línea 102:		Línea 102:
	Archivo:Diferenciasfiinitaseuler.jpg\|Aproximación numérica de la cuación de ondas		Archivo:Diferenciasfiinitaseuler.jpg\|Aproximación numérica de la cuación de ondas
	</gallery>		</gallery>
		+
		+	[[Categoría:Grado en Ingeniería Civil y Territorial]]
		+	[[Categoría:Ecuaciones Diferenciales]]

← Revisión anterior		Revisión del 09:12 28 may 2013
Línea 96:		Línea 96:
	surf(X,T,sol)		surf(X,T,sol)
	}}		}}
		+
		+	== Ejemplo ==
		+
		+	<gallery>
		+	Archivo:Diferenciasfiinitaseuler.jpg\|Aproximación numérica de la cuación de ondas
		+	</gallery>

← Revisión anterior		Revisión del 09:03 28 may 2013
Línea 59:		Línea 59:
	o escrito en otra forma		o escrito en otra forma
	<math>\begin{equation}\label{diferencia8}\left\{ \begin{array}{lll} U_{n+1}=U_n+dtV_n \\ V_{n+1} = V_n+dt(-KU_n+F(t_n)) \\ U_0=\left( \begin{array}{llll} x_1^2(1-x_1) \\ x_2^2(1-x_2) \\ x_3^2(1-x_3) \\ x_4^2(1-x_4) \end{array} \right) , \;\; V_0=\left( \begin{array}{llll} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array} \right) \end{array} \right.\end{equation}</math>		<math>\begin{equation}\label{diferencia8}\left\{ \begin{array}{lll} U_{n+1}=U_n+dtV_n \\ V_{n+1} = V_n+dt(-KU_n+F(t_n)) \\ U_0=\left( \begin{array}{llll} x_1^2(1-x_1) \\ x_2^2(1-x_2) \\ x_3^2(1-x_3) \\ x_4^2(1-x_4) \end{array} \right) , \;\; V_0=\left( \begin{array}{llll} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array} \right) \end{array} \right.\end{equation}</math>
		+
		+	== Programa para OCTAVE o MATLAB ==
		+
		+	{{matlab\|codigo=
		+	L1=0; L2=1; % extremos del intervalo
		+	% Discretizacion espacialN=10; % número de subintervalosdx=(L2-L1)/N; % paso de discretizacionx=L1+dx:dx:L2-dx; % malla espacialU=x.^2.*(1-x); % posición inicialV=x; % velocidad inicial
		+	% Discretización en tiempoTF=3; % tiempo finaldt=dx/10; % paso de tiempot=0:dt:TF; % malla temporal
		+	% Matriz discretización de -u_xx K=(1/(dx^2))(2diag(ones(1,N-1)) -diag(ones(1,N-2),1)-diag(ones(1,N-2),-1));K(1,1)=1/dx^2;
		+	g=4./(1+x.^2); % solucion en x=0
		+	% Guardamos la solucion en la matriz solsol(1,:)=[U(1),U,0];f=zeros(1,N-1);% Bucle temporalfor k=1:length(t)-1 f(1)=-(2t(k)exp(-t(k)^2)/dx); f(N-1)=(t(k)/((1+t(k)^2)dx^2)); Z=U'+dtV'; W=V'+dt(-KU'+g'exp(-t(k))+f'); %+(4./(1+X.^2)).exp(-T)); U=Z'; V=W'; sol(k+1,:)=[U(1)-2t(k)exp(-t(k))*dx , U, t(k)/(1+t(k)^2)];end
		+	% Dibujo de la solucionx1=L1:dx:L2[X,T]=meshgrid(x1,t)surf(X,T,sol)
		+	}}

@@ Línea 49: / Línea 49: @@
 El siguente paso es resolver este sistema de 4 ecuaciones diferenciales de segundo orden.Empezamos primeramente utilizando el m\'etodo de Euler. Esto implica pasar de este sistema de 4 ecuaciones diferenciales de segundo orden a un sistema de 8 ecuaciones diferenciales de primer orden . Introducimos las funciones $u'_1(t)=v_1(t), \; u'_2(t)=v_2(t), \; u'_3(t)=v_3(t), \; u'_4(t)=v_(t),$, las nuevas icognitas son $8 $ funciones que las designamos por $u_1(t)$, $u_2(t)$, $u_3(t)$, $u_4(t)$, $v_1(t)$, $v_2(t)$, $v_3(t)$, $v_4(t)$ que satisfacen el sistema de 8 ecuaciones diferenciales de primer orden
-<math>\begin{equation}\label{diferencia6}\left\{ \begin{array}{llllllll} u_1'(t)=v_1(t) \\ u_2'(t)=v_2(t) \\ u_3'(t)=v_3(t) \\ u_4'(t)=v_4(t) \\  v'_1(t)=-\frac{-u_1(t)+2u_1(t)-u_2(t)}{dx^2}+\frac{4}{1+x_1^2}e^{-t} -\frac{2te^{-t^2}}{dx}\vspace{0.2cm} \\  v'_2(t)=-\frac{-u_1(t)+2u_2(t)-u_3(t)}{dx^2}+\frac{4}{1+x_2^2}e^{-t} \vspace{0.2cm} \\  v'_3(t)=-\frac{-u_2(t)+2u_3(t)-u_4(t)}{dx^2}+\frac{4}{1+x_3^2}e^{-t} \vspace{0.2cm} \\  v'_4(t)=\frac{-u_3(t)+2u_4(t)}{dx^2}+\frac{4}{1+x_1^2}e^{-t}+\frac{t}{1+t^2} \frac{1}{dx^2}  \end{array} \right.\end{equation}
+<math>\begin{equation}\label{diferencia6}\left\{ \begin{array}{llllllll} u_1'(t)=v_1(t) \\ u_2'(t)=v_2(t) \\ u_3'(t)=v_3(t) \\ u_4'(t)=v_4(t) \\  v'_1(t)=-\frac{-u_1(t)+2u_1(t)-u_2(t)}{dx^2}+\frac{4}{1+x_1^2}e^{-t} -\frac{2te^{-t^2}}{dx} \\  v'_2(t)=-\frac{-u_1(t)+2u_2(t)-u_3(t)}{dx^2}+\frac{4}{1+x_2^2}e^{-t}  \\  v'_3(t)=-\frac{-u_2(t)+2u_3(t)-u_4(t)}{dx^2}+\frac{4}{1+x_3^2}e^{-t}  \\  v'_4(t)=\frac{-u_3(t)+2u_4(t)}{dx^2}+\frac{4}{1+x_1^2}e^{-t}+\frac{t}{1+t^2} \frac{1}{dx^2}  \end{array} \right.\end{equation}
 </math>
 LLamando $V(t)= \left(   \begin{array}{llll} v_1(t) \\ v_2(t) \\ v_3(t)  \\ v_4(t) \end{array}  \right)$, \ref{diferencia5} es equivalente a
-<math>\begin{equation}\label{diferencia7}\left\{ \begin{array}{ll}U'(t)=V(t) \vspace{0.2cm} \\V'(t)=-KU(t)+F(t)  \vspace{0.2cm}  \\ U(0)=\left( \begin{array}{llll} x_1^2(1-x_1) \\ x_2^2(1-x_2) \\ x_3^2(1-x_3) \\ x_4^2(1-x_4) \end{array}      \right) , \;\; V(0)=\left( \begin{array}{llll} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array}      \right)  \end{array}  \right.\end{equation}</math>
+<math>\begin{equation}\label{diferencia7}\left\{ \begin{array}{ll}U'(t)=V(t)  \\V'(t)=-KU(t)+F(t)  \\ U(0)=\left( \begin{array}{llll} x_1^2(1-x_1) \\ x_2^2(1-x_2) \\ x_3^2(1-x_3) \\ x_4^2(1-x_4) \end{array}      \right) , \;\; V(0)=\left( \begin{array}{llll} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array}      \right)  \end{array}  \right.\end{equation}</math>
 Si tomamos un paso temporal $dt$, definimos $t_j=0+jdt$, $U_n$ y $V_n$ designa la soluci\'on aproximada de $U(t_n)$ y $V(t_n)$ respectivamente,   el m\'etodo de Euler aplicado a \ref{diferencia7} queda en la forma
 <math>\left(  \begin{array}{ll} U_{n+1}  \\ V_{n+1}  \end{array} \right) =\left(  \begin{array}{ll} U_{n}  \\ V_{n}  \end{array} \right)+dt \left(  \begin{array}{ll} U'_{n}(t_n)  \\ V'_{n}(t_n)  \end{array} \right)=\left(  \begin{array}{ll} U_{n}  \\ V_{n}  \end{array} \right)+dt\left(  \begin{array}{ll} V_{n}  \\ -KU_n+F(t)  \end{array} \right)</math>
@@ Línea 58: / Línea 58: @@
 o escrito en otra forma
-<math>\begin{equation}\label{diferencia8}\left\{ \begin{array}{lll} U_{n+1}=U_n+dtV_n   \vspace{0.2cm}\\ V_{n+1} = V_n+dt(-KU_n+F(t_n))   \vspace{0.2cm} \\ U_0=\left( \begin{array}{llll} x_1^2(1-x_1) \\ x_2^2(1-x_2) \\ x_3^2(1-x_3) \\ x_4^2(1-x_4) \end{array}      \right) , \;\; V_0=\left( \begin{array}{llll} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array}      \right)  \end{array}  \right.\end{equation}</math>
+<math>\begin{equation}\label{diferencia8}\left\{ \begin{array}{lll} U_{n+1}=U_n+dtV_n   \\ V_{n+1} = V_n+dt(-KU_n+F(t_n))    \\ U_0=\left( \begin{array}{llll} x_1^2(1-x_1) \\ x_2^2(1-x_2) \\ x_3^2(1-x_3) \\ x_4^2(1-x_4) \end{array}      \right) , \;\; V_0=\left( \begin{array}{llll} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array}      \right)  \end{array}  \right.\end{equation}</math>

@@ Línea 39: / Línea 39: @@
 Si llamamos
-<math> U(t)= \left(   \begin{array}{llll} u_1(t) \\ u_2(t) \\ u_3(t)  \\ u_4(t) \end{array}  \right), \qquad  K= \left(  \begin{array}{cccc} 1 & -1 & 0 & 0 \\-1 & 2 & -1 & 0 \\    0 & -1 & 2 & -1 \\   0 & 0 & -1 & 2 \; \end{array} \right), </math>
+<math> U(t)= \left(   \begin{array}{llll} u_1(t) \\ u_2(t) \\ u_3(t)  \\ u_4(t) \end{array}  \right), \qquad  K= \left(  \begin{array}{cccc} 1 & -1 & 0 & 0 \\-1 & 2 & -1 & 0 \\    0 & -1 & 2 & -1 \\   0 & 0 & -1 & 2 \; \end{array} \right), \qquad  F(t)=\left( \begin{array}{llll}  \frac{4}{1+x_1^2}e^{-t} -\frac{2te^{-t^2}}{dx} \\   \frac{4}{1+x_2^2}e^{-t} \\ \frac{4}{1+x_3^2}e^{-t} \\  \frac{4}{1+x_4^2}e^{-t}  +\frac{t}{1+t^2} \frac{1}{dx^2}  \end{array}  \right), </math>
-−
-−
-<math> F(t)=\left( \begin{array}{llll}  \frac{4}{1+x_1^2}e^{-t} -\frac{2te^{-t^2}}{dx} \\   \frac{4}{1+x_2^2}e^{-t} \\ \frac{4}{1+x_3^2}e^{-t} \\  \frac{4}{1+x_4^2}e^{-t}  +\frac{t}{1+t^2} \frac{1}{dx^2}  \end{array}  \right), </math>
-−
 el sistema \ref{diferencia4}, junto con la condici\'on inicial, puede escribirse en la forma
@@ Línea 51: / Línea 47: @@
 == Discretización temporal ==