El tanque de tormentas de Arroyofresno (Grupo 71)

Trabajo realizado por estudiantes
Título	El tanque de tormentas de Arroyofresno. Grupo 71
Asignatura	Teoría de Campos
Curso	2025-26
Autores	David Santafé Palacios Pedro Suñé Pérez Beatriz Bernal Castañeda Nicole Filip Palacios
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

El Tanque de Tormentas de Arroyofresno es una gran infraestructura subterránea del sistema de saneamiento de Madrid diseñada para almacenar temporalmente el exceso de agua durante lluvias intensas. Su función principal es evitar inundaciones, proteger el medio ambiente y reducir la carga sobre la red de alcantarillado y las depuradoras, reteniendo el agua hasta que pueda ser tratada de forma segura.

1 Funcionamiento del Tanque de Arroyofresno:

TANQUE TORMENTAS

El Tanque de Tormentas de Arroyofresno forma parte del Plan de Saneamiento de Madrid. Su principal función es reducir el riesgo de inundaciones y proteger el medio ambiente durante lluvias intensas. El tanque se activa automáticamente cuando la red de alcantarillado unitaria, que transporta aguas pluviales y residuales conjuntamente, supera su capacidad máxima. Los contaminantes arrastrados al tanque de tormentas provienen principalmente de: sólidos en suspensión, materia orgánica; metales pesados (zinc, plomo, cobre), nutrientes (nitrógeno y fósforo) y agentes específicos críticos para el hormigón (cloruros y sulfatos).

El agua contaminada permanece en el tanque de tormentas de forma temporal, hasta que los contaminantes decantan, así como, la Estación Depuradora de Aguas Residuales pueda asumir el caudal restante para su posterior bombeo y tratamiento final. Los tanques de tormentas son comunes en grandes ciudades como: Barcelona (España): Cuenta con numerosos tanques de tormentas como el tanque de la Zona Franca, con objetivo de evitar inundaciones en las ciudades. Berlín (Alemania): Emplea un sistema mixto de tanques y canales de almacenamiento subterráneo. París (Francia): Posee grandes infraestructuras subterráneas "Grand Collecteur" de similar funcionamiento.

El almacenamiento de estas aguas sucias provoca que las columnas de hormigón estén expuestas a los agentes agresivos. Los principales mecanismos de degradación son: Corrosión de la armadura por cloruros que penetran el hormigón permitiendo la presencia de oxígeno y humedad provocando la oxidación. Ataque por sulfatos que reaccionan con los componentes del cemento dando lugar a agrietamientos internos y desintegrando la matriz del hormigón.

Para garantizar la vida útil de estas estructuras críticas, se adoptan rigurosas medidas de protección: Aumento de espesor del recubrimiento de hormigón para retrasar la llegada de contaminantes; hormigón de baja porosidad para reducir la permeabilidad; empleo de aceros inoxidables en la armadura; recubrimientos superficiales en las superficies del hormigón para crear una barrera física contra el agua contaminada.

2 Presión y fuerzas estructurales sobre columnas

2.1 Representación del campo de presiones:

Cuando el tanque está lleno, el campo escalar de presión en el agua viene dado por:

[math]P(z) = P_0 - \rho_{\text{agua}} \, g \, z,\qquad z \in [-H,0][/math]

donde \(P_0\) es la presión atmosférica, \(\rho_{\text{agua}}\) la densidad del agua y \(g\) la aceleración de la gravedad. Como \(P\) sólo depende de la coordenada vertical \(z\), en la superficie cilíndrica de una columna todos los puntos de una misma cota tienen la misma presión. A continuación representamos este campo como un mapa de colores sobre un cilindro de radio [math]r_c = 0{,}75 m[/math].

CAMPO DE PRESIONES

% Parámetros
rcolumna = 0.75;      % en m
H = 22;               % en m 
Patm = 101325;        % Pa
rhoagua = 1000;       % kg/m^3
g = 9.81;             % m/s^2

% Malla en coordenadas cilíndricas
ptostheta = 200;
ptosz = 200;
theta = linspace(0,2*pi,ptostheta);
z = linspace(-H,0,ptosz);
[Theta,Z] = meshgrid(theta,z);

% Coordenadas cartesianas de la superficie cilíndrica
R = rcolumna;
X = R*cos(Theta);
Y = R*sin(Theta);

% Campo de presión P(z) 
P = Patm - rhoagua * g .* Z;                

% Superficie con mapa de color según P
figure('Color','w','Position',[200 200 700 600]);
h = surf(X, Y, Z, P, 'EdgeColor','none');  
colorbar_handle = colorbar;
colorbar_handle.Label.String = 'Presión (Pa)';

% Etiquetas
axis equal
xlabel('x (m)')
ylabel('y (m)')
zlabel('z (m)')
title('Campo de presión P(z) sobre la superficie cilíndrica de una columna')
view(40,20)

2.2 Gradiente de presión y representación como campo vectorial:

El gradiente del campo de presión se define como:

[math]\nabla P = \left( \frac{\partial P}{\partial x},\frac{\partial P}{\partial y},\frac{\partial P}{\partial z} \right)[/math]

Para el campo hidrostático \(P(z) = P_0 - \rho_{\text{agua}} g z\) sólo hay dependencia en \(z\), de modo que:

[math]\nabla P = (0,\,0,\,-\rho_{\text{agua}} g)[/math]

Es un campo vectorial constante dirigido hacia abajo. A continuación lo representamos sobre la superficie cilíndrica de una columna.

REPRESENTACIÓN GRADIENTE DE PRESIONES

% CAMPO DEL GRADIENTE DE PRESIÓN SOBRE UNA COLUMNA

% Datos
H = 22;                % [m] 
rcolumna = 0.75;      % [m] 
rhoagua = 1000;         % [kg/m^3] 
g = 9.81;         % [m/s^2] 

% Mallado cilíndrico
ptosTheta = 20;
ptosZ = 20;
theta = linspace(0, 2*pi, ptosTheta);
z = linspace(-H, 0, ptosZ);
[Theta, Z] = meshgrid(theta, z);
X = rcolumna*cos(Theta);
Y = rcolumna*sin(Theta);

% Gradiente de presión 
Gx = zeros(size(X));
Gy = zeros(size(Y));
Gz = -ones(size(Z));   % en -z

% Representación
figure;
surf(X, Y, Z, 'FaceAlpha', 0.2, 'EdgeColor', 'none');  
hold on;
quiver3(X, Y, Z, Gx, Gy, Gz, 0.7, 'k');                
axis equal;
xlabel('x [m]');
ylabel('y [m]');
zlabel('z [m]');
title('Gradiente de presión sobre la superficie de una columna');
grid on;
hold off;

2.3 Fuerza total agua sobre columna:

Para la fuerza total trabajamos con presiones manométricas:

[math] p(z) = P(z) - P_0 = \rho_{\text{agua}} g |z| = -\rho_{\text{agua}} g z,\qquad z\in[-H,0]. [/math]

En la superficie lateral de la columna de radio \(r_c\) el elemento de área es \(dA = r_c\,d\theta\,dz\). La fuerza total normal (como magnitud escalar) es

[math] F = \int_{\text{superficie}} p(z)\,dA = 2\pi r_c \rho_{\text{agua}} g \int_{-H}^{0}(-z)\,dz = \pi r_c \rho_{\text{agua}} g H^2. [/math]

Con \(r_c = 0{,}75\ \text{m}\), \(H = 22\ \text{m}\), \(\rho_{\text{agua}} = 1000\ \text{kg/m}^3\) y \(g = 9{,}81\ \text{m/s}^2\):

[math] F \approx 1{,}12\cdot 10^7\ \text{N} = 1{,}12\cdot 10^4\ \text{kN} \approx 1{,}14\cdot 10^3\ \text{toneladas-fuerza}. [/math]

3 Infiltración y corrosión en columnas

3.1 Campo de concentración infiltrada:

El campo de concentración infiltrada describe cómo los contaminantes presentes en el agua del tanque penetran a través del hormigón de las columnas y cómo su concentración varía dentro del material, desde la superficie externa hasta llegar a la armadura de acero situada en el interior. Es, por tanto, un campo escalar que asigna a cada punto del hormigón de la columna (definido por sus coordenadas) un valor de concentración de contaminantes.

Su formula es:

[math] C_{\text{col}}(r,z) = C_{0}\left(1 + \alpha \frac{|z|}{H}\right)\exp\!\left[-\lambda(\rho_{c}-r)\right] [/math]

SECCION VERTICAL

SECCION TRANSVERSAL

% Parámetros del problema
H    = 22;                
radiocolumna = 0.75;             
radioarmadura = 0.65;             
concentracionsup = 0.3;                 
alpha = 3;                
lambda = 10;             

% Mallado para la sección vertical (r,z)
nr = 300; 
nz = 300;
r = linspace(0, radiocolumna, nr);
z = linspace(-H, 0, nz);     % fondo a superficie
[R, Z] = meshgrid(r, z);

% Campo de concentración en agua
Cconagua = concentracionsup * (1 + alpha * abs(Z) / H);

% Campo infiltrado en columna
Cinfcol = Cconagua .* exp(-lambda * (radiocolumna - R));

%  FIGURA 1: SECCIÓN VERTICAL
figure;
imagesc(r, z, Cinfcol); %muestra matriz como imagen
colorbar;
xlabel('r (m)');
ylabel('z (m)');
title('Campo de concentración infiltrada C_{col}(r,z) - Sección vertical');
% marcar la armadura
hold on;
plot([radioarmadura radioarmadura], [-H 0], 'k--', 'LineWidth', 2);
text(radioarmadura+0.02, -H+2, 'Armadura', 'Color','k');


%  FIGURA 2: SECCIÓN TRANSVERSAL (x–y)
% malla en el plano transversal
nx = 300; 
ny = 300;
x = linspace(-radiocolumna, radiocolumna, nx);
y = linspace(-radiocolumna, radiocolumna, ny);
[X, Y] = meshgrid(x, y);
R2 = sqrt(X.^2 + Y.^2); %distancia radial al eje de la columna

% seleccionamos una profundidad
profundidad = -H/2;
Conanguaenz = concentracionsup * (1 + alpha * abs(profundidad)/H); %concentracion del agua en z= -H/2 (profundidad)
Cinfcol_xy = Conanguaenz .* exp(-lambda * (radiocolumna - R2)); %concentración infiltrada dentro del hormigón en cada punto

% enmascarar puntos fuera del cilindro
Cinfcol_xy(R2 > radiocolumna) = NaN;

% figura
figure;
imagesc(x, y, Cinfcol_xy);
set(gca,'YDir','normal'); %en direccion normal
axis equal;
colorbar;
xlabel('x (m)');
ylabel('y (m)');
title(['C_{col}(x,y) a z = ' num2str(profundidad) ' m']);

% dibujar límite de la armadura
hold on;
th = linspace(0, 2*pi, 300);
plot(radioarmadura*cos(th), radioarmadura*sin(th), 'k--', 'LineWidth',2);
text(radioarmadura*1.1, 0, 'Armadura', 'Color','k');

3.2 Campo del gradiente de concentración infiltrada:

La concentración de contaminantes en el agua varía con la profundidad según

[math]C_{\text{agua}}(z) = C_0\left(1 + \alpha \frac{|z|}{H}\right),\qquad z\in[-H,0][/math]

donde [math]C₀ = 0{,}3 kg/m³[/math] y [math]α = 3[/math]. En el interior del hormigón de la columna, la concentración infiltrada se modela como

[math]C_{\text{col}}(r,z) = C_{\text{agua}}(z)\,\mathrm{e}^{-\lambda(\rho_a-r)},\qquad 0\le r\le \rho_a[/math]

con [math]ρₐ = 0{,}65 m[/math] (posición de la armadura) y [math]λ = 10 m⁻¹[/math] (coeficiente de atenuación en el hormigón).

El gradiente de esta concentración en coordenadas cilíndricas \((r,z)\) es

[math]\nabla C_{\text{col}}(r,z) = \left(\frac{\partial C_{\text{col}}}{\partial r},\,\frac{\partial C_{\text{col}}}{\partial z}\right) = \big(\lambda\,C_{\text{col}}(r,z),\, -\frac{C_0\alpha}{H}\,\mathrm{e}^{-\lambda(\rho_a-r)}\big)[/math]

En una sección vertical que pasa por el eje de la columna (plano \(y=0\)), tomamos la coordenada horizontal \(x\) de manera que \(r = |x|\) y representamos el campo vectorial \(\nabla C_{\text{col}}(x,z)\).

CAMPO GRADIENTE CONCENTRACION INFILTRADA

% GRADIENTE DE LA CONCENTRACIÓN INFILTRADA EN UNA SECCIÓN VERTICAL

% Parámetros en SI
H     = 22;        
concentracionsup = 0.3;       
alpha = 3;         
lambda = 10;       
rarmadura = 0.65;     

% Mallado en la sección vertical (plano y = 0)
nx = 41;
nz = 41;
x  = linspace(-rarmadura, rarmadura, nx);   
z  = linspace(-H, 0, nz);           
[X, Z] = meshgrid(x, z);
R = abs(X);  

% Máscara (Todo lo que esté fuera del radio de la columna se marcará como NaN)
mask = (R <= rarmadura);

% Concentración en el agua
Concentracionagua = concentracionsup*(1 + alpha*abs(Z)/H);   % [kg/m^3]

% Concentración infiltrada en el hormigón
concentracionhorm = Concentracionagua .* exp(-lambda*(rarmadura - R));   % [kg/m^3]

% Derivadas analíticas en (r,z)
dC_dr = lambda .* concentracionhorm;      %derivada parcial respecto de r                 
dC_dz = -(concentracionsup*alpha/H) .* exp(-lambda*(rarmadura - R));     %derivada parcial respecto de z

% Paso de componente radial a componente en x (en y=0)(ahora esta en
% direccion radial, para poder usar quiver necesitamos que este en
% coordenadas cartesianas)
signX = sign(X);
signX(signX == 0) = 1;   
dC_dx = dC_dr .* signX;

% Aplicar máscara 
dC_dx(~mask) = NaN;
dC_dz(~mask) = NaN;

% Representación del campo vectorial en la sección (x,z)
figure;
quiver(X, Z, dC_dx, dC_dz, 0.5, 'filled');
axis equal;
xlabel('x [m]');
ylabel('z [m]');
title('Gradiente de concentración \nabla C_{col}(x,z) en una sección vertical');
set(gca,'YDir','normal');  %coloca el z=0 mas arriba que z=-H 
grid on;

3.3 Función de concentración en la armadura en función de la profundidad:

En la armadura de la columna (radio \(r = \rho_a\)), la concentración infiltrada coincide con la concentración del agua, ya que el factor exponencial vale 1:

[math]C_{\text{acero}}(z) = C_{\text{col}}(\rho_a,z) = C_{\text{agua}}(z) = C_0\left(1 + \alpha \frac{|z|}{H}\right),\qquad z\in[-H,0][/math]

Tomamos [math]C₀ = 0{,}3 kg/m³[/math], [math]α = 3[/math] y [math]H = 22 m[/math]. El umbral crítico para el inicio de la corrosión es

[math]C_{\text{crit}} = 0{,}05\ \text{kg/m}^3[/math]

Dibujamos \(C_{\text{acero}}(z)\) en función de la profundidad \(|z|\) y añadimos la línea horizontal \(C_{\text{crit}}\). En este caso, como \(C_{\text{acero}}(z) > C_{\text{crit}}\) para toda la profundidad, toda la armadura se encuentra en zona de riesgo de corrosión.

CONCENTRACION ARMADURA Y UMBRAL CRITICO

% CONCENTRACIÓN EN LA ARMADURA Y UMBRAL CRÍTICO
 % Parámetros (SI)
H     = 22;       
concentracionsup    = 0.3;      
alpha = 3;       
umbralcritico= 0.05;    

% Profundidad positiva desde la superficie hasta el fondo
nPts   = 200;
p  = linspace(0, H, nPts);                    
cacero = concentracionsup*(1 + alpha*p/H);                 

% Representación
figure;
plot(p, cacero, 'LineWidth', 2);
hold on;
yline(umbralcritico, 'r--', 'C_{crit}', 'LineWidth', 1.5);

xlabel('Profundidad  [m]');
ylabel('C_{acero} [kg/m^3]');
title('Concentración en la armadura C_{acero}(z) y umbral crítico C_{crit}');
grid on;

% Comentario en la consola sobre la zona de riesgo
if all(cacero > umbralcritico)
    disp('Toda la armadura está en zona de riesgo de corrosión (C_acero > C_crit para toda la profundidad).');
else
    disp('Solo la parte de la armadura donde C_acero > C_crit está en zona de riesgo de corrosión.');
end
hold off;

3.4 Superficies de isoconcentración dentro de una columna:

Las isosuperficies de concentración son superficies tridimensionales dentro del agua o del hormigón de las columnas que unen todos los puntos donde la concentración del contaminante tiene un mismo valor constante. En el Tanque de Tormentas de Arroyofresno describen las zonas del volumen (tanto del agua como del hormigón) donde la concentración es igual a un cierto nivel, lo que permite identificar regiones con igual grado de contaminación.

ISOSUPERFICIE DE CONCENTRACION DE UNA COLUMNA

% Isosuperficies 3D de C_col dentro de la columna
% Parámetros del enunciado
rcol = 0.75;    % m (radio de la columna)
rarmd = 0.65;    % m (radio de la armadura)
H     = 22;      % m (profundidad)
ConSup    = 0.3;     % kg/m3
alpha = 3;
lambda = 10;     % 1/m

% Resolución del mallado
Nx = 80; Ny = 80; Nz = 120;

% Mallado cartesiano
x = linspace(-rcol, rcol, Nx);
y = linspace(-rcol, rcol, Ny);
z = linspace(-H, 0, Nz);

[X, Y, Z] = meshgrid(x, y, z);

% Coordenada radial
R = sqrt(X.^2 + Y.^2);

% Concentración en el agua: Cagua(z)
Cagua = ConSup * (1 + alpha * abs(Z) / H);

% Concentración infiltrada en columna: Ccol(r,z)
Ccol = Cagua .* exp(-lambda * (rcol - R));

% Fuera del cilindro no se representa
Ccol(R > rcol) = NaN;

% Valores de isosuperficies a mostrar
isoSup = [0.1 0.2 0.3 0.4];

% Colores (RGBA)
colores = [ ...
    0.2 0.5 1.0 0.35;  % azul
    0.3 0.9 0.3 0.35;  % verde
    1.0 0.6 0.1 0.35;  % naranja
    0.7 0.2 0.7 0.35]; % morado

figure; hold on;

% Crear cada isosuperficie
for i = 1:length(isoSup)
    val = isoSup(i);

    % Extraer la isosuperficie
    p = patch(isosurface(X, Y, Z, Ccol, val));

    % Estética de la superficie
    p.FaceColor = colores(i,1:3);
    p.FaceAlpha = colores(i,4); %controla la opacidad 
    p.EdgeColor = 'none';
end

% Dibujar el cilindro de la armadura (rho = rho_a)
theta = linspace(0, 2*pi, 120);
z_cyl = linspace(-H, 0, 80);
[TH, ZZ] = meshgrid(theta, z_cyl);

Xc = rarmd * cos(TH);
Yc = rarmd * sin(TH);
Zc = ZZ;

surf(Xc, Yc, Zc, 'FaceColor', 'y', 'FaceAlpha', 0.95, 'EdgeColor', 'none');

% Ajustes de gráfico
xlabel('x (m)');
ylabel('y (m)');
zlabel('z (m)');
title('Isosuperficies 3D de C_{col} dentro de la columna');
axis equal;
view(3);
grid on;

3.5 Vida útil de la columna en la zona crítica:

Se tiene que el diámetro de las barras de armadura es [math] d_0 = 20 mm [/math] . Se considera como pérdida crítica cuando la barra pierde el 25% de su diámetro original. Calculamos la pérdida crítica debida a la corrosión:

[math]Δd_0 = 0.25 * 20 = 5 mm [/math]

Estos 5 mm son material total perdido respectó al diámetro de la barra y como la corrosión actua del exterior hacia el interior, se obtiene un radio de material corroído de 2.5 mm. La velocidad de corrosión dada es [math] v_c = 0.1 \frac{mm}{año} [/math]. Por último, calculamos el tiempo (t) que tarda en alcanzar esa perdida crítica del 25%.

[math]tiempo = t = \frac{Δd_0}{v_c} = \frac{5}{0.1} = 50 años [/math]

3.6 Masa total de contaminantes infiltrados en una columna:

Para calcular la masa de contaminantes infiltrados en una columna se ha usado código matlab debido a la ayuda que ofrece para cálculos mas complejos como en este caso.

Tras ejecutar el código, se obtiene una masa de contaminantes infiltrados por columna de 6.74 kg (6.739290 kg). Así mismo obtenemos que la masa infiltrada para el mismo volumen de la columna pero en el agua es de 29.16 kg (29.157907 kg). Si se comparan se observa que es bastante mayor la masa infiltrada en agua que en la columna. [math] \frac{MasaColumna}{MasaAguaEquiv}= \frac{6.74}{29.16} [/math]

% Datos
rcolumna = 0.75;    % m
H     = 22.0;    % m
C0    = 0.3;     % kg/m^3
alpha = 3.0;
lambda = 10.0;   % 1/m

% I_rho = rho_c/lambda - (1 - exp(-lambda*rho_c))/lambda^2
I_rho = rcolumna./lambda - (1 - exp(-lambda.*rcolumna))./(lambda.^2);

% J_z = C0 * H * (1 + alpha/2)
J_z = C0 .* H .* (1 + alpha./2);

% Masa contaminantes infiltrados en la columna
M_col = 2 * pi .* I_rho .* J_z;   % [kg]

% Masa contaminantes en volumen equivalente de agua
M_agua_eq = pi .* rcolumna.^2 .* J_z; % [kg]

% Relacion contaminantes columna y contaminantes vol agua equiv
fraction = M_col ./ M_agua_eq;
percentage = fraction * 100;

% Resultados
fprintf('I_rho = %.12g   (unidad: m^2)\n', I_rho);
fprintf('J_z   = %.12g   (unidad: kg/m^2 - resultado de int_z Cagua dz)\n', J_z);
fprintf('Masa infiltrada (M_col) = %.6f kg\n', M_col);
fprintf('Masa en agua equivalente (M_agua_eq) = %.6f kg\n', M_agua_eq);
fprintf('Fracción M_col / M_agua_eq = %.6f (%.3f%%)\n', fraction, percentage);

4 Decantación y distribución de contaminantes

4.1 Campo de concentración de contaminantes en el agua:

La concentración de contaminantes en el agua está gobernada por el proceso de decantación: las partículas más densas descienden con el tiempo y se acumulan en el fondo del tanque. Como consecuencia, la concentración solo depende de la profundidad z, y viene dada por:

[math]C_{\text{agua}}(z) = C_0\left(1 + \alpha \frac{|z|}{H}\right),\qquad z\in[-H,0][/math]

donde z=0 corresponde a la superficie y z=-H al fondo del tanque.

Dado que en esta región z≤0, se cumple ∣z∣=-z, por lo que la expresión se convierte en:

[math]C_{\text{agua}}(z) = C_0\left(1 - \alpha \frac{z}{H}\right)[/math]

Esta relación es lineal y creciente con la profundidad. El gráfico de C_{\text{agua}}(z) muestra una recta que aumenta desde la superficie hasta el fondo, reflejando el mecanismo físico de sedimentación: los contaminantes más pesados tienden a concentrarse en la parte inferior. En una representación espacial, esto se observa como:

              colores más intensos o densos cerca del fondo (mayor concentración),
              colores más suaves cerca de la superficie (menor concentración).

Por tanto, la distribución forma capas horizontales de igual concentración, coherentes con un proceso de sedimentación estratificada dentro del tanque.

PARED DEL TANQUE

CONCENTRACIÓN CONTAMINANTE DEL AGUA

% Concentración de contaminantes en el agua

% Parámetros del enunciado
Csup = 0.3;       % kg/m^3
alpha = 3;
H = 22;         % m
L = 290;        % m

% (a) Representación de la función Cagua(z)
z = linspace(-H, 0, 200);
Cz = Csup * (1 + alpha * abs(z)/H);  
figure;
plot(z, Cz,'r','LineWidth', 2);
xlabel('z (m)');
ylabel('C_{agua}(z)  (kg/m^3)');
title('Concentración de contaminantes en el agua');
grid on;

% (b) Mapa de colores en una pared del tanque

% Mallado para pared vertical
nx = 200;
nz = 200;
x = linspace(0, L, nx);   
z2 = linspace(0, -H, nz);  
[X, Z] = meshgrid(x, z2);

% Concentración (solo depende de z)
C = Csup * (1 + alpha * abs(Z)/H);

figure;
pcolor(X, Z, C);
shading interp;      
colorbar;
xlabel('x (m)');
ylabel('z (m)');
title('Mapa de colores de la concentración C_{agua}(z)');
set(gca,'YDir','normal');   % para que 0 esté arriba

4.2 Gradiente de concentración de contaminante en el agua:

El gradiente de un campo escalar indica la dirección en la que la magnitud crece más rápidamente. Para esta situación: [math]C_{\text{agua}}(z) = C_0\left(1 - \alpha \frac{z}{H}\right)[/math]

solo depende de la coordenada vertical z, por lo que el gradiente tendrá únicamente una componente en esa dirección:

[math]\nabla C_{\text{agua}} = \left(0, 0, \frac{dC}{dz}\right)[/math]

Al derivar obtenemos:

[math]\frac{dC}{dz} = -\frac{C_0 \alpha}{H}[/math]

Este valor es constante y negativo, lo que significa: La concentración aumenta hacia abajo. El gradiente apunta hacia la dirección donde la concentración crece, es decir, hacia el fondo del tanque. Físicamente, esto representa el flujo sedimentario de los contaminantes y la acumulación en zonas profundas.

GRADIENTE CONCENTRACIÓN CONTAMINANTES EN AGUA

%– Gradiente del campo de concentración

% Parámetros del problema
C0 = 0.3;      % kg/m^3
alpha = 3;
H = 22;        % m
L = 290;       % longitud de la pared

% Gradiente analítico
dCdz = -C0*alpha/H;   % constante

% Malla 2D para la pared vertical
nx = 60; 
nz = 50;
x = linspace(0, L, nx);
z = linspace(0, -H, nz);   % superficie 0 → fondo -H

[X, Z] = meshgrid(x, z);

% Campo de concentración (solo depende de z)
C = C0 * (1 + alpha*abs(Z)/H);

% Componentes del gradiente
U = zeros(size(X));           % dC/dx = 0
V = dCdz * ones(size(X));     % dC/dz = constante

% Representación gráfica
figure;
pcolor(X, Z, C);
shading interp;
colorbar;
hold on;

% Reducir densidad de flechas
skip = 3;

quiver(X(1:skip:end,1:skip:end), ...
       Z(1:skip:end,1:skip:end), ...
       U(1:skip:end,1:skip:end), ...
       V(1:skip:end,1:skip:end), ...
       'k', 'AutoScale', 'on');

xlabel('x (m)');
ylabel('z (m)');
title('Gradiente del campo de concentración \n \nabla C_{agua}');
set(gca,'YDir','normal'); % para que 0 esté arriba

4.3 Superficies de isoconcentración en el agua:

Concentración en función de la profundidad:[math]C_{\text{agua}}(z) = C_0\left(1 + \alpha \frac{|z|}{H}\right),\qquad z\in[-H,0][/math]

α=3. Despejando z de la ecuación de concentración:

[math] \begin{aligned} C_{\text{agua}} &= C_0\left(1 + \alpha \frac{|z|}{H}\right) \\ \frac{C_{\text{agua}}}{C_0} - 1 &= \alpha \frac{|z|}{H} \\ |z| &= \frac{H}{\alpha}\left(\frac{C_{\text{agua}}}{C_0} - 1\right) \\ z &= -\frac{H}{\alpha}\left(\frac{C_{\text{agua}}}{C_0} - 1\right) \end{aligned} [/math]

Eso significa que: C=0.4 aparecerá cerca de la superficie, C=1.0 aparecerá más cerca del fondo.

SUPERFICIES DE ISOCONCENTRACIÓN EN EL AGUA

% Superficies de isoconcentración en el agua del tanque

% Dimensiones del tanque
L = 290;   % longitud (m)
W = 140;   % anchura  (m)
H = 22;    % profundidad (m)

% Parámetros de la concentración
C0 = 0.3;      % concentración superficial (kg/m^3)
alpha = 3;     % factor de incremento hacia el fondo

% Mallado 3D del tanque
Nx = 40; Ny = 40; Nz = 80;
x = linspace(-L/2, L/2, Nx);
y = linspace(-W/2, W/2, Ny);
z = linspace(-H, 0, Nz);  % z va desde el fondo (-H) hasta superficie (0)

[X, Y, Z] = meshgrid(x, y, z);

% Campo de concentración Cagua(z)
Cagua = C0 .* (1 + alpha .* abs(Z) ./ H);

% Valores de isosuperficie
isoVals = [0.4, 0.6, 0.8, 1.0];

% Figura 3D
figure;
hold on;

colors = lines(length(isoVals));

for k = 1:length(isoVals)
    p = patch(isosurface(X, Y, Z, Cagua, isoVals(k)));
    isonormals(X, Y, Z, Cagua, p);
    set(p, 'FaceColor', colors(k,:), 'EdgeColor', 'none', ...
           'FaceAlpha', 0.35);
end

% Formato
xlabel('x (m)');
ylabel('y (m)');
zlabel('z (m)');
title('Superficies de isoconcentración en el tanque');
colormap turbo;
colorbar;
grid on;
view(45, 25);
axis tight
axis vis3d;
daspect([1 1 0.1])

hold off;

4.4 Masa total de contaminantes en el agua del tanque:

Consideramos que el fondo del tanque está en la cota \(z = -22\ \text{m}\) y que el nivel del agua alcanza la cota \(z = -7\ \text{m}\). Por tanto, el agua ocupa la región

[math] -22 \le z \le -7, [/math]

es decir, una columna de agua de espesor \(15\ \text{m}\).

La concentración de contaminantes en el agua viene dada por

[math] C_{\text{agua}}(z) = C_0\left(1 + \alpha \frac{|z|}{H}\right),\qquad C_0 = 0{,}3\ \text{kg/m}^3,\quad \alpha = 3,\quad H = 22\ \text{m}. [/math]

Como en todo el agua se cumple \(z \le 0\), tenemos \(|z| = -z\) y, por tanto,

[math] C_{\text{agua}}(z) = C_0\left(1 - \alpha \frac{z}{H}\right),\qquad z\in[-22,-7]. [/math]

En planta, el tanque tiene dimensiones \(L = 290\ \text{m}\), \(W = 140\ \text{m}\). Hay aproximadamente

[math] N \approx \frac{L}{15}\,\frac{W}{15} [/math]

columnas cilíndricas de radio \(r_c = 0{,}75\ \text{m}\), que restan volumen al agua. El área efectiva ocupada por agua en planta es

[math] A_{\text{fondo}} = LW - N\,\pi r_c^2. [/math]

La masa total de contaminantes en el agua se obtiene integrando la concentración en el volumen ocupado por el agua:

[math] M_{\text{tot}} = \iiint C_{\text{agua}}(z)\,dV = A_{\text{fondo}} \int_{-22}^{-7} C_0\left(1 - \alpha \frac{z}{H}\right)\,dz. [/math]

Evaluando la integral:

[math] \int_{-22}^{-7} \left(1 - \alpha \frac{z}{H}\right)\,dz = \left[\, z - \frac{\alpha}{2H}z^2 \right]_{-22}^{-7} \approx 44{,}66, [/math]

de modo que

[math] M_{\text{tot}} = A_{\text{fondo}}\,C_0\cdot 44{,}66. [/math]

Sustituyendo los valores numéricos \(L = 290\ \text{m}\), \(W = 140\ \text{m}\), \(r_c = 0{,}75\ \text{m}\), \(C_0 = 0{,}3\ \text{kg/m}^3\), \(\alpha = 3\), \(H = 22\ \text{m}\), se obtiene:

[math] A_{\text{fondo}} \approx 4{,}03\cdot 10^4\ \text{m}^2, [/math]
[math] M_{\text{tot}} \approx 5{,}40\cdot 10^5\ \text{kg}. [/math]

Por tanto, la masa total de contaminantes en el agua del tanque es aproximadamente

[math] M_{\text{tot}} \approx 5{,}40\cdot 10^5\ \text{kg} \;\approx\; 5{,}4\cdot 10^2\ \text{toneladas}. [/math]

4.5 Caso pasadas 24 horas:

Después de 24 horas de decantación, el 60% de la masa total de contaminantes se deposita en el fondo formando una capa de lodo. El lodo tiene densidad aparente

[math] \rho_{\text{lodo}} = 1400\ \text{kg/m}^3. [/math]

Partimos de la masa total de contaminantes calculada en el apartado anterior:

[math] M_{\text{tot}} \approx 5{,}40\cdot 10^5\ \text{kg}. [/math]

(a) Masa de contaminantes depositados

Se deposita el 60% de la masa total:

[math] M_{\text{dep}} = 0{,}60\,M_{\text{tot}} \approx 3{,}24\cdot 10^5\ \text{kg}. [/math]

(b) Volumen de lodo

El volumen de lodo se obtiene mediante

[math] V_{\text{lodo}} = \frac{M_{\text{dep}}}{\rho_{\text{lodo}}} \approx \frac{3{,}24\cdot 10^5}{1400} \approx 2{,}31\cdot 10^2\ \text{m}^3. [/math]

(c) Espesor promedio de la capa de lodo (teniendo en cuenta las columnas)

El área efectiva del fondo donde se deposita el lodo es la misma que en el apartado (14):

[math] A_{\text{fondo}} = LW - N\,\pi r_c^2 \approx 4{,}03\cdot 10^4\ \text{m}^2. [/math]

El espesor medio de lodo por evento de lluvia es

[math] e_1 = \frac{V_{\text{lodo}}}{A_{\text{fondo}}} \approx \frac{2{,}31\cdot 10^2}{4{,}03\cdot 10^4} \approx 5{,}7\cdot 10^{-3}\ \text{m} \;\approx\; 5{,}7\ \text{mm}. [/math]

(d) Espesor acumulado en un año y frecuencia de limpieza

Si el tanque se usa de media 8 veces al año y no se limpia entre eventos:

[math] e_{\text{año}} = 8\,e_1 \approx 8\cdot 5{,}7\cdot 10^{-3} \approx 4{,}6\cdot 10^{-2}\ \text{m} \;\approx\; 4{,}6\ \text{cm}. [/math]

Como

[math] e_{\text{año}} \approx 4{,}6\ \text{cm} \lt 30\ \text{cm}, [/math]

una limpieza anual del tanque sería suficiente con este modelo.

% EJERCICIO (15): LODO EN EL FONDO DEL TANQUE (USANDO M_tot DEL APARTADO 14)

% Datos geométricos
L   = 290;      % [m]
W   = 140;      % [m]
rc  = 0.75;     % [m]
sep = 15;       % [m]

% Área efectiva del fondo (igual que en el apartado 14)
Ncol   = (L/sep)*(W/sep);
A_cols = Ncol * pi*rc^2;
A_fondo = L*W - A_cols;        % [m^2]

% Masa total de contaminantes (resultado del 14)
M_tot = 5.3968e5;              % [kg] ≈ 5,40·10^5 kg

% (a) Masa depositada
M_dep = 0.60 * M_tot;          % [kg]

% (b) Volumen de lodo
rho_lodo = 1400;               % [kg/m^3]
V_lodo = M_dep / rho_lodo;     % [m^3]

% (c) Espesor por evento
e1 = V_lodo / A_fondo;         % [m]

% (d) Espesor anual (8 eventos)
nEventos = 8;
e_anual  = nEventos * e1;      % [m]
e_limite = 0.30;               % [m] umbral de limpieza

fprintf('Masa total (14):          M_tot   = %.2f kg\n', M_tot);
fprintf('(a) Masa depositada:      M_dep   = %.2f kg\n', M_dep);
fprintf('(b) Volumen de lodo:      V_lodo  = %.2f m^3\n', V_lodo);
fprintf('(c) Espesor por evento:   e1      = %.4f m (%.2f mm)\n', ...
        e1, e1*1000);
fprintf('(d) Espesor anual (8 ev): e_anual = %.4f m (%.2f cm)\n', ...
        e_anual, e_anual*100);
if e_anual > e_limite
    fprintf('-> No es suficiente una limpieza anual.\n');
else
    fprintf('-> Una limpieza anual es suficiente.\n');
end