Visualizacion de campos escalares y vectoriales (grupo 3A)

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.

Marta de Castro Pérez/ Alejandra García-Page Acevedo/ Silvia Pinedo Gil/ Ana María Ragolta Villarroya/

1 Enunciado

Consideramos una placa rectangular 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 transversales. Supondremos lo siguiente: [math] \vec a=\frac{\vec i}{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 i[/math].

2 Mallado

x=-0.5:0.1:0.5;       % Creamos los intervalos x= [-1/2,1/2] y=[0,2]
y=0:0.1:2;             
[xx,yy]=meshgrid(x,y); %Creamos la matriz a partir del intervalo 
mesh(xx,yy,0*xx)       % Dibujamos malla
axis([-2,2,-1,3])      % seleccionamos región

Figura 1 Mallado de la placa

3 Temperatura del sólido

Como ya hemos definido previamente las variables x e y creamos directamente el campo temperatura ([math]T(\rho,\theta)=-\log(\rho+0.1)[/math] ) a partir de esas variables. Sin embargo esta función viene dada en coordenadas polares por lo que realizamos el cambio de variable : [math]\rho = x^2 + y^2[/math]

T=-log(sqrt(xx.^2+yy.^2)+0.1); 
surf(xx,yy,f)          % Dibujamos la superficie 
axis([-2,2,-1,3])      % Seleccionamos la región

En el gráfico las altas temperaturas aparecen en rojo y las más bajas en azul. El foco está situado en el origen y por eso es la zona más caliente.

Figura 2 Distribución de la temperatura en la malla

Figura 3 Distribución de la temperatura en la malla vista desde arriba

4 Gradiente y curvas del campo

Calculamos el gradiente del campo derivando respecto de las dos variables de las que depende

%Como ya tenemos definida la malla y el campo, calculamos directamente el gradiente
Tx=-xx./((sqrt(xx.^2+yy.^2)).*(sqrt(xx.^2+yy.^2)+0.1));
Ty=-yy./((sqrt(xx.^2+yy.^2)).*(sqrt(xx.^2+yy.^2)+0.1));
%Dibujamos el Vector gradiente (figura 3)
quiver(xx,yy,Tx,Ty)
hold on 
contour(xx,yy,T)       % lineas de nivel

plot(x,x-x,'m','linewidth',1);
plot(x,2+x-x,'m','linewidth',1);
plot(-0.5+y-y,y,'m','linewidth',1);
plot(0.5+y-y,y,'m','linewidth',1);
axis([-1,1,-0.5,2.5]) %Numeramos los ejes

Figura 4 Representación del gradiente y las curvas de nivel del campo

5 Campo de vectores adicional

Consideremos ahora el campo de vectores [math]\vec u(x,y)=\frac{\sin(\pi y)}{10}\vec i[/math] y lo dibujamos en los puntos del mallado del sólido para determinar el desplazamiento de los puntos del mallado.

%U(xx,yy)=(0.1*sin(pi.*yy))i/10 + 0.*uu j
%definimos las componentes del campo U
Ux=(0.1*sin(pi.*yy))/10;
Uy=0.*uu;
quiver(xx,yy,Ux,Uy)
axis([-1,1,-0.5,2.5])
view([0,0,1])

Figura 4 Campo de vectores en los puntos del mallado del sólido

6 Desplazamiento sobre el mallado

A continuación consideramos el vector [math]\vec u(x,y)=\frac{\sin(\pi y)}{10}\vec i[/math] como un vector de desplazamiento sobre la malla y comparamos con la figura 1

[xx,yy]=meshgrid(x,y); %dibujamos la malla inicial 
subplot(1,2,1)
mesh(xx,yy,0*xx)
axis([-1,1,-0.5,2.5]) %malla con el vector de desplazamiento aplicado
subplot(1,2,2)
Ux=0.1*sin(pi*yy);
Uy=0*uu;
mesh(xx+Ux,yy+Uy,0*xx)
axis([-1,1,-0.5,2.5]) %Numeración de ejes

Figura 5 Comparación de la malla sin y con el campo de desplazamiento

7 Cálculo de la divergencia y el rotacional

La divergencia resulta nula puesto que el vector u sólo tiene componente en x, y esa componente depende de la variable y. La divergencia representa en este caso el incremento de volumen en cada punto del sólido bajo la acción del campo vectorial [math]\vec u(x,y)=\frac{\sin(\pi y)}{10}\vec i[/math] .

En nuestra malla la divergencia es nula porque al tratarse de un plano todos los puntos se mueven de la misma forma,lo que hace que el cambio (desplazamiento) sea igual en todos los puntos.

El caso del rotacional es distinto. El rotacional de un vector representa la tendencia de un campo vectorial a rotar alrededor de un punto o eje.

R=abs((0.1*pi).*cos(pi.*yy)); %función rotacional de U
surf(xx,yy,R)
axis([-1,1,-0.5,2.5]) %numeramos los ejes
view([0,0,1])
colorbar %añadimos la leyenda de colores

Figura 6 Representación del rotacional

8 Tensor de deformaciones

Definimos la parte simétrica del tensor gradiente de [math]\vec u[/math], que se denomina tensor de defrmaciones : [math]\epsilon(\vec u)=(\nabla \vec u + \nabla \vec u^t)/2[/math] Todo esto en un medio elástico lineal, isótropo y homogéneo donde 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], dibujamos 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].

Las tensiones sobre un sólido vienen definidas por un tensor que depende de los coeficientes de Lamé (que son dos constantes elásticas que caracterizan por completo el comportamiento elástico de un sólido isótropo para pequeñas deformaciones). Aplicado a la placa rectangular puesto que tiene una divergencia nula, por lo explicado anteriormente, dependerá solo de Eij siendo ésta la parte antisimétrica del tensor gradiente del desplazamiento. Se puede observar analíticamente que las componentes de sigma son las mismas que las de E. Para hallar las tensiones normales en la dirección del eje i ,aplicamos i*E*i resultando cero pues es la proyección del vector E*i al solo depender de la coordenada j sobre la dirección i es 0. Para hallar las tensiones normales en la dirección del eje j aplicamos j*E*j también es cero, puesto que el vector E*j solo depende de la componente i por lo que la proyección sobre el eje j es nula. Las tensiones tangenciales (componentes del vector tensión rasantes al plano de corte) vienen definidas por la ecuación |σ⋅i→−(i→⋅σ⋅i→)i→| si son respecto al plano ortogonal a i y por la ecuación |σ⋅j→−(j→⋅σ⋅j→)j→| respecto al plano ortogonal a j. Dado que en nuestro caso |σ⋅i→−(i→⋅σ⋅i→)i→| el segundo sumando es nulo ,las tensiones tangenciales quedan de la forma σ*i originando la siguiente gráfica:

subplot(1,2,1)
quiver(xx,yy,xx*0,Uy)
axis([-1,1,-0.5,2.5])
view([0,0,1])
subplot(1,2,2) 
R=abs((0.1*pi).*cos(pi*yy));
surf(xx,yy,R)
axis([-1,1,-0.5,2.5])
colorbar