Visualización de campos escalares y vectoriales en elasticidad (Grupo C3)
| |
La finalidad del trabajo es el estudio, mediante el uso de operadores diferenciales tales como la divergencia y el rotacional, de los cambios producidos por la temperatura [math]T(x,y,z)[/math], asi como los desplazamientos [math]\vec u(\rho,\theta)[/math], sobre una placa plana (en dimensión 2) que ocupa el anillo comprendido entre las circunferencias centradas en el origen de radios 2 y 1.
De esta forma:
- La temperatura viene definida por la función:
[math]T(x, y, x)=(8−y^2+2y)e^{−x^2+2x−1}[/math]
- Vamos a suponer que la fuerza aplicada sobre la placa ha provocado un desplazamiento de los puntos
de la misma dado por el vector de desplazamientos vienen dados por:
[math]u(ρ, θ) = (1 − ρ)^2\vec g_{\rho}.[/math]
Si definimos [math]r_0(\rho,\theta)[/math] el vector de posición de los puntos de la placa inicial, la posición de cada punto [math](\rho,\theta)[/math] de la placa después de haber aplicado el desplazamiento está determinada por: [math] \vec r (x,y)= \vec r_{0}(x,y)+\vec u(x,y). [/math]
1 Situación inicial
Se toma una placa plana que ocupa el anillo circular centrado en el origen y comprendido entre los radios 1 y 2. Para definir gráficamente dicho sólido, se diseñó el siguiente programa, haciendo uso de las coordenadas polares.
h=0.1; %Paso del muestreo
u=1:h:2; %u [1,2]
v=0:h:2*pi+h; %v [0,2*pi]
[RO,THETA]=meshgrid(u,v); %Se generan las matrices de u y v
figure(1)
X=RO.*cos(THETA); %Parametrizamos el anillo
Y=RO.*sin(THETA);
mesh(X,Y,0*X) %Dibujamos la malla
axis([-3,3,-3,3]) %Definimos los ejes con su dominio correspondiente
axis equal
view(2) %Visualizamos el mallado
2 Comportamiento de la temperatura
h=0.1; %Paso del muestreo
u=1:h:2; %u [1,2]
v=0:h:2*pi+h; %v [0,2*pi]
[RO,THETA]=meshgrid(u,v); %Se generan las matrices de u y v
figure(1)
X=RO.*cos(THETA); %Parametrizamos el anillo
Y=RO.*sin(THETA);
T=((8-Y.^2+2.*Y).*exp(-X.^2+2.*X-1)); %Definimos el campo escalar
surf(X,Y,T)
axis([-3,3,-3,3]) ; %Definimos los ejes con su dominio correspondiente
axis equal
view(2) %Visualizamos el mallado
El resultado obtenido es el gráfico que se puede observar a continuación:
Podemos ver que las temperaturas más elevadas se sitúan en la zona del primer cuadrante, prácticamente entorno al punto (1,1). Sin embargo, este aspecto se puede tratar con más precisión estudiando el gradiente de temperaturas, que trataremos a continuación.
- Variación de la temperatura
El gradiente de temperaturas indica la dirección de máxima variación de la misma en cada punto, es decir, en que dirección (y sentido) aumenta en mayor proporción. La expresión del gradiente queda como:
[math]∇T= \frac{∂T}{∂x} i ⃗+ \frac{∂T}{∂y}j ⃗+ \frac{∂T}{∂z} \vec k=(8-y^2+ 2y)(-2x+2) e^{-x^2+2x-1}\vec i+ (-2y+2)e^{-x^2+2x-1}j ⃗ [/math]
Con la ayuda de MATLAB, podemos visualizar gráficamente el gradiente:
h=0.1; % sampling step
u=1:h:2; % sampling of the interval [1,2]
v=0:h:2*pi+h; % sampling of the interval [0,2*pi]
[RO,THETA]=meshgrid(u,v); % matrixes of ro and teta coordinates
figure(1)
X=RO.*cos(THETA); % parametrization
Y=RO.*sin(THETA);
T=((8-Y.^2+2.*Y).*exp(-X.^2+2.*X-1)); % The scalar field
contour(X,Y,T) % Draw the level sets
hold on
% El gradiente de T viene dado por el campo -exp(-y) según el versor j
Tx=(8-Y.^2+2.*Y).*(-2.*X+2).*exp(-X.^2+2.*X-1); % x-component of the vector field
Ty=(-2.*Y+2).*exp(-X.^2+2.*X-1); % y-component of the vector field
quiver(X,Y,Tx,Ty) % Draw the vector field
axis([-3,3,-3,3]) % select region for drawing
view(2) % See the picture from the top
colorbar
Tal y como vemos en la gráfica, los vectores gradiente asociados a cada punto del sólido tienen una dirección ortogonal a las líneas equipotenciales, es decir, aquellas que unen puntos con igual temperatura. Por otro lado, en las proximidades de la zona más caliente del anillo convergen los vectores, ya que el gradiente se anula. Esto se debe a que en el máximo, el valor de la temperatura no crece.
3 Deformaciones
- Campo desplazamientos \vec u
[math]u ⃗= (1- ρ)^2 g ⃗_ρ [/math]
Recordemos que:[math] g ⃗_ρ=cos(θ)•i ⃗+sin(θ)•j ⃗+0•k ⃗ [/math]
h=0.1; % paso de muestreo
u=1:h:2; % u [1,2]
v=0:h:2*pi+h; % v [0,2*pi]
[RO,THETA]=meshgrid(u,v); % matrices de las coordenadas RO y THETA
figure(1)
X=RO.*cos(THETA); % parametrización en polares
Y=RO.*sin(THETA);
hold on
mesh(X,Y,0*X) %Dibujamos la malla
axis([-3,3,-3,3]) %Definimos los ejes que se observarán en la gráfica
axis equal
view(2) % Vista en planta
fx=((1-RO).^2.*cos(THETA));
fy=((1-RO).^2.*sin(THETA));
quiver(X,Y,fx,fy) % Dibujamos el campo vectorial
axis([-3,3,-3,3]) % Seleccionamos el dominio del dibujo
view(2)
hold off
El campo u ⃗ actúa en dirección radial, ya que sólo tiene componentes según g ⃗_ρ. Este comportamiento se estudiará posteriormente con la divergencia y el rotacional.
Si utilizamos el comando subplot, se puede comparar la figura antes y después de aplicarle el campo u ⃗:
h=0.1; %Paso del muestreo
u=1:h:2; %u [1,2]
v=0:h:2*pi+h; %v [0,2*pi]
[RO,THETA]=meshgrid(u,v); %Se generan las matrices de u y v
figure(1)
X=RO.*cos(THETA); %Parametrizamos el anillo
Y=RO.*sin(THETA);
subplot(1,3,1);
mesh(X,Y,0*X) %Dibujamos la malla del anillo sin deformar
axis([-3,3,-3,3]) %Definimos los ejes con su dominio correspondiente
axis equal
view(2)
subplot(1,3,2); %se abre el subplot
fx=(RO.^2-RO+1).*cos(THETA); %tras calcular a mano el r desplazado, se introduce
fy=(RO.^2-RO+1).*sin(THETA);
mesh(fx,fy,0*fx);% Dibujo de malla desplazada
axis([-3,3,-3,3])% Definimos los ejes con su dominio correspondiente
axis equal
view(2);
subplot(1,3,3);
hold on %abrimos hold on para superponer las dos graficos siguientes
mesh(fx,fy,0*fx) %Dibujamos la malla deformada
axis([-3,3,-3,3]) %Definimos los ejes con su dominio correspondiente
axis equal
view(2)
quiver(X,Y,fx,fy) % Draw the vector field
axis([-3,3,-3,3]) % select region for drawing
axis equal
view(2)
hold off
- Divergencia
La divergencia indica la variación de volumen local en el sólido debido al desplazamiento. Ésta viene dada por
[math]∇• u ⃗= \frac{1}{√g} \frac{∂}{∂x^i} (√g u^i )=\frac{1}{p}\frac{∂}{∂ρ}[ρ(1-ρ)^2]= \frac{1}{ρ}+3ρ-4[/math]
h=0.1; %Paso del muestreo
u=1:h:2; %u [1,2]
v=0:h:2*pi+h; %v [0,2*pi]
[RO,THETA]=meshgrid(u,v); %Se generan las matrices de u y v
figure(1)
X=RO.*cos(THETA); %Parametrizamos el anillo
Y=RO.*sin(THETA);
mesh(X,Y,0*X) %Dibujamos la malla
axis equal %Definimos los ejes con su dominio correspondiente
f=(1./RO)-4+3.*RO;% Módulo de las tensiones en la dirección radial
surf(X,Y,f)
view(3)
Cuando el valor de ρ es igual a 1, la divergencia de u ⃗ se anula (y también el campo), tal y como se puede ver en ambos gráficos. Por otro lado, cuando ρ es igual a 2, la divergencia adquiere su valor máximo. Esto indica que el campo de deformación tiene más incidencia en la zona exterior del anillo que no en la interior.
- Rotacional
El rotacional indica la cantidad de giro sobre una superficie y su expresión es la siguiente:
[math]\nabla \times \vec u=\begin{bmatrix} \vec g_{\rho} & \vec g_{\theta} & \vec g_{z} \\ \frac{\delta}{\delta \rho} & \frac{\delta}{\delta \theta} & \frac{\delta}{\delta z} \\ (1- ρ)^2 & \ 0 & 0 \end{bmatrix}=\vec {0}\ [/math]:
La interpretación que podemos dar es que no se produce giro entorno al vector normal de la superficie. Y es que es lógico, ya que, como indicamos anteriormente, el campo u ⃗ viene expresado según g ⃗_p
4 Tensiones
- Tensor de tensiones σ
[math]σ= λ∇•u ⃗ 1+2μϵ[/math]
| |