Onda transversal plana (Grupo 54)

Trabajo realizado por estudiantes
Título	Derformación plana. Grupo
Asignatura	Teoría de Campos
Curso	2025-26
Autores	Jorge Muñoz Jimenez Daniel Galarza Polo Armando de Tomás
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

En este documento representaremos la sección longitudinal de un arco comprendido entre los radios [math]1\lt\sqrt{x^2+y^2}\lt2[/math] Supondremos que están definidas dos cantidades físicas.

• La Temperatura
• Los Desplazamientos

La temperatura [math]T(x, y)[/math] viene dada por la ecuación:

[math]T(x, y)=(x-y)^2[/math]

Los desplazamientos [math]u(x, y)[/math] producidos por la acción de una fuerza determinada.

Al definir el vector de posición de los puntos de la placa antes de que se produzca cualquier deformación [math]\vec{r_{0}}(x, y)=x\vec{i}+y\vec{j}[/math] , la posición de cada punto [math](x, y)[/math] de la placa después de la deformación vendrá dada mediante la ecuación [math]\vec{r}(x, y)=\vec{r_{0}}(x, y)+\vec{u}(x, y,t)[/math]

Conocemos [math]\vec{u}(\rho, \theta )=\frac{1}{5}(\rho - 1) \rho\vec{e_{\rho}}[/math]

1 Mallado

A continuación podemos ver la representación del mallado que contiene a los puntos interiores del sólido realizado con MatLab.

Hacemos uso de un paso de muestreo, el intervalo entre punto y punto es [math]h=\frac{1}{10}[/math] para las variables [math]x[/math] e [math]y[/math].

% Radios del arco
r1 = 1;
r2 = 2;
%Divisores
Nr = 10;
Nt = 40;
%Crear vectores
r = linspace (r1, r2, Nr+1);
theta = linspace(0, pi, Nt+1);
%Mallado
figure;
hold on;
axis equal; 
title('Mallado');
xlabel('x');
ylabel('y');
%Lineas radiales
for i = 1:length(r)
     x = r(i) * cos(theta);
     y = r(i) * sin(theta);
     plot(x, y, 'k'); 
end 
for j = 1:length(r)
     x = r(i) * cos(theta);
     y = r(i) * sin(theta);
     plot(x, y, 'g'); 
end 

grid on;

2 Temperatura

La temperatura viene dada por la siguiente expresión [math] T(x,y) = (x - y)^2 [/math], que depende únicamente de x e y.

r1 = 1;           
r2 = 2; 
% divisiones radiales        
Nr = 10; 
% divisiones angulares         
Nt = 60;          
r = linspace(r1, r2, Nr+1);
theta = linspace(0, pi, Nt+1);
[R, TH] = meshgrid(r, theta);
% Convertir a coordenadas cartesianas
X = R .* cos(TH);
Y = R .* sin(TH);
%Calcular la temperatura T = (x - y)^2
T = (X - Y).^2;
%Grafica de la temperatura
figure;
pcolor(X, Y, T);
shading interp;
colorbar;
axis equal;
title('Temperatura T = (x - y)^2 en el semiarco');
xlabel('x');
ylabel('y');

3 [math]\nabla T[/math]

Para poder calcular el Gradiente de T, primero debemos obtener las derivadas parciales. Derivadas parciales:

[math] \frac{\partial T}{\partial x} = T_x = 2(x - y), \qquad \frac{\partial T}{\partial y} = T_y = -2(x - y). [/math]

Por tanto el gradiente es

[math] \nabla T(x,y) = (T_x, T_y) = \bigl(2(x - y),\, -2(x - y)\bigr) = 2(x - y)\,(1, -1). [/math]

Tras realizar los cálculos obtenemos que [math]\nabla T[/math] = 0 , siendo este un punto crítico. El Gradiente de T podemos confirmar que es ortogonal a la dirección de las curvas de nivel en todos los puntos.

Como el producto escalar es cero, [math]\nabla T[/math] es ortogonal a la dirección de las curvas de nivel en todos los puntos.

Calcularemos [math]\nabla T[/math] y se dibujará como un campo vectorial, representado con Matlab:

r1 = 1;          
r2 = 2;           
Nr = 10;           
Nt = 60;        
r = linspace(r1, r2, Nr+1);
theta = linspace(0, pi, Nt+1);
[R, TH] = meshgrid(r, theta);
%Coordenadas cartesianas
X = R .* cos(TH);
Y = R .* sin(TH);

T = (X - Y).^2;
figure; 
pcolor(X, Y, T);
shading interp;
colorbar;
axis equal;
title('Gradiente');
xlabel('x');
ylabel('y');

%Gradiente en coord enadas cartesianas
dTdx = 2 * (X - Y);     
dTdy = -2 * (X - Y);    

hold on;
% Curvas de nive
contour(X, Y, T, 15, 'k', 'LineWidth', 0.8);
quiver(X(1:step:end,1:step:end), Y(1:step:end,1:step:end), ...
       dTdx(1:step:end,1:step:end), dTdy(1:step:end,1:step:end), ...
       0.6, 'r', 'LineWidth', 1);

legend({'Mapa T (pcolor)','Curvas de nivel','Gradiente \nabla T'}, 'Location','best');
hold off;

4 Campo de vectores

Consideramos ahora el campo de vectores [math]\vec{u}(\rho, \theta )=\frac{1}{5}(\rho - 1) \rho\vec{e_{\rho}}[/math] , el nuevo mallado del sólido será:

r1 = 1;
r2 = 2;
Nr = 10;
Nt = 40;
r = linspace(r1, r2, Nr+1);
theta = linspace(0, pi, Nt+1);

figure;
hold on;
axis equal;
title('Campo de vectores');
xlabel('x');
ylabel('y');

% Dibujar líneas radiales
for i = 1:length(r)
    x = r(i) * cos(theta);
    y = r(i) * sin(theta);
    plot(x, y, 'k');
end

for j = 1:length(theta)
    x = r .* cos(theta(j));
    y = r .* sin(theta(j));
    plot(x, y, 'g'); 
end

xlim([-3 3]);
ylim([-1 2.5]);
grid on;

% Mallado
[TH, R] = meshgrid(theta, r);
X = R' .* cos(TH');
Y = R' .* sin(TH');
u_r = (1/5) * (R' - 1) .* R';
u_x = u_r .* cos(TH');
u_y = u_r .* sin(TH');

step = 2;
Xq = X(1:step:end, 1:step:end);
Yq = Y(1:step:end, 1:step:end);
Uq = u_x(1:step:end, 1:step:end);
Vq = u_y(1:step:end, 1:step:end);

quiver(Xq, Yq, Uq, Vq, 'r');

hold off;

5 Desplazamiento antes y después

Si [math]\vec{u}[/math] determina el desplazamiento que sufre cada punto del sólido, dibujar el sólido antes y después del desplazamiento.

Cada punto del sólido inicialmente está en:

[math] \vec{r}_0 = (x,y) = (\rho \cos\theta,\; \rho \sin\theta). [/math]

Tras deformarse, su nueva posición es:

[math] \vec{r} = \vec{r}_0 + \vec{u}. [/math]

El desplazamiento es:

[math] \vec{r} = (\rho + u_\rho)\,\vec{e}_\rho, [/math]

Por tanto, el nuevo radio es:

[math] \rho_{\text{nuevo}} = \rho + \tfrac{1}{5}(\rho - 1)\rho. [/math]

h = 0.1;
r = 1:h:2;
npuntos = round(pi/h)+1;
ang = linspace(0,pi,npuntos);
[rho,theta] = meshgrid(r,ang);
x = rho.*cos(theta);
y = rho.*sin(theta);
desprad = (1/5).*(rho-1).*rho;

% Proyección del desplazamiento en cartesianas.
despx = desprad.*cos(theta);
despy = desprad.*sin(theta);

X = x+despx;
Y = y+despy;
figure('Color','w');
limitesejes=[-3 3 -1 3];

subplot(1,2,1)
mesh(x,y,0*x); 
view(2)
axis equal
axis(limitesejes)
set(gca,'DataAspectRatio',[1 1 1]);
title('Antes de la Deformación');
set(findobj(gca,'Type','Surface'),'EdgeColor','k');

subplot(1,2,2)
mesh(X,Y,0*X);
view(2)
axis equal
axis(limitesejes)
set(gca,'DataAspectRatio',[1 1 1]);
title('Después de la Deformación');
set(findobj(gca,'Type','Surface'),'EdgeColor','b');

6 Divergencia

Calcularemos la [math]\nabla \cdot \vec{u}[/math] en todos los puntos del sólido.

Fórmula de la divergencia en cilíndricas:

Para un campo [math]\vec{u} = u_r(r,\varphi,z)\,\hat{e}_r + u_\varphi(r,\varphi,z)\,\hat{e}_\varphi + u_z(r,\varphi,z)\,\hat{e}_z[/math], la divergencia viene dada por

[math] \nabla \cdot \vec{u} = \frac{1}{r}\frac{\partial}{\partial r}\!\left(r\,u_r\right) + \frac{1}{r}\frac{\partial u_\varphi}{\partial \varphi} + \frac{\partial u_z}{\partial z}. [/math]

En nuestro caso [math]u_\varphi = 0[/math] y [math]u_z = 0[/math], por tanto

[math] \nabla \cdot \vec{u} = \frac{1}{r}\frac{d}{dr}\!\left(r\,u_r(r)\right). [/math]

Sustitución de [math]u_r[/math]

Tenemos [math]u_r(r) = \tfrac{1}{5}(r - 1)r = \tfrac{1}{5}(r^2 - r)[/math].

Calculemos [math]r\,u_r[/math]:

[math] r\,u_r = r \cdot \tfrac{1}{5}(r^2 - r) = \tfrac{1}{5}(r^3 - r^2). [/math]

Derivando respecto a [math]r[/math]:

[math] \frac{d}{dr}(r\,u_r) = \tfrac{1}{5}(3r^2 - 2r). [/math]

Dividiendo por [math]r[/math]:

[math] \nabla \cdot \vec{u} = \frac{1}{r} \cdot \tfrac{1}{5}(3r^2 - 2r) = \tfrac{1}{5}(3r - 2). [/math]

Así que la expresión cerrada de la divergencia es

[math] \nabla \cdot \vec{u}(r) = \tfrac{1}{5}(3r - 2), [/math]

válida para todo punto del dominio (independiente de [math]\varphi[/math] y [math]z[/math]).

Análisis en el dominio [math]r \in [1,2][/math]

• Es una función lineal creciente en [math]r[/math].
• En [math]r = 1[/math]:

[math] \nabla \cdot \vec{u}(1) = \tfrac{1}{5}(3 - 2) = \tfrac{1}{5} = 0.2. [/math]

• En [math]r = 2[/math]:

[math] \nabla \cdot \vec{u}(2) = \tfrac{1}{5}(6 - 2) = \tfrac{4}{5} = 0.8. [/math]

Por tanto, los puntos con mayor divergencia son los del borde exterior [math]r = 2[/math] .

r1 = 1;
r2 = 2;
Nr = 10;      
Nt = 40;     
% Crear vectores
r = linspace(r1, r2, Nr+1);
theta = linspace(0, pi, Nt+1);
% Creamos malla R,TH y convertimos a X,Y (mismo formato que en tus scripts anteriores)
[TH, R] = meshgrid(theta, r);  % TH,R de tamaño (Nr+1) x (Nt+1)
X = R' .* cos(TH');            % X,Y tamaño (Nt+1) x (Nr+1)
Y = R' .* sin(TH');
% Calcular rho en cada punto (alternativa segura)
rho = sqrt(X.^2 + Y.^2);
% Divergencia 
div_u = (1/5) * (3 * rho - 2);
figure;
pcolor(X, Y, div_u);   
shading interp;
colorbar;
axis equal;
title('Divergencia');
xlabel('x'); ylabel('y');

hold on;
plot(X, Y, 'k:', 'LineWidth', 0.5);    
plot(X', Y', 'g:', 'LineWidth', 0.5);  
hold off;

Los puntos con mayor divergencia, se pueden observar en el limite superior del arco de radio dos. Esta representados con círculos rojos a lo largo de toda la figura.

.

7 Rotacional

Calcular [math]\lvert \nabla \times \vec{u} \rvert[/math] en todos los puntos del sólido. ¿Qué puntos sufren un mayor rotacional?

Rotacional en coordenadas cilíndricas

[math] \nabla \times \vec{u} = \left( \frac{1}{r}\frac{\partial u_z}{\partial \varphi} - \frac{\partial u_\varphi}{\partial z} \right)\hat{e}_r \;+\; \left( \frac{\partial u_r}{\partial z} - \frac{\partial u_z}{\partial r} \right)\hat{e}_\varphi \;+\; \left( \frac{1}{r}\frac{\partial}{\partial r}(r\,u_\varphi) - \frac{1}{r}\frac{\partial u_r}{\partial \varphi} \right)\hat{e}_z. [/math]

Sustitución de nuestro campo

En nuestro caso [math]u_\varphi \equiv 0[/math], [math]u_z \equiv 0[/math], y [math]u_r[/math] depende sólo de [math]r[/math] (no depende de [math]\varphi[/math] ni de [math]z[/math]). Por tanto:

• [math](\nabla \times \vec{u})_r = \frac{1}{r}\frac{\partial 0}{\partial \varphi} - \frac{\partial 0}{\partial z} = 0.[/math]

• [math](\nabla \times \vec{u})_\varphi = \frac{\partial u_r}{\partial z} - \frac{\partial 0}{\partial r} = 0 - 0 = 0.[/math] (porque [math]u_r[/math] no depende de [math]z[/math])

• [math](\nabla \times \vec{u})_z = \frac{1}{r}\left(\frac{\partial}{\partial r}(r \cdot 0) - \frac{\partial u_r}{\partial \varphi}\right) = \frac{1}{r}(0 - 0) = 0.[/math] (porque [math]u_\varphi = 0[/math] y [math]u_r[/math] no depende de [math]\varphi[/math])

Por tanto:

[math] \nabla \times \vec{u}(r,\varphi) = \vec{0} [/math]

en todo el dominio.

Módulo del rotacional y puntos con mayor rotacional

El módulo es

[math]\lvert \nabla \times \vec{u} \rvert = 0[/math]

en todos los puntos. No hay puntos con mayor rotacion

figure;
pcolor(X, Y, rot_z); 
shading interp;

colormap(parula);          
caxis([-0.1 0.1]);           
colorbar;
axis equal;
title('Rotacional (∇×u)_z — Campo irrotacional');
xlabel('x'); ylabel('y');

8 Tensiones normales en la dirección [math]\vec{e}_\rho \; y \; \frac{1}{\rho}\,\vec{e}_\theta[/math]

Divergencia en coordenadas cilíndricas:

[math] \operatorname{div}\vec{u} = \frac{1}{\rho}\frac{\partial}{\partial \rho}(\rho u_\rho) + \frac{1}{\rho}\frac{\partial u_\theta}{\partial \theta}. [/math]

Como [math]u_\theta = 0[/math] y [math]u_\rho[/math] depende sólo de [math]\rho[/math],

[math] \operatorname{div}\vec{u} = \frac{1}{\rho}\frac{d}{d\rho}(\rho u_\rho). [/math]

Con [math]u_\rho = \tfrac{1}{5}(\rho^2 - \rho)[/math] se obtiene

[math] \rho u_\rho = \tfrac{1}{5}(\rho^3 - \rho^2), \qquad \frac{d}{d\rho}(\rho u_\rho) = \tfrac{1}{5}(3\rho^2 - 2\rho), [/math]

por tanto:

[math] \operatorname{div}\vec{u} = \tfrac{1}{5}(3\rho - 2). [/math]

Componentes del tensor de deformación:

• [math]\varepsilon_{\rho\rho} = \dfrac{\partial u_\rho}{\partial \rho}[/math].
• [math]\varepsilon_{\theta\theta} = \dfrac{u_\rho}{\rho}[/math] (porque [math]u_\theta = 0[/math]).
• [math]\varepsilon_{\rho\theta} = \varepsilon_{\theta\rho} = \tfrac{1}{2}\left( \frac{1}{\rho}\frac{\partial u_\rho}{\partial \theta} + \frac{\partial u_\theta}{\partial \rho} - \frac{u_\theta}{\rho} \right) = 0[/math]

Nulo por no haber dependencia en [math]\theta[/math] y [math]u_\theta = 0[/math]

Calculamos:

[math] \varepsilon_{\rho\rho} = \frac{d}{d\rho}\left[\tfrac{1}{5}(\rho^2 - \rho)\right] = \tfrac{1}{5}(2\rho - 1), [/math]

[math] \varepsilon_{\theta\theta} = \frac{u_\rho}{\rho} = \tfrac{1}{5}(\rho - 1). [/math]

Componentes del tensor de tensiones:

[math] \sigma_{\rho\rho} = (\operatorname{div}\vec{u}) + 2\varepsilon_{\rho\rho}, \qquad \sigma_{\theta\theta} = (\operatorname{div}\vec{u}) + 2\varepsilon_{\theta\theta}. [/math]

Sustituyendo los valores:

[math] \sigma_{\rho\rho} = \tfrac{1}{5}(3\rho - 2) + 2\cdot\tfrac{1}{5}(2\rho - 1) = \tfrac{1}{5}(3\rho - 2 + 4\rho - 2) = \tfrac{1}{5}(7\rho - 4). [/math]

[math] \sigma_{\theta\theta} = \tfrac{1}{5}(3\rho - 2) + 2\cdot\tfrac{1}{5}(\rho - 1) = \tfrac{1}{5}(3\rho - 2 + 2\rho - 2) = \tfrac{1}{5}(5\rho - 4). [/math]

Por tanto las tensiones son:

• Tensión normal en dirección [math]\vec{e}_\rho[/math]:

[math] \vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho = \sigma_{\rho\rho} = \tfrac{1}{5}(7\rho - 4). [/math]

• Tensión normal en dirección [math]\tfrac{1}{\rho}\vec{e}_\theta[/math]:

[math] \left(\frac{1}{\rho}\vec{e}_\theta\right)\cdot \sigma \cdot \left(\frac{1}{\rho}\vec{e}_\theta\right) = \frac{1}{\rho^2}\,\sigma_{\theta\theta} = \frac{1}{\rho^2}\cdot \frac{1}{5}(5\rho - 4) = \frac{5\rho - 4}{5\rho^2}. [/math]

r = linspace(1,2,11);
theta = linspace(0,pi,41);
[R,TH] = meshgrid(r,theta);
R = R'; TH = TH';
X = R .* cos(TH);
Y = R .* sin(TH);
u_r = (1/5)*(R.^2 - R);

eps_rr = (1/5)*(2*R - 1);
eps_tt = u_r ./ R;


div_u   = eps_rr + eps_tt;
sigma_rr = div_u + 2*eps_rr;
sigma_tt = div_u + 2*eps_tt;

T_er  = sigma_rr;
T_eth = sigma_tt ./ (R.^2);
[max_er, i1]  = max(T_er(:));
[max_eth,i2]  = max(T_eth(:));
x1 = X(i1); y1 = Y(i1);
x2 = X(i2); y2 = Y(i2);

figure;

subplot(1,3,1); hold on; axis equal; grid on;
for i = 1:length(r), plot(r(i)*cos(theta), r(i)*sin(theta),'k'); end
for j = 1:length(theta), plot(r*cos(theta(j)), r*sin(theta(j)),'g'); end
title('Mallado'); xlim([-3 3]); ylim([-1 2.5]);

subplot(1,3,2);
scatter(X(:),Y(:),80,T_er(:),'filled'); axis equal; colorbar; grid on;
hold on; plot(x1,y1,'kx','LineWidth',2);
title('\sigma_{rr}');

subplot(1,3,3);
scatter(X(:),Y(:),80,T_eth(:),'filled'); axis equal; colorbar; grid on;
hold on; plot(x2,y2,'kx','LineWidth',2);
title('\sigma_{\theta\theta}/\rho^2');

fprintf('Max sigma_rr = %.5f en (x,y) = (%.4f, %.4f)\n', max_er, x1, y1);
fprintf('Max sigma_tt/rho^2 = %.5f en (x,y) = (%.4f, %.4f)\n', max_eth, x2, y2);

9 Tensiones tangenciales respecto al plano ortogonal a [math]\vec{e}_\rho[/math]

Dibujar las tensiones tangenciales respecto al plano ortogonal a [math]\vec{e}_\rho[/math], es decir

[math] \left|\, \sigma \cdot \vec{e}_\rho \;-\; (\vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho)\,\vec{e}_\rho \,\right| [/math].

¿Dónde son mayores? Compáralas con los puntos de mayor deformación de la malla.

La cantidad pedida es la magnitud del componente tangencial del vector tensión sobre la superficie cuya normal es [math]\vec{e}_\rho[/math].

El vector sobre la superficie normal [math]\vec{n} = \vec{e}_\rho[/math] es

[math] \vec{t} = \sigma \cdot \vec{e}_\rho . [/math]

Su componente normal es [math](\vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho)\,\vec{e}_\rho[/math]. Por tanto, el componente tangencial es

[math] \vec{t}_{\text{tan}} = \sigma \cdot \vec{e}_\rho \;-\; (\vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho)\,\vec{e}_\rho , [/math]

y la magnitud pedida es [math]\lVert \vec{t}_{\text{tan}} \rVert[/math].

Cálculo en coordenadas cilíndricas

En la resolución anterior obtuvimos el tensor de tensiones en la base polar ([math]\vec{e}_\rho[/math], [math]\vec{e}_\theta[/math]) :

[math] \sigma_{ij} = \begin{pmatrix} \sigma_{\rho\rho} & \sigma_{\rho\theta} \\ \sigma_{\theta\rho} & \sigma_{\theta\theta} \end{pmatrix} = \begin{pmatrix} \dfrac{7\rho - 4}{5} & 0 \\ 0 & \dfrac{5\rho - 4}{5} \end{pmatrix}, [/math]

es decir [math]\sigma_{\rho\theta} = \sigma_{\theta\rho} = 0[/math].

Entonces

[math] \sigma \cdot \vec{e}_\rho = \sigma_{\rho\rho}\,\vec{e}_\rho + \sigma_{\theta\rho}\,\vec{e}_\theta = \sigma_{\rho\rho}\,\vec{e}_\rho + 0 \cdot \vec{e}_\theta, [/math]

y

[math] (\vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho)\,\vec{e}_\rho = \sigma_{\rho\rho}\,\vec{e}_\rho. [/math]

Por tanto el vector tangencial es

[math] \vec{t}_{\text{tan}} = \sigma \cdot \vec{e}_\rho - \sigma_{\rho\rho}\,\vec{e}_\rho = \vec{0}, [/math]

y su magnitud es

[math] \left\|\, \sigma \cdot \vec{e}_\rho - (\vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho)\,\vec{e}_\rho \,\right\| = 0 \quad \text{en todo el dominio}. [/math]

Conclusión: las tensiones tangenciales sobre las superficies normales a [math]\vec{e}_\rho[/math] son nulas en todos los puntos del dominio (para el campo [math]\vec{u}[/math] dado).

¿Dónde son mayores?

Dado que la magnitud del componente tangencial es exactamente cero en todo el dominio, no hay puntos donde sean mayores: son nulas en todo el arco.

Por consiguiente, no hay correlación con los puntos de mayor deformación de la malla — las mayores deformaciones radiales existen en [math]\rho = 2[/math] (o crecen con [math]\rho[/math]), pero las tensiones tangenciales pedidas son cero en todas partes.

Comparación

• La magnitud de la tensión tangencial pedida,

[math] \left\| \sigma \cdot \vec{e}_\rho - (\vec{e}_\rho \cdot \sigma \cdot \vec{e}_\rho)\,\vec{e}_\rho \right\| [/math]

resulta exactamente cero en todo el dominio.

⇒ No hay puntos donde las tensiones tangenciales sean mayores; son nulas en todas partes.

• Las deformaciones (componentes del tensor [math]\varepsilon[/math]) son no nulas y aumentan con [math]\rho[/math]. En particular,

[math] \varepsilon_{\rho\rho} = \frac{2\rho - 1}{5}, \qquad \varepsilon_{\theta\theta} = \frac{\rho - 1}{5}. [/math]

Las mayores deformaciones ocurren en el borde exterior [math]\rho = 2[/math]

r = linspace(1,2,40);
theta = linspace(0,pi,80);

[R,TH] = meshgrid(r,theta);
R = R'; TH = TH';
X = R .* cos(TH);
Y = R .* sin(TH);

eps_rr = (1/5)*(2*R - 1);   
eps_tt = (1/5)*(R - 1);   
div_u = eps_rr + eps_tt;


sigma_rr = div_u + 2*eps_rr;      
sigma_tt = div_u + 2*eps_tt;   
sigma_rt = zeros(size(R));       

T_er = sigma_rr;                   
T_eth = (1./R.^2) .* sigma_tt;     
T_tang = zeros(size(R));       
figure;

subplot(1,3,1);
scatter(X(:),Y(:),40,T_er(:),'filled');
title('\sigma_{rr}');
axis equal; colorbar;
xlabel('x'); ylabel('y');


subplot(1,3,2);
scatter(X(:),Y(:),40,T_eth(:),'filled');
title('(1/\rho e_\theta)^T \sigma (1/\rho e_\theta)');
axis equal; colorbar;
xlabel('x'); ylabel('y');

subplot(1,3,3);
scatter(X(:),Y(:),40,T_tang(:),'filled');
title('Tensión tangencial respecto al plano normal a e_\rho');
axis equal; 
colorbar;
xlabel('x'); 
ylabel('y');

10 Tensiones tangenciales respecto al plano ortogonal a [math]\frac{1}{\rho}\,\vec{e}_\theta[/math]

Dibujar las tensiones tangenciales respecto al plano ortogonal a [math]\tfrac{1}{\rho}\,\vec{e}_\theta[/math], es decir

[math] \left| \,\sigma \cdot \left(\tfrac{1}{\rho}\,\vec{e}_\theta\right) \;-\; \left( \left(\tfrac{1}{\rho}\,\vec{e}_\theta\right) \cdot \sigma \cdot \left(\tfrac{1}{\rho}\,\vec{e}_\theta\right) \right) \left(\tfrac{1}{\rho}\,\vec{e}_\theta\right) \right| [/math].

¿Dónde son mayores? Compáralas con los puntos de mayor deformación de la malla.

Del apartado anterior tenemos, en la base ([math]\vec{e}_\rho[/math], [math]\vec{e}_\theta[/math]),

[math] \sigma = \begin{pmatrix} \sigma_{\rho\rho} & 0 \\ 0 & \sigma_{\theta\theta} \end{pmatrix}, \qquad \sigma_{\rho\rho} = \frac{7\rho - 4}{5}, \qquad \sigma_{\theta\theta} = \frac{5\rho - 4}{5}. [/math]

El vector unitario tangencial escalado que aparece en el enunciado es

[math] \vec{m} = \frac{1}{\rho}\,\vec{e}_\theta . [/math]

Cálculo del vector tangencial

Primero calculamos [math]\sigma \cdot \vec{m}[/math] , [math]\vec{m} = (0,\;1/\rho)[/math], por lo que

[math] \sigma \cdot \vec{m} = \sigma_{\rho\rho} \cdot 0 \,\vec{e}_\rho + \sigma_{\theta\theta} \cdot \frac{1}{\rho}\,\vec{e}_\theta = \frac{\sigma_{\theta\theta}}{\rho}\,\vec{e}_\theta. [/math]

El escalar [math]\vec{m}\cdot\sigma\cdot\vec{m}[/math] es

[math] \vec{m}\cdot\sigma\cdot\vec{m} = \frac{1}{\rho}\vec{e}_\theta \cdot \left(\frac{\sigma_{\theta\theta}}{\rho}\vec{e}_\theta\right) = \frac{\sigma_{\theta\theta}}{\rho^2}. [/math]

Por tanto

[math] (\vec{m}\cdot\sigma\cdot\vec{m})\,\vec{m} = \frac{\sigma_{\theta\theta}}{\rho^2} \cdot \frac{1}{\rho}\,\vec{e}_\theta = \frac{\sigma_{\theta\theta}}{\rho^3}\,\vec{e}_\theta. [/math]

La diferencia vectorial es:

[math] \sigma \cdot \vec{m} - (\vec{m}\cdot\sigma\cdot\vec{m})\,\vec{m} = \frac{\sigma_{\theta\theta}}{\rho}\vec{e}_\theta - \frac{\sigma_{\theta\theta}}{\rho^3}\vec{e}_\theta = \sigma_{\theta\theta}\left(\frac{1}{\rho} - \frac{1}{\rho^3}\right)\vec{e}_\theta . [/math]

La magnitud pedida en valor absoluto es,

[math] \left\|\sigma \cdot \vec{m} - (\vec{m}\cdot\sigma\cdot\vec{m})\,\vec{m}\right\| = \left|\sigma_{\theta\theta}\right|\, \frac{\rho^2 - 1}{\rho^3}. [/math]

Sustituyendo [math]\sigma_{\theta\theta} = (5\rho - 4)/5[/math],

[math] T(\rho) = \frac{5\rho - 4}{5} \cdot \frac{\rho^2 - 1}{\rho^3} = \frac{(5\rho - 4)(\rho^2 - 1)}{5\rho^3}. [/math]

Observa que [math]T(\rho)[/math] depende solo de [math]\rho[/math] .

¿Dónde son mayores?

Para [math]\rho \in [1,2][/math], la función [math]T(\rho)[/math] se anula en [math]\rho = 1[/math] y es positiva para [math]\rho \gt 1[/math]. Evaluando en el extremo exterior [math]\rho = 2[/math]:

[math] T(2) = \frac{(5\cdot 2 - 4)(2^2 - 1)}{5\cdot 2^3} = \frac{(10 - 4)\cdot 3}{40} = \frac{18}{40} = 0.45. [/math]

Un análisis numérico muestra que [math]T(\rho)[/math] crece monótonamente en [math][1,2][/math] y alcanza su máximo en [math]\rho = 2[/math].

Las deformaciones calculadas eran:

[math] \varepsilon_{\rho\rho} = \frac{2\rho - 1}{5}, \qquad \varepsilon_{\theta\theta} = \frac{\rho - 1}{5}, [/math]

ambas crecen con [math]\rho[/math] en [math][1,2][/math] y alcanzan sus valores máximos en [math]\rho = 2[/math].

Comparación

Las tensiones tangenciales respecto al plano ortogonal a [math]\tfrac{1}{\rho}\vec{e}_\theta[/math] son máximas en el borde exterior [math]\rho = 2[/math], coincidiendo con los puntos de mayor deformación . Es decir, en este caso los máximos de [math]T(\rho)[/math] y de las deformaciones ocurren ambos en [math]\rho = 2[/math].

r = linspace(1,2,100);         
theta = linspace(0,pi,200);
[R,TH] = meshgrid(r,theta);
R = R';
TH = TH';
X = R .* cos(TH);
Y = R .* sin(TH);

sigma_tt = (1/5) .* (5.*R - 4);       
T_tang = (sigma_tt ./ R) .* (1 - 1./(R.^2));   

eps_rr = (1/5) .* (2.*R - 1);
eps_tt = (1/5) .* (R - 1);
trace_u = eps_rr + eps_tt;   % = 1/5 (3 rho - 2);
subplot(1,2,1);
contourf(X, Y, T_tang, 30, 'LineColor','none');
axis equal; colorbar;
title('Magnitud de la tensión tangencial respecto a (1/\rho e_\theta)');
xlabel('x');
ylabel('y');

hold on;
plot([1*cos(theta); 2*cos(theta)], [1*sin(theta); 2*sin(theta)], 'k', 'LineWidth', 0.6);
hold off;


subplot(1,2,2);
contourf(X, Y, trace_u, 30, 'LineColor','none');
axis equal; colorbar;
title('Deformación (traza = \nabla\cdot u)');
xlabel('x'); ylabel('y');
hold on;
plot([1*cos(theta); 2*cos(theta)], [1*sin(theta); 2*sin(theta)], 'k', 'LineWidth', 0.6);
hold off;

rho_vals = [1, 1.2, 1.4, 1.6, 1.8, 2.0];
fprintf('rho

11 Densidad

Si la densidad de la placa viene dada por

[math] d(\rho,\theta) = 1 + e^{\rho^2}\cos\theta, [/math]

calcular la masa aproximando la integral numéricamente.

El dominio es [math]\rho \in [1,2][/math], [math]\theta \in [0,\pi][/math] .

La masa [math]M[/math] se obtiene integrando la densidad sobre el área (usar coordenadas polares, elemento de área [math]dA = \rho\, d\rho\, d\theta[/math]):

[math] M = \int_{\theta=0}^{\pi} \int_{\rho=1}^{2} \left( 1 + e^{\rho^2}\cos\theta \right)\rho \, d\rho \, d\theta. [/math]

Separamos la integral en dos términos:

[math] M = \int_{0}^{\pi} \int_{1}^{2} \rho \, d\rho \, d\theta \;+\; \int_{0}^{\pi} \int_{1}^{2} e^{\rho^2}\cos\theta \, \rho \, d\rho \, d\theta. [/math]

Observa que [math]\int_{0}^{\pi} \cos\theta \, d\theta = 0[/math]. Por tanto el segundo término se anula y

[math] M = \int_{0}^{\pi} d\theta \int_{1}^{2} \rho \, d\rho = \pi \cdot \left[ \frac{\rho^2}{2} \right]_{\rho=1}^{2} = \pi \cdot \frac{4 - 1}{2} = \frac{3\pi}{2}. [/math]

Resultado

[math] M = \frac{3\pi}{2} \approx 4.71238898038469. [/math]

Resolución con Matlab de la integral:

f = @(th,r) (1 + exp(r.^2).*cos(th)).*r;

M = integral2(f, 0, pi, 1, 2);

12 Aplicación del trabajo

Para dar una aplicación concreta, interpretaremos el dominio del problema (arco circular definido entre los radios 1 y 2) como una sección de la corteza terrestre, y consideraremos que el campo de desplazamientos 𝑢 representa el movimiento inducido por ondas sísmicas, en particular las ondas S, que son ondas transversales que se propagan durante un terremoto.

Aplicado a lo anteriormente calculado:

La función temperatura:

Podría representar, en este contexto, un gradiente térmico asociado a zonas profundas de la corteza donde existe actividad volcánica o magmática

El campo de desplazamientos:

Puede interpretarse como un desplazamiento radial provocado por el paso de ondas sísmicas S

La divergencia:

Representa cambios locales de volumen en el material

El rotacional:

Representa la tendencia del material a rotar, torcerse o cizallarse, algo que es característico de las ondas S Si es elevado --> deformación intensa

Cálculo del tensor de deformaciones y tensor de tensiones:

Nos permite determinar cómo se distribuyen las tensiones internas en el material

Masa de la placa (densidad variable):

Representa una distribución no uniforme del material

13 Conclusión

Este trabajo no solo nos permite practicar conceptos de cálculo vectorial, tensores y campos, sino que además representa de manera simplificada una herramienta real para estudiar cómo se comportan los materiales bajo deformaciones, incluyendo uno de los fenómenos naturales más relevantes: los terremotos.

14 Póster

Descargar PDF

15 Bibliografía

https://simula.industriales.upm.es/apuntes/mcd.pdf

https://riunet.upv.es/bitstreams/baf345fa-326d-4f11-aa52-7655bbb6f4b5/download

https://books.google.com/books/about/Introduction_to_Continuum_Mechanics.html?id=lEhh-hjG6EgC

https://ceae.colorado.edu/~amadei/CVEN5768/PDF/NOTES3.pdf

https://www.civil.northwestern.edu/people/rudnicki/Continuum/cmbook_11_03_2011.pdf