Visualización de campos escalares y vectoriales en elasticidad G8
 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.
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:
- 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].
- La temperatura del sólido viene dada por [math]T(\rho,\theta)=e^{-y}[/math]. Dibujarla.
- 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.
- 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.
- 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.
- 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.
- Calcular [math]|\nabla \times \vec u|[/math] en todos los puntos del sólido y dibujarlo. ¿Qué puntos sufren un mayor rotacional?
- 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].
- 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.
- 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.
Contenido
1 Inicio
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
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á cte para cada valor de yo.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
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.
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
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
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
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 este caso son iguales a 1, y por tanto el tensor de tensiones es:[math] \sigma_{ij}=\frac{\pi}{10} cos (\pi y) \vec i [/math]
