Visualizacion de campos escalares y vectoriales en elasticidad (Grupo A - 12)
Para el correcto estudio de la visualización de campos escalares y vectoriales, debemos suponer una placa rectangular plana que, en este caso, ocupará la región [-1/2, 1/2]x[0,2]. Así mismo, definimos sobre ella dos magnitudes físicas, la temperatura,T(x,y,t) un campo escalar; y los desplazamientos,[math]\vec u(x,y,t)[/math],un campo vectorial.
Partiendo de [math]\vec r_{0}(x,y)[/math], vector de posición de los puntos de la placa en reposo, la posición de cada punto de la placa en un instante t vendrá dado 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 transversal:
[math]\vec u(x,y)=1/10 sin(πy)\vec i [/math]
Contenido
- 1 Representación del mallado del sólido
- 2 Representación de un campo escalar
- 3 Cálculo y representación del gradiente y curvas de nivel de una función escalar.
- 4 Representación de campos vectoriales
- 5 Representación de cada punto del sólido antes y despúes del desplazamiento determinado por el campo vectorial [math]\vec u[/math].
- 6 Calculo y representación de la divergencia en todos los puntos del sólido
- 7 Cálculo del rotacional
- 8 Tensor de tensiones
- 9 Tensiones tangenciales respecto al plano ortogonal a [math]\vec i[/math].
- 10 Tensiones tangenciales respecto al plano ortogonal a [math]\vec j[/math].
1 Representación del mallado del sólido
Para representar la placa antes mencionada, creamos dos vectores con el comando linspace, que representan los intervalos de las variables x e y donde la placa está definida, en este caso, tomaremos un muestreo 1/10. Creamos la malla donde la placa está definida, con el comando meshgrid, y la representamos mediante el comando surf tomando como Z una matriz de ceros, ya que la placa es plana.
x=linspace(-1/2,1/2,10); %Límites de la placa en el eje X
y=linspace(0,2,20); %Límites de la placa en el eje Y
[X,Y]=meshgrid(x,y); %Mallado de la placa
Z=zeros(20,10); %Altura de cada punto de la placa
surf(X,Y,Z) %Representación de la placa
axis([-1,1,-1,3]) %Límites de representación
view(2) %Punto de vista vertical
2 Representación de un campo escalar
En este caso, un campo escalar adecuado es la temperatura de la placa, la cual puede ser definida por la función:
[math] T(ρ; θ) = -log(ρ + 0.1)[/math]
Podremos traducirla a Matlab mediante la función “inline”. Como la función depende de variables cilíndricas, la podemos escribir directamente en coordenadas cartesianas con los cambios de variables adecuados. Aplicando esta función a las matrices de coordenadas de la placa, obtenemos Z1, la matriz imagen de la función T. Por último, representamos la superficie con el comando “surf” y la visualizamos.
f = inline(' -log( sqrt( x.^2 + y.^2 ) + 0.1)','x','y'); %función de temperatura
Z1 = f( X , Y ); %Matriz dla temperatura de los puntos de la placa
surf( X , Y , Z1) %Representación de la temperatura en la placa
axis( [ -1,1 , -1,3 ]) %Límites de representación
view(2) %Punto de vista vertical
3 Cálculo y representación del gradiente y curvas de nivel de una función escalar.
Una curva de nivel es aquella línea que en una figura, une todos los puntos que tienen igualdad de altura, para representarlas, se usará el comando “contour”, en el cuál se usarán las variables X, Y, Z1 definidas anteriormente, y por último se añade el número de línes de nivel que se quieren representar.
Contour( X , Y , Z , 20 ) %Representación de las líneas de nivel
axis([ -1,1 , -1,3 ]) %Límites de representación
view(2) %Punto de vista vertical
Procedemos pues al cálculo y representación del gradiente, para ello se han de calcular las derivadas parciales con respecto a X e Y de la función T, en este caso, serán las funciones gx, gy, las cuales son traducidas a matlab con “inline”. Aplicando las matrices del mallado a estas funciones, se obtienen las matrices de las coordenadas de el campo vectorial gradiente de T, las cuales pueden ser representadas en el recinto de la placa, (limitado por las matrices del mallado) con el comando “quiver3”.
gx = inline('-(x ./ (x.^2 + y.^2 + 0.1*sqrt(x.^2 + y.^2)))','x','y'); %Derivada parcial de X de la función temperatura
gx = inline('-(x ./ (x.^2 + y.^2 + 0.1*sqrt(x.^2 + y.^2)))','x','y'); %Derivada parcial de Y de la función temperatura
U = gx( X , Y ); %Componentes X del gradiente de la temperatura
V = gy( X , Y ); %Componentes Y del gradiente de la temperatura
W = zeros( size( V )); %Componentes Z del gradiente de la temperatura
quiver3( X , Y , Z , U , V , W ) %Representación del gradiente de la temperatura
axis([ -1,1 , -1,3 ]) %Límites de representación
view(2) %Punto de vista vertical
Podemos observar finalmente que las curvas de nivel son ortogonales al campo vectorial gradiente.
4 Representación de campos vectoriales
Ahora, se va 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 transversal:
[math]\vec u(x,y)=1/10 sin(πy)\vec i [/math]
En este caso, el campo vectorial está solo definido en la componente [math]\vec i[/math]. Por tanto se pasa la función de esta componente a matlab, nombrándola, en este caso como fx. Así mismo fy será la componente [math]\vec j[/math], en este caso no está definida, luego será 0. Damos valores a las funciones dentro de la región de la placa, obteniendo U1 y V1 y finalmente representamos con “quiver” en la region que contiene la placa.
fx = inline(' sin( pi*y ) / 10' , 'x' , 'y'); %Función de la componente X del campo U
fy = inline('0' , 'x' , 'y'); %Función de la componente Y del campo U
U1 = fx( X , Y ); %Componentes X del campo vectorial U
V1 = fy( X , Y ); %Componentes Y del campo vectorial U
quiver( X , Y , U1 , V1 ) %Representación del campo vectorial U
axis([ -1,1 , -1,3 ]) %Límites de representación
view(2) %Punto de vista vertical
5 Representación de cada punto del sólido antes y despúes del desplazamiento determinado por el campo vectorial [math]\vec u[/math].
Al aplicar el campo vectorial U, antes definido, a la placa ejemplo, se obtiene una segunda placa, la cual es modificada por el campo vectorial. Dado que el campo vectorial ejemplo U sólo tiene componente en x, la unica coordenada modificada de los puntos de la placa será la componente en x, por ello, llamaremos X1 a la componente en x después del desplazamiento. La obtenemos sumando la X inicial y el desplazamiento. La Y y la Z son las mismas, usamos el comando surf con las variables X1, Y, Z.
X1=X+U1; %Componente X modificada
surf(X1,Y,Z) %Representación de la placa modificada
axis([-1,1,-1,3]) %Límites de representación
view(2) %Punto de vista vertical
6 Calculo y representación de la divergencia en todos los puntos del sólido
En coordenadas cartesianas, podemos definir la divergencia como:
[math]\nabla\cdot\vec u= \frac{ \partial \vec u_1 }{ \partial x} + \frac{ \partial \vec u_2 }{ \partial y}+ \frac{ \partial \vec u_3 }{ \partial z}=0
[/math]
La divergencia de un campo vectorial es una medida del cambio de volumen local debido al desplazamiento, en este caso, la divergencia es 0, por tanto la representación es 0 en todos los puntos.
7 Cálculo del rotacional
En coordenadas cartesianas, se puede definir el rotacional como:
[math]\nabla\times\vec{u}= \begin{pmatrix} \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{pmatrix}= \begin{pmatrix} \vec{i} & \vec{j} & \vec{k} \\ \frac{ \partial}{\partial x} & \frac{\partial }{\partial y} & \frac{\partial }{\partial z}\\ \frac{sen(πy)}{10}& 0 & 0 \end{pmatrix}=- \frac{π}{10} cos(πy) \vec k [/math]
Para representar el rotacional se ha de hallar su módulo, en este caso el módulo del campo vectorial U sólo tiene componente en z, luego se aplica el valor absoluto:
[math]|\nabla\times\vec{u}|=| \frac{π}{10} cos(πy) |[/math]
Introducimos con el comando “inline” la componente Z del rotacional, y la aplicamos a los puntos del sólido. Los puntos que maximizan una función son los que hacen máxima la función coseno, es decir, puntos de los ejes y iguales a los números enteros en este caso, como se puede comprobar en el dibujo.
hz=inline('-pi*cos(pi*y)/10','x','y','z'); %Función del rotacional en cada punto de la placa
Z2=abs(hz(X,Y,Z)); %Módulo del rotacional
surf(X,Y,Z2); %Representación del módulo del rotacional en cada punto de la placa
axis([-1,1,-1,3]) %Límites de representación
view(2) %Punto de vista vertical
8 Tensor de tensiones
Definimos la función épsilon como la parte simétrica del tensor gradiente de [math]\vec u[/math]. El tensor de tensiones se define como [math]σ = λ* \nabla.U*δ+2*μ*ε[/math] ( tomando λ y μ = 1 )
9 Tensiones tangenciales respecto al plano ortogonal a [math]\vec i[/math].
Para hallar las tensiones tangenciales al plano ortogonal a [math]\vec i[/math], usamos la fórmula [math]|\sigma \cdot \vec i-(\vec i \cdot \sigma \cdot \vec i) \vec i|[/math] para calcularlas, saliendo con solo componente en y .Construimos las matrices Kx Ky y Kz con las componentes de las tensiones, y después las dibujamos usando el comando quiver3 para representar campos vectoriales en 3 dimensiones.
Kx=(pi/10)*cos(pi*Y)/10; %Componente Y de las tensiones tangenciales del eje i
Ky=zeros(size(Kx)); %Componente X de las tensiones tangenciales del eje i
Kz=zeros(size(Kx)); %Componente Z de las tensiones tangenciales del eje i
quiver3(X,Y,Z,Kx,Ky,Kz); %Representación del campo de las tensiones tangenciales del eje i
axis([-1,1,-1,3]) %Límites de representación
view(2) %Punto de vista vertical
10 Tensiones tangenciales respecto al plano ortogonal a [math]\vec j[/math].
Para obtener y dibujar estas tensiones, las tangenciales respecto al plano ortogonal a [math]\vec j[/math], el método es igual al de las tensiones tangenciales respeto al plano ortognal [math]\vec i[/math], solo que usaremos para hallarlas la fórmula [math]|\sigma \cdot \vec j-(\vec j \cdot \sigma \cdot \vec j) \vec j|[/math]. Luego para representarlas lo haremos igual también, construyendo matrices con las componentes de las tensiones y representándolas sobre los punos de coordenadas X Y Z con el comando quiver3.
Ky=(pi/10)*cos(pi*Y)/10; %Componente Y de las tensiones tangenciales del eje j
Kx=zeros(size(Ky)); %Componente X de las tensiones tangenciales del eje j
Kz=zeros(size(Ky)); %Componente Z de las tensiones tangenciales del eje j
quiver3(X,Y,Z,Kx,Ky,Kz); %Representación del campo de las tensiones tangenciales del eje j
axis([-1,1,-1,3]) %Límites de representación
view(2) %Punto de vista vertical
Comparando las dos figuras anteriores se puede observar que hay una deformación máxima común en y=2; sin embargo, las deformaciones mínimas son diferentes: en las tensiones del eje x están en y=0; y=1,5 mientras que en las tensiones del eje y, el mínimo esta en y=1.
