Revisión del 15:02 1 mar 2015

1 Introducción

1.1 Objetivos y metodología

El objetivo de este trabajo es estudiar las concentraciones de los reactivos y productos de una reacción química a lo largo de la evolución en el tiempo de dicha reacción. Para ello modelizaremos el proceso empleando ecuaciones diferenciales y teniendo en cuenta la Ley de Acción de Masas y el Principio de Conservación de la materia.

Al ser las ecuaciones y sistemas que se nos presentan de difícil resolución analítica, utilizaremos el programa informático Matlab para resolverlos mediante métodos numéricos como el de Euler, el del Trapecio, el de Runge-Kutta y el de Heun.

Plasmaremos las soluciones en gráficos que nos permitirán dar una interpretación físico química a los resultados.

1.2 Ley de acción de masas

Enunciada por Guldberg y Waage en 1864, establece que dada una reacción química reversible en equilibrio, a temperatura constante, la relación de concentraciones de los reactivos y productos tiene un valor constante.

1.3 Principio de conservación de la masa

En una reacción química ordinaria la masa permanece constante, siendo la masa consumida de los reactivos igual a la masa obtenida de los productos.

2 Reacción I

2.1 Interpretación del problema y deducción de las ecuaciones diferenciales

Deducir a partir de la ley de accion de masas y del principio de conservacion de la masa que las concentraciones de A y B deben satisfacer las ecuaciones

[math]x'(t)+y'(t)=0[/math]
[math]y'(t)=k_1x(t)y(t)[/math]

.....................

2.2 Obtención del problema de valor inicial

Una vez demostrado el origen de las ecuaciones (1) y (2) en el apartado anterior, obtendremos el problema de valor inicial (PVI) mediante las manipulaciones detalladas a continuación:

[math]x'(t)+y'(t)=0\quad\quad(1)[/math]

[math]y'(t)=k_1x(t)y(t)\quad\quad(2)[/math]

Integrando (1) obtenemos (3):

[math]\int_\mathbb{D} (\frac{dx(t)}{dt} + \frac{dy(t)}{dt})\,dt =\int_\mathbb{D} 0\,dt \quad \Rightarrow \quad x(t)+y(t)=k_2\quad \Rightarrow \quad x(t)=k_2-y(t)\quad\quad(3)[/math]

Sustituyendo (3) en (2), obtenemos (4):

[math]y'(t)=k_1x(t)y(t)\quad \Rightarrow \quad y'(t)=k_1(k_2-y(t))y(t)\quad \Rightarrow \quad y'(t)=k_1k_2y(t)-k_1y^2(t)\quad\quad(4)[/math]

Conocemos las concentraciones iniciales de A y de B, por lo que podemos obtener el valor de la constante k₂:

[math][A]:\quad x(0)=1 mol/l;\quad\quad [B]:\quad y(0)=0,01 mol/l[/math]

[math]x(t)+y(t)=k_2\quad \Rightarrow \quad x(0)+y(0)=k_2\quad \Rightarrow \quad 1+0,01=k_2\quad \Rightarrow \quad k_2=1,01[/math]

Por último, conocemos también el valor de k₁, por lo que podremos plantear a continuación el correspondiente PVI:

[math]k_1=1[/math]

[math](PVI)\begin{cases}y'(t)=1,01y(t)-y^2(t) \\ \\ y(0) = 0,01 \end{cases}\quad (t,y)\in[0,10]\times\mathbb{R}[/math]

2.3 Unicidad de la solución

Intersección B con dominio

Aplicamos a continuación el Teorema de Existencia y Unicidad:

Sea [math]f(t,y)=y'(t) \Rightarrow f(t,y)=1.01y(t)-y^2(t)\quad[/math], y [math]\quad\frac{\partial }{\partial y} f(t,y) = 1.01-2y(t)[/math].

Y sea [math]B[(t_0,y_0),r],\quad r\gt0,\quad en (t_0,y_0)=(0,0.01)[/math].

[math]f(t,y)[/math] es continua en [math] \mathbb{B}\cap \mathbb{D}\quad \Rightarrow \quad \exists y(t)\quad /\quad y'(t)=1,01y(t)-y^2(t)[/math]

[math]\frac{\partial }{\partial y} f(t,y)[/math] es continua en [math] \mathbb{B}\cap \mathbb{D}\quad \Rightarrow \quad \exists ! y(t)\quad /\quad y'(t)=1,01y(t)-y^2(t)[/math]

Con lo que queda demostrado que el (PVI) tiene solución, y además esta es única.

2.4 Resolución numérica del problema de valor inicial

2.4.1 Método de Euler

2.4.2 Método del trapecio

2.4.3 Método de Runge-Kutta

Representación gráfica y(t). Método Runge-Kutta

%Método de Runge-Kutta
clc

%Función f(t)=y'(t)
fyt=inline('1.01*y-y^2','y');

%Datos del problema
y0=0.01;
t0=0;
tN=10;
h=0.1;

%Variable independiente (t)
t=t0:h:tN;

%Variable dependiente (y)
y=zeros(1,length(t));
y(1)=y0;

for i=1:length(t)-1
    k1=fyt(y(i));
    k2=fyt(y(i)+k1/2*h);
    k3=fyt(y(i)+k2/2*h);
    k4=fyt(y(i)+k3*h);
    y(i+1)=y(i)+h/6*(k1+2*k2+2*k3+k4);
end

%Representación gráfica
hold on
plot(t,y)
legend('Sustancia B: y(t)')
xlabel('Tiempo (s)')
ylabel('Concentración (mol/l)')
grid on

3 Reacción II