Deformaciones de una columna.Grupo 28

De MateWiki
Revisión del 22:48 9 dic 2024 de Rodrigo Prado Fornos (Discusión | contribuciones) (Variación de temperatura)

Saltar a: navegación, buscar
Trabajo realizado por estudiantes
Título Deformaciones de una columna . Grupo 28
Asignatura Teoría de Campos
Curso 2024-25
Autores
  • Mateo Navarro Díaz
  • María Victoria González Junco
  • Fernando Benítez Pérez
  • Rodrigo Prado Fornos
  • Claudia Elimar Manrique
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura


Para el siguiente articulo, consideraremos una placa rectangular plana ocupando la región en el espacio plano
[math](x, y) ∈ [-1, 1]×[0, 10][/math]. A continuación definimos dos cantidades físicas; por un lado la temperatura dada
como

[math] T(\rho, \theta)=\sin(2π\rho)[/math]

en cartesianas como

[math] T(x, y)=\sin(2π\sqrt{x^2 + y^2})[/math]
También, definimos la densidad de la placa como
[math] d(\rho, \theta)=(2-\rho)(4-\cos(4(\rho+π/2)))[/math]

Por otro lado, hemos de tener en cuenta los desplazamientos [math]\vec{u}(\rho, \theta)[/math] producidos por la acción de una fuerza
determinada. Estos desplazamientos se corresponden con el campo [math]\vec{u}(\rho, \theta)= \rho sin(\theta)sin(2π\rho/50)e_\vec{\theta}[/math] .




1 Representación del mallado del sólido

Mallado de la placa

Definimos el mallado que representa el interior del sólido mediante dos vectores x e y los cuales representan el intervalo en el que está definida nuestra placa, tomando un muestreo h=1/10. Al crear la malla, al representarla con el comando surf, tenemos que tener en cuenta que es plana dando valor 0 a la altura de los puntos de la misma.

h=0.1;                       %muestreo
x=-1:h:1;                    %eje x de la placa
y=0:h:10;                     %eje y de la placa
[xx,yy]=meshgrid(x,y);       %mallado de la placa
figure(1)
surf(xx,yy,0*xx);            %representación de la placa
axis([-2,2,0,10])             %ejes del rectángulo
view(2)


2 Gradiente de la temperatura

El cálculo del gradiente de un campo escalar se expresa como la derivada parcial respecto de cada coordenada de dicho campo en función de la base ortonormal orientada positiva [math]\vec{i},\vec{j},\vec{k}[/math] (COORDENADAS CARTESIANAS). Es decir:

[math] \nabla T=\frac{\partial{T}}{\partial{x}}\vec{i}+\frac{\partial{T}}{\partial{y}}\vec{j}[/math]


[math] T=\sin(2π\sqrt{x^2 + y^2})[/math]


[math] \nabla T=\frac{2πx cos(2π\sqrt{x^2 + y^2})}{\sqrt{x^2 + y^2}}\vec{i}+\frac{2πy cos(2π\sqrt{x^2 + y^2})}{\sqrt{x^2 + y^2}}\vec{j}[/math]

A continuación se muestra el código completo de matlab, del cual resultan las siguientes imágenes en las que podemos ver las curvas de nivel y donde se encuentra su máximo y mínimo, a partir de los colores de las gráficas y a partir de el propio matlab. Este nos índica al ejecutar el programa que el máximo es: 0.99996

C.Temperatura C.Nivel
clc
clear
h = 0.2; %MUESTREO
x =-1.0:h:1.0; %DOMINIO DE X
y =0.0:h:10.0; %DOMINIO DE Y
[X,Y]= meshgrid(x,y); % CREACIÓN DEL MALLADO
T=sin(2.*pi.*(sqrt(X.^2+Y.^2))); %TEMPERATURA EN CADA PUNTO DEL MALLADO
mesh(X,Y,0*X) %CREACIÓN DE LA MALLA
subplot(1,3,1) %GRÁFICO SUPERIOR
surf(X,Y,T); %DEGRADACIÓN DE COLORES
axis ([-1,1,-0.5,10.5]) 
view(2)
colorbar %BARRA DE COLORES
%TÍTULO PRIMERA GRÁFICA Y EJES
title('CAMPO DE TEMPERATURAS')
xlabel('eje X')
ylabel('eje Y')
subplot(1,3,2) %GRÁFICO INFERIOR
contour(X,Y,T,11);
axis ([-1.0,1.0,-0.5,10.5]);
%TÍTULO SEGUNDA GRÁFICA Y EJES 
title('CURVAS DE NIVEL')
xlabel('eje X')
ylabel('eje Y')
colorbar %BARRA DE COLORES
Tmax=max(max(T)) %PUNTO MÁXIMO
subplot(1,3,3)
axis ([-1.25,1.25,-0.5,10.5]) 
view(2)
%CURVAS DE NIVEL ENUMERADAS
[M,c]=contour(X,Y,T,[0:0.75:10],'ShowText','on') 
%TÍTULO TERCERA GRÁFICA Y EJES
title('CURVAS DE NIVEL NUMERADAS')
xlabel('eje X')
ylabel('eje Y')
colorbar %BARRA DE COLORES


A continuación, se va a demostrar como el gradiente, previamente calculado, es perpendicular a las líneas de nivel anteriores. Al no verse claro en la imagen original, con la cantidad de líneas de flujo, se opta por aplicar un zoom para una correcta visualización.


Mallado de la placa
clear;clc;
h = 2/10; %MUESTREO
x = [-1:h:1]; %DOMINIO DE X
y = [0:h:12]; %DOMINIO DE Y
[X,Y]= meshgrid(x,y); %CREACIÓN DEL MALLADO
T=sin(2.*pi.*(sqrt(X.^2+Y.^2))); %TEMPERATURA EN CADA PUNTO DEL MALLADO
contour(X,Y,T,11); %CURVAS DE NIVEL
dx=2*pi*X.*cos(2.*pi.*sqrt(X.^2+Y.*2))./sqrt(X.^2+Y.*2); %PARCIAL DE X
dy=2*pi*Y.*cos(2.*pi.*sqrt(X.^2+Y.*2))./sqrt(X.^2+Y.*2); %PARCIAL DE Y
%TÍTULO Y EJES
title('Gradiente de temperatura');
xlabel('Eje X');
ylabel('Eje Y');
% Representación de la temperatura y las curvas de nivel
hold on
quiver(x,y,dx,dy);
axis equal
colorbar


3 Ley de Fourier

Gracias al segundo principio de la termodinámica conocemos el calor emitido por dos focos, uno cálido y otro frío, solo puede ir en un sentido, de mayor a menor calor. Y la ley expresada por el matemático y físico Jean-Baptiste Joseph Fourier trata de establecer una relación entre la distancia, tiempo y el flujo de calor entre focos. Dicha ley viene expresada [math]\vec{Q}=−K∇T[/math], donde K es la constante de conductividad térmica que dependerá del material a estudio. En nuestro caso, la sección del solido tiene por K el valor asociado 1.

Se observa que las flechas van en sentido contrario al gradiente.

Figura 3
clear;clc;
%Definimos las variables x e y usando como paso de muestreo h
 h=2/10;
 x=[-1:h:1];
 y=[0:h:10];
%Creación del mallado
 [X,Y]=meshgrid(x,y);
 mesh(X,Y,0.*X);
%Definimos la función T
 T=sin(2.*pi.*(sqrt(X.^2+Y.^2)));
%Curvas de nivel
 contour(X,Y,T,11);
hold on
%Calculo del gradiente de T
 dx=2*pi*X.*cos(2.*pi.*sqrt(X.^2+Y.*2))./sqrt(X.^2+Y.*2);
 dy=2*pi*X.*cos(2.*pi.*sqrt(X.^2+Y.*2))./sqrt(X.^2+Y.*2);
%Constante de conductividad térmica y ley de Fourier 
 k=1;
  Q1=-k.*(dx);
  Q2=-k.*(dy);
%Representación del campo vectorial
 quiver(x,y,Q1,Q2,'r');
 axis equal;
hold off
%Título y nombre de los ejes
 title('Energía calorífica');
  xlabel('Eje X');
  ylabel('Eje Y');


4 Variación de temperatura

La temperatura [math]T(\rho) = \sin(2\pi \rho)[/math] varía cuando cambia [math]T(\rho)[/math], y su mayor cambio ocurre cuando la derivada [math]\frac{\partial T}{\partial \rho} = 2\pi \cos(2\pi \rho)[/math] alcanza valores máximos o mínimos. Esto sucede cuando [math]\rho[/math] es un número entero, ya que en esos puntos el coseno toma sus valores máximos. La variación es radial porque [math]T[/math] solo depende de [math]\rho[/math]. El gradiente de temperatura apunta hacia afuera o hacia adentro, dependiendo del signo de [math]\cos(2\pi \rho)[/math], lo que indica la dirección de mayor cambio. Por lo tanto, la dirección de máximo crecimiento es perpendicular a las circunferencias de [math]\rho=n[/math], donde 𝑛 es un número entero. Esto es lo que se observa en las flechas ubicadas en la circunferencia de radio [math]\rho = 1[/math].

Variación de temperatura.
% Parámetros en coordenadas cartesianas
x = linspace(-1, 1, 200);
y = linspace(0, 10, 200);
% Crear una malla para la sección
[X, Y] = meshgrid(x, y);
% Conversión a coordenadas cilindricas
R = sqrt(X.^2 + Y.^2); % rho
Theta = atan2(Y, X); % theta (arctan(y/x))
% Temperatura en coordenadas polares
T = sin(2 * pi * R);
% Gradiente de T con respecto a rho
dT_drho = 2 * pi * cos(2 * pi * R);
% Dibujar el mapa de temperatura en coordenadas cartesianas
figure;
contourf(X, Y, T, 20, 'LineColor', 'none'); % Mapa de temperatura
colorbar;
hold on;
% Definir las circunferencias de máxima variación
max_rhos = 1:10;
for rho_val = max_rhos
    num_puntos = 300; 
    theta_circle = linspace(0, 2*pi, num_puntos);
    x_circle = rho_val * cos(theta_circle); 
    y_circle = rho_val * sin(theta_circle);
    valid_idx = (x_circle >= -1 & x_circle <= 1 & y_circle >= 0 & y_circle <= 10);
    x_circle = x_circle(valid_idx);
    y_circle = y_circle(valid_idx);
    
    % Marcar los puntos de máxima variación con rojo
    plot(x_circle, y_circle, 'ro', 'MarkerFaceColor', 'r', 'MarkerSize', 2);
    
    
    if rho_val == 1
        num_flechas = 20; 
        theta_flechas = linspace(0, 2*pi, num_flechas + 1); 
        
        for t = theta_flechas(1:end-1)
            x_start = rho_val * cos(t); 
            y_start = rho_val * sin(t); 
            dx = 0.4 * cos(t);
            dy = 0.4 * sin(t); 
            quiver(x_start, y_start, dx, dy, 'r', 'LineWidth', 1, 'MaxHeadSize', 1.5);
        end
    end
end
title('Temperatura y Dirección de Máxima Variación (Coordenadas Cartesianas)');
xlabel('x');
ylabel('y');
axis equal;
xlim([-1.5, 1.5]);
ylim([0, 10]);
hold off;


5 Campo de desplazamientos

Campo de desplazamientos y puntos fijos.
% Definir la malla
h = 0.3; 
x = -1:h:1;  
y = 0:h:10;  

% Crear la malla en coordenadas cartesianas
[Mx, My] = meshgrid(x, y);

% Convertir a coordenadas cilíndricas
rho = sqrt(Mx.^2 + My.^2);  
theta = atan2(My, Mx);  

% Desplazamiento en coordenadas cilíndricas
u_theta = rho .* sin(theta) .* sin((2 * pi * rho) / 50);  % Componente en theta

% Convertir a coordenadas cartesianas
up = u_theta .* (-sin(theta));  % Componente en x
uth = u_theta .* cos(theta);  % Componente en y

% Graficar el campo de desplazamientos (flechas)
figure;
h_quiver = quiver(Mx, My, up, uth, 1.5, 'Color', [0, 0, 1], 'AutoScale', 'on', 'LineWidth', 0.5);  
h_quiver.AutoScaleFactor = 2.5; 

hold on;

% Graficar la región [−1, 1] × [0, 10]
plot([-1 1], [0 0], 'k', 'LineWidth', 2); 
plot([-1 1], [10 10], 'k', 'LineWidth', 2);  
plot([-1 -1], [0 10], 'k', 'LineWidth', 2);  
plot([1 1], [0 10], 'k', 'LineWidth', 2);  
% Etiquetas y título
axis equal;
axis([-2, 2, -0.5, 12.5]);
xlabel('Eje de las X');
ylabel('Eje de las Y');
title('Campo de Desplazamientos de la Columna');

% Detectar puntos fijos (donde los desplazamientos son cero)
fixed_points = (u_theta == 0);  % Puntos fijos cuando u_theta = 0

% Graficar los puntos fijos (solo aquellos donde u_theta es cero)
h_points = plot(Mx(fixed_points), My(fixed_points), 'ro', 'MarkerFaceColor', 'r', 'MarkerSize', 4); 
legend([h_quiver, h_points], {'Campo de desplazamientos', 'Puntos fijos'}, 'Location', 'best');
grid on;
hold off;


6 Sólido antes y después del desplazamiento

Campo de desplazamientos y puntos fijos.

En esta ocasión representaremos el desplazamiento de la columna definido por el campo vectorial 𝑢(𝜌,𝜃)=𝜌sin(𝜃)sin(2𝜋𝜌/50)𝑒𝜃

Hemos ejecutado el mallado para generar los puntos de dominio, aplicamos la fórmula del campo vectorial a dicho mallado para obtener la nueva posición del sólido,y por último, nos hemos apoyado en subgráficos para comparar el sólido antes y después del desplazamiento.

% Parametrización 
h = 0.1; 
rho = -1:h:1; 
tt = 0:h:10;

% Mallado
[Mrho, Mtt] = meshgrid(rho, tt);

% Campo de desplazamientos 
u_theta = 5 * (Mrho .* sin(Mtt) .* sin(Mrho .* 2 * pi ./ 50)); 
rdrho = Mrho; 
rdtt = Mtt + u_theta;

% Grafica sólido antes del desplazamiento
subplot(1, 2, 1);
mesh(Mrho, Mtt, zeros(size(Mrho))); 
axis equal;
axis([-1.5, 1.5, 0, 10]); 
title('Sólido antes del desplazamiento');
xlabel('(Eje X)');
ylabel('(Eje Y)');
view(2); 
grid on;

% Grafica sólido después del desplazamiento
subplot(1, 2, 2);
mesh(rdrho, rdtt, zeros(size(rdrho))); 
axis equal;
axis([-1.5, 1.5, 0, 10]); 
title('Sólido después del desplazamiento');
xlabel('(Eje X)');
ylabel('(Eje Y)');
view(2); 
grid on;


7 Divergencia del campo de desplazamientos

La divergencia de un campo vectorial es una medida de la cantidad de flujo que sale o entra en un punto de la región definida.


En este caso la divergencia resultante es [math] \nabla\cdot\vec u = cosθsin( \frac{2πρ}{50})\ [/math]

Puntos mínimo y máximo de la divergencia : [math] \frac{\partial \nabla\cdot\vec u}{\partial θ}=0 [/math] ; [math] \frac{\partial \nabla\cdot\vec u}{\partial ρ }=0 [/math]

Mínimo(θ,ρ)=(0,-1) ; Máximo(θ,ρ)=(0,1)

Punto nulo de la divergencia : [math] \nabla\cdot\vec u = 0 [/math] ; si ρ=0 el [math] sin( \frac{2πρ}{50})\ = 0 [/math] lo que anularía el valor de la divergencia.


Divergencia
Divergencia
% Definición de la malla en coordenadas cartesianas
x = -1:0.2:1;         % Coordenadas x
y = 0:0.1:10;         % Coordenadas y
[xx, yy] = meshgrid(x, y);

% Conversión de la malla cartesiana (x, y) a coordenadas cilíndricas (r, theta)
r = sqrt(xx.^2 + yy.^2);         % Radio en coordenadas cilíndricas
theta = atan2(yy, xx);           % Ángulo en coordenadas cilíndricas

% Definir la divergencia en coordenadas cilíndricas
diver = cos(theta) .* sin(2 * pi * r / 50);

% Encontrar los puntos de divergencia máxima, mínima y nula
[minVal, minIdx] = min(diver(:));       % Valor mínimo y su índice lineal
[maxVal, maxIdx] = max(diver(:));       % Valor máximo y su índice lineal
nullIdx = find(abs(diver) < 1e-6);      % Índices donde la divergencia es nula (aproximación)

% Convertir los índices lineales a índices de matriz (para la malla)
[rowMin, colMin] = ind2sub(size(diver), minIdx); % Índices del mínimo
[rowMax, colMax] = ind2sub(size(diver), maxIdx); % Índices del máximo
[rowNull, colNull] = ind2sub(size(diver), nullIdx); % Índices donde la divergencia es nula

% Coordenadas correspondientes en (x, y)
minPoint = [xx(rowMin, colMin), yy(rowMin, colMin)];
maxPoint = [xx(rowMax, colMax), yy(rowMax, colMax)];
nullPoints = [xx(rowNull, colNull), yy(rowNull, colNull)];

% Mostrar los resultados
fprintf('Divergencia mínima:\n Valor: %.3f, Punto: (%.3f, %.3f)\n', minVal, minPoint(1), minPoint(2));
fprintf('Divergencia máxima:\n Valor: %.3f, Punto: (%.3f, %.3f)\n', maxVal, maxPoint(1), maxPoint(2));
fprintf('Divergencia nula en %d puntos:\n', size(nullPoints, 1));
disp(nullPoints);

% Representación de la divergencia
figure;
surf(xx, yy, diver);             % Malla original, pero con valores de divergencia corregidos
shading interp;                  % Suavizado
colorbar;                        % Barra de colores
hold on;

% Marcar los puntos de interés
plot3(minPoint(1), minPoint(2), minVal, 'ro', 'MarkerSize', 8, 'MarkerFaceColor', 'r'); % Mínimo
plot3(maxPoint(1), maxPoint(2), maxVal, 'go', 'MarkerSize', 8, 'MarkerFaceColor', 'g'); % Máximo
for i = 1:size(nullPoints, 1)  % Marcar puntos nulos
    plot3(nullPoints(i, 1), nullPoints(i, 2), 0, 'bo', 'MarkerSize', 6, 'MarkerFaceColor', 'b');
end

% Etiquetas del gráfico
title('Divergencia en coordenadas cilíndricas con puntos de interés');
xlabel('x');
ylabel('y');
zlabel('Divergencia');
legend('Divergencia', 'Mínimo', 'Máximo', 'Nulo');


8 Rotacional del campo

  • [math]\triangledown \times \vec{u}(\rho,\theta,z)=\frac{1}{\rho}\begin{vmatrix} \vec{e}_{\rho} & \rho \vec{e}_{\theta } & \vec{e}_{z}\\ \frac{\partial }{\partial \rho } & \frac{\partial }{\partial \theta } & \frac{\partial }{\partial z}\\ u_{\rho } & \rho u_{\theta } & u_{z} \end{vmatrix} [/math]

Aplicando nuestro campo escalar, obtenemos:

  • [math]\triangledown \times \vec{u}(\rho,\theta,z)=\frac{1}{\rho}\begin{vmatrix} \vec{e}_{\rho} & \rho \vec{e}_{\theta } & \vec{e}_{z}\\ \frac{\partial }{\partial \rho } & \frac{\partial }{\partial \theta } & \frac{\partial }{\partial z}\\ 0 & p^2 sin(\theta)sin( \frac{2πρ}{50}) & 0 \end{vmatrix} = (2sin(\theta)sin( \frac{2πρ}{50}) + psin(\theta)cos( \frac{2πρ}{50})*\frac{2πρ}{50})\vec{e}_{z}\\ [/math]


  • [math]\left \|\triangledown \times \vec{u}(\rho,\theta,z)\right \| = sin(\theta)\sqrt{4sin^2(\frac{2πρ}{50})+p^2cos^2(\frac{2πρ}{50})*\frac{4π^2}{50^2}+4psin(\frac{2πρ}{50})cos(\frac{2πρ}{50})*\frac{2π}{50}} [/math]


Rotacional
% Dimensiones de la sección transversal
x = -1:0.1:1;  % Coordenadas x
y = 0:0.1:10;  % Coordenadas y
[xx, yy] = meshgrid(x, y);

% Convertir coordenadas cartesianas a cilíndricas
r = sqrt(xx.^2 + yy.^2);         % Radio en coordenadas cilíndricas
theta = atan2(yy, xx);           % Ángulo en radianes (atan2 para cuadrante correcto)

% Definir los límites de r y theta
r(r == 0) = eps;  % Evitar división por cero en r

% Cálculo del rotacional
% (nabla x u)_z = 2 * sin(theta) * sin(2*pi*r/50) + (2*pi*r/50) * sin(theta) * cos(2*pi*r/50)
rot_z = 2 .* sin(theta) .* sin(2 * pi * r / 50) + (2 * pi * r / 50) .* sin(theta) .* cos(2 * pi * r / 50);

% Calcular el valor absoluto del rotacional
abs_rot_z = abs(rot_z);

% Representación gráfica del valor absoluto del rotacional
figure;
contourf(xx, yy, abs_rot_z, 20, 'LineColor', 'none');  % Gráfico de contornos rellenos
colorbar;                                             % Barra de colores
title('Valor absoluto del rotacional (∇ × ⃗u_z) en la sección transversal');
xlabel('x (m)');
ylabel('y (m)');

% Añadir información sobre los ejes
hold on;
quiver(xx, yy, -yy ./ r, xx ./ r, 0.5, 'k');  % Direcciones del campo vectorial (rotacional)
hold off;


9 Tensiones Normales

El tensor de tensiones [math]σ_{ij}[/math] se puede escribir a través de la fórmula:
[math]σ = λ∇ · \vec{u} 1 + 2µϵ[/math]

Donde [math]ϵ(\vec{u})= (∇\vec{u} + ∇\vec{u}^t)/2 [/math], es la parte simétrica del tensor gradiente de [math] \vec{u}[/math].


Empezaremos calculando
[math]ε(\vec u)= \begin{pmatrix} 0 & 0 & 0\\ 0 & cos(\theta)sin(\frac{2\pi*p}{50}) & 0 \\ 0 & 0 & 0 \end{pmatrix}[/math]


[math]σ=λ\nabla·\vec u 1+2με = λ*cos(\theta)sin(\frac{2\pi*p}{50})·\begin{pmatrix} 1 & 0 & 0\\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix} + 2μ\begin{pmatrix}0 & 0 & 0 \\ 0 & cos(\theta)sin(\frac{2\pi*p}{50}) & 0 \\ 0 & 0 & 0 \end{pmatrix} [/math]


[math]= cos(\theta)sin(\frac{2\pi*p}{50}) (\begin{pmatrix} λ & 0 & 0\\ 0 & λ & 0 \\ 0 & 0 & λ \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0\\ 0 & 2μ & 0 \\ 0 & 0 & 0 \end{pmatrix}) = cos(\theta)sin(\frac{2\pi*p}{50})\begin{pmatrix} λ & 0 & 0\\ 0 & 2μ+λ & 0 \\ 0 & 0 & λ \end{pmatrix}→{λ=μ=1}→cos(\theta)sin(\frac{2\pi*p}{50})\begin{pmatrix} 1 & 0 & 0\\ 0 & 3 & 0 \\ 0 & 0 & 1 \end{pmatrix}[/math]


Las tensiones normales en la dirección de los ejes [math](\vec e_ρ,\vec e_\theta,\vec e_z)[/math] son:

  • [math]\vec e_ρ · σ ·\vec e_\rho = \begin{pmatrix} 1 & 0 & 0 \end{pmatrix}·\begin{pmatrix} cos(\theta)sin(\frac{2\pi*p}{50}) & 0 & 0\\ 0 & 3*cos(\theta)sin(\frac{2\pi*p}{50}) & 0 \\ 0 & 0 & cos(\theta)sin(\frac{2\pi*p}{50}) \end{pmatrix}·\begin{pmatrix} 1 \\0 \\ 0 \end{pmatrix} = cos(\theta)sin(\frac{2\pi*\rho}{50}) [/math]


  • [math]\vec e_\theta · σ ·\vec e_\theta = \begin{pmatrix} 0 & 1 & 0 \end{pmatrix}·\begin{pmatrix} cos(\theta)sin(\frac{2\pi*p}{50}) & 0 & 0\\ 0 & 3*cos(\theta)sin(\frac{2\pi*p}{50}) & 0 \\ 0 & 0 & cos(\theta)sin(\frac{2\pi*p}{50}) \end{pmatrix}·\begin{pmatrix} 0 \\1 \\ 0 \end{pmatrix} = 3*cos(\theta)sin(\frac{2\pi*\rho}{50}) [/math]


  • [math]\vec e_z · σ ·\vec e_z = \begin{pmatrix} 0 & 0 & 1 \end{pmatrix}·\begin{pmatrix} cos(\theta)sin(\frac{2\pi*p}{50}) & 0 & 0\\ 0 & 3*cos(\theta)sin(\frac{2\pi*p}{50}) & 0 \\ 0 & 0 & cos(\theta)sin(\frac{2\pi*p}{50}) \end{pmatrix}·\begin{pmatrix} 0 \\0 \\ 1 \end{pmatrix} = cos(\theta)sin(\frac{2\pi*\rho}{50}) [/math]

10 Tensiones tangenciales al plano ortogonal al \( \vec{e_\rho} \)

Para calcular las tensiones tangenciales respecto al plano ortogonal a \( \vec{e_\rho} \), utilizaremos la fórmula:[math]\ \left| \sigma ·\vec i -(\vec i·\sigma·\vec i)\vec i \right| [/math].
en la que sustituiremos \(\vec{i}\) por \( \vec{e_\rho} \)

Entonces para comenzar a operar sabemos que el tensor tensiones es el: σ=\begin{pmatrix}\cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) & 0 & 0 \\ 0 & 3\cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) & 0 \\ 0 & 0 & \cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) \end{pmatrix}.

Primero, calculamos el producto \( \sigma \cdot \vec{e_\rho} \):=\begin{pmatrix}\cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) \\ 0 \\ 0\end{pmatrix} y luego el segundo sumando \((\vec{e_\rho} \cdot \sigma \cdot \vec{e_\rho}) \vec{e_\rho}\)= \begin{pmatrix} \cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) \\ 0 \\ 0 \end{pmatrix}


Finalmente, restamos para obtener las tensiones por medio de la fórmula \( \vec{\tau} \)= \begin{pmatrix} \cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) \\ 0 \\ 0 \end{pmatrix} - \begin{pmatrix} \cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right) \\ 0 \\ 0 \end{pmatrix} = \begin{pmatrix} 0 \\ 0 \\ 0 \end{pmatrix}.


Por lo tanto, no existen tensiones tangenciales no nulas en este caso.


11 Tensión de Von Mises

La tensión de Von Mises se define mediante la fórmula:
[math] \sigma_{VM} = \sqrt{ \frac { ( \sigma_1 - \sigma_2 )^2 + ( \sigma_2 - \sigma_3 )^2 + ( \sigma_3 - \sigma_1 )^2 }{2} }[/math];
donde \( \sigma_1 \), \( \sigma_2 \), y \( \sigma_3 \) son los autovalores del tensor de tensiones \( \sigma \).

Este resultado permite calcular la tensión de Von Mises indicando los puntos donde la tensión alcanza valores máximos y así podemos determinar si un material soportará las cargas aplicadas o si llegará a fallar.


Tensión Von Mises
% Definir la sección transversal
x = linspace(-1, 1, 200);
y = linspace(0, 10, 500);
[X, Y] = meshgrid(x, y);
rho = sqrt(X.^2 + Y.^2);
theta = atan2(Y, X);

% Componentes del tensor de tensiones en coordenadas polares
T_rr = rho .* cos(theta) .* sin(2 * pi * rho / 50);
T_rtheta = 3 * cos(theta) .* sin(2 * pi * rho / 50);
T_zz = T_rr;

% Inicializar los autovalores
sigma_1 = zeros(size(X));
sigma_2 = zeros(size(X));
sigma_3 = zeros(size(X));

% Calcular los autovalores en cada punto
for i = 1:size(X, 1)
    for j = 1:size(X, 2)
        % Construir el tensor de tensiones
        stress_tensor = [T_rr(i, j), T_rtheta(i, j), 0;
                         T_rtheta(i, j), T_rr(i, j), 0;
                         0, 0, T_zz(i, j)];
                     
        % Calcular los autovalores
        eigenvalues = eig(stress_tensor);
        sigma_1(i, j) = eigenvalues(1);
        sigma_2(i, j) = eigenvalues(2);
        sigma_3(i, j) = eigenvalues(3);
    end
end

% Calcular la tensión de Von Mises
VM = sqrt(0.5 * ((sigma_1 - sigma_2).^2 + (sigma_2 - sigma_3).^2 + (sigma_3 - sigma_1).^2));

% Encontrar el valor máximo y su posición
[VM_max, idx] = max(VM(:));
[max_row, max_col] = ind2sub(size(VM), idx);
max_x = X(max_row, max_col);
max_y = Y(max_row, max_col);

% Graficar en 3D
figure;
surf(X, Y, VM, 'EdgeColor', 'none'); % Gráfico de superficie
colormap('parula');
colorbar;
hold on;

% Asegurar que el punto rojo aparece en la gráfica
z_offset = VM_max + 0.01 * (max(VM(:)) - min(VM(:))); % Evitar problemas de superposición
plot3(max_x, max_y, z_offset, 'ro', 'MarkerSize', 10, 'MarkerFaceColor', 'r');

% Añadir texto para mayor claridad
text(max_x, max_y, z_offset, sprintf('VM Máx: %.2f', VM_max), 'Color', 'red', 'FontSize', 10);

% Personalizar el gráfico
xlabel('x');
ylabel('y');
zlabel('Tensión de Von Mises (VM)');
title('Distribución de Tensión de Von Mises en 3D');
legend(sprintf('Máximo VM: %.2f en (%.2f, %.2f)', VM_max, max_x, max_y), 'Location', 'best');
grid on;
view(45, 30); % Ajustar ángulo de la vista


12 Campo de fuerzas

El desplazamiento de la placa viene dado por un vector fuerza [math]\ \vec F [/math] que actua en el interior de la placa, El cálculo de dicho vector de fuerzas se realiza haciendo una aproximación de la ecuación de la elasticidad lineal, que es : [math]\ \vec F = \nabla ·\sigma_{ij} [/math], es decir que habrá que aplicar la divergencia a la matriz del antes mencionado tensor tensiones [math] \ \vec{F}_\rho = -\frac{\partial \sigma_{\rho \rho}}{\partial \rho} \cdot \vec{e}_\rho - \frac{1}{\rho}\frac{\partial \sigma_{\rho \theta}}{\partial \theta} \cdot \vec{e}_\rho - \frac{\sigma_{\theta \theta}}{\rho} \cdot \vec{e}_\rho = -\left( \frac{\partial \cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right)}{\partial \rho} + \frac{1}{\rho} \frac{\partial 0}{\partial \theta} + \frac{3\cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right)}{\rho} \right) \cdot \vec{e}_\rho = -\left(\frac{2\pi}{50}\cos(\theta) \cos\left(\frac{2\pi \rho}{50}\right) + \frac{3\cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right)}{\rho}\right) \cdot \vec{e}_\rho [/math]


[math] \ \vec{F}_\theta = -\frac{\partial \sigma_{\theta \rho}}{\partial \rho} \cdot \vec{e}_\theta - \frac{1}{\rho}\frac{\partial \sigma_{\theta \theta}}{\partial \theta} \cdot \vec{e}_\theta - \frac{\sigma_{\rho \theta}}{\rho} \cdot \vec{e}_\theta = -\left( \frac{\partial 0}{\partial \rho} + \frac{1}{\rho} \frac{\partial 3\cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right)}{\partial \theta} + \frac{0}{\rho} \right) \cdot \vec{e}_\theta = -\left(\frac{-3\sin(\theta) \sin\left(\frac{2\pi \rho}{50}\right)}{\rho}\right) \cdot \vec{e}_\theta [/math]

[math] \ \vec{F}_z = -\frac{\partial \sigma_{z \rho}}{\partial \rho} \cdot \vec{e}_z - \frac{1}{\rho}\frac{\partial \sigma_{z \theta}}{\partial \theta} \cdot \vec{e}_z - \frac{\partial \sigma_{z z}}{\partial z} \cdot \vec{e}_z = -\left( \frac{\partial 0}{\partial \rho} + \frac{1}{\rho} \frac{\partial 0}{\partial \theta} + \frac{\partial \cos(\theta) \sin\left(\frac{2\pi \rho}{50}\right)}{\partial z} \right) \cdot \vec{e}_z = 0 \cdot \vec{e}_z [/math]


Después de hacer la divergencia vemos que la fuerza no tiene componente \( \vec{F}_z \) ya que es nula es por eso que solo actuará en \( \vec{F}_\rho \) y \( \vec{F}_\theta \)

13 Masa Total de la Placa

Queremos calcular la masa de la placa rectangular [−1,1]×[0,10]. La masa se obtiene integrando el campo escalar: [math] d(\rho, \theta)=(2-\rho)(4-\cos(4(\theta+π/2)))[/math] (densidad de la placa) a lo largo de la superficie. Como la placa es rectangular, la podemos parametrizar utilizando coordenadas cartesianas.[math] d(x, y)=(2-\sqrt{x^2+y^2})(4-\cos(4(arctan(\frac{x}{y})+π/2)))[/math]


[math] Masa=\int_{-1}^{1} \int_{0}^{10} ​ (2-\sqrt{x^2+y^2})(4-\cos(4(arctan(\frac{x}{y})+π/2)))dxdy = 184.771 [/math]