Visualización de campos escalares y vectoriales sobre semianillo- grupo B4

De MateWiki
Saltar a: navegación, buscar
Trabajo realizado por estudiantes
Título Visualización de campos escalares y vectoriales sobre un semianillo-Grupo 4
Asignatura Teoría de Campos
Curso 2021-22
Autores

Sandra González De Mendoza 3471, Gema González Esteban 3472, Lucia Monge Mansilla 3539

Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

En este artículo veremos cómo distintos campos escalares y vectoriales afectan a una placa plana con forma de semianillo. Dicha placa ocupa la mitad de un anillo circular centrado en el origen y comprendido entre los radios 1 y 2, en el plano y≥0. En ella tenemos definidas dos cantidades físicas: la temperatura [math]T(x,y)[/math], que viene dada por:

[math]T(x,y)=log((x-3)^2+2)[/math]

y los desplazamientos [math]\vec{u}(x, y)[/math] producidos por la acción de un campo ya conocido en nuestro caso:

[math]\vec{u}(ρ,θ)=\frac{ρ-1}{5}sen(θ)\vec{e_{θ}}[/math]

De esta forma, si definimos [math]\vec{r_{0}}(x,y)[/math] el vector de posición de los puntos de la placa antes de la deformación, la posición de cada punto [math](x,y)[/math] de la placa después de la deformación viene dada por:

[math]\vec{r}(x,y)=\vec{r_{0}}(x, y)+\vec{u}(x, y)[/math]

Analizaremos las curvas de nivel de la temperatura, las deformaciones provocadas por el campo de desplazamientos, las tensiones causadas por un tensor de deformaciones, etc. Para poder realizar todo lo comentado anteriormente, debemos llevar a cabo diferentes cálculos respecto al campo vectorial [math]\vec{u}(x, y)[/math] previamente definido; como gradientes, rotacionales y divergencias.

1 Mallado. Visualización de la placa.

Representaremos los puntos del interior del sólido con un mallado, siendo el resultado un semianillo. Tomar los ejes en el rectángulo (x, y) ∈ [−3; 3] × [−1; 3].

%Definimos los intervalos de rho y teta
u=1:0.1:2;
v=0:pi/50:pi;
%Creamos el mallado
[U,V]=meshgrid(u,v);
figure(1);
%Las funciones de X e Y en coordenadas cilíndricas
X=U.*cos(V);
Y=U.*sin(V);
%El mallado queda representado dentro de los intervalos [-3,3]x[-1,3]
mesh(X,Y,0*X)
axis([-3,3,-1,3])
view(2)
mallado

2 Curvas de nivel de la temperatura.

Una vez definida la temperatura,

[math]T(x,y)=log((x-3)^2+2)[/math]

dibujamos con el comando "contour" sus curvas de nivel. Observando la gráfica, podemos comprobar que la temperatura máxima se alcanza en el punto (-2,0) con un valor de 3.1416 ºC.

%Definimos los intervalos de rho y teta
r=1:0.1:2;
t=0:pi/50:pi;
%Creamos el mallado
[R,T]=meshgrid(r,t);
figure(1);
%Las funciones de X e Y en coordenadas cilíndricas
X=R.*cos(T);
Y=R.*sin(T);
%El mallado queda representado dentro de los intervalos [-3,3]x[-1,3]
mesh(X,Y,0*X)
axis([-3,3,-1,3])
view(2)
%Definimos la función de temperatura
f=log10((X-3).^2+2);
subplot(1,2,1)
surf(X,Y,f)
view(2)
axis([-3,3,-1,3])
colorbar
subplot(1,2,2)
contour(X,Y,f,30)
colorbar
axis([-3,3,-1,3])
Tmax=max(max(T));
fprintf('La temperatura máxima alcanzada es de: %1.4f ºC\n' , Tmax)
temperatura

3 Gradiente de la temperatura.

Calculamos el gradiente de la función de temperatura y observamos donde crece más rápido. Dibujaremos el gradiente como campo vectorial, obteniendo un campo perpendicular a las curvas de nivel del apartado anterior. El gradiente es la suma de las derivadas parciales:
[math]Gradiente(T)={\frac{dT}{dX}+\frac{dT}{dY}}[/math]

Calculamos las deriavdas parciales obteniendo:

[math]{\frac{dT}{dX}}={\frac{2(x-3)}{log[(x-3)^2+2]}}[/math]; [math]{\frac{dT}{dY}}=0[/math]
%creamos los intervalos de las variables
u=1:0.1:2;
v=0:pi/50:pi;
%creamos mallado y gráfica
[U,V]=meshgrid(u,v);
figure(1);
X=U.*cos(V);
Y=U.*sin(V);
mesh(U,V,0*U)
axis([-3,3,-1,3])
%Campo vectorial (gradiente): 2(x-3)ln10 * ((x-3)2+2)
%introducimos el código de la función
f=log10((X-3).^2+2);
%derivadas parciales de la función
Tx=((X-3).*2)./(log(10)*(((X-3).^2)+2));
Ty=0;
%dibujamos la gráfica
contour(X,Y,f,30);
view(2)
axis equal
axis ([-3,3,-1,3])
hold on
%representamos el campo vectorial
quiver(X,Y,Tx,Ty);
axis equal
axis ([-3,3,-1,3])
view(2)
hold off
gradiente

4 Campo de desplazamientos en los puntos del sólido.

Este apartado no ha sido realizado debido a que en los datos ya nos dan el campo vectorial [math]\vec{u}(ρ,θ)=\frac{ρ-1}{5}sen(θ)\vec{e_{θ}}[/math] con lo que no es necesario calcularlo de nuevo.

5 Visualización campo de vectores.

A continuación se dibujará el campo [math]\vec{u}(ρ,θ)[/math] sobre los puntos del mallado del sólido. Siendo el campo vectorial:

[math]\vec{u}(ρ,θ)=\frac{ρ-1}{5}sen(θ)(-sen(\vec{i})+cos(\vec{j}))[/math]
%creamos los intervalos de las variables
u=1:0.1:2
v=0:pi/50:pi;
%mallado
[U,V]=meshgrid(u,v); 
%parametrización
X=U.*cos(V);        
Y=U.*sin(V);
FX=((-U+1)./5).*(sin(V).^2); 
FY=((U-1)./5).*sin(V).*cos(V);
%mostramos los vectores en la gráfica
quiver(X,Y,FX,FY)
axis([-3,3,-1,3]);
view(2)
xlabel('Eje x')
ylabel('Eje y')
title('Campo de desplazamientos sobre la placa')
campo vectorial

6 Sólido previo y posterior a los desplazamientos.

Ahora observaremos las deformaciones producidas en nuestra placa debido a la actuación del campo vectorial [math]\vec{u}(ρ,θ)[/math]. Calcularemos las coordenadas de X y de Y, y a su vez, las coordenadas después de la deformación a las que llamaremos UX y UY. En la gráfica representaremos las placas antes y después de la deformación además de compararlas, haciendo uso del comando subplot y el comenaod hold on.

%creamos los intervalos de las variables
h=50; 
u=linspace(1,2,h); 
v=linspace(0,pi,h);  
%mallado
[U,V] = meshgrid(u,v);
x = U.*cos(V);
y = U.*sin(V);
%función
Ux=(-(U+1)./5).*(sin(V).^2); 
Uy=((U-1)./5).*sin(V).*cos(V);
%representación gráfica
%antes de la deformación
subplot(1,3,1)
surf(x,y,0*x);
view(2)
axis equal
axis([-3,3,-1,3])
title('placa antes de la deformación');
%después de la deformación
subplot(1,3,2)
UX=x+Ux;
UY=y+Uy;
surf(UX,UY,0*UX);
view(2)
axis equal
axis([-3,3,-1,3])
title('placa después de la deformación');
%comparación
subplot(1,3,3)
hold on
plot3(x,y,x*0);
plot3(UX,UY,0*UX);
view(2)
axis equal
axis([-3,3,-1,3])
title('comparación');
hold off
desplazamientos

7 Divergencia campo vectorial sobre placa.

La divergencia define la diferencia de volumen debido al desplazamiento, en nuestro caso, al estar trabajando en 2D, nos muestra el incremento de área. La divergencia se define como:

[math]Divergencia(T)={\frac{dTx}{dX}+\frac{dTy}{dY}+\frac{dTz}{dz}}[/math]
%creamos los intervalos de las variables
h=50; 
u=linspace(1,2,h); 
v=linspace(0,pi,h);  
%mallado
[U,V] = meshgrid(u,v);
x = U.*cos(V);
y = U.*sin(V);
%divergencia
DIVu=(1./U).*(((U-1)./5).*cos(V));
%gráfica
surf(x,y,DIVu)
view(2)
axis equal
axis([-3,3,-1,3])
colorbar
%divergencias máximas y mínimas
DIVuMAX=max(DIVu);
DIVuMIN=min(DIVu);
title('Divergencia del campo vectorial')
divergencia

8 Rotacional .

El operador rotacional muestra la tendencia del campo vectorial a inducir rotación de un punto sobre un vector. Como en dos dimensiones no se aprecia este cambio, hemos decidido realizarlo también en 3-D.

[math]|\nabla×\vec u(ρ,θ)| =\frac{(2ρ-1)sen(θ)}{5-ρ}[/math]
%creamos los intervalos de las variables
h=50; 
u=linspace(1,2,h); 
v=linspace(0,pi,h);  
%mallado
[U,V] = meshgrid(u,v);
x = U.*cos(V);
y = U.*sin(V);
%añadimos el rot
ROTu=abs((sin(V)./(5-U)).*(2*U-1));
%gráfica en 2-D
subplot(1,2,1);
surf(x,y,ROTu);
axis equal
view(2);
colorbar;
title('rotacional en 2-D');
%gráfica en 3-D
subplot(1,2,2);
surf(x,y,ROTu);
colorbar
title('rotacional en 3-D');
Maximo=max(max(ROTu));
rotacional

9 Tensor de tensiones. Tensiones normales.

Ahora hallaremos las tensiones a las que se ve sometida la placa en un medio elástico lineal, isótropo y homogéneo donde los desplazamientos permiten escribir el tensor de la siguiente forma:

[math]σ(x,y)=λ∇.(\vec u)1 + 2με[/math]

Sabiendo que ε es la parte simétrica del ∇\(\vec u\) y conociendo λ=μ=1 (debido a las propiedades del material), calculamos las tensiones en las direcciones de los ejes y las representamos en el dibujo.

%Definimos las variables 
u=1:0.1:2;
v=0:pi/50:pi;
%mallado
[U,V]=meshgrid(u,v); 
X=U.*cos(V);        
Y=U.*sin(V);
Z=zeros(size(X)); 
clf
%calculamos las componentes de la matriz sigma 
a=(((U-1).*cos(V))./(5.*U)); %Elemento (1,1,1) de la matriz sigma
b=((((U-1)./5.*U)+((2.*U-2)./5)).*cos(V)); %Elemento (2,2,2) de la matriz sigma
c=((U-1).*cos(V))./(5.*U); %Elemento (3,3,3) de la matriz sigma
subplot(1,3,1);
surf(X,Y,a)
axis equal
view(2)
colorbar
subplot(1,3,2);  
surf(X,Y,b)
axis equal
view(2)
colorbar
subplot(1,3,3)
surf(X,Y,c)
axis equal
view(2)
colorbar
tensiones normales

10 Tensión tangencial respecto al plano ortogonal X=0

Calculamos las tensiones tangenciales, que vienen dadas por:

[math]|σ(\vec e_{ρ})-((\vec e_{ρ})σ(\vec e_{ρ}))(\vec e_{ρ})|; [/math] siendo esto igual a 0.

11 Tensión de Von Mises.

Conociendo que la tensión de Von Mises se define como:

Tensión Von Mises

donde σ1, σ2 y σ3 son los autovalores de σ (también conocidos como tensiones principales). Se trata de una magnitud escalar que se suele usar como indicador para saber cuando un material inicia un comportamiento plástico (y no elástico puro). Pintar la tensión de Von Mises y señalar en qué punto se alcanza el mayor valor.

h=0.1;
u=[1:h:2];
v=[0:h:pi];
[U,V]=meshgrid(U,V);
%Creación de la matriz sigma
s=zeros(3,3);
%Creamos la matriz de von mises 
vm=zeros(32,11); %Tiene un tamaño 32x11 porque es el tamaño del mallado
for i=1:32*11
UU=U(i)';
VV=V(i)';
s(1,1)=(((UU-1).*cos(VV))./(5.*UU));
s(1,2)=(((-UU+2)./5).*sin(VV));
s(2,1)=s(1,2);
s(2,2)=((((UU-1)./(5.*UU))+((2.*UU-2)./5)).*cos(VV));
s(3,3)=((UU-1).*cos(VV))./(5.*UU);
  [v,d]=eig(s);
  vm(i)=sqrt(((d(1,1)-d(2,2))^2+(d(2,2)-d(3,3))^2+(d(3,3)-d(1,1)).^2)/2);
endfor
xx=U.*cos(V);
yy=U.*sin(V);
subplot(1,2,1)
surf(xx,yy,vm)
colorbar
subplot(1,2,2)
surf(xx,yy,vm)
colorbar
view(2)
Tensión VM en placa

12 Campo de fuerzas.

Sobre nuestra placa actúa un campo de fuerzas [math]\vec{F}[/math], que equivale aproximadamente a la ecuación de elasticidad, lo que produce un desplazamiento, como podemos observar en la gráfica. Definimos este campo de fuerzas como [math]\vec{F}=-▽σ[/math]. Siendo σ los vectores que componen las filas de la matriz sigma.

%Creamos los intervalos de las variables
r=1:0.1:2;
t=0:pi/50:pi;
%hacemos el mallado
[R,T]=meshgrid(r,t); 
X=R.*cos(T);        
Y=R.*sin(T);
%Definimos las funciones x e y del campo F
mx=((-R+3)./(5.*R)).*cos(T);
my=(((-4.*R)+4)./(5.*R)).*sin(T)-((R-1)./(5.*(R.^2))).*sin(T);
%dibujamos la gráfica de vectores
quiver(X,Y,mx,my);
axis([-3,3,-1,3]);
view(2);
title('Campo de fuerzas F')
Campo de fuerzas