Usuario:Albertodiaz.r

Trabajo realizado por estudiantes
Título	Campos escalares y vectoriales en elasticidad (Grupo 14)
Asignatura	Teoría de Campos
Curso	2023-24
Autores	ALBERTO DÍAZ REÑONES, MIGUEL ÁNGEL DE GREGORIO ÁVILA, ABDALLAH ATTAR ALTARAZI, EDUARDO CASADO PÉREZ, ALEJANDRO GARCÍA FERNÁNDEZ
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

Estudio de Campos Escalares y Vectoriales en Elasticidad

Para este artículo, vamos a ver las distintas aplicaciones matemáticas que simbolizan fenómenos físicos, sobre una placa rectangular plana en dimensión 2, comprendida entre los intervalos [-1,1] x [0,12] para las variables x e y. Procedemos a definir dos cantidades físicas que van a ayudar a su modelización:

• Por un lado la temperatura T: [math] T(x,y) = log(1+x^2) + log(1+(y-4) ^2) [/math];
• Por otro lado, la posición de cada punto tras la deformación, definida como [math] \vec r_d(x,y) [/math] = [math] \vec r_0(x,y) [/math] + [math] \vec u(x,y) [/math], siendo [math] \vec r_0(x,y) [/math] el vector de posición antes de la deformación, y [math] \vec u(x,y) [/math] el desplazamiento por la acción de una fuerza externa.

Adicionalmente, vamos a suponer que esa fuerza externa aplicada ha provocado un desplazamiento ondulatorio de los puntos representados por el vector: [math] \vec u(x,y,t) [/math] = [math] \vec a [/math] * sin(Π*k*([math] \vec d [/math] * [math] \vec r_0(x,y) [/math] - vt), donde [math] \vec a [/math] es la amplitud, k > 0 es el nº de onda y [math] \vec d [/math] es un vector unitario que marca a dirección de propagación y v es la velocidad de propagación. La variable t es usada para definitir el tiempo, pero usaremos t=0 en los apartados iniciales, luego la fuerza aplicada es [math] \vec u(x,y,t) [/math] = [math] \vec a [/math] * sin(Π*k*([math] \vec d [/math] * [math] \vec r_0(x,y) [/math]). Como supondremos que es una onda transversal en la que la dirección de propagación es ortogonal a la amplitud, [math] \vec a [/math] tomará el valor de [math] \vec 1/3i [/math], k tendrá valor unitario (k=1) y [math] \vec d [/math] será igual a [math] \vec 1/3j [/math]

1 MALLADO DE LOS PUNTOS INTERIORES DEL SÓLIDO

El primer paso que hemos de hacer es condicionar un paso de muestreo, en este caso h=2/10. Una vez decidido el valor, procedemos a representar nuestro sólido (la placa rectangular plana). Definimos los ejes [-1,1] x [0,12]. El paso de muestro es el nº de divisiones por eje que se hacen, luego el sólido queda dividido en 5 partes. Finalmente, diseñamos una malla que contiene una superficie con todas las combinaciones de los vectores x e y, similar una matriz con las dimensiones las de los vectores de las variables.

MALLADO DE LA PLACA

%Definimos el paso de muestreo y las variables
h=2/10;
x=[-1:h:1];
y=[0:h:12];
% Se realiza el mallado.
[Mx,My]=meshgrid(x,y);
%figure(1)
% Se representa el mallado con tres ejes, de ellos uno nulo.
mesh(Mx,My,Mx*0);
% Se pone título a la gráfica.
title('Mallado de la placa');
% Se da nombre  a los ejes.
xlabel('Eje de las X');
ylabel('Eje de las Y');
% Se da equidistancia a los ejes.
axis equal;
% Se define el rango de visión de la gráfica.
axis([-1.5,1.5,-0.5,12.5]);
% Se cambia a vista cenital
view(2)

2 CURVAS DE NIVEL DE LA TEMPERATURA. GRADIENTE DE T [math] ∇(T) [/math]

En este apartado, primeramente, vamos a calcular y representar las curvas de nivel de la temperatura (usando el comando contour entre otros). Veremos en que punto o puntos es máxima, y a continuación, vamos a representar el gradiente de T [math]∇(t)[/math]. En la imagen de la izquierda, hemos situado las curvas de nivel de la Temperatura, y en la parte derecha, el gradiente de T.

CURVAS DE NIVEL DE LA Tª

%Definimos las variables 
h = 2/10;
x=(-1:h:1);
y=(0:h:12);
% Creación del mallado
[X,Y]= meshgrid(x,y);
T= log(1 + X.^2) + log(1 + (Y - 4).^2);
% Escribimos el título del gráfico
title('Curvas de nivel de la temperatura');
% Escribimos los nombres de los ejes
xlabel('Eje X');
ylabel('Eje Y');
% Establecemos el límite de los ejes
axis([-1,1,0,12]); 
% Representación de la temperatura y las curvas de nivel
subplot(1,2,1);
mesh(X,Y,T);
subplot(1,2,2); 
contour(X,Y,T,50);
%añadimos la barra de colores
colorbar
axis equal

3 LEY DE FOURIER

De acuerdo con la Ley de Fourier, la energía calorífica [math]\vec{Q}[/math] viaja de acuerdo con la formula [math]\vec Q [/math]=−k[math]∇T[/math] . Donde k es la constante de conductividad térmica de la placa, que supondremos que es k=1.

Para el desarrollo de este apartado, utilizaremos nuevamente el código del programa anterior, conservando las variables h, x e y, junto con las matrices que hemos generado con la función meshgrid(), X e Y, además de la matriz T que representa la temperatura.

Después, procederemos a calcular el gradiente en coordenadas cartesianas, siguiendo la siguiente expresión.

[math]\frac {\partial T} {\partial x} \vec{i} + \frac {\partial T} {\partial y} \vec{j} = \frac{2x} {1+x^2} \vec{i}+\frac{2(y-4)} {1+(y-4) ^2} \vec{j}[/math]

A continuación, crearemos las variables [math]dx[/math] y [math]dy[/math] en función de las matrices X e Y. Y por último, utilizaremos el comando quiver() que nos permitirá representar el gradiente de la temperatura de una manera detallada.

Para comprobar la ortogonalidad entre las curvas de nivel del campo de temperatura y el gradiente, usaremos el comando hold on para poder representar las curvas de nivel y el gradiente en un mismo gráfico.

GRADIENTE

%Se definen las variables.
h=2/10;
x=[-1:h:1];
y=[0:h:12];
%Definimos las matrices con las que vamos a trabajar
[X,Y]=meshgrid(x,y);
T=log(1+X.^2)+log(1+(Y-4).^2);
%Derivada Parcial respecto X
dx=-(2.*X)./(X.^2+1);
%Derivada Parcial respecto Y
dy=-(2.*(Y-4))./(1+(Y-4).^2);
figure(3)
contour(X,Y,T,15)
hold on
quiver(X,Y,dx,dy)
hold off
title('Gradiente de la Placa')
xlabel('EJE X')
ylabel('EJE Y')
zlabel('EJE Z')
% Se define el rango de visión de la gráfica.
axis ([-1.5,1.5,-0.5,12.5]);
axis equal
% Se le aplica una barra con colores asignados a cada temperatura.
colorbar;

4 CAMPO DE VECTORES DE DESPLAZAMIENTO

Dibujar el campo de vectores en los puntos del mallado del sólido, en t=0.

Iniciaremos mediante el cálculo de [math]\vec{u}[/math]. Luego, procederemos a la representación de las variables h, x e y, emplearemos el comando meshgrid() para obtener las matrices X e Y. Después estableceremos el campo de desplazamientos utilizando las variables ux y uy en relación con las matrices X e Y. A continuación, emplearemos el comando quiver() para visualizar el campo de vectores de desplazamiento. Para concluir, ajustaremos los ejes según los límites especificados y les asignaremos nombres, al igual que al gráfico.

CAMPO DE VECTORES

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[Mx,My]=meshgrid(x,y);
ux=(1/3).*sin(pi/3.*My);
uy=0.*My;
figure(4)
% Campo de vectores.
quiver(Mx,My,ux,uy);
% Se pone título a la gráfica.
title('CAMPO DE VECTORES EN T=0');
xlabel('EJE X')
ylabel('EJE Y')
axis equal;
axis([-1.5,1.5,-0.5,12.5]);

5 SÓLIDO ANTES Y DESPUÉS DEL DESPLAZAMIENTO DADO POR EL CAMPO DE VECTORES

Dibujar el sólido antes y después del desplazamiento dado por el campo de vectores [math]\vec{u}[/math] (en t=0). Dibujar ambos en la misma figura usando el comando subplot.

Empezamos definiendo las regiones y los parámetros del sólido, antes y después del desplazamiento del sólido. Con el comando subplot representamos ambas situaciones y una representación de las dos juntas.

DESPLAZAMIENTO Y COMPARACIÓN

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[Mx,My]=meshgrid(x,y);
ux=(1/3).*sin(pi/3.*My);
uy=0.*My;
figure(5)
%Antes de la deformación
subplot(2,2,1)
mesh(Mx,My,Mx*0);
title('ANTES DEL DESPLAZAMIENTO');
xlabel('EJE X')
ylabel('EJE Y')
axis equal;
axis([-1.5,1.5,-0.5,12.5]);
view(2)
%Después del desplazamiento
subplot(2,2,2)
mesh(Mx+ux,My+uy,0*My);
title('DESPUÉS DEL DESPLAZAMIENTO');
xlabel('EJE X')
ylabel('EJE Y')
axis equal;
axis([-1.5,1.5,-0.5,12.5]);
view(2)
%Deformación
subplot(2,2,3)
mesh(Mx,My,Mx*0);
hold on
mesh(Mx+ux,My+uy,0*My,'LineWidth',1);
hold off
title('COMPARACIÓN');
xlabel('EJE X')
ylabel('EJE Y')
axis equal;
axis([-1.5,1.5,-0.5,12.5]);
view(2)

6 DIVERGENCIA DE U ( [math] ∇·\vec{u} [/math] )

Dibujar [math]∇[/math]·[math]\vec{u}[/math] en t=0. Determinar analíticamente los puntos en los que la divergencia de [math]\vec{u}[/math] es máxima, mínima y nula. La divergencia es una medida del cambio de volumen local debido al desplazamiento.

¿Se puede apreciar esto en la gráfica?

Calcularemos la divergencia a través de la siguiente expresión: [math]∇·\vec{u} = \frac{∂}{∂x}(\frac{1}{3}·sen(\frac{Π}{3}y))+\frac{∂}{∂y}(0) = 0[/math]

Obtenemos un color uniforme en la gráfica debido a que el laplaciano nos ha dado 0. Por lo tanto la placa tendrá un color uniforme y no habrá cambio de volumen debido al desplazamiento en ningún punto.

DIVERGENCIA

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[Mx,My]=meshgrid(x,y);
% Calculamos la divergencia y da cero.
D=0.*Mx+0.*My;
%Divergencia
surf(Mx,My,D)
%El shading flat hace que se vea uniforme
shading flat
title('DIVERGENCIA DE LA PLACA');
xlabel('EJE X');
ylabel('EJE Y');
axis equal;
axis([-1.5,1.5,-0.5,12.5]);
colorbar;
view(2)

7 CÁLCULO DEL ROTACIONAL EN TODOS LOS PUNTOS DEL SÓLIDO

En este apartado buscamos calcular el rotacional de los puntos del sólido en t=0. Definimos primeramente la función [math]\vec{u}(x,y)=\frac{1}{3}sen(\frac{Π}{3}y)\vec{i}[/math]

Para calcular el rotacional de un campo de desplazamientos, aplicamos la fórmula:

[math]∇ × \vec{u}= \begin{vmatrix}\vec{i} & \vec{j}& \vec{k}\\ \frac{\partial }{\partial x} & \frac{\partial }{\partial y}& \frac{\partial }{\partial z}\\ u_{1}& u_{2}& u_{3} \end{vmatrix} = \begin{vmatrix} \vec{i} & \vec{j} & \vec{k}\\ \frac{\partial }{\partial x} & \frac{\partial }{\partial y} &\frac{\partial }{\partial z} \\ \frac{1}{3}sen(\frac{Π}{3}y) & 0 & 0\end{vmatrix} = \frac{1}{3}sen(\frac{Π}{3}y) \vec{k}[/math]

Siendo el módulo: [math]|∇ × \vec{u}|= \frac{1}{3}sen(\frac{Π}{3}y)[/math]

ROTACIONAL

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[Mx,My]=meshgrid(x,y);
R=(pi/9).*cos((pi/3).*My);
%Calculamos el rotacional R:
surf(Mx,My,R);
shading flat
title('ROTACIONAL');
xlabel('EJE X');
ylabel('EJE Y');
axis equal;
axis([-1.25,1.25,-0.25,12.5]);
colorbar;
view(2)

8 TENSIONES NORMALES

En un medio elástico lineal, isótropo y homogéneo los desplazamientos permiten escribir el tensor de tensiones σij a través de la fórmula [math]σ=λ∇·\vec{u}1 + 2μ∈[/math], donde 1 es el tensor identidad, λ y µ son los coeficientes de Lamé.

Calculamos ∈ como la parte simétrica del tensor gradiente de [math]\vec u[/math]. Definimos la parte simétrica del tensor gradiente de [math]\vec u [/math], conocido como tensor de deformaciones: ∈([math]\vec u [/math])=(∇[math]\vec u [/math]+∇[math]\vec u^t [/math])/2. Procedo primero a calcular el gradiente de T [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {\partial u_1 \over \partial x} & {\partial u_1 \over \partial y} & {\partial u_1 \over \partial z}\\ {\partial u_2 \over \partial x} & {\partial u_2 \over \partial y} & {\partial u_2 \over \partial z}\\ {\partial u_3 \over \partial x} & {\partial u_3 \over \partial y} & {\partial u_3 \over \partial z}\\ \end{pmatrix} [/math] Debido a que nuestras 2ª y 3ª fila son 0 (ya que no tenemos componente [math]\vec u_2 ni u_3 [/math], la matriz resultante es: [math]\vec u [/math] : [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0}\\ {0} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Haciendo la matriz traspuesta queda también: [math]\vec ∇u^t [/math] : [math] \mathbb (∇\vec u^t) = \; \begin{pmatrix} {0} & {0} & {0}\\ {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Con estas 2 matrices, defino la parte simétrica del tensor gradiente de [math]\vec u [/math], conocido como tensor de deformaciones: ∈([math]\vec u [/math])=(∇[math]\vec u [/math]+∇[math]\vec u^t [/math])/2 = [math] \mathbb = \; \begin{pmatrix} {0} & {\frac{Π}{18}cos(\frac{Π}{3}y)} & {0}\\ {\frac{Π}{18}cos(\frac{Π}{3}y)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Uso ahora la fórmula del tensor de tensiones [math]σ=λ∇·\vec{u}1 + 2μ∈[/math]. Calculando la divergencia de u ([math] ∇· \vec u [/math]) resulta 0, ya que nuestra única componente de [math]\vec{u}(x,y)=\frac{1}{3}·sen(\frac{Π}{3}y-v·t)\vec{i}[/math] es [math]\vec u_1 [/math], y la fórmula de la divergencia es [math] ∇·\vec{u}= \frac{\partial u_{1} }{\partial\ x}+\frac{\partial u_{2}}{\partial y}+\frac{\partial u_{3}}{\partial z} [/math], luego el tensor de tensiones finalmente es igual a [math]σ=2∈ [/math] ya que μ=1 [math] \mathbb (σ) = \; \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0}\\ {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Dibujamos ahora las tensiones normales en la dirección que marca el eje [math] \vec i [/math],[math] \vec j [/math], [math] \vec k [/math]

[math]\vec({i}·σ·\vec{i})=\begin{pmatrix} 1 & 0 & 0 \end {pmatrix} · \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0}\\ {\frac{Π}{9}·(\frac{Π}{3}y)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math] · [math] \begin{pmatrix} 1\\ 0\\ 0\\ \end {pmatrix} [/math] = 0

[math]\vec({j}·σ·\vec{j})=\begin{pmatrix} 0 & 1 & 0 \end {pmatrix} · \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0}\\ {\frac{Π}{9}·(\frac{Π}{3}y)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math] · [math] \begin{pmatrix} 0\\ 1\\ 0\\ \end {pmatrix} [/math] = 0

[math]\vec({k}·σ·\vec{k})=\begin{pmatrix} 0 & 0 & 1 \end {pmatrix} · \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y)} & {0}\\ {\frac{Π}{9}·(\frac{Π}{3}y)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math] · [math] \begin{pmatrix} 0\\ 0\\ 1\\ \end {pmatrix} [/math] = 0

Ya que todas son nulas, no podemos representarlo gráficamente en Matlab.

9 TENSIONES TANGENCIALES

Este apartado consiste en calcular las tensiones generadas respecto al plano ortogonal a [math]\vec{i}[/math] y representarlas en el caso de que existiesen. Conseguimos hallarlas mediante la siguiente fórmula [math]|σ·\vec{i}-(\vec{i}·σ·\vec{i})·\vec{i}|[/math]. En nuestro caso el sustraendo va a ser igual a cero por lo cual el único término que queda es el minuendo, esto quiere decir que no hay ningún tipo de tensión en la dirección del eje [math]\vec{i}[/math], lo que implica que el único componente de tensión es tangencial al plano ortogonal al eje x. [math]|σ·\vec{i}-(\vec{i}·σ·\vec{i})·\vec{i}|[/math] = [math]|\frac{π}{9}·cos(π·y/3)\vec{j}-0|[/math] = [math]|\frac{π}{9}·cos(π·y/3)\vec{j}|[/math]. Vamos a mostrar a continuación la representación gráfica de las tensiones respecto a los diferentes puntos de la placa:

TENSIONES TANGENCIALES

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[Mx,My]=meshgrid(x,y);
Ttg=(pi/9).*cos((pi/3).*My); %Tensióntangencial respecto al plano ortogonal a i
quiver(Mx,My,Ttg.*0,Ttg);
title('Tensión respecto al plano ortogonal a i')
xlabel
ylabel
axis equal
axis([-1.5,1.5,-0.5,12.5])
view(2)

10 TENSIÓN DE VON MISSES

Para este apartado vamos a trabajar con la fórmula de Von Mises

[math]σ_{VM}=\sqrt{\frac{(σ_{1}-σ_{2})^2+(σ_{2}-σ_{3})^2+(σ_{3}-σ_{1})^2}2}[/math]

que nos indica cuando un material inicia un comportamiento plástico en lugar de tener las cualidades de un material elástico puro. La fórmula de Von Mises depende [math]σ_{1}[/math],[math]σ_{2}[/math] y [math]σ_{3}[/math] que son las tensiones principales en el material que a su vez son los autovalores asociados a la matriz de σ:

[math]\mathbb σ=\; \begin{pmatrix} {0}&{\frac{π}{9}·cos(\frac{π·y}{3})}&{0}\\ {\frac{π}{9}·cos(\frac{π·y}{3})}&{0}&{0}\\ {0}&{0}&{0}\\ \end{pmatrix}[/math], siendo los autovalores [math]σ_{1}=0,σ_{2}=\frac{π}{9}·cos(\frac{π·y}{3}),σ_{3}=-\frac{π}{9}·cos(\frac{π·y}{3})[/math]

Si sustituimos los valores que hemos hallado en las tensiones principales en la fórmula de Von Mises nos da el siguiente resultado:

[math]σ_{VM}=\sqrt{\frac{(0-\frac{π}{9}·cos(\frac{π·y}{3}))^2+(\frac{π}{9}·cos(\frac{π·y}{3})-(-\frac{π}{9}·cos(\frac{π·y}{3}))^2+(-\frac{π}{9}·cos(\frac{π·y}{3})-0)^2}2}=\sqrt{3}·\frac{π}{9}·cos(\frac{π·y}{3})[/math]

VON MISES

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[Mx,My]=meshgrid(x,y);
M_VonMises=0.*My;
VonMises=inline('(((T1-T2)^2+(T2-T3)^2+(T3-T1)^2)/2)^(1/2)','T1','T2','T3'); %Definimos la
función de Von Mises siendo T1,T2,T3 las tensiones elementales
[a,b]=size(My);
for i=1:a
for j=1:b
Deformaciones=[[0;(pi/9).*cos((pi/3).*My(i,j));0],[(pi/9).*cos((pi/3).*My(i,j));0;0],[0;0;0]];
Lamb=eig(Deformaciones);
T1=Lamb(1,1);
T2=Lamb(2,1);
T3=Lamb(3,1);
M_VonMises(i,j)=VonMises(T1,T2,T3);
end
end
surf(Mx,My,M_VonMises)
shading flat
title('Von Mises')
xlabel=('Eje x')
ylabel=('Eje y')
axis equal
axis([-1.5,1.5,-0.5,12.5])
colorbar
view(2)
max(max(M_VonMises))

Vemos que son varios los puntos de la placa donde la función de Von Mises toma su valor máximo y esto es cada vez que el [math]cos(\frac{π·y}{3})=1[/math] siendo su valor [math]σ=0.6046[/math]

11 CAMPO DE FUERZAS [math] \vec F [/math] CAUSANTES DEL DESPLAZAMIENTO USANDO LA ECUACIÓN DE LA ELASTICIDAD LINEAL

Para este apartado, vamos a usar la ecuación de la elasticidad lineal, definida por [math]\vec{F}=\frac{\partial^2\vec{u}}{\partial t^2}-∇·σ[/math]

Definimos primeramente la función [math]\vec{u}(x,y)=\frac{1}{3}sen(\frac{Π}{3}-vt)·\vec{i}[/math]. Lo primero que hemos de hacer, es calcular la segunda derivada, luego:

• Primera derivada: [math]\frac{\partial\vec{u}}{\partial t}=\frac{v}{3}cos(\frac{Π}{3}y - vt)\vec{i}[/math]
• Segunda derivada: [math]\frac{\partial^2\vec{u}}{\partial t^2}=\frac{-v^2}{3}sen(\frac{Π}{3}y - vt)\vec{i}[/math]

Ahora que tenemos calculado la segunda derivada, tenemos que calcular el ∇·σ, que se obtiene al hacer la divergencia de los vectores cuyas componentes son las filas de la matriz σ. Luego, tenemos que calcular el tensor de tensiones del apartado 8. Para ello:

• 1º: Cálculo del gradiente de [math]\vec u [/math] : [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {\partial u_1 \over \partial x} & {\partial u_1 \over \partial y} & {\partial u_1 \over \partial z}\\ {\partial u_2 \over \partial x} & {\partial u_2 \over \partial y} & {\partial u_2 \over \partial z}\\ {\partial u_3 \over \partial x} & {\partial u_3 \over \partial y} & {\partial u_3 \over \partial z}\\ \end{pmatrix} [/math]

Como la 2ª y 3ª fila son 0, debido a que nuestro vector [math] \vec u [/math] solo tiene componente [math] \vec u_1 [/math] ([math] \vec u_2 [/math] y [math] \vec u_3 [/math] son 0), la matriz queda así: Gradiente de [math]\vec u [/math] : [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0}\\ {0} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Hacemos la matriz traspuesta de [math]\vec ∇u [/math] : [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {0} & {0} & {0}\\ {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Ahora, con estas 2 matrices, puedo definir la parte simétrica del tensor gradiente de [math]\vec u [/math], conocido como tensor de deformaciones: ∈([math]\vec u [/math])=(∇[math]\vec u [/math]+∇[math]\vec u^t [/math])/2 = [math] \mathbb = \; \begin{pmatrix} {0} & {\frac{Π}{18}cos(\frac{Π}{3}y-vt)} & {0}\\ {\frac{Π}{18}cos(\frac{Π}{3}y-vt)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Ahora usaremos la fórmula del tensor de tensiones [math]σ=λ∇·\vec{u}1 + 2μ∈[/math]. Si calculamos la divergencia de u ([math] ∇· \vec u [/math]) nos resulta 0, ya que nuestra única componente de [math]\vec{u}(x,y)=\frac{1}{3}·sen(\frac{Π}{3}y-vt)\vec{i}[/math] es [math]\vec u_1 [/math], y la fórmula de la divergencia es [math] ∇·\vec{u}= \frac{\partial u_{1} }{\partial\ x}+\frac{\partial u_{2}}{\partial y}+\frac{\partial u_{3}}{\partial z} [/math], luego el tensor de tensiones finalmente es igual a [math]σ=μ·2∈[/math] = [math] \mathbb (σ) = μ· \; \begin{pmatrix} {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0}\\ {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

• 2º: Continuo calculando la divergencia de los vectores cuyas componentes son las filas de la matriz σ, siendo ∇·σ:

[math]∇· \begin{pmatrix}0 & \frac{Π}{9}cos(\frac{Π}{3}y-vt)) & 0\end{pmatrix}=-μ\frac{Π^2}{27}sen(\frac{Π}{3}y-vt)[/math]

[math]∇· \begin{pmatrix}\frac{Π}{9}cos(\frac{Π}{3}y-vt) & 0 & 0\end{pmatrix}=0 [/math]

[math]∇· \begin{pmatrix}0 & 0 & 0\end{pmatrix}=0 [/math]

• 3º: Aplico la fórmula de la elasticidad lineal: [math]\vec{F}=\frac{\partial^2\vec{u}}{\partial t^2}-∇·σ[/math]

[math]\vec{F}=\frac{-v^2}{3}sen(\frac{Π}{3}y - vt)+μ·\frac{Π^2}{27}sen(\frac{Π}{3}y-vt)[/math]

Suponemos ahora que la [math]\vec{F}=0[/math], vamos a poner la velocidad en términos de λ y μ, resultando finalmente: v = [math]\frac{Π}{3}·\sqrtμ[/math]

Vamos a suponer ahora que la onda fuese longitudinal, para ello tomamos [math]\vec a = \frac{1}{3}j [/math] y calcularemos ahora la velocidad también en términos de λ y μ. Esta vez hemos de calcular nuestra nueva [math]\vec u [/math] que será igual a [math]\vec{u}(x,y)=\frac{1}{3}sen(\frac{Π}{3}-vt)\vec{j}[/math] Repito los mismos pasos que anteriormente:

• Primera derivada: [math]\frac{\partial\vec{u}}{\partial t}=\frac{v}{3}cos(\frac{Π}{3}y - vt)\vec{j}[/math]
• Segunda derivada: [math]\frac{\partial^2\vec{u}}{\partial t^2}=\frac{-v^2}{3}sen(\frac{Π}{3}y - vt)\vec{j}[/math]

Ahora que tenemos calculado la segunda derivada, tenemos que calcular el ∇·σ, que se obtiene al hacer la divergencia de los vectores cuyas componentes son las filas de la matriz σ. Luego, tenemos que calcular el tensor de tensiones del apartado 8. Como la 1ª y 3ª fila son 0, debido a que nuestro vector [math] \vec u [/math] solo tiene componente [math] \vec u_2 [/math] ([math] \vec u_1 [/math] y [math] \vec u_3 [/math] son 0), la matriz queda así:

• 1º: Cálculo del Gradiente de [math]\vec u [/math] : [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {0} & {0} & {0}\\ {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Hacemos la matriz traspuesta de [math]\vec ∇u [/math] : [math] \mathbb (∇\vec u) = \; \begin{pmatrix} {0} & {0} & {0}\\ {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

En este caso, como las matrices coinciden, el tensor deformaciones es finalmente igual a: ∈ [math] \mathbb = \; \begin{pmatrix} {0} & {0} & {0}\\ {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Ahora usaremos la fórmula del tensor de tensiones [math]σ=λ∇·\vec{u}1 + 2μ∈[/math]. Si calculamos la divergencia de u (∇⋅σ) vemos que ya no es nula, sino que la componente [math]\vec u_2 [/math] se puede derivar [math]\partial y [/math], luego la divergencia de [math]\vec u [/math] es ∇·[math]\vec u [/math] = [math]\frac{Π}{9}cos(\frac{Π}{3}y-vt)[/math]

El tensor de tensiones finalmente es igual a [math]σ=λ∇·\vec{u}1 + 2μ∈[/math] [math] \mathbb = \; \begin{pmatrix} {\frac{Π}{9}cos(\frac{Π}{3}y-vt)·λ} & {0} & {0}\\ {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)·(2μ+λ)} & {0}\\ {0} & {0} & {\frac{Π}{9}cos(\frac{Π}{3}y-vt)·λ}\\ \end{pmatrix} [/math]

• 2º: Continuo calculando la divergencia de los vectores cuyas componentes son las filas de la matriz σ, siendo ∇·σ:

[math]∇· \begin{pmatrix}\frac{Π}{9}cos(\frac{Π}{3}y-vt)λ & 0 & 0\end{pmatrix}=0 [/math]

[math]∇· \begin{pmatrix}0 & \frac{Π}{9}cos(\frac{Π}{3}y-vt)(2μ+λ)) & 0\end{pmatrix}=-\frac{Π^2}{27}(2μ+λ)sen(\frac{Π}{3}y-vt)[/math]

[math]∇· \begin{pmatrix}0 & 0 & \frac{Π}{9}cos(\frac{Π}{3}y-vt)λ\end{pmatrix}=0 [/math]

• 3º: Aplico la fórmula de la elasticidad lineal: [math]\vec{F}=\frac{\partial^2\vec{u}}{\partial t^2}-∇·σ[/math]

[math]\vec{F}=\frac{-v^2}{3}sen(\frac{Π}{3}y - vt)\vec{j}[/math] + [math]\frac{Π^2}{27}(2μ+λ)sen(\frac{Π}{3}y-vt)[/math]

Suponemos ahora que la [math]\vec{F}=0[/math], vamos a poner la velocidad en términos de λ y μ, resultando finalmente: v = [math]\frac{Π}{3}·\sqrt{2μ+λ}[/math]

12 MÓDULO DE DESPLAZAMIENTO TRANSVERSAL VARIABLE SEGÚN EL TIEMPO

Para este apartado, con los cálculos realizados del apartado 11, vamos a calcular el módulo de desplazamiento transversal, (en este caso, según la dirección [math] \vec i [/math] ), según el intervalo de tiempo t ∈ [0,10]. Vamos también a representarlo gráficamente en Matlab, viendo los distintos valores numéricos (así como también analíticamente). Hemos adjuntado 3 imágenes, ya que hemos dado 3 valores distintos al vector linspace: 5, 10 y 100 puntos; y en cada una de las imágenes, podemos ver como a medida que aumentan los puntos, tenemos una gráfica totalmente distinta, más nítida y más similar a una función senoidal.

%Defino la velocidad=1, para tener solo una variable que es el tiempo t
v = 1;
%Definir la función en función del tiempo
t = linspace(0, 10, 100); %Hemos usado también 5 y 10 puntos, así como 100 puntos
u = (1/3) * sin(pi/3 - v * t); % Nuestra función u
%Creamos el gráfico
figure;
%Trazar la función en función del tiempo
plot(t,u); 
title('Función u en función del tiempo');
xlabel('Tiempo (s)');
ylabel('Valor de la función');
grid on; % Mostrar la cuadrícula