Placa Plana (Grupo 09)

Trabajo realizado por estudiantes
Título	Placa plana. Grupo 09
Asignatura	Teoría de Campos
Curso	2025-26
Autores	Rafael Gonzalez Gomez Mario Belinchón Buendía Alejandro Morales Tari Pablo Márquez Blanco
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

1 Introducción

Consideramos la sección longitudinal de una viga en voladizo de geometría rectangular (en dimensión 2). Esta viga ocupa la región mostrada en la Figura 1, y permanece anclada en la pared vertical situada en el lado izquierdo. La zona que vamos a estudiar se describe mediante:

[math] (x,y) \in [0,4] \times [f(x),, g(x)], \qquad f(x)=\frac{x}{8}, \qquad g(x)= 2-\frac{x}{8} [/math]

En esta placa se consideran dos magnitudes físicas relevantes. La primera es el campo de temperaturas, modelado por la función:

[math] T(x,y)=\big(1+(y-1)^2\big),(4-x) [/math]

La segunda magnitud corresponde al campo de desplazamientos [math]\vec u(x,y)[/math], producido por una fuerza externa —cuyo valor exacto no conocemos— y que deforma la placa.

Para describir la deformación, sea [math]\vec r_0(x,y)=x,\vec i + y,\vec j[/math] el vector de posición inicial de cada punto. Entonces, tras aplicarse la fuerza, la nueva posición de un punto queda dada por:

[math] \vec r_d(x,y)=\vec r_0(x,y)+\vec u(x,y) [/math]

Supondremos que la fuerza aplicada induce un desplazamiento expresado, en coordenadas polares [math](\rho,\theta)[/math], mediante:

[math] \vec u(\rho,\theta)= -\frac{1}{20},\rho^2\cos\theta ; \vec e_\theta [/math]

Este modelo permite estudiar tanto la distribución térmica en la viga como la deformación resultante, proporcionando una base adecuada para analizar su comportamiento estructural.

2 Mallado

Para comenzar, construiremos un mallado de la sección longitudinal que nos permita representar los puntos interiores del sólido. Para ello, hemos empleado un paso de muestreo de [math]h = \tfrac{1}{10}[/math] y hemos definido la malla sobre el rectángulo que delimita los ejes coordenados:

[math] (x, y) \in [-1,, 5] \times [-1,, 3] [/math]

Este mallado servirá como base para discretizar la región de interés y facilitar los cálculos posteriores.

h = 1/10;
color_malla = [0, 0.6, 0.6]; % Color 

% Coordenadas
u = 0 : h : 4;
v = 0 : h : 2;
[U, V] = meshgrid(u, v);

% Calculamos los bordes 
y_abajo = U ./ 8;
y_arriba = 2 - (U ./ 8);

% Interpolamos
factor_altura = V ./ 2; % Va de 0 (abajo) a 1 (arriba)

X = U;
Y = y_abajo + factor_altura .* (y_arriba - y_abajo);

% Visualización
figure('Color', 'w'); hold on;

% Dibujar las líneas VERTICALES
plot(X, Y, 'Color', color_malla, 'LineWidth', 0.5);

% Dibujar las líneas HORIZONTALES
plot(X', Y', 'Color', color_malla, 'LineWidth', 0.5);

% Dibujar el contorno 
plot(u, u./8, 'k-', 'LineWidth', 2); % Borde Abajo
plot(u, 2 - u./8, 'k-', 'LineWidth', 2); % Borde Arriba
plot([0 0], [0 2], 'k-', 'LineWidth', 2); % Borde Izquierdo
plot([4 4], [0.5 1.5], 'k-', 'LineWidth', 2); % Borde Derecho

% Configuración final
axis([-1 5 -1 3]); % Zoom/Encuadre exacto
xlabel('x'); ylabel('y');
title('Mallado de la placa plana');
grid on;
box on;
hold off;

3 Curvas de nivel

A continuación se representan las curvas de nivel de la temperatura con el fin de identificar visualmente la zona donde esta alcanza su valor máximo. A partir de la gráfica puede apreciarse que los puntos donde la temperatura es mayor son:

[math](x,y) = (0,0) \quad \text{y} \quad (0,2)[/math]

Además, calculamos el gradiente de la temperatura para analizar la dirección de mayor variación de esta magnitud. El gradiente viene dado por:

[math] \nabla T(x,y) = \left( \frac{\partial T}{\partial x},, \frac{\partial T}{\partial y} \right) = \Big( -\big(1+(y-1)^2\big),; (4-x)\cdot 2(y-1) \Big). [/math]

En la representación gráfica se observa que el vector gradiente es ortogonal a las curvas de nivel, tal como establece la teoría: apunta siempre en la dirección de crecimiento más rápido de la temperatura.

Figura 2

clc; clear; close all;
% Definición 
h=0.1;               
x=0:h:4;                  
y=0:h:2;           
[X,Y] = meshgrid(x, y);

% Funciones  
f=@(x)x./8;  
g=@(x)2-x./8;  


EnPlaca=Y>=f(X)&Y< g(X);

% Definición de la Temperatura y el Gradiente

T=(1+(Y-1).^2).*(4-X);


T_plot=T;
T_plot(~EnPlaca)=NaN;

% Cálculo numérico del gradiente para los vectores 

U=-(1+(Y-1).^2);          
V=2.*(Y-1).*(4-X);   


U(~EnPlaca)=NaN;
V(~EnPlaca)=NaN;

% Gráfica
figure(1);
clf; 
hold on;
axis equal; 
axis([-0.5 4.5 -0.5 2.5]);
grid on;
title({'Isotermas y Gradiente \nabla T', 'T(x,y) = (1+(y-1)^2)(4-x)'});
xlabel('x'); ylabel('y');

% Mapa de calor 
[C_fill, h_fill] = contourf(X, Y, T_plot, 20, 'LineStyle', 'none');
colormap(jet); 
cb=colorbar;
ylabel(cb,'Temperatura T');

% Curvas de nivel
[C, h_cont]=contour(X,Y,T_plot,10,'k','LineWidth',1.5);
clabel(C,h_cont,'FontSize',9,'Color','k','FontWeight','bold');

% Campo Vectorial 

step=3; 
quiver(X(1:step:end,1:step:end),Y(1:step:end, 1:step:end), ...
U(1:step:end,1:step:end),V(1:step:end,1:step:end), ...
1.2,'k','LineWidth',1); % '1.2' es el factor de escala de las flechas

% Dibujar bordes de la placa
plot(x,f(x),'k-','LineWidth',2);
plot(x,g(x),'k-','LineWidth',2); 
plot([0 0],[f(0) g(0)],'k-','LineWidth',2); 
plot([4 4],[f(4) g(4)],'k-','LineWidth',2); 

% Señalar Máximo
plot(0,0,'rp','MarkerSize',12,'MarkerFaceColor','r');
plot(0,2,'rp','MarkerSize',12,'MarkerFaceColor','r');
text(0.1,0,'Max T','Color','white','FontWeight','bold');

hold off;

4 Ley de Fourier

Según la ley de Fourier, el flujo de energía calorífica [math]\vec{Q}[/math] viene dado por:

[math]\vec{Q} = -\kappa , \nabla T[/math]

donde [math]\kappa[/math] representa la conductividad térmica del material. En nuestro caso, consideramos [math]\kappa = 1[/math], de modo que la expresión se simplifica a:

[math]\vec{Q} = - \nabla T[/math]

Bajo esta suposición, en la representación gráfica puede apreciarse que el flujo calorífico es simplemente el opuesto del gradiente de la temperatura, señalando siempre la dirección de máxima disminución térmica.

% Configuración del Mallado
h = 0.2;
x = 0 : h : 4;
y = 0 : h : 2;
[X, Y] = meshgrid(x, y);

% Definir la Geometría 
y_abajo = X ./ 8;
y_arriba = 2 - (X ./ 8);


% Definir la Temperatura 
T_func = @(x,y) (1 + (y - 1).^2) .* (4 - x);
T = T_func(X, Y);

% Cálculos Físicos 
% A) Calculamos el Gradiente (nabla T) mediante derivadas numéricas
[dTdx, dTdy] = gradient(T, h, h);

% Calculamos el Flujo de Calor 
Qx = -dTdx;
Qy = -dTdy;

% Limpieza 

Qx(~dentro) = NaN;
Qy(~dentro) = NaN;
T(~dentro) = NaN;

% Visualización
figure('Color', 'w', 'Name', 'Flujo de Calor Fourier');
hold on;

% FONDO: Mapa de Calor 

contourf(X, Y, T, 20, 'LineStyle', 'none');
colormap('jet'); % Paleta de colores: Azul (frío) -> Rojo (caliente)

% Barra de color explicativa
c = colorbar;
c.Label.String = 'Temperatura T(x,y)';
c.Label.FontSize = 10;

% El Campo Vectorial.
quiver(X, Y, Qx, Qy, 'k', 'LineWidth', 1.2, 'AutoScaleFactor', 1.5);

% Dibujo del contorno.
plot([0 4], [0 0.5], 'k-', 'LineWidth', 2); 
plot([0 4], [2 1.5], 'k-', 'LineWidth', 2); 
plot([0 0], [0 2], 'k-', 'LineWidth', 2); 
plot([4 4], [0.5 1.5], 'k-', 'LineWidth', 2); 

% Configuración Final
axis equal;
axis([-0.5 4.5 -0.5 2.5]); 
grid on;
xlabel('Eje x'); ylabel('Eje y');
title({'Vector Flujo de Calor', ...
'El calor fluye desde las capas más calientes hacia las más frías'});

hold off;

5 Gradiente Térmico

El punto de máxima variación térmica corresponde al lugar donde el gradiente de la temperatura alcanza su valor máximo. En la práctica, esto se obtiene igualando la derivada del gradiente a cero.

Mientras que la variación térmica sigue la dirección del gradiente, el flujo de calor más intenso se orienta en sentido contrario, es decir, hacia donde disminuye la temperatura.

Paso 1: Definición de la función para maximizar

Llamamos a la función:

[math]G(x,y) = |\nabla T|^2[/math]

Desarrollando los términos:

[math]G(x,y) = \big(1+(y-1)^2\big)^2 + \big(2(y-1)(4-x)\big)^2[/math]

[math]G(x,y) = (y^2 - 2y + 2)^2 + 4(y-1)^2(4-x)^2[/math]

Paso 2: Cálculo de los puntos críticos

Para maximizar la función, derivamos y resolvemos el sistema:

[math]\frac{\partial G}{\partial x} = -8 (y-1)^2 (4-x) = 0[/math] [math]\frac{\partial G}{\partial y} = 4 (y-1) \big(1 + (y-1)^2 + 2 (4-x)^2 \big) = 0[/math]

Al resolver este sistema, encontramos que los mínimos se sitúan sobre la recta:

[math]y = 1[/math]

Como no hay máximos aislados dentro del dominio, los máximos deben encontrarse en la frontera:

Borde izquierdo: [math]x=0[/math]

Borde derecho: [math]x=4[/math]

Curva inferior: [math]y = x/8[/math]

Curva superior: [math]y = 2 - x/8[/math]

Tras analizar los límites del recinto, los máximos de la función del gradiente térmico se localizan en el borde izquierdo, [math]x=0[/math]. Las coordenadas de los máximos son:

P(1): x=0.00, y=0.00, Dirección del vector: [-2.00, -8.00] P(2): x=0.00, y=2.00, Dirección del vector: [-2.00, 8.00]

Paso 3: Cálculo de la máxima variación térmica

El valor de la norma del gradiente es igual para P(1) y P(2). Calculamos, por ejemplo, para P(1) = (0,0):

[math]\nabla T(0,0) = \big(-(1 + (0-1)^2), 2(0-1)(4-0)\big) = (-2,-8)[/math]

[math]|\nabla T(0,0)| = \sqrt{(-2)^2 + (-8)^2} = \sqrt{4 + 64} = 2\sqrt{17} \approx 8.246[/math]

Por lo tanto, la máxima variación térmica en la placa se produce en los puntos P(1) y P(2) y tiene un valor aproximado de 8.246.

Figura 4:Máxima Variación de T

clc; clear; close all;

h=0.05;
[X,Y]=meshgrid(0:h:4,0:h:2);
f=X./8;g=2-X./8;
EnPlaca=Y>=f&Y<=g;

% Gradiente y Norma
Tx=-(1+(Y-1).^2); 
Ty=2.*(Y-1).*(4-X);

% Gradiente
NormaGrad=sqrt(Tx.^2+Ty.^2);
NormaGrad(~EnPlaca)=0; 

% Máximo
max_val=max(NormaGrad(:));
[filas,cols]=find(NormaGrad==max_val);

% puntos máximos
x_max=X(filas,cols);
y_max=Y(filas,cols);

% Componentes del vector dirección 
u_dir=Tx(filas,cols);
v_dir=Ty(filas,cols);

% Visualización
figure('Name','Máxima Variación de T');
hold on;axis equal;axis([-0.5 4.5 -0.5 2.5]);grid on;
title({'Puntos de Máxima Variación de Temperatura', '(Color = Magnitud del Gradiente T)'});
xlabel('x'); ylabel('y');

% Mapa de calor
NormaGrad(~EnPlaca)=NaN;
contourf(X,Y,NormaGrad,20,'LineStyle','none');
colormap(jet);colorbar;

% Contorno
plot([0 4 4 0 0],[0 0.5 1.5 2 0],'k','LineWidth',1.5);

% Puntos rojos y dirección
plot(x_max, y_max,'ro','MarkerSize',7.5,'MarkerFaceColor','r');

% Flechas de dirección
quiver(x_max,y_max,u_dir,v_dir,0.5,'k','LineWidth',2,'MaxHeadSize',0.5);

% Resultados
fprintf('La variación máxima es %.2f\n', max_val);
fprintf('Ocurre en los puntos:\n');
for i=1:length(x_max)
    fprintf('P(%d): x=%.2f, y=%.2f. Dirección del vector: [%.2f, %.2f]\n',i,x_max(i),y_max(i),u_dir(i),v_dir(i));
end

6 Campo de Desplazamientos

El campo de desplazamientos de un sólido se representa mediante un campo vectorial, que indica la dirección y sentido del desplazamiento de cada punto mediante vectores (representados como flechas rojas).

Como se puede observar, los puntos que se encuentran más próximos al punto de aplicación de la fuerza experimentan los desplazamientos más grandes, mientras que los situados más alejados apenas se mueven.

h = 0.1; 
x = 0:h:4;
y = 0:h:2;
[X, Y] = meshgrid(x, y);

% Definición de la región 
Condicion = (Y >= X./8) & (Y <= 2 - X./8);

% Aplicamos Máscara
X(~Condicion) = NaN;
Y(~Condicion) = NaN;

% Cálculo del Campo de Desplazamientos 

Rho = sqrt(X.^2 + Y.^2);
Theta = atan2(Y, X);


U_theta_val = -1/20 .* Rho.^2 .* cos(Theta);


Ux = -U_theta_val .* sin(Theta);
Uy = U_theta_val .* cos(Theta);


% Calcular la posición deformada

% Posición final = Posición inicial + Desplazamiento
X_final = X + Ux;
Y_final = Y + Uy;


% Comparativa

figure('Name', 'Deformación de la Placa', 'Color', 'w');

% PLACA ORIGINAL
subplot(1, 2, 1);
% Usamos mesh con Z=0 para dibujar la 'rejilla' plana
mesh(X, Y, zeros(size(X)), 'EdgeColor', 'b');
view(2); % Vista superior (2D)
axis equal; grid on;
title('Placa antes del desplazamiento');
xlabel('Eje X'); ylabel('Eje Y');
axis([-0.5 4.5 -0.5 2.5]);

% PLACA DEFORMADA 
subplot(1, 2, 2);
mesh(X_final, Y_final, zeros(size(X)), 'EdgeColor', 'b');
view(2);
axis equal; grid on;
title('Placa después del desplazamiento');
xlabel('Eje X'); ylabel('Eje Y');
axis([-0.5 4.5 -0.5 2.5]);


 Apartado 5 

figure('Name', 'Campo Vectorial de Desplazamientos');
hold on; axis equal; grid on;
title('Campo de Desplazamientos');

% Dibujamos contorno original de referencia
plot(X(:), Y(:), '.k', 'MarkerSize', 1);

% Dibujamos los vectores 
quiver(X, Y, Ux, Uy, 0, 'r', 'LineWidth', 1.5);

xlabel('Eje X'); ylabel('Eje Y');
xlim([-0.5 4.5]); ylim([-0.5 2.5]);

7 Desplazamiento del Sólido

Debido a la fuerza de magnitud desconocida, el sólido experimenta un desplazamiento descrito por el campo vectorial:

[math]\vec{u}(\rho, \theta) = -\frac{1}{20},\rho^2 \cos \theta ; \vec{e}_\theta[/math]

El sólido se deforma siguiendo este campo de desplazamientos. Como era de esperar, los puntos más próximos al punto de aplicación de la fuerza son los que sufren los mayores desplazamientos, mientras que los situados más lejos apenas se mueven.

Se adjuntan las gráficas que muestran el desplazamiento sufrido. (Nota: el desplazamiento inferior no se aprecia en la visualización actual; sería necesario ajustar los ejes en el código para visualizarlo correctamente).

% Definición de Geometría y Desplazamientos
h = 0.1; % Paso del mallado
x = 0:h:4;
y = 0:h:2;
[X, Y] = meshgrid(x, y);

% Definición de la Región 

% Aplicamos máscara 
X(~Condicion) = NaN;
Y(~Condicion) = NaN;

% Cálculo de Desplazamientos

Rho = sqrt(X.^2 + Y.^2);
Theta = atan2(Y, X);

% Valor escalar de u_theta
U_theta_val = -1/20 .* Rho.^2 .* cos(Theta);

% Proyección a cartesianas 
Ux = -U_theta_val .* sin(Theta);
Uy = U_theta_val .* cos(Theta);

% --- Posición Final Deformada ---
X_final = X + Ux;
Y_final = Y + Uy;

% Gráfica Comparativa de Sólidos Mallados 
figure('Name', 'Comparación Antes/Después', 'Color', 'w', 'Position', [100, 100, 1000, 500]);

% Placa Original
subplot(1, 2, 1);
hold on;

mesh(X, Y, zeros(size(X)), 'EdgeColor', [0 0.7 0.9], 'FaceColor', 'none');

view(2); 
axis equal; 
grid on;
box on; 

title('Placa antes del desplazamiento', 'FontWeight', 'bold');
xlabel('Eje X');
ylabel('Eje Y');
axis([-1.5 4.5 0 2.5]);

% Placa Deformada
subplot(1, 2, 2);
hold on;
% Dibujamos la malla usando las coordenadas finales (X_final, Y_final)
mesh(X_final, Y_final, zeros(size(X)), 'EdgeColor', [0 0.7 0.9], 'FaceColor', 'none');

view(2);
axis equal;
grid on;
box on;

title('Placa después del desplazamiento', 'FontWeight', 'bold');
xlabel('Eje X');
ylabel('Eje Y');
axis([-1.5 4.5 0 2.5]); 

sgtitle('Visualización del Sólido Mallado bajo Deformación'); % Título

8 Divergencia

La divergencia [math]\nabla \cdot \vec{u}[/math] mide la variación de volumen local en el sólido. En nuestro caso, el cálculo da como resultado:

[math]\nabla \cdot \vec{u} = \frac{y}{20}[/math]

Este resultado indica que la placa experimenta una expansión positiva, que aumenta linealmente con la altura, debido al estiramiento de las fibras superiores.

En la representación gráfica se observa la distribución de la divergencia mediante franjas de color:

En azul se encuentra la base, donde el volumen permanece prácticamente constante.

Hacia la parte superior, el color se degrada progresivamente hasta un rojo intenso, indicando el punto de máxima expansión, independiente de la posición horizontal 𝑥 x.

Figura 7

clc; clear; close all;
% Definir Malla y Placa
h=0.05;
[X,Y]=meshgrid(0:h:4, 0:h:2);
f=X./8;g=2-X./8; % Límites
EnPlaca=Y>=f&Y<=g; % Máscara lógica

% Divergencia
Div=Y./20; 
Div(~EnPlaca)=NaN;

% Visualización
figure('Name','Divergencia');clf;hold on;
axis equal; axis([-0.5 4.5 -0.5 2.5]);box on;
title('Divergencia: Cambio de Volumen');
xlabel('X'); ylabel('Y');
[C, h_cont]=contourf(X, Y, Div, 20,'LineStyle','none');
colormap(jet);
c=colorbar;ylabel(c,'Valor de Divergencia');

% Dibujar bordes negros
x_b = linspace(0, 4, 200);
plot(x_b, x_b./8,'k','LineWidth',1);
plot(x_b, 2 - x_b./8,'k','LineWidth',1);
plot([0 0],[0 2],'k','LineWidth',1);
plot([4 4],[0.5 1.5],'k','LineWidth',1);

% Señalar puntos clave
plot(0,2,'kp','MarkerSize',12,'MarkerFaceColor','r');
text(0.1, 2.1,'Max Expansión','FontWeight','bold');


plot(0,0,'kp','MarkerSize',12,'MarkerFaceColor','b');
text(0.1,-0.1,'Divergencia Nula','FontWeight','bold');
hold off;

9 Rotacional

Se calcula el rotacional [math]|\nabla \times \vec{u}|[/math] en todos los puntos del sólido.

Los resultados muestran que los puntos más alejados del origen y alineados con el eje 𝑥 x experimentan los mayores valores de rotacional, indicando que estas zonas sufren la mayor rotación local debido al campo de desplazamientos aplicado.

% Mallado y Región
h = 0.05;
x = 0:h:4;
y = 0:h:2;
[X, Y] = meshgrid(x, y);

% Definición de la placa
Condition = (Y >= X./8) & (Y <= 2 - X./8);
X(~Condition) = NaN;
Y(~Condition) = NaN;

% Definición del Mallado y Región
h = 0.05;
x = 0:h:4;
y = 0:h:2;
[X, Y] = meshgrid(x, y);

% Definición de la placa 
Condition = (Y >= X./8) & (Y <= 2 - X./8);
X(~Condition) = NaN;
Y(~Condition) = NaN;

% 2. Cálculo del Rotacional (Magnitud)
%

10 Tensores Deformacionales

Se define el tensor de deformaciones como la parte simétrica del gradiente del campo de desplazamientos [math]\vec{u}[/math]:

[math]\varepsilon(\vec{u}) = \frac{1}{2} \big( \nabla \vec{u} + (\nabla \vec{u})^T \big)[/math]

En un medio elástico lineal, isótropo y homogéneo, los desplazamientos permiten expresar el tensor de tensiones [math]\sigma[/math] mediante la relación de Hooke generalizada:

[math]\sigma = \lambda , (\nabla \cdot \vec{u}) , I + 2 \mu , \varepsilon(\vec{u})[/math]

donde:

[math]I[/math] es el tensor identidad en [math]\mathbb{R}^3[/math],

[math]\lambda[/math] y [math]\mu[/math] son los coeficientes de Lamé, que dependen de las propiedades elásticas del material.

Tomando [math]\lambda = \mu = 1[/math], se representan las tensiones normales en las direcciones de los ejes coordenados:

[math]\vec{i} \cdot \sigma \cdot \vec{i}[/math] — tensión normal en la dirección del eje 𝑥 x,

[math]\vec{j} \cdot \sigma \cdot \vec{j}[/math] — tensión normal en la dirección del eje 𝑦 y,

[math]\vec{k} \cdot \sigma \cdot \vec{k}[/math] — tensión normal en la dirección del eje 𝑧 z.

Solo se representan las componentes no nulas, mostrando cómo se distribuyen las tensiones normales en el sólido bajo el campo de desplazamientos aplicado.

Figura 9

clc; clear; close all;

% Geometría de la placa
h = 0.1;                        
x = 0:h:4;
num_puntos_y = round(2/h);
s = linspace(0,1,num_puntos_y+1);% parámetro vertical [0,1]

[U,S] = meshgrid(x,s);
f = U./8;                        
g = 2 - U./8;                   

X = U;                           
Y = f + S.*(g - f);           

% Campo de desplazamientos aplicado

u = (X.*Y)/20;                   % componente en x
v = -(X.^2)/20;                  % componente en y

% Derivadas analíticas
ux_x = Y/20;                    
ux_y = X/20;                   
vy_x = -(X/10);                 
vy_y = zeros(size(X));           

div_u = ux_x + vy_y; % divergencia = y/20

% Tensor de deformaciones
eps_xx = ux_x;                  
eps_yy = vy_y;                  
eps_xy = 0.5*(ux_y + vy_x);      

% Tensiones 
lambda = 1; mu = 1;
sigma_xx = lambda*div_u + 2*mu*eps_xx;   % = 3y/20
sigma_yy = lambda*div_u + 2*mu*eps_yy;   % = y/20
sigma_zz = lambda*div_u;                 % = y/20 (eps_zz=0)
sigma_xy = 2*mu*eps_xy;                  % = -x/20
sigma_yx = sigma_xy;                     % simétrico
sigma_xz = zeros(size(X));               % no hay variación en z
sigma_zx = zeros(size(X));               % idem

% Tensión tangencial respecto al plano ortogonal a i 

tau_i = abs(sigma_yx);                   % = | -x/20 |

% Gráficos

cm = turbo;
lw = 1.2;

% Mallado y contornos de la placa
figure('Name','Mallado de la placa','Color','w');
mesh(X,Y,zeros(size(X)),'EdgeColor',[0,0.6,0.7],'FaceColor','none'); hold on;
plot(X(1,:),Y(1,:),'k','LineWidth',2);       % borde inferior
plot(X(end,:),Y(end,:),'k','LineWidth',2);    % borde superior
plot(X(:,1),Y(:,1),'k','LineWidth',2);        % borde izquierdo (pared)
plot(X(:,end),Y(:,end),'k','LineWidth',2);    % borde derecho
axis([-0.5 4.5 -0.5 2.5]); axis equal; grid on; box on;
xlabel('x'); ylabel('y'); title('Mallado y contorno de la placa');

% Tensiones normales
figure('Name','Tensiones normales','Color','w');
tiledlayout(1,3,'Padding','compact','TileSpacing','compact');

nexttile;
surf(X,Y,sigma_xx); shading interp; colormap(cm); colorbar;
title('\sigma_{xx}'); xlabel('x'); ylabel('y'); zlabel('\sigma_{xx}');
axis tight; daspect([1 1 0.35]); view(45,30);

nexttile;
surf(X,Y,sigma_yy); shading interp; colormap(cm); colorbar;
title('\sigma_{yy}'); xlabel('x'); ylabel('y'); zlabel('\sigma_{yy}');
axis tight; daspect([1 1 0.35]); view(45,30);

nexttile;
surf(X,Y,sigma_zz); shading interp; colormap(cm); colorbar;
title('\sigma_{zz}'); xlabel('x'); ylabel('y'); zlabel('\sigma_{zz}');
axis tight; daspect([1 1 0.35]); view(45,30);
sgtitle('Distribución de tensiones normales (\lambda=\mu=1)','FontWeight','bold');

% Tensión tangencial respecto a i 
figure('Name','Tension tangencial respecto a i','Color','w');
surf(X,Y,tau_i); shading interp; colormap(cm); colorbar;
title('\tau_{i} = |\sigma_{xy}|'); xlabel('x'); ylabel('y'); zlabel('\tau_{i}');
axis tight; daspect([1 1 0.35]); view(45,30); grid on;


figure('Name','Contornos 2D de tensiones','Color','w');
tiledlayout(1,4,'Padding','compact','TileSpacing','compact');

nexttile; contourf(X,Y,sigma_xx,20,'LineColor','none'); axis equal tight; colorbar;
title('\sigma_{xx}'); xlabel('x'); ylabel('y');

nexttile; contourf(X,Y,sigma_yy,20,'LineColor','none'); axis equal tight; colorbar;
title('\sigma_{yy}'); xlabel('x'); ylabel('y');

nexttile; contourf(X,Y,sigma_zz,20,'LineColor','none'); axis equal tight; colorbar;
title('\sigma_{zz}'); xlabel('x'); ylabel('y');

nexttile; contourf(X,Y,tau_i,20,'LineColor','none'); axis equal tight; colorbar;
title('\tau_{i}'); xlabel('x'); ylabel('y');

11 Tensiones tangenciales

Se representan las tensiones tangenciales respecto al plano ortogonal al eje [math]\vec{i}[/math], definidas como:

[math]|\sigma \cdot \vec{i} - (\vec{i} \cdot \sigma \cdot \vec{i}) , \vec{i}|[/math]

En la gráfica se muestran únicamente las componentes no nulas, permitiendo visualizar las tensiones tangenciales que actúan sobre los planos perpendiculares al eje 𝑥 x.

% Definición del Mallado
h = 0.15; 
x = 0:h:4;
y = 0:h:2;
[X, Y] = meshgrid(x, y);

% Definición de la placa 
Condition = (Y >= X./8) & (Y <= 2 - X./8);
X(~Condition) = NaN;
Y(~Condition) = NaN;

% Cálculo de la Tensión Tangencial 

Sig_xy = -X ./ 20;

% Definición del Vector Tangencial

U_vec = zeros(size(X)); % Componente X es 0
V_vec = Sig_xy; % Componente Y es el valor de la tensión

% Gráfica
figure('Name', 'Apartado 10: Vector Tensión Tangencial', 'Color', 'w');

% Usamos quiver para dibujar las flechitas azules

quiver(X, Y, U_vec, V_vec, 'b', 'LineWidth', 1);

% Configuración de la vista
view(2); 
axis equal; 
grid on; 

% Títulos 
title('Tensiones tangenciales (Plano ortogonal a i)');
xlabel('Eje X');
ylabel('Eje Y');

% Ajuste de límites 
axis([-0.5 4.5 -0.5 2.5]);

12 Tensión de Von Mises

La tensión de Von Mises se define mediante la expresión:

[math] \sigma_{VM} = \sqrt{\frac{(\sigma_1 - \sigma_2)^2 + (\sigma_2 - \sigma_3)^2 + (\sigma_3 - \sigma_1)^2}{2}} [/math]

Esta fórmula combina las diferencias entre las tensiones principales (𝜎1,𝜎2,𝜎3σ1,σ2,σ3) y permite calcular el esfuerzo equivalente que experimenta el material, útil para evaluar la resistencia frente a la plasticidad o el fallo.

Figura 11

clc; clear; close all;

% Malla
h = 0.1;
[x, y] = meshgrid(0:h:4, 0:h:2);   % ojo: x -> columnas, y -> filas
X = x; Y = y;

% Bordes de la placa
f = X./8;
g = 2 - X./8;
EnPlaca = (Y >= f) & (Y <= g);

% Definición del campo de desplazamientos en cartesianas
ux=(1/20).*X.*Y;
uy=-(1/20).*X.^2;

% Derivadas analíticas:
% ux_x=d(u_x)/dx=y/20
% ux_y=d(u_x)/dy=x/20
% uy_x=d(u_y)/dx=-x/10
% uy_y=d(u_y)/dy 0
ux_x=Y./20;
ux_y=X./20;
uy_x=-X./10;
uy_y=zeros(size(X));

% Divergencia (útil para sigma)
divU=ux_x+uy_y;

% Tensor de deformaciones simétrico 
% epsilon = 1/2(grad u + grad u^T)
% epsilon_xx = ux_x
% epsilon_yy = uy_y
% epsilon_xy = 1/2*(ux_y + uy_x)
eps_xx = ux_x;
eps_yy = uy_y;
eps_xy = 0.5*(ux_y + uy_x);

% Para el caso 3D consideramos eps_zz = 0
eps_zz = zeros(size(X));

% Coeficientes de Lamé
lambda=1;
mu=1;

% Construcción del tensor de tensiones 
%sigma = lambda*div(u)*I + 2*mu*eps
sigma_xx = lambda.*divU + 2*mu.*eps_xx;
sigma_yy = lambda.*divU + 2*mu.*eps_yy;
sigma_zz = lambda.*divU + 2*mu.*eps_zz;
sigma_xy = 2*mu.*eps_xy;    % = sigma_yx

% Von Mises y matriz para autovalores
VM = NaN(size(X));
S1 = NaN(size(X)); S2 = NaN(size(X)); S3 = NaN(size(X));

% Loop sobre puntos dentro de la placa para calcular autovalores
[idx_i,idx_j]=find(EnPlaca);
npts=numel(idx_i);

for k = 1:npts
    i = idx_i(k); j = idx_j(k);   % i: fila (y index), j: columna (x index)
    % construir sigma 3x3 en ese punto
    S = [sigma_xx(i,j), sigma_xy(i,j), 0;sigma_xy(i,j), sigma_yy(i,j),0;0,0,sigma_zz(i,j)];
    e = eig(S);          
    % ordenar autovalores (no estrictamente necesario, pero ordeno descendente)
    e = sort(e,'descend');
    s1 = e(1); s2 = e(2); s3 = e(3);
    % tensión de Von Mises (fórmula)
    sigmaVM = sqrt( ((s1-s2).^2 + (s2-s3).^2 + (s3-s1).^2) / 2 );
    VM(i,j) = sigmaVM;
    S1(i,j) = s1; S2(i,j) = s2; S3(i,j) = s3;
end

VM(~EnPlaca) = NaN;
% Buscar máximo de Von Mises dentro de la placa
[maxVM, indMax] = max(VM(:));
[rowMax, colMax] = ind2sub(size(VM), indMax);
xMax = X(rowMax, colMax);
yMax = Y(rowMax, colMax);

% Visualización
figure('Name','Tensión de Von Mises','Color','w'); clf; hold on;
axis equal; axis([-0.5 4.5 -0.5 2.5]); box on;
title('Tensión de Von Mises \sigma_{VM}');
xlabel('x'); ylabel('y');

% Contourf del campo VM
numLevels = 30;
[C, hC] = contourf(X, Y, VM, numLevels, 'LineStyle','none');
colormap(jet);
c = colorbar; ylabel(c, '\sigma_{VM}');
% Bordes de la placa
xb = linspace(0,4,200);
plot(xb, xb./8, 'k-', 'LineWidth', 1.5);         
plot(xb, 2 - xb./8, 'k-', 'LineWidth', 1.5);    
plot([0 0], [0 2], 'k-', 'LineWidth', 1.5);  
plot([4 4], [0.5 1.5], 'k-', 'LineWidth', 1.5); 

% Punto de máxima tensión de Von Mises
plot(xMax, yMax, 'rp', 'MarkerSize', 14, 'MarkerFaceColor', 'r');
text(xMax+0.08, yMax+0.08, sprintf('Máx VM = %.3g', maxVM), 'FontWeight','bold','Color','r');

% Mejora estética
set(gcf,'Color','w');
shading interp;

hold off;

% Información
fprintf('Máximo Von Mises = %.6g en (x,y) = (%.4g, %.4g)\n', maxVM, xMax, yMax);

13 Campo de Fuerzas

El campo de fuerzas [math]\vec{F}[/math] que actúa sobre la placa, y que es responsable del desplazamiento observado, se puede aproximar mediante la ecuación de la elasticidad lineal:

[math]\vec{F} = - \operatorname{div} \sigma[/math]

donde [math]\operatorname{div} \sigma[/math] es la divergencia de un campo tensorial. En coordenadas cartesianas, se expresa como:

[math]\operatorname{div} \sigma = (\operatorname{div} \sigma)_1 , \vec{i} + (\operatorname{div} \sigma)_2 , \vec{j} + (\operatorname{div} \sigma)3 , \vec{k}, \quad (\operatorname{div} \sigma)i = \sum{j=1}^{3} \frac{\partial \sigma{ij}}{\partial x_j}[/math]

Esta formulación permite calcular, a partir del tensor de tensiones, el campo de fuerzas internas que genera los desplazamientos de la placa.

%% Apartado 12: Campo de Fuerzas Volumétricas (F = -div Sigma)
clc; clear; close all;

% 1. Definición del Mallado
h = 0.1;
x = 0:h:4;
y = 0:h:2;
[X, Y] = meshgrid(x, y);

% Definición de la placa (Trapecio)
Condition = (Y >= X./8) & (Y <= 2 - X./8);
X(~Condition) = NaN;
Y(~Condition) = NaN;

lambda = 1; mu = 1;

% 2. Recalculamos las Tensiones (Necesarias para derivar)
% Usamos las fórmulas analíticas que ya conocemos
Div_u = Y ./ 20;
Eps_xx = Y ./ 20;
Eps_yy = zeros(size(X));

Sig_xx = lambda .* Div_u + 2 .* mu .* Eps_xx; % = 3y/20
Sig_yy = lambda .* Div_u + 2 .* mu .* Eps_yy; % = y/20
Sig_xy = -X ./ 20; % = -x/20 (cizalladura)

% 3. Cálculo de la Divergencia del Tensor (Numéricamente)
% F = -div(sigma)
% Fx = -( d(Sxx)/dx + d(Sxy)/dy )
% Fy = -( d(Sxy)/dx + d(Syy)/dy )

% La función gradient calcula derivadas centrales.
% Sintaxis: gradient(Matriz, paso_x, paso_y)
[dSxx_dx, dSxx_dy] = gradient(Sig_xx, h, h);
[dSxy_dx, dSxy_dy] = gradient(Sig_xy, h, h);
[dSyy_dx, dSyy_dy] = gradient(Sig_yy, h, h);

Fx = - (dSxx_dx + dSxy_dy);
Fy = - (dSxy_dx + dSyy_dy);

% 4. Comprobación de Magnitud (Para ver si es cero)
MagnitudF = sqrt(Fx.^2 + Fy.^2);
max_F = max(MagnitudF(:));

fprintf('--- Resultado Apartado 12 ---\n');
fprintf('La fuerza máxima calculada es: %e\n', max_F);

if max_F < 1e-10
fprintf('CONCLUSIÓN: La fuerza es prácticamente NULA (Equilibrio).\n');
subtitle_text = 'Fuerza \approx 0 (Equilibrio Estático)';
else
subtitle_text = 'Existen fuerzas volumétricas';
end

% 5. Gráfica
figure('Name', 'Apartado 12: Campo de Fuerzas', 'Color', 'w');
hold on;

% Dibujamos los vectores.
% Como son casi cero, MATLAB pintará puntos. Para ver "algo",
% a veces se normalizan, pero lo correcto es mostrar que no hay fuerza.
quiver(X, Y, Fx, Fy, 'r');

% Dibujar bordes
plot([0 4 4 0 0], [0 0.5 1.5 2 0], 'k--', 'LineWidth', 0.5); % Caja aproximada

title('Campo de Fuerzas Volumétricas');
subtitle(subtitle_text);
xlabel('Eje X'); ylabel('Eje Y');
axis equal; axis([-0.5 4.5 -0.5 2.5]);
grid on;

% Añadimos una nota en la gráfica explicándolo
text(1, 1, {'Fuerza Nula', '(puntos indican magnitud ~0)'}, ...
'Color', 'r', 'FontWeight', 'bold', 'HorizontalAlignment', 'center');

hold off;

14 Densidad

Figura 10

clc; clear; close all;
% Mallado
h=0.05; 
[X,Y]=meshgrid(0:h:4,0:h:2);

% Forma de la placa
Condition=(Y>=X./8)&(Y<=2-X./8);

% Definición de la Densidad : d(x,y)=(4-x)*|y|
Densidad=(4-X).*abs(Y);

% Filtrado para la Gráfica
Densidad(~Condition)=NaN;

% Visualización 3D
figure('Name','Densidad', 'Color','w');

% Superficie 3D
surf(X,Y,Densidad);

% Configuración
shading interp; % Elimina las líneas negras para que se vea suave
colormap(jet);  % Mapa de colores (De azul a rojo)
colorbar;       % Barra de referencia
caxis([0 8]);   % Fija la escala de colores

% Etiquetas y Título
title({'Densidad'},'FontSize',13);
xlabel('Eje X','FontWeight','bold');
ylabel('Eje Y','FontWeight','bold');
zlabel('Densidad (kg/m^2)','FontWeight','bold');

% Ajustes
axis tight;      % Ajusta los ejes al contenido
view(-45, 30);   % Ángulo de vista 
grid on;
box on;

% Proyección de sombra en el suelo
hold on;
% Dibujamos el contorno plano abajo del todo (en z=0) para referencia
contour(X,Y,Densidad,20,'z','Offset',0); 
hold off;

15 Ejemplos de uso en la ingeniería

Las vigas en voladizo destacan por su forma y por la versatilidad que ofrecen, lo que las convierte en un elemento muy empleado en distintos ámbitos. Algunos de sus usos más frecuentes son:

• Arquitectura: se utilizan en balcones, terrazas, marquesinas o cubiertas de acceso. • Maquinaria: forman parte de brazos de grúas y de las palas de los aerogeneradores. • Ingeniería: se aplican en la construcción de puentes en voladizo y en diferentes componentes estructurales.