Diferencia entre revisiones de «Difusión de una sustancia contaminante (Grupo 10)»

De MateWiki
Saltar a: navegación, buscar
Línea 133: Línea 133:
 
Con <math>N+1</math> ecuaciones.
 
Con <math>N+1</math> ecuaciones.
 
==== Resolución por método del trapecio ====
 
==== Resolución por método del trapecio ====
{{matlab|codigo=
 
clear all
 
%%Problema de difusión de barrera
 
%Datos
 
L=5;
 
T=6;
 
%Discretización temporal y espacial
 
dx=0.1;
 
N=L/dx;
 
dt=dx^2/2;
 
x=0:dx:L;
 
t=0:dt:T;
 
%Matriz de coeficientes
 
K=2*diag(ones(1,N+1))-diag(ones(1,N),1)-diag(ones(1,N),-1);
 
K(N,N-1)=-2;
 
K(1,2)=-2;
 
K=(1/dx^2)*K;
 
%Condición inicial
 
U0=(3.*(x>3))';
 
sol(1,:)=[U0'];
 
U=U0;
 
%Bucle iterativo aplicando método de diferencias finitas con método del trapecio
 
for j=length(t)-1
 
U=(eye(N+1)+(dt/2)*K)\(U-(dt/2)*K*U);
 
sol(j+1,:)=[U'];
 
end
 
[xx,tt]=meshgrid(x,t);
 
figure(1)
 
mesh(xx,tt,sol);
 
xlabel('x');
 
ylabel('y');
 
grid
 
 
}}
 
 
==== Resolución por Euler explícito ====
 
 
 
{{matlab|codigo=
 
{{matlab|codigo=
 
clear all
 
clear all

Revisión del 16:53 16 may 2014

La difusión implica un movimiento espontáneo y desordenado de moléculas individuales y, especialmente, nos interesa la difusión de moléculas de una sustancia disuelta en otra.

La difusión de una sustancia puede considerarse análoga al flujo de calor, y la ley de Fick establece que el ritmo de difusión por unidad de superficie, en dirección perpendicular a ésta, es proporcional al gradiente de la concentración de esa sustancia en esa dirección. La concentración es la masa de la sustancia por unidad de volumen, y el gradiente de concentración es la variación de concentración por unidad de distancia.

Para determinar como se distribuye la concentración de una sustancia contaminante tomaremos como ejemplo un tubo largo compuesto por dos sustancias de las cuáles una de ellas es contaminante.


1 Simplificación del problema

Para que este problema resulte más sencillo aplicaremos una serie de condiciones, facilitando el cálculo del mismo.

En primer lugar denotaremos su concentración [math] u [/math], concentración de contaminante en cada posición del tubo y constante en cada sección, medida en [math] \frac{mol}{m^2*s} [/math].

Supondremos que el tubo mide L=5 m, ocupando el intervalo [math] x \in (0,L) [/math], y que la concentración es la misma en cualquier punto de la sección transversal del tubo. Por tanto, la concentración dependerá unicamente de dos variables [math] u = u(x,t) [/math]. En los extremos se ha colocado un aislante que hace que el flujo de contaminante al exterior del tubo sea nulo. Además, consideraremos la simplificación de que se trata de un tubo aislado por sus laterales, pero abierto por las bases del cilindro.


1.1 Ley de conservación de masa

Tras una serie de experimentos realizados por Lavoisier, se pusieron de manifiesto que si tenemos en cuenta todas las sustancias que forman parte en una reacción química y todos los productos formados, nunca varía la masa. Esta ley exige que no haya una generación espontánea de materia. Es decir, en la región del espacio estudiada, no existe entrada o salida por un punto cualquiera (frontera o interior), la masa en esta región se mantiene constante. Esta es la ley de la conservación de la masa, que podemos enunciarla, pues, de la siguiente manera:
"En toda reacción química la masa se conserva, esto es, la masa total de los reactivos es igual a la masa total de los productos"
La ley de conservación de masas exige que no haya una generación espontánea de materia. Es decir, en la región del espacio estudiada, no existe entrada o salida por un punto cualquiera (frontera o interior), la masa en esta región se mantiene constante.

Esta suposición, en términos de variación de masa, puede quedar determinada como:

Variación de la masa en el tiempo = flujo de masa a través de la frontera por unidad de tiempo + flujo de masa generada en el interior por unidad de tiempo.

1.2 Ley de Fick

Se trata de una ley cuantitativa en forma diferencial que describe la difusión de materia o energía en un medio en el que inicialmente existe equilibrio químico o físico. Esta ley se trata de una generalización de los fenómenos descritos por la Ley de Fourier, ya que esta es una consideración concreta del calor.

La ley de Fick anuncia que en el caso de existir diferencias de concentración en cualquier especie, el paso de estas especies se llevará a cabo desde las regiones con mayor concentración a las de menor. Formulado matemáticamente:

[math]F(x,t)=-D\frac{\partial u}{\partial x}[/math]

Donde [math]F[/math] es el flujo de la sustancia y [math]D[/math] el coeficiente de difusión, que representa la facilidad en que un soluto particular se mueve por un disolvente determinado.

La ley de Fick establece que el flujo de difusión del contaminante es proporcional a la variación de concentración: [math] F(x,t)=-D\frac{\partial u}{\partial x} [/math] donde D es el coeficiente de difusión (medido en [math] \frac{m^2}{s} [/math]), que depende de las propiedades físicas de los compuestos, y que se supone constante D=1.

2 Modelización de la ecuación diferencial

Aplicando el principio de conservación de la masa y la ley de Fourier, definimos la variable [math] u(x,t)[/math] en unidades de sustancia por área, es decir, concentración por superficie del tubo.

Consideramos [math] A [/math] como variable superficie, y tomamos un pequeño trozo, definido por [math]\Delta x[/math]. Entonces, para obtener la variación de la masa en el tiempo:

[math] A u(x,t) \Delta x [/math]

Y derivando en el tiempo:

[math] A u_t(x,t) \Delta x [/math]

obtenemos el primer término de la igualdad anteriormente enunciada.

En el segundo, consideramos nulo el flujo de masa en los puntos interiores por unidad de tiempo, pues suponemos la tubería perfecta y sin sumideros.

Como la superficie lateral del tubo está aislada, solamente habrá flujo de sustancia en la dirección del eje [math]x[/math]. Si [math]F(x,t)\gt 0[/math] pensamos que el flujo va hacia la derecha, si [math]F(x,t)\lt 0[/math], este va a hacia la izquierda. El flujo de calor en un intervalo [math][x, x + ∆x ] [/math] viene dado dado por:

[math]F(x,t)A-F(x+\Delta x, t)A[/math]

Igualando los resultados anteriores:

[math]A u_t(x,t) \Delta x = F(x,t)A-F(x+\Delta x, t)A[/math]

y seguidamente pasamos dividiendo [math]A \Delta x[/math]:

[math]u_t(x,t)= \frac {F(x,t)-F(x+\Delta x, t)}{\Delta x}[/math]

Si hacemos [math]\Delta x \rightarrow 0 [/math]:

[math]u_t(x,t)= -F_x(x,t)[/math]

Aplicando la ley de Fick:

[math]u_t(x,t)=-\frac{\partial }{\partial x}(-D u_x(x,t)) = D u_xx(x,t) [/math]

Y de esta forma podemos obtener la ecuación en derivadas parciales que nos define la concentración de sustancia:

[math]u_t(x,t)- D u_xx(x,t)= 0 [/math]

2.1 Planteamiento del problema bien propuesto

Para la resolución de ésta ecuación necesitaremos mínimo tres condiciones, dos de las cuales serán de contorno, y la otra será la condición inicial.Si suponemos [math]u_0(x,0)[/math], entonces:

[math] u_0(x,0)=\begin{cases} 0, si x\leq 3\\ 3, si x\gt3 \end{cases} [/math]

Suponemos que el flujo será nula al haber una tapa en los extremos del tubo, con lo que el problema bien propuesto quedaría::

[math]\begin{cases} u_t(x,t)- D u_xx(x,t)= 0 \\ u_0(x,0), x\in (0,L) \\ u_x(0,t)=0, u_x(L,t)=0, t\geq 0 \end{cases} [/math]

3 Resolución numérica

3.1 Caso particular 1: fronteras aisladas y concentraciones iniciales conocidas

En el primer caso particular que estudiaremos, trataremos de obtener una resolución numérica del problema anteriormente planteado:

[math]\begin{cases} u_t(x,t)- D u_xx(x,t)= 0 \\ u_0(x,0), x\in (0,L) \\ u_x(0,t)=0, u_x(L,t)=0, t\geq 0 \end{cases} [/math]

3.1.1 Resolución por el método de las diferencias finitas

La resolución por el método de las diferencias finitas se basa en expresar el problema continuo a un sistema discreto de ecuaciones diferenciales de primer orden.

Particularizando para nuestro caso con las condiciones dadas, el problema planteado sería:

[math]\begin{cases} u'_n(t)+\frac{-u_{n-1}(t)+2u_n(t)+u_{n+1}(t)}{h^2} \\ \frac{u_1(t)-u_-1(t)}{2h} \rightarrow u_1(t)=u_-1(t) \\ \frac{u_{N+1}(t)-u_{N-1}(t)}{2h} \rightarrow u_{N+1}(t)=u_{N-1}\\ u_n(0)=u^0(x) \end{cases} [/math]

Por lo que tenemos un sistema de la forma:

[math]\begin{cases} U'(t)+KU(t)=F(t)=0\\ U(0)=U^0 \end{cases} [/math]

Con [math]N+1[/math] ecuaciones.

3.1.1.1 Resolución por método del trapecio

clear all
%%Problema de difusión de barrera
%Datos
L=5;
T=6;
%Discretización temporal y espacial
dx=0.1;
N=L/dx;
dt=dx^2/2;
x=0:dx:L;
t=0:dt:T;
%Matriz de coeficientes
K=2*diag(ones(1,N+1))-diag(ones(1,N),1)-diag(ones(1,N),-1);
K(N,N-1)=-2; 
K(1,2)=-2;
K=(1/dx^2)*K;
%Condición inicial
U0=(3.*(x>3))';
sol(1,:)=[U0'];
U=U0;
%Bucle iterativo aplicando método de diferencias finitas con método del trapecio
for j=length(t)-1
U=(eye(N+1)+(dt/2)*K)\(U-(dt/2)*K*U);
sol(j+1,:)=[U'];
end
[xx,tt]=meshgrid(x,t);
figure(1)
mesh(xx,tt,sol);
xlabel('x');
ylabel('y');
grid


3.2 Conservación de la masa total de contaminante

[math]\frac d {dt} \int_0^5 u(x,t)\,dx= \int_0^5 u_{xx}(x,t)\,dx= u_x(x,t)|_0^5=0[/math]