Diferencia entre revisiones de «La presa de El Atazar (Grupo 1)»

r0 = 200;  % Radio base de la presa (aproximado)
b = 35;    % Curvatura del arco parabólico
H = 134;   % Altura de la presa
theta = linspace(3*pi/4, 5*pi/4, 100);  % Ángulo en radianes
z = linspace(0, H, 100);  % Altura en metros

% Creación del mallado en 2D para los ángulos y las alturas
[Z, Theta] = meshgrid(z, theta);

% Cálculo del radio para cada punto en la superficie
R = r0 + b * (1 - (Z.^2) / H^2);

% Conversión de coordenadas cilíndricas a cartesianas
X = R .* cos(Theta);
Y = R .* sin(Theta);

% Crear la figura
figure;
surf(X, Y, Z, 'FaceColor', 'b', 'EdgeColor', 'none');  % Gráfico 3D de la superficie
xlabel('Eje X (m)');
ylabel('Eje Y (m)');
zlabel('Eje Z (m)');
title('Superficie parametrizada de la presa');

% Configuración visual
axis equal;                      % Asegura la misma escala en todos los ejes
view(3);                          % Vista en 3D
grid on;                          % Muestra la cuadrícula

% Parámetros físicos
rho = 1000;  % Densidad del agua (kg/m^3)
g = 9.81;    % Aceleración de la gravedad (m/s^2)

% Cálculo de la presión en función de la altura
P = rho * g * (H - Z);  % Presión hidrostática en la superficie

% Visualización de la superficie con colores representando la presión
figure;
surf(X, Y, Z, P, 'EdgeColor', 'none');  % Superficie con mapa de colores basado en P
xlabel('Eje X (m)');
ylabel('Eje Y (m)');
zlabel('Eje Z (m)');
title('Presión sobre la superficie de la presa');

% Configuración de colores y visualización
colormap(jet);     % Paleta de colores para la presión
colorbar;          % Barra de colores para interpretar la presión
axis equal;        % Igualar las escalas de los ejes
view(3);           % Vista tridimensional
grid on;           % Mostrar la cuadrícula

% Derivadas numéricas para obtener vectores tangentes
[Rx_theta, Rx_z] = gradient(X);
[Ry_theta, Ry_z] = gradient(Y);
[Rz_theta, Rz_z] = gradient(Z);
% Producto cruzado de vectores tangentes para obtener el vector normal
Nx = Ry_theta .* Rz_z - Rz_theta .* Ry_z;
Ny = Rz_theta .* Rx_z - Rx_theta .* Rz_z;
Nz = Rx_theta .* Ry_z - Ry_theta .* Rx_z;
% Magnitud del vector normal
N_magnitude = sqrt(Nx.^2 + Ny.^2 + Nz.^2);
% Vector normal unitario
Nx_unit = Nx ./ N_magnitude;
Ny_unit = Ny ./ N_magnitude;
Nz_unit = Nz ./ N_magnitude;
% Campo de fuerza de presión
Fx = -P .* Nx_unit;
Fy = -P .* Ny_unit;
Fz = -P .* Nz_unit;
% Representación del campo vectorial sobre la superficie
figure;
quiver3(X, Y, Z, Fx, Fy, Fz, 0.5, 'Color', 'yellow'); % Escala ajustable con el factor 2
hold on;
surf(X, Y, Z, P, 'EdgeColor', 'none', 'FaceAlpha',1); % Superficie semitransparente
colorbar;
xlabel('X (m)');
ylabel('Y (m)');
zlabel('Z (m)');
title('Campo vectorial de la fuerza de presión sobre la presa');
view(3);
grid on;
% Fijamos theta = 0 (corte longitudinal)
theta_cut = 0; % Corte en el plano x-z
X_cut = X(:, 51); % Extraemos el corte longitudinal (theta cerca de 0)
Z_cut = Z(:, 51); % Correspondiente en altura
P_cut = P(:, 51); % Presión en el plano
% Vectores normales en el corte
Nx_cut = Nx(:, 51); % Componente x del vector normal
Nz_cut = Nz(:, 51); % Componente z del vector normal
N_magnitude_cut = sqrt(Nx_cut.^2 + Nz_cut.^2);
% Normalización del vector normal
Nx_unit_cut = Nx_cut ./ N_magnitude_cut;
Nz_unit_cut = Nz_cut ./ N_magnitude_cut;
% Campo de fuerza de presión en el corte
Fx_cut = -P_cut .* Nx_unit_cut;
Fz_cut = -P_cut .* Nz_unit_cut;
% Gráfica del campo vectorial en el corte
figure;
quiver(X_cut, Z_cut, Fx_cut, Fz_cut, 'b', 'LineWidth', 1.5); % Vectores de fuerza
hold on;
plot(X_cut, Z_cut, 'k-', 'LineWidth', 2); % Contorno de la presa en el plano
xlabel('X (m)');
ylabel('Z (m)');
title('Campo de fuerza de presión en el plano longitudinal');
grid on;