Explotación Minera (G15-C)

Trabajo realizado por estudiantes
Título	Explotación minera. Grupo 15-C
Asignatura	Ecuaciones Diferenciales
Curso	Curso 2014-15
Autores	Belén Salamanca, M.Rosario Ruiz Serrano , Almudena Román, Carmén Rocio LLanes, Elena Suta, Sergio Fernández
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

1 Introducción

1.1 Relación entre Producción-Material extraido

Dada la ecuación diferencial proporcionada en el problema, donde se tiene la derivada de Q con respecto al tiempo, la función de producción es la propia derivada, ya que Q es la función de distribución y P es la densidad, con lo cual la relación entre ambas es la derivada, esto es:

dx/dt=rQ log(K/Q)=P

1.2 Valor de la tasa intrínseca de crecimiento (r)

Como ya se ha comentado anteriormente la función P(producción) es la derivada de la cantidad extraída de material,

dx/dt=rQ log(K/Q)=P

sabemos que la producción máxima es de 240, por tanto, derivando la función P e igualando a 0 obtenemos donde se produce el máximo

dP/dt=0    t= c/r

que introducido en la ecuación P e igualando a 240 se obtiene el valor de la tase de crecimiento r pedido

r= 240e/k =0.06

1.3 Relación entre la producción y el volumen de toneladas extraidas según el modelo de Gompertz

Mediante el modelo de Gompertz, podemos expresar la producción de minerales (P en función de Q), como hemos visto en el apartado anterior. En la gráfica se observa cómo dicha producción aumenta hasta llegar a su máximo y a continuación disminuye. Para representar la curva necesitamos el siguiente código matlab:

b=240
k=10875;
r=b*exp(1)/k;
q=0:1:10875;
N=length(q);
P=zeros(1,N);
for i=1:N
    P(i)=r*q(i)*log(k/q(i));
end
plot(q,P)

1.4 Relación entre la producción y el volumen de toneladas extraídas según el modelo de Verhulst

Volvemos a dibujar la curva del modelo de Gompertz descrita en el apartado anterior, y en una segunda ventana la curva de Verhulst. Para ello, necesitamos definir una nueva función, así como una nueva tasa intrínseca de crecimiento. [math] Q'=rQ(1-\frac{Q}{k}) [/math]

Añadimos una tercera ventana en la que se dibujan las gráficas de ambos modelos, para poder compararlas. Necesitamos el siguiente código matlab para obtener los resultados:

b=240;
k=10875;
r1=b*exp(1)/k;%Tasa intrínseca de crecimiento según el modelo de Gompertz 
r2=(240*4)/10875; %Tasa intrinseca de creciemiento según el modelo de Verhulst
q=0:1:10875;
N=length(q);
P=zeros(1,N);
for i=1:N
P(i)=r1*q(i)*log(k/q(i)); %Gompertz
PV(i)=r2*q(i)*(1-q(i)/k);  %Verhulst
end
subplot(1,3,1)
plot(q,P,'k') %Curva para la función de Gompertz
xlabel('cantidad')            
ylabel('produccion')         
subplot(1,3,2)
plot(q,PV,'r') %Curva para la función de Verhulst
xlabel('cantidad (tn)')      
ylabel('produccion (tn/año)') 
subplot(1,3,3) %Comparación entre las curvas de ambos modelos
plot(q,P,'k')                   
xlabel('cantidad')             
ylabel('produccion')         
hold on %Para superponer gráficas                   
plot(q,PV,'r')              
xlabel('cantidad (tn)')        
ylabel('produccion (tn/año)')  
legend('Modelo Gompertz','Modelo Verhulst','Location','best') %Leyenda
hold off

1.5 Problema de valor inicial y aproximación de la cantidad de material extraído por el método de Euler

Sea el problema de valor inicial dQ/dt=r*Q*log(k/Q) tal que Q(0)=0.1, la gráfica representa la variación de cantidad de mineral extraible a lo largo del tiempo. Como se puede observar, en un primer momento podemos extraer mucha cantidad de mineral hasta que, al llegar a un tiempo t y debido a que los recursos minerales limitados, la tasa de rendimiento de extracción del mineral es inferior a 25 toneladas/año (estamos alcanzando la máxima cantidad extraible). Por lo tanto, ya no es rentable seguir con la extracción.

t0=0;
h=1/12;
k=10875;
b=240;
r=b*exp(1)/k;
t(1)=t0;
q(1)=0.1;
c=log(log(k/q(1)));
Q=k/(exp(exp(-r*t+c)));
i=1;
while 1
    q(i+1)=q(i)+h*r*q(i)*log(k/q(i));
    t(i+1)=t(i)+h;
    if i>1&&abs((r*q(i)*log(k/q(i)))-25)<0.1&&abs((r*q(i-1)*log(k/q(i-1))))>abs((r*q(i)*log(k/q(i))));
      
        break
    end
    i=i+1;
end
[t',q']
plot(t,q)

1.6 Aproximación de la cantidad de material extraído con los métodos de Runge kutta y Heun

clear all
%Condiciones iniciales
t0=0;
q0=0.1;
h=1/12;
k=10875;
b=240;
r=b*exp(1)/k;
%Variable dependiente
t=t0;
q(1)=q0;
z(1)=q0;
i=1;
while 1
K1=r*q(i)*(log(k)-log(q(i)));
K2=r*(q(i)+1/2*K1*h)*(log(k)-log(q(i)+1/2*K1*h));
K3=r*(q(i)+1/2*K2*h)*(log(k)-log(q(i)+1/2*K2*h));
K4=r*(q(i)+K3*h)*(log(k)-log(q(i)+K3*h));
q(i+1)=q(i)+h/6*(K1+2*K2+2*K3+K4);
t(i+1)=t(i)+h;
    if i>1&&abs(r*q(i)*(log(k)-log(q(i)))-25)<0.1&&r*q(i-1)*(log(k)-log(q(i-1)))>r*q(i)*(log(k)-log(q(i)))
        break
    end
   i=i+1;
end
i=1;
while 1
    K1=r*z(i)*(log(k)-log(z(i)));
    K2=r*(z(i)+K1*h)*(log(k)-log(z(i)+K1*h));
    z(i+1)=z(i)+h/2*(K1+K2);
    t(i+1)=t(i)+h;
    if i>1&&abs(r*z(i)*(log(k)-log(z(i)))-25)<0.1&&r*z(i-1)*(log(k)-log(z(i-1)))>r*z(i)*(log(k)-log(z(i)))
        break
    end
i=i+1;
end
hold on
plot(t,q,'b*')
plot(t,z,'r+')
legend('RK4','HEUN','location','best')
hold off

1.7 Cantidad de mineral cuando el tiempo tiende a infinito

tN=input('Introduce un valor del tiempo muy grande:');
t=0:1:tN;                             
N=length(t);                             
Q=zeros(1,N);                            
r=240*exp(1)/10875;                      
Q0=0.1;                                  
K=10875;                                 
Q(1)=Q0;                                 
Q=K*exp(exp(-r*t)*(log(Q0/K)));          
for i=1:N                                
    Q(i)=K*exp(exp(-r*t(i))*(log(Q0/K)));
end                                      
Q(tN)                                                
plot(t,Q)

1.8 Representación de la función P(t)

Se adjunta el código matlab que da la función de producción en función del tiempo.

clear all
t=0:1:200;
n=length(t);
Q0=0.1;
K=10875;
r=240*exp(1)/K;
C=log(log(K/Q0));
for i=1:n
    P(i)=(K*r*exp(-r*t(i)+C))/(exp(exp(-r*t(i)+C)));
    
end
plot(t,P)

La función resultante es: El punto de máxima producción se encuentra en P=239.9906, en un tiempo de aproximadamente 40 años, que se obtiene en matlab mediante el comando 'max'.

1.9 Cantidad de mineral sin extraer una vez transcurrida la vida útil

Utilizando el programa del apartado 1.5 podemos obtener el valor de la cantidad de material extraído en t=25 que nos da un valor de Q(25)=820ton. Sabemos que tenemos K= 1087ton. en el total del yacimiento, por tanto restando Q(25) menos el total del yacimiento deducimos que quedarán aproximadamente 10056ton. de material sin extraer ==1.11 Reajuste a datos reales