Visualización de campos escalares y vectoriales en elasticidad G8

De MateWiki
Saltar a: navegación, buscar
Warning.png Este artículo está en versión beta. El autor de este artículo no lo ha terminado todavía, por favor no lo edites hasta que elimine este mensaje.


Visualización de campos escalares y vectoriales en elasticidad.


Mallado de una placa plana rectangular

Consideramos una placa plana (en dimensión 2) que ocupa la región [math] [-1/2,1/2] \times [0,2][/math]. En ella vamos a suponer que tenemos definidas dos cantidades físicas: la temperatura [math]T(x,y,t)[/math], que depende de las dos variables espaciales [math](x,y)[/math] y el tiempo [math]t[/math], y los desplazamientos [math]\vec u(x,y,t)[/math]. De esta forma, si definimos [math]r_0(x,y)[/math] el vector de posición de los puntos de la placa en reposo, la posición de cada punto [math](x,y)[/math] de la placa en un instante de tiempo [math]t[/math] viene dada por: [math] \vec r (x,y,t)= \vec r_{0}(x,y)+\vec u(x,y,t). [/math] Vamos a suponer que sobre la placa se ha aplicado una fuerza que ha provocado una vibración de manera que los desplazamientos vienen dados por la onda: [math] \vec u(x,y,t) = \vec a \sin(\vec b \cdot \vec r_0-ct), [/math] donde [math]\vec a[/math] se conoce como amplitud, [math]\vec b[/math] es la fase que indica la dirección de propagación y [math]c/|\vec b|[/math] es la velocidad de propagación.

Si [math]\vec a [/math] es paralelo a [math]\vec b[/math] diremos que la onda es longitudinal mientras que si es perpendicular hablaremos de onda transversal.

En este trabajo vamos a centrarnos en las ondas longitudinales. Supondremos lo siguiente: [math] \vec a=\frac{\vec j}{10}, \qquad \vec b= \pi \vec j, \qquad t=0. [/math] En este caso, [math]\vec u(x,y)=\frac{\sin(\pi y)}{10}\vec j[/math].

Se pide:

  1. Dibujar un mallado que represente los puntos interiores del sólido (ver figura). Tomar como paso de muestreo [math]h=1/10[/math] para las variables [math]x[/math] e [math]y[/math].
  2. La temperatura del sólido viene dada por [math]T(\rho,\theta)=e^{-y}[/math]. Dibujarla.
  3. Calcular [math]\nabla T[/math] y dibujarlo como campo vectorial. Calcular las curvas de nivel de T y observar gráficamente que [math]\nabla T[/math] es ortogonal a dichas curvas.
  4. Consideramos ahora el campo de vectores [math]\vec u=\frac{\sin(y)}{10}\vec j[/math]. Dibujar el campo de vectores en los puntos del mallado del sólido.
  5. 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. Dibujar ambos en la misma figura usando el comando subplot.
  6. Calcular la [math]\nabla \cdot \vec u[/math] en todos los puntos del sólido y dibujarla. ¿Qué puntos tienen mayor divergencia? La divergencia es una medida del cambio de volumen local debido al desplazamiento.
  7. Calcular [math]|\nabla \times \vec u|[/math] en todos los puntos del sólido y dibujarlo. ¿Qué puntos sufren un mayor rotacional?
  8. Definamos [math]\epsilon(\vec u)=(\nabla \vec u + \nabla \vec u^t)/2[/math], la parte simétrica del tensor gradiente de [math]\vec u[/math], que se denomina tensor de deformaciones. En un medio elástico lineal, isótropo y homogéneo los desplazamientos permiten escribir el tensor de tensiones [math]\sigma_{ij}[/math] a través de la fórmula:

[math] \sigma_{ij}=\lambda \nabla \cdot \vec u \delta_{ij} + 2\mu \epsilon_{ij}, [/math] donde [math]\lambda[/math] y [math]\mu[/math] son los conocidos como coeficientes de Lamé que dependen de las propiedades elásticas de cada material. Tomando [math]\lambda=\mu=1[/math], dibujar las tensiones normales en la dirección que marca el eje [math]\vec i[/math], es decir [math]\vec i \cdot \sigma \cdot \vec i[/math] y las tensiones normales en la dirección que marca el eje [math]\vec j[/math], es decir [math]\vec j \cdot \sigma \cdot \vec j[/math].

  1. Dibujar las tensiones tangenciales respecto al plano ortogonal a [math]\vec i[/math], es decir [math]|\sigma \cdot \vec i-(\vec i \cdot \sigma \cdot \vec i) \vec i|[/math]. ¿Donde son mayores? Compararlas con los puntos de mayor deformación de la malla.
  2. Dibujar las tensiones tangenciales respecto al plano ortogonal a [math]\vec j[/math], es decir [math]|\sigma \cdot \vec j-(\vec j \cdot \sigma \cdot \vec j) \vec j|[/math]. ¿Donde son mayores? Compararla con los puntos de mayor deformación de la malla.

1 Inicio

Mallado de una placa rectangular de región [math] [-1/2,1/2] \times [0,2][/math]

El objetivo de este trabajo es poder visualizar tanto campos escalares como vectoriales en elasticidad. Pero antes de nada, se debe dejar dos conceptos claros: Un campo escalar es la representación de una distribución espacial de una magnitud escalar, es decir, de un número, asociando un valor a cada punto del espacio. Un campo vectorial es la representación de una distribución espacial de una magnitud vectorial, asociando un vector a cada punto del espacio.

Para poder llevar acabo el objetivo de este trabajo, primero hemos de visualizar la placa rectangular dada, ya que en ella es donde vamos a dibujar los distintos tipos de campos. La mejor forma de visualizar una placa, en nuestro caso rectangular, es a partir del dibujo de un mallado. Esta mallado viene determinado por una región, [math] [-1/2,1/2] \times [0,2][/math], la cual vamos a dibujar introduciendo el siguiente código en Matlab:

x=-0.5:0.1:0.5;
y=0:0.1:2;
[xx,yy]=meshgrid(x,y);
figure(1)
mesh(xx,yy,0*xx);
axis([-2,2,-1,3]);
view(2)


Al ser una placa en 2D debemos abatirla, pues si no la veríamos en 3D.

2 Efecto de la Temperatura

Campo de Temperaturas sobre la placa

La temperatura del sólido viene definida por el campo de temperaturas [math]T(x,y)=e^{-y}[/math]. que es un campo escalar que depende sólo de la variable y, por lo tanto será constante para cada valor de y.Observamos que cuando la variable y disminuye, el valor de la temperatura aumenta, esto se debe a que la función es una exponencial negativa en y .

x=-0.5:0.1:0.5;              % Intervalo [-0.5,0.5]
y=0:0.1:2;                   % Intervalo[0,2]
[xx,yy]=meshgrid(x,y)        % matrices del mallado
figure(1)
f=exp(-yy);                  % Campo escalar
surf(xx,yy,f)                % visualización 
axis([-2,2,-1,3])            % región de visualización
view(2)


3 Divergencia de la función Temperatura

Representación del campo vectorial

Se define la divergencia de la función T, como la suma de las derivadas parciales de dicha función, :[math]\nabla\cdot\vec T = \frac{\partial T_x}{\partial x}+ \frac{\partial T_y}{\partial y}+\frac{\partial T_z}{\partial z}[/math] Por tanto, al calcularla para este caso será: [math]\nabla\cdot\vec T = -e^{-y} [/math] La representación de dicha función como campo vectorial, vendría dada por:

x=-0.5:0.1:0.5;                %Intervalo [-0.5,0.5]
y=0:0.1:2;                     %Intervalo [0,2]
[xx,yy]=meshgrid(x,y);         %Mallado
figure(1)       
fy=-exp(-yy)                   %derivada con respecto de y
fx=0*xx                        %derivdad con respecto de x
quiver(xx,yy,fx,fy)            %visualización del campo en 2D
axis([-2,2,-1,3])              %región de visualización
view(2)


Sean el conjunto de puntos (x,y) donde una función toma un determinado valor; este conjunto es llamado conjunto de nivel, de tal forma que si se representan en el plano OXY se denominan curvas de nivel.

Curvas de nivel
x=-0.5:0.1:0.5;          %Intervalo [-0.5,0.5]
y=0:0.1:2;               %Intervalo [0,2]
[xx,yy]=meshgrid(x,y);   %Mallado
figure(1)
f=-exp(-yy);             %derivada con respecto de y
contour(xx,yy,f)         %representación de las curvas de nivel 
hold on
plot(xx,xx-xx,'k','linewidth',1);
plot(xx,2+xx-xx,'k','linewidth',1);
plot(-0.5+yy-yy,yy,'k','linewidth',1);
plot(0.5+yy-yy,yy,'k','linewidth',1);
axis([-2,2,-1,3])
view(2)


En vista de las figuras adjuntas, se observa que la representación vectorial de [math]\nabla T[/math] es perpendicular a las curvas de nivel.

4 Influencia del desplazamiento

Desplazamiento de los puntos

Ahora consideramos [math]\vec u(x,y)=\frac{\sin(\pi y)}{10}\vec j[/math] como un nuevo campo de vectores para determinar el desplazamiento de los puntos del mallado.

Al ser el vector unitario [math]\vec j [/math] el desplazamiento será en dirección vertical. Por ello se mantiene la coordenada [math] "x" [/math] del punto constante.

Para obtener la gráfica utilizamos el siguiente código en Matlab:

x=-0.5:0.1:0.5;      %genera vectores
 y=0:0.1:2;
 [xx,yy]=meshgrid(x,y);    %genera malla
 figure(3)
 fx=0*xx;
 fy=(sin(pi.*yy))/10
 quiver(xx,yy,fx,fy)
 axis([-2,2,-1,3])
 view(4)

5 Comparación de la malla sin y con desplazamiento

Al introducir los comandos adecuados en Matlab se obtienen dos gráficas continuas, una corresponde a la malla sin desplazamiento y la otra a la malla con desplazamiento. Podemos observar que las diferencias son poco apreciables, pero existen; las diferencias más significativas se encuentran en el centro de la malla donde los puntos han sufrido una máxima compresión, por el contrario los puntos más alejados del centro son los que han sufrido el mayor estiramiento. Por otra parte existen algunos puntos intermedios donde la deformación es nula o muy difícil de apreciar.

comparación de la malla sin y con desplazamiento
x=-0.5:0.1:0.5;
y=0:0.1:2;
[X,Y]=meshgrid(x,y)
subplot(1,2,1)        %subventana 1
mesh(X,Y,0*X);
view(2)
xlabel('eje x')
ylabel('eje y')
title('MALLA SIN DESPLAZAMIENTO')
axis image

DX=0*X;
DY=0.1.*sin(pi*Y);
subplot(1,2,2)        %subventana 2
mesh(X+DX, Y+DY,0*X);
view(2)
xlabel('eje x')
ylabel('eje y')
title('MALLA CON DESPLAZAMIENTO')
axis image


6 Divergencia del Campo

Divergencia 3-D
Divergencia 2-D

En primer lugar, definimos la divergencia de un campo vectorial en un punto como el escalar que resulta de derivar parcialmente cada una de las componentes del campo: [math]\nabla\cdot\vec U = \frac{\partial U_x}{\partial x}+ \frac{\partial U_y}{\partial y}+\frac{\partial U_z}{\partial z}[/math]

Físicamente, la divergencia representa la cantidad de flujo neto que fluye a través de una superficie cerrada; es decir, la divergencia no es más que la diferencia entre el flujo saliente y el entrante en dicha superficie.

En nuestro caso, el valor de la divergencia de nuestro campo es : [math]\nabla\cdot\vec U = \frac{\pi cos (\pi y)}{10}\vec j [/math]

Los comandos de MATLAB necesarios para obtener la representación gráfica de la divergencia de los puntos del sólido tratado son los siguientes:

>> x=-0.5:0.1:0.5;
>> y=0:0.1:2;
>> [xx,yy]=meshgrid(x,y);
>> figure(1)
>> f=(cos(pi.*yy).*pi)/10;
>> surf(xx,yy,f)
>> axis([-2,2,-1,3])
>> view(2)


Además, una correcta representación gráfica de la divergencia nos mostrará la expansión del campo en el recinto escogido. En nuestro caso, obtenemos las figuras, en 3-D y 2-D respectivamente, que se muestran a la derecha.


Partiendo de todo lo anterior, y de que la divergencia nos muestra el cambio de volumen local en los puntos del sólido debido al desplazamiento, podemos concluir que el conjunto de puntos con mayor divergencia será el que haya sufrido mayores variaciones volumétricas, que en nuestro caso se corresponde con el intervalo en Y [0,9;1,2].

7 Rotacional

El rotacional de un campo vectorial, [math]|\nabla \times \vec U|[/math], es un campo vectorial que se define y calcula: [math] \nabla\times \vec U=\frac{1}{\sqrt{g}} \left| \begin{matrix} \vec g_1 & \vec g_2 & \vec g_3 \\ & & \\ \frac{\partial}{\partial u} & \frac{\partial}{\partial v} & \frac{\partial}{\partial w} \\ & & \\ U_1 & U_2 & U_3 \end{matrix}\right| [/math] En nuestro caso, el rotacional sería la siguiente expresión: [math] \nabla\times \vec U=\left| \begin{matrix} \vec i & \vec j & \vec k \\ & & \\ \frac{\partial}{\partial x} & \frac{\partial}{\partial y} & \frac{\partial}{\partial z} \\ & & \\ U_x & U_y & U_z \end{matrix}\right| [/math] Siendo [math] U_x=0 [/math] [math] U_y=e^{-y} [/math] [math]U_z=0[/math] y por tanto el resultado del rotacional es: [math]|\nabla \times \vec U|=0[/math]

El rotacional define una circulación,por tanto, si es nulo significa que el campo vectorial no se mueve, es decir, permanece estático.

8 Tensiones normales y tangenciales

El tensor de deformaciones es aquel definido como [math]\epsilon(\vec u)=(\nabla \vec u + \nabla \vec u^t)/2[/math], es decir, la parte simétrica del tensor gradiente de [math] \vec U [/math].

Como el campo vectorial estudiado se encuentra en un medio elástico, lineal, isótropo y homogéneo, esto permite la descripción de un desplazamiento definido por la fórmula: [math] \sigma_{ij}=\lambda \nabla \cdot \vec u \delta_{ij} + 2\mu \epsilon_{ij}, [/math] donde [math]\lambda[/math] y [math]\mu[/math], conocidos como los coeficientes de Lamé, dependen de las propiedades elásticas de cada material. En nuestro caso, ambos son iguales a 1 (λ=μ=1). A continuación, procedemos al cálculo de las tensiones normales y tangenciales:

Tensores normales

[math] \sigma_{x}=\frac{\pi}{10} cos (\pi y) [/math]
[math] \sigma_{y}=\frac{3 \pi}{10} cos (\pi y) [/math]

Tensores tangenciales

[math]|\sigma \cdot \vec i-(\vec i \cdot \sigma \cdot \vec i) \vec i|=\frac{\pi}{10} cos (\pi y) - \frac{\pi}{10} cos (\pi y) = 0 [/math]
[math]|\sigma \cdot \vec j-(\vec j \cdot \sigma \cdot \vec j) \vec j| =\frac{3 \pi}{10} cos (\pi y) - \frac{3 \pi}{10} cos (\pi y) = 0 [/math]

Con lo anterior, observamos que los tensores tangenciales en las dos direcciones propuestas (respecto al plano ortogonal a OX y respecto al plano ortogonal a OY) son nulos, lo que implica que no existen tensiones tangenciales al campo en dichas direcciones.