Visualización de campos escalares y vectoriales en elasticidad. (Grupo 42)

Trabajo realizado por estudiantes
Título	Visualización de Campos Escalares y Vectoriales en elasticidad, Grupo(42)
Asignatura	Teoría de Campos
Curso	2023-24
Autores	Andrés Ruiz, Jorge Martin, Nadir Ahnihan
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

1 Introducción

La redacción de este artículo tiene por objetivo el estudio del desplazamiento que adquiere un sólido al someterse a la acción de un campo de fuerzas externas. Para ello, se tienen dos cantidades físicas dependientes de las variables x e y:

La temperatura

[math] T(x,y)=3log(1+(x-1)^2)+log(1+(y-8)^2)[/math]

El campo de desplazamientos [math]\vec{u}(x,y)=\vec{a}·sin(\pi\cdot k(d·\vec{r0}(x, y)−vt))[/math] .

Tomando las variable del campo de desplazamientos los siguiente valores:

[math]\vec{a}=x/3 \vec{i} [/math];

[math]\vec{d}=1/12 \vec{j}[/math];

[math]k=1[/math];

Emplearemos el software de programación y cálculo numérico Matlab/Octave, que nos simplifica las operaciones complejas y además nos periten generar figuras en 2D y 3D que ilustran a la perfección las cálculos analíticos.

2 Representación de la placa rectangular plana.

La figura adjunta ilustra un mallado que representa los puntos interiores del sólido. Tomando como superficie del mallado el rectángulo (x, y) ∈ [−1; 1] × [0; 12] y como paso de muestreo h = 2/10 para las variables x e y.

Figura 1. Representación del mallado

%Se define el paso de muestreo para las variables x e y
h=2/10;
%Se definen los parámetros que representan la superficie de la placa
x=[-1:h:1];
y=[0:h:12];
%Se crea el mallado
[X,Y]=meshgrid(x,y);
%Dibujo de la malla, al ser mesh() un comando que requiere tres elementos
%de entrada, se toma una de las creadas y se multiplica por 0
mesh(X,Y,Y*0)
%Se definen los ejes, el título y el visualizado en dos dimensiones
axis([-5,5,0,14]);
title('Mallado de la placa')
xlabel('Eje X')
ylabel('Eje Y')
view(2);

3 Representaciones asociadas la temperatura

Asociada a la ecuación de la temperatura [math] T(x,y)=3log(1+(x-1)^2)+log(1+(y-8)^2)[/math]

empelada para este estudio, vamos a visualizar dos representaciones asociadas:

3.1 Curva de nivel

A continuación, vamos a dibujar las curvas de nivel asociadas a la temperatura,

[math] T(x,y)=3log(1+(x-1)^2)+log(1+(y-8)^2)[/math]

Distinguiendo con una escala de colores los diferentes valores que toma, siendo los más cálidos los puntos en los que mayor temperatura alcanza, y los más fríos donde es menor la temperatura.

El valor máximo de la temperatura es 9.0027 y se obtiene en el punto (-1,0).

Figura 2: Campo de temperaturas y sus curvas de nivel

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se define la función temperatura, T
T=3*log(1+(X-1).^2)+log(1+(Y-8).^2);
%Se representan las curvas de nivel de la temperatura
subplot(1,2,1)
surf(X,Y,T)
colorbar
view(2)
shading flat
subplot(1,2,2)
contour(X,Y,T,30)
title('Curvas de nivel')
%Tmax, valor máximo de la temperatura
Tmax=max(max(T))

3.2 Gradiente

El gradiente es una campo vectorial, que aplicado a la temperatura de nuestro estudio indica la razón del cambio de la temperatura del por unidad de distancia. El gradiente en coordenadas cartesianas se obtiene de la siguiente manera:

[math]\nabla T(x,y,z) =\frac{d∂}{dx} +\frac{d∂}{dy} + \frac{d∂}{dz}[/math]

Si lo aplicamos a nuestro estudio, obtenemos que:

[math]\nabla T(x,y) =\frac{d∂}{dx} +\frac{d∂}{dy}=\gt3 \cdot \frac{2(x-1)}{1+(x-1)^2} \vec i + \frac{2(y-8)}{1+(y-8)^2} \vec j [/math]

Observamos gráficamente que ∇T es ortogonal a dichas curvas.

Figura 3: Gradiente de la temperatura

Figura 4: Gradiente de la temperatura aumentado

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
T=3*log(1+(X-1).^2)+log(1+(Y-8).^2);
%Se calcula el gradiente de la temperatura
[TX,TY]=gradient(T);
%Se representa el gradiente sobre las curvas de nivel de la temperatura
hold on
quiver(X,Y,TX,TY)
contour(X,Y,T,30)
title('Gradiente de la temperatura')
hold off

4 Ley de Fourier

Jean-Baptiste Joseph Fourier, fue un matemático y físico francés conocido por sus trabajos sobre la descomposición de funciones periódicas en series trigonométricas convergentes llamadas Series de Fourier, método con el cual consiguió resolver la ecuación del calor. De acuerdo a la Ley de Fourier la energía calorífica [math]\vec{Q}[/math] viaja de acuerdo a la fórmula

[math]\vec Q = -k *\nabla T[/math]
(siendo k una constante con valor 1)

Hemos calculamos esta energía calorífica aplicada a los datos de nuestro estudio y la hemos representado como un campo vectorial.

Figura 5: Campo vectorial de la energía calorífica

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
T=3*log(1+(X-1).^2)+log(1+(Y-8).^2);
[TX,TY]=gradient(T);
%Se calcula la energía calorífica, Q, por la Ley de Fourier
k=-1;
QX=k*TX;
QY=k*TY;
%Se dibuja Q como campo vectorial
quiver(X,Y,QX,QY)
title('Campo vectorial de la energía calorífica')

5 Campo de vectores

A continuación, dibujamos el campo de vectores en los puntos del mallado del sólido, cuando el tiempo es nulo, es decir, en t = 0.

Figura 6: Campo de vectores en t=0

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se definen la variable a y las constantes d y k
a=X/3;
d=1/12;
k=1;
%Se definen los desplazamientos, u, para t=0;
ui=X/3.*sin(k*(pi*d).*Y);
uj=0.*Y;
%Se dibuja el campo de vectores en los puntos del mallado del sólido
hold on
axis([-5,5,0,14]);
title('Campo de vectores en los puntos del mallado, en t=0')
view(2);
mesh(X,Y,Y*0)
quiver(X,Y,ui,uj,'Linewidth',1.5)
hold off

6 Desplazamiento del sólido

Para este apartado, hemos representado en primer lugar el sólido en su posición inicial, es decir, antes de que el campo de fuerzas externo actúa sobre él, y en segunda estancia el solido desplazado debido a la acción de las fuerzas, en t=0 amabas figuras.

Figura 7: Comparación del sólido antes y después del desplazamiento

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
ui=X/3.*sin((pi/12).*Y);
uj=0.*Y;
%El vector posición en el instante inical, r0, es igual a Xi+Yj
%Se define el vector posición después del desplazamiento, rd
Xd=X+ui;
Yd=Y+uj;
%Se dibuja el mallado antes del desplazamiento
subplot(1,2,1)
mesh(X,Y,Y*0)
view(2);
axis([-5,5,0,14]);
title('Antes del desplazamiento')
%Se dibuja el mallado después del desplazamiento
subplot(1,2,2)
mesh(Xd,Yd,Yd*0)
view(2);
axis([-5,5,0,14]);
title('Después del desplazamiento')

7 Estudio de la divergencia

Para ver la divergencia, vamos a dibujar [math]\nabla\cdot\vec u[/math] en t = 0. Pudiendo apreciar en la grafica como la divergencia es una medida del cambio de volumen local debido al desplazamiento, que adopta la forma de una parábola con vértice en Y=6, donde la divergencia es máxima, y en la figura se corresponde con la tonalidad mas cálida mientras que el mínimo punto de divergencia de [math] \vec u [/math] esta situado en los extremos de la placa(x=0 y x=12), y su valor corresponde a 0, correspondiendose con los colores mas fríos de la figura. . Para el desarrollo de este apartado hemos utilizado que:

la divergencia: [math]\nabla\cdot\vec u= \frac{\partial}{\partial{x}}({\frac{x}{3} ({sen(\frac{πy}{12}})})[/math], es decir, [math]\nabla\cdot\vec u= \frac{1}{3} ({sen(\frac{πy}{12}}) [/math]

Figura 8: Divergencia

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se define la divergencia del vector u en t=0
div=(1/3)*sin((pi/12)*Y);
%Se dibuja la divergencia
subplot(1,2,1)
surf(X,Y,div)
shading flat
title('Divergencia')
%Se dibuja la divergencia en 2D
subplot(1,2,2)
surf(X,Y,div)
view(2)
shading flat
colorbar
%Se calculan los valores máximos y mínimos y los puntos nulos de la divergencia
DivMax=max(max(div))
DivMin=min(min(div))
%La divergencia es nula en los puntos mínimos

8 Rotacional de [math] \vec u [/math]

El rotacional de un campo vectorial indica la tendencia de un campo a inducir rotación alrededor de un punto. Viene determinad por la siguiente expresión:

[math]\nabla×\vec u(x,y) = {\begin{vmatrix} \vec i & \vec j & \vec k \\ \frac{\partial}{\partial x} & \frac{\partial}{\partial y} & \frac{\partial}{\partial z} \\ \vec u_x & \vec u_y & \vec u_z \end{vmatrix}} [/math]

Si aplicamos esto a nuestro campo vectorial [math] \vec u [/math], en t=0, obtenemos:

[math] \nabla×\vec u(x,y) = -\frac{\pi \cdot x}{36} \cdot cos(\frac{\pi\cdot y}{12}) \vec k[/math]

El valor máximo del rotacional es 0.0873 y se alcanza en los puntos (-1, 12) y (1, 0)

Figura 9: Módulo del rotacional

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se define el rotacional de u en t=0
rot=((pi/36)*X).*cos((pi/12)*Y);
%Se dibuja el rotacional
subplot(1,2,1)
surf(X,Y,rot)
shading flat
%Se dibuja el rotacional en 2D
subplot(1,2,2)
surf(X,Y,rot)
shading flat
view(2)
colorbar
title('Rotacional')
%Se calcula el valor máximo del rotacional
Rotmax=max(max(rot))

9 Tensiones tangenciales

Definamos [math] є(\vec{u}) = (∇\vec{u} + \triangledown \vec{u}^{t})/2[/math], la parte simétrica del tensor gradiente de ū conocido como tensor de deformaciones. 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]σ=λ·\nabla·\vec u I+2με[/math], donde 1 es el tensor identidad en el conjunto de vectores libres del espacio R 3 y λ, µ son los conocidos como coeficientes de Lamé que dependen de las propiedades elásticas de cada material. A pesar de que los desplazamientos son planos (es decir ū no tiene componente en la dirección de k)las tensiones no tienen por qué ser planas y puede haber tensiones en la dirección ortogonal al plano de la placa. Tomando λ = µ = 1, dibujar las tensiones normales en la dirección que marca el eje [math] \vec i [/math], es decir [math] \vec i · \sigma · \vec i [/math], las tensiones normales en la dirección que marca el eje [math] \vec j [/math], es decir [math] \vec j · \sigma · \vec j [/math] y las correspondientes al eje [math] \vec k [/math] , es decir [math] \vec k · \sigma · \vec k [/math] (dibujar las que no son nulas).

Se calcula el gradiente de [math]e(\vec{u})[/math] y su transpuesto. [math] \nabla\vec u= \begin{pmatrix}\frac {\partial u_1}{\partial x} & \frac{\partial u_1}{\partial y} & \frac{\partial u_1}{\partial z} \\ \frac{\partial u_2}{\partial x} & \frac{\partial u_2}{\partial y} & \frac{\partial u_2}{\partial z} \\ \frac{\partial u_3}{\partial x} & \frac{\partial u_3}{\partial y} & \frac{\partial u_3}{\partial z}\end{pmatrix} = \begin{pmatrix}\frac{1}{3}sin(\frac{π}{12}y) & \frac{πx}{36}cos(\frac{π}{12}y) & 0 \\ 0 & 0 & 0\\ 0 & 0 & 0\end{pmatrix} [/math]

[math]\nabla \vec u ^t= \begin{pmatrix}\frac{1}{3}sin(\frac{π}{12}y) & 0 & 0 \\ \frac{πx}{36}cos(\frac{π}{12}y) & 0 & 0 \\ 0 & 0 & 0\end{pmatrix} [/math]

Por lo que el tensor deformación es: [math] є(\vec{u}) = \begin{pmatrix}\frac{1}{3}sin(\frac{π}{12}y) & \frac{πx}{72}cos(\frac{π}{12}y) & 0 \\ \frac{πx}{72}cos(\frac{π}{12}y) & 0 & 0 \\ 0 & 0 & 0\end{pmatrix} [/math]

Aplicando la fórmula del tensor de tensiones, [math]\sigma = \left(\begin{matrix} \frac{1}{3}sin(\frac{πy}{12}) & 0 & 0 \\ & \\0 & \frac{1}{3}sin(\frac{πy}{12}) & 0 \\ & \\0 & 0 & \frac{1}{3}sin(\frac{πy}{12})\end{matrix} \right) + \left (\begin{matrix} \frac{2}{3}sin(\frac{πy}{12}) & \frac{πx}{36}cos(\frac{πy}{12}) & 0 \\&\\ \frac{πx}{36}cos(\frac{πy}{12}) & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right) [/math]

Sumando sale que [math]σ[/math] = [math] \left (\begin{matrix} (sin(\frac{πy}{12}) & (\frac{πx}{36})cos(\frac{πy}{12}) & 0 \\&\\ (\frac{πx}{36})cos(\frac{πy}{12}) & (\frac{1}{3})sin(\frac{πy}{12}) & 0 \\&\\ 0 & 0 & (\frac{1}{3})sin(\frac{πy}{12})\end{matrix} \right) [/math]

El módulo de la tensión normal según la dirección [math]\vec i[/math]:

[math] \vec i · \sigma · \vec i = \sin (\frac{y·\pi}{12}) [/math]

El módulo de la tensión normal según la dirección [math]\vec j[/math]:

[math] \vec j · \sigma · \vec j = 1/3\sin (\frac{y·\pi}{12}) [/math]

El módulo de la tensión normal según la dirección [math]\vec k[/math]:

[math] \vec k · \sigma · \vec k = 1/3\sin (\frac{y·\pi}{12}) [/math]

Figura 10: tension normal en la dirección i

Figura 11: tension normal en la dirección j

Figura 12: tension normal en la dirección k

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se definen las tensiones normales en la dirección que marcan los ejes i, j
%y k
Ti=sin((pi/12)*Y);
Tj=(1/3)*sin((pi/12)*Y);
Tk=(1/3)*sin((pi/12)*Y);
%Se dibuja la tensión normal en la dirección que marca el eje i
subplot(1,2,1)
surf(X,Y,Ti)
shading flat
colorbar
title('Tensiones normales dirección i')
subplot(1,2,2)
surf(X,Y,Ti)
shading flat
colorbar
view(2)
%Se dibuja la tensión normal en la dirección que marca el eje j
figure
subplot(1,2,1)
surf(X,Y,Tj)
shading flat
colorbar
title('Tensiones normales dirección j')
subplot(1,2,2)
surf(X,Y,Tj)
shading flat
colorbar
view(2)
%Se dibuja la tensión normal en la dirección que marca el eje k
figure
subplot(1,2,1)
surf(X,Y,Tk)
shading flat
colorbar
title('Tensiones normales dirección k')
subplot(1,2,2)
surf(X,Y,Tk)
shading flat
colorbar
view(2)

10 Representación de las tensiones normales

En primer lugar calculamos las tensiones tangenciales respecto al plano ortogonal a [math]\vec{i}[/math], es decir [math]|σ ·\vec{i} − (\vec{i} · σ ·\vec{i})\vec{i}|[/math], en t = 0. Después dibujamos sólo las que no son nulas.

[math] \left| \sigma \cdot \vec{i} - (\vec{i} \cdot \sigma \cdot \vec{i}) \vec{i} \right|=\left| \left (\begin{matrix} sin(\frac{πy}{12}) & \frac{πx}{36}cos(\frac{πy}{12}) & 0 \\&\\ \frac{πx}{36}cos(\frac{πy}{12}) & \frac{1}{3}sin(\frac{πy}{12}) & 0 \\&\\ 0 & 0 & \frac{1}{3}sin(\frac{πy}{12})\end{matrix} \right) \cdot \begin{pmatrix}1\\0 \\ 0\end{pmatrix} - sen(\frac{\pi y}{12})\cdot \begin{pmatrix}1\\0 \\ 0\end{pmatrix} \right |= \frac{πx}{36}cos(\frac{πy}{12}) [/math]

Figura 13: tensiones tangenciales en la dirección i en t=0, no nulas

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se define la tensión tangencial respecto al plano ortogonal a i en t=0
Tti=(pi/36)*X.*cos((pi/12)*Y);
%Se dibujan las tensiones tangeciales no nulas
subplot(1,2,1)
surf(X,Y,Tti)
shading flat
subplot(1,2,2)
surf(X,Y,Tti)
view(2)
colorbar
shading flat
title('Tensión tangencial dirección i')

11 Tensión de Von Mises

La tensión de Von Mises se define por la fórmula [math] \sigma_V =\sqrt{\frac{(\sigma_1-\sigma_2)^2+(\sigma_2-\sigma_3)^2+(\sigma_3-\sigma_1)^2}{2}}[/math], donde σ1, σ2 y σ3 son los autovalores de σ (también conocidos como tensiones principales). Se trata de una magnitud escalar que se suele usar como indicador para saber cuando un material inicia un comportamiento plástico (y no elástico puro).En la figura adjunta, se puede diferenciar la distribución de la Tensión de Von Mises tanto en 3D como en 2D, ambas figuras coloreadas cada parte con el color correspondiente asociado a la tensión en ese punto, siendo los colores mas cálidos los que indican una mayor tensión de Von Mises y los mas fríos una menor Tensión de Von Mises, como indica la barra de color adjunta.

Siendo [math]σ=\begin{pmatrix} \sin(\frac{\pi*y}{12}) & 0 & 0\\ 0 & 1/3\sin(\frac{\pi*y}{12}) & 0 \\ 0 & 0 & 1/3\sin(\frac{\pi*y}{12}) \end{pmatrix}[/math]

Esta tensión alcanza su máximo en el punto 0.666667.

Figura 14: Tensión de Von Mises

clc
clear all
close all
h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
Ti=sin((pi/12).*Y);
Tj=(1/3)*sin((pi/12).*Y);
Tk=(1/3)*sin((pi/12).*Y);
%Se crean matrices de unos y ceros
TVM=ones(size(y,2),size(x,2));
M=zeros(3);
%Se nombran los valores de cada componente de la tensión de Von Mises
for i=1:length(y)
for j=1:length(x)
u=Ti(i,j);
v=Tj(i,j);
w=Tk(i,j);
M=[u 0 0; 0 v 0; 0 0 w];
%Se obtienen los autovalores
[p,e]=eig(M);
a1=e(1,1);
a2=e(2,2);
a3=e(3,3);
%Se calcula la tensión de Von Mises
VonMises=sqrt( ((a1-a2)^2+(a2-a3)^2+(a3-a1)^2) * 1/2 );
TVM(i,j)=VonMises;
end
end
%Se dibuja la tensión de Von Mises
subplot(1,2,1)
surf(X,Y,TVM)
shading flat
title('Tensión de Von Mises')
subplot(1,2,2)
surf(X,Y,TVM)
shading flat
view(2)
colorbar
%Se calcula la tensión de Von Mises máxima
TVMmax=max(max(TVM))

12 Campo de fuerzas que actúa sobre la placa

En la placa, actúa un campo de fuerzas [math]\vec{F}[/math], que se aproxima por la ecuación de la elasticidad lineal

[math]\vec{F}=\frac{∂^2\vec{u}}{∂t^2}-∇· σ[/math]

donde [math]∇ · σ[/math] es el campo vectorial que se obtiene al hacer la divergencia del tensor de tensores, [math]σ[/math].

Siendo el vector [math]\vec{u}=\vec{a}·sin(π\cdot k(d·\vec{r0}(x, y)−vt))[/math] y los desplazamientos: [math]\vec{a}=x/3 \vec{i} [/math] , [math]\vec{d}=1/12 \vec{j}[/math] y [math]k=1[/math]; Sustituyendo en el vector [math]\vec u[/math]:

[math]\vec{u}=\frac{x}{3}sin (\frac{π \cdot y}{12}-vt)\vec{i}.[/math]

Como resultado, el tensor deformaciones es:

[math] Ԑ = \; \begin{pmatrix} {\frac{2}{3}sin(\frac{πy}{12}-vt)} & {\frac{πx}{36}cos(\frac{πy}{12}-vt)} & {0}\\ {\frac{πx}{36}cos(\frac{πy}{12}-vt)} & {0} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Por lo tanto el tensor de tensiones y su divergencia son:

[math]σ=λ∇·\vec{u}1 + 2μԐ = \; \begin{pmatrix} {sin(\frac{πy}{12}-vt)} & {\frac{πx}{36}cos(\frac{πy}{12}-vt)} & {0}\\ {\frac{πx}{36}cos(\frac{πy}{12}-vt)} & {\frac{1}{3}sin(\frac{πy}{12}-vt)} & {0}\\ {0} & {0} & {\frac{1}{3}sin(\frac{πy}{12}-vt)}\\ \end{pmatrix} [/math]

[math] ∇· σ= \begin{pmatrix} \ \frac{-xπ^2}{432}sin(\frac{πy}{12}-vt)\vec{i} & \frac{π}{18} \cdot cos(\frac{πy}{12}-vt)\vec{j} & 0\vec{k} \end{pmatrix} [/math]

A continuación se calcula la derivada de segundo grado del vector [math]\vec u[/math] en función de t:

[math]\frac{∂\vec{u}}{∂t}= \frac{-vx}{3}cos (\frac{π \cdot y}{12}-vt)\vec{i}.[/math]

[math]\frac{∂^2\vec{u}}{∂t^2}= \frac{-v^2x}{3} (sin(\frac{π \cdot y}{12}-vt))\vec{i}.[/math]

Se sustituyen los términos calculados en la ecuación de elasticidad lineal con la condición que [math]\vec F[/math]=0;

[math]\vec{F}=\frac{∂^2\vec{u}}{∂t^2}-∇· σ = \frac{-v^2x}{3} (sin(\frac{πy}{12}-vt))\vec{i}. - \begin{pmatrix} \ \frac{-xπ^2}{432}sin(\frac{πy}{12}-vt)\vec{i} & \frac{π}{18}cos(\frac{πy}{12}-vt)\vec{j} & 0\vec{k} \end{pmatrix} =0 [/math]

[math]\frac{π}{18}cos(\frac{πy}{12}-vt)\vec{j}.[/math] =0

Sí y solo sí [math]cos(\frac{πy}{12}-vt) [/math]= 0

Como resultado da que la velocidad de propagación es [math]V= \frac{π}{12}\vec{j}.[/math]

Si tomamos la onda longitudinal en lugar de transversal, es decir, con [math]\vec{a}=1/3 \vec{j} [/math]

El vector [math]\vec{u}[/math] es igual a:

[math]\vec{u}=\frac{1}{3}sin (\frac{π \cdot y}{12}-vt)\vec{j}.[/math]

Siendo, por lo tanto, sus derivadas de primer y segundo orden:

[math]\frac{∂\vec{u}}{∂t}=\frac{-v}{3}cos (\frac{π \cdot y}{12}-vt)\vec{j}.[/math]

[math]\frac{∂^2\vec{u}}{∂t^2}=\frac{-v^2}{3}(sin(\frac{π \cdot y}{12}-vt))\vec{j}.[/math]

Siendo el nuevo tensor deformación: [math] Ԑ = \; \begin{pmatrix} {0} & {0} & {0}\\ {0} & {\frac{π}{18}cos(\frac{πy}{12}-vt)} & {0}\\ {0} & {0} & {0}\\ \end{pmatrix} [/math]

Siendo, en este caso, el tensor de tensiones: [math]σ=λ∇·\vec{u}1 + 2μԐ[/math] [math] \mathbb = \; \begin{pmatrix} {\frac{π}{36}sin(\frac{πy}{12}-vt)} & {0} & {0}\\ {0} & {\frac{5π}{36}cos(\frac{πy}{12}-vt)} & {0}\\ {0} & {0} & {\frac{π}{36}sin(\frac{πy}{12}-vt)}\\ \end{pmatrix} [/math]

Y su divergencia: [math] ∇· σ= \begin{pmatrix} \ 0\vec{i} & \frac{-5π^2}{432}sen(\frac{π*y}{12}-vt)\vec{j} & 0\vec{k} \end{pmatrix} [/math]

Por último al usar la ecuación de la elasticidad lineal de nuevo: [math]\vec{F}=\frac{∂^2\vec{u}}{∂t^2}-∇· σ = \frac{-v^2}{3}(sin(\frac{\pi \cdot y}{12}-vt))\vec{j}. - \begin{pmatrix} \ 0\vec{i} & \frac{-5π^2}{432}sen(\frac{π*y}{12}-vt)\vec{j} & 0\vec{k} \end{pmatrix} [/math]

Se llega a que el resultado de la velocidad de propagación no varía: [math]V= \frac{π}{12}\vec{j}.[/math]

13 Desplazamiento transversal

Fijado ahora el punto (1/2, 1), hemos calculado en primera estancia el módulo del desplazamiento transversal (en la dirección [math]\vec{i}[/math] ) a lo largo del tiempo (correspondiente al intervalo [math]t∈[0, 10][/math] y hemos representado la función que a cada t le asocia dicho desplazamiento. Sustituyendo la velocidad de propagación el vector de desplazamientos en el punto (1/2,1) es :

[math]\vec{u}(x,y)=\frac{x}{3}sen(\frac{π}{12}y-πvt)·\vec{i}[/math] .

Figura 15: Gráfica tiempo-desplazamiento

h=2/10;
x=[-1:h:1];
y=[0:h:12];
[X,Y]=meshgrid(x,y);
%Se define t
t=[0:h:10];
%Se define el desplazamiento en función de t
d=1/6*sin((pi/12)-pi.*t*(pi/12));
%Se dibuja la gráfica Tiempo-Desplazamiento
plot(t,d);
title('Desplazamiento en función del tiempo');
xlabel('Tiempo');
ylabel('Desplazamiento');