Diferencia entre revisiones de «Movimiento de partículas. Grupo 22C»

Revisión del 16:26 5 dic 2014

Consideramos un conjunto de 20 partículas que inicialmente se encuentran en los puntos de coordenadas (x_i,y_i,z_i) = (cos(2πi/10),sin(2πi/10), i/10), i=1, 2, ..., 20 respecto a la base ortonormal {e₁,e₂,e₃}. Supondremos que las partículas están unidas por alambres de masa despreciable de manera que su posición relativa no cambia.

1 VISUALIZACIÓN DE UN SISTEMA DE PARTÍCULAS

Creamos un programa para obtener representación de la posición de las partículas Los ejes de la siguiente figura están fijados en la región [-2,2]x[-2,2]x[0,2]. construimos los vectores con las coordenadas X, Y, Z de los 20 puntos y los representamos.

Visualizacion de un sistema de particulas

%la primera parte del programa calcula las coordenadas de los puntos
%materiales.
%crea una matriz de 3 filas y 20 columnas, donde cada columna son las
%tres coordenadas de los puntos.
coordenadas=zeros(3,20);
for i=1:20
coordenadas(1,i)=(cos((2*pi*i)/10));
coordenadas(2,i)=(sin((2*pi*i)/10));
coordenadas(3,i)=(i/10);
end
clear('i');
%dibujamos la helice en 3d con los ejes pedidos, 
%creamos los vectores de coordeandas X, Y, Z
coorX=zeros(1,20);
coorY=zeros(1,20);
coorZ=zeros(1,20);
for i=1:20
coorX(1,i)=coordenadas(1,i);
coorY(1,i)=coordenadas(2,i);
coorZ(1,i)=coordenadas(3,i);
end
clear('i');
%pintamos esos vectores
figure(1);
plot3(coorX,coorY,coorZ, '--. b');
axis([-2,2,-2,2,0,2]);
hold on

2 CENTRO DE MASAS DE UN SISTEMA DE PARTÍCULAS

Sabiendo que las partículas tienen masa creciente m_i=10+i/10, para calcular el centro de masas hacemos uso de la siguiente fórmula: (∑_iⁿr_im_i)/M donde M es la masa total y ri= la distancia del punto al eje correspondiente x,y ó z. El código de Matlab será el siguiente:

Visualización del centro de masas

%calculamos una matriz con las masas
masas=zeros(1,20);
for i=1:20
masas(1,i)=(10+i/10);
end
clear('i');
%calculamos cuanto aporta cada particula al centro de masas
%y el centro de masas
auxXG=zeros(1,20);
auxYG=zeros(1,20);
auxZG=zeros(1,20);
for i=1:20
auxXG(1,i)=(masas(1,i)*coorX(1,i));
auxYG(1,i)=(masas(1,i)*coorY(1,i));
auxZG(1,i)=(masas(1,i)*coorZ(1,i));
end
clear('i');
masatotal=sum(masas);
auxXG=sum(auxXG);
auxYG=sum(auxYG);
auxZG=sum(auxZG);
%calculamos las coordenads del centro de masas y lo pintamos
XG=((1/masatotal)*auxXG);
YG=((1/masatotal)*auxYG);
ZG=((1/masatotal)*auxZG);
plot3(XG,YG,ZG,'. g ','linewidth',10)
hold off
figure(2);

3 ROTACIÓN DE UN SISTEMA DE PARTÍCULAS

3.1 Matriz rotación con eje ω=e₃ y angulo θ= π/16

Calculamos la matriz de rotación a partir de un eje ω=e₃ y un ángulo θ= π/16. El código que nos permite obtener dicha matriz de rotación es esta:

Visualización de un sistema de puntos rotados con eje ω=e₃ y angulo θ= π/16

%esta parte del programa calcula la primera matriz de rotación
%calculamos la matriz de la rotación a partir de un eje de vector w
%y un ángulo theta
w=[0, 0 , 1];
theta=pi/16;
matriz1=eye(3);
matriz2=zeros(3,3);
matriz3=zeros(3,3);
for i=1:3
matriz2(1,i)=(w(1)*w(i));
matriz2(2,i)=(w(2)*w(i));
matriz2(3,i)=(w(3)*w(i));
end
clear('i');
matriz3(1,1)=0;
matriz3(1,2)=(-1*w(3));
matriz3(1,3)=(w(2));
matriz3(2,1)=(w(3));
matriz3(2,2)=0;
matriz3(2,3)=(-1*w(1));
matriz3(3,1)=(-1*w(2));
matriz3(3,2)=(w(1));
matriz3(3,3)=0;
%construimos la matriz de rotación
rotacion1=zeros(3,3);
matrizaux1=(cos(theta)*matriz1);
matrizaux2=((1-cos(theta))*matriz2);
matrizaux3=(sin(theta)*matriz3);
rotacion1=matrizaux1+matrizaux2+matrizaux3;
%separamos los vectores directores de cada particula y les hacemos la rotación
coorxR1=zeros(1,20);
cooryR1=zeros(1,20);
coorzR1=zeros(1,20);
vaux=zeros(3,1);
for i=1:20
vaux(1,1)=coorX(1,i);
vaux(2,1)=coorY(1,i);
vaux(3,1)=coorZ(1,i);
rotpunto=rotacion1*vaux;
coorxR1(1,i)=rotpunto(1);
cooryR1(1,i)=rotpunto(2);
coorzR1(1,i)=rotpunto(3);
clear('rotpunto');
end
clear('i');
plot3(coorxR1,cooryR1, coorzR1, '--. g')
axis([-2,2,-2,2,0,2])
hold on
plot3(coorX, coorY, coorZ, '--p b')
hold off

La matriz de componentes del tensor de la primera rotación es:

   0.9808   -0.1951         0
   0.1951    0.9808         0
        0         0    1.0000

3.2 Matriz rotación con eje ω=e₁ y angulo θ= π/16

Código:

Visualización de un sistema de puntos rotados con eje ω=e₁ y angulo θ= π/16

%segunda matriz de rotacion, eje e1
w=[1, 0 , 0];
theta=pi/16;
matriz1=eye(3);
matriz2=zeros(3,3);
matriz3=zeros(3,3);
for i=1:3
matriz2(1,i)=(w(1)*w(i));
matriz2(2,i)=(w(2)*w(i));
matriz2(3,i)=(w(3)*w(i));
end
clear('i');
matriz3(1,1)=0;
matriz3(1,2)=(-1*w(3));
matriz3(1,3)=(w(2));
matriz3(2,1)=(w(3));
matriz3(2,2)=0;
matriz3(2,3)=(-1*w(1));
matriz3(3,1)=(-1*w(2));
matriz3(3,2)=(w(1));
matriz3(3,3)=0;
%construimos la matriz de rotación
rotacion2=zeros(3,3);
matrizaux1=(cos(theta)*matriz1);
matrizaux2=((1-cos(theta))*matriz2);
matrizaux3=(sin(theta)*matriz3);
rotacion2=matrizaux1+matrizaux2+matrizaux3;
clear ('vaux');
%separamos los vectores directores de cada particula y les hacemos la rotación
coorxR2=zeros(1,20);
cooryR2=zeros(1,20);
coorzR2=zeros(1,20);
vaux=zeros(3,1);
for i=1:20
vaux(1,1)=coorX(1,i);
vaux(2,1)=coorY(1,i);
vaux(3,1)=coorZ(1,i);
rotpunto=rotacion2*vaux;
coorxR2(1,i)=rotpunto(1);
cooryR2(1,i)=rotpunto(2);
coorzR2(1,i)=rotpunto(3);
clear('rotpunto');
end
clear('i');
plot3(coorxR2,cooryR2, coorzR2, '--p r')
axis([-2,2,-2,2,0,2])
hold on
plot3(coorX, coorY, coorZ, '--. b')
hold off

La matriz de componentes del tensor de la segunda rotación es:

   1.0000         0         0
        0    0.9808   -0.1951
        0    0.1951    0.9808

3.3 Matriz rotación con eje ω=e₂ y angulo θ= π/16

Calculamos la matriz de rotación a partir de un eje ω=e₂ y un ángulo θ= π/16. El código que nos permite obtener dicha matriz de rotación es esta:

Visualización de un sistema de puntos rotados con eje ω=e₂ y ángulo θ= π/16

%calculamos la matriz de la tercera rotación
%a partir de un vector w y un eje angulo theta
w=[0, 1, 0];
theta=pi/16;
matriz1=eye(3);
matriz2=zeros(3,3);
matriz3=zeros(3,3);
for i=1:3
matriz2(1,i)=(w(1)*w(i));
matriz2(2,i)=(w(2)*w(i));
matriz2(3,i)=(w(3)*w(i));
end
clear('i');
matriz3(1,1)=0;
matriz3(1,2)=(-1*w(3));
matriz3(1,3)=(w(2));
matriz3(2,1)=(w(3));
matriz3(2,2)=0;
matriz3(2,3)=(-1*w(1));
matriz3(3,1)=(-1*w(2));
matriz3(3,2)=(w(1));
matriz3(3,3)=0;
%construimos la matriz de rotación
rotacion3=zeros(3,3);
matrizaux1=(cos(theta)*matriz1);
matrizaux2=((1-cos(theta))*matriz2);
matrizaux3=(sin(theta)*matriz3);
rotacion3=matrizaux1+matrizaux2+matrizaux3;
clear('vaux');
%separamos los vectores directores de cada particula y les hacemos la rotación
coorxR3=zeros(1,20);
cooryR3=zeros(1,20);
coorzR3=zeros(1,20);
vaux=zeros(3,1);
for i=1:20
vaux(1,1)=coorX(1,i);
vaux(2,1)=coorY(1,i);
vaux(3,1)=coorZ(1,i);
rotpunto=rotacion3*vaux;
coorxR3(1,i)=rotpunto(1);
cooryR3(1,i)=rotpunto(2);
coorzR3(1,i)=rotpunto(3);
clear('rotpunto');
end
clear('i');
plot3(coorxR3,cooryR3, coorzR3, '--p k')
axis([-2,2,-2,2,0,2])
hold on
plot3(coorX, coorY, coorZ, '--. b')
hold off

La matriz de componentes del tensor de la tercera rotación es:

   0.9808         0    0.1951
        0    1.0000         0
  -0.1951         0    0.9808

3.4 Matriz rotación con eje ω=e₁+e₂+e₃ y ángulo θ= π/16

Visualizacion de un sistema de puntos rotados con eje ω=e₁+e₂+e₃ y ángulo θ= π/16

%calculamos la matriz de la cuarta rotación
%a partir de un vector w y un eje angulo theta
a=1/sqrt(3);
w=[a, a, a];clc
theta=pi/16;
matriz1=eye(3);
matriz2=zeros(3,3);
matriz3=zeros(3,3);
for i=1:3
matriz2(1,i)=(w(1)*w(i));
matriz2(2,i)=(w(2)*w(i));
matriz2(3,i)=(w(3)*w(i));
end
clear('i');
matriz3(1,1)=0;
matriz3(1,2)=(-1*w(3));
matriz3(1,3)=(w(2));
matriz3(2,1)=(w(3));
matriz3(2,2)=0;
matriz3(2,3)=(-1*w(1));
matriz3(3,1)=(-1*w(2));
matriz3(3,2)=(w(1));
matriz3(3,3)=0;
%construimos la matriz de rotación
rotacion4=zeros(3,3);
matrizaux1=(cos(theta)*matriz1);
matrizaux2=((1-cos(theta))*matriz2);
matrizaux3=(sin(theta)*matriz3);
rotacion4=matrizaux1+matrizaux2+matrizaux3;
clear('vaux');
%separamos los vectores directores de cada partícula y les hacemos la rotación
coorxR4=zeros(1,20);
cooryR4=zeros(1,20);
coorzR4=zeros(1,20);
vaux=zeros(3,1);
for i=1:20
vaux(1,1)=coorX(1,i);
vaux(2,1)=coorY(1,i);
vaux(3,1)=coorZ(1,i);
rotpunto=rotacion4*vaux;
coorxR4(1,i)=rotpunto(1);
cooryR4(1,i)=rotpunto(2);
coorzR4(1,i)=rotpunto(3);
clear('rotpunto');
end
clear('i');
plot3(coorxR4,cooryR4, coorzR4, '--p m')
axis([-2,2,-2,2,0,2])
hold on
plot3(coorX, coorY, coorZ, '--. b')
clear ('matriz1','matriz2','matriz3','w','matrizaux1','matrizaux2','matrizaux3');
hold off
%fin del apartado 3

La matriz de componentes del tensor de la tercera rotación es:

   0.9872   -0.1062    0.1190
   0.1190    0.9872   -0.1062
  -0.1062    0.1190    0.9872

4 DEFINICIÓN DE LA VELOCIDAD COMO UN TENSOR ANTISIMÉTRICO

Dado un sistema que gira alrededor de un eje genérico con vector de dirección unitario [math]\overrightarrow{\omega}[/math] cuya variación angular viene dada por la función θ(t) donde t ∈ [0, π] es el tiempo, comprobamos analíticamente que la velocidad de los puntos del sistema está relacionada con su posición mediante un tensor antisimétrico A de la siguiente forma [math]v_{i}(t)=\frac{dr_{i}(t)}{dt}=A\cdot r_{i}(t)[/math] para todo i = 1, 2, ..., 20.

Para poder expresar la velocidad angular [math]\Omega [/math] en forma vectorial, primero multiplicamos su módulo [math]\dot{\Theta}(t)[/math] por su vector director unitario [math]\overrightarrow{\omega}[/math] [math]\rightarrow \Omega =\left | \dot{\Theta}(t) \right |\cdot \overrightarrow{\omega}[/math].

Por el campo de velocidades sabemos que [math]\overrightarrow{V_{i}}=\overrightarrow{\Omega}\times \overrightarrow{OP}=\overrightarrow{\Omega}\times \overrightarrow{r_{i}}[/math].

En forma tensorial [math]\overrightarrow{V_{i}}=(\overrightarrow{\Omega }\times )\cdot \overrightarrow{r_{i}}[/math].

A continuación se demuestra que [math]\overrightarrow{\Omega }\times [/math] es un tensor antisimétrico donde [math]\Omega =\overrightarrow{\omega}\times [/math] :

[math]\Omega _{ij}=\overrightarrow{g_{i}}\cdot [\Omega \cdot \overrightarrow{g_{j}}]=\overrightarrow{g_{i}}\cdot [(\overrightarrow{\omega}\times )\cdot \overrightarrow{g_{j}}]=\overrightarrow{g_{i}}\cdot (\overrightarrow{\omega}\times \overrightarrow{g_{j}})=[\overrightarrow{g_{i}},\overrightarrow{\omega},\overrightarrow{g_{j}}]=-[\overrightarrow{g_{i}},\overrightarrow{g_{j}},\overrightarrow{\omega}]=-\omega ^{k}[\overrightarrow{g_{i}},\overrightarrow{g_{j}},\overrightarrow{g_{k}}]=-\omega ^{k}\pm \sqrt{g}\varepsilon _{ijk}=-\omega ^{1}\cdot \varepsilon _{ij1}-\omega ^{2}\cdot \varepsilon _{ij2}-\omega ^{3}\cdot \varepsilon _{ik3}[/math] [math]-\omega^{1}\cdot \bigl(\begin{smallmatrix} 0 & 0 & 0\\ 0 & 0 & 1\\ 0 & -1 & 0 \end{smallmatrix}\bigr)-\omega ^{2}\cdot \bigl(\begin{smallmatrix} 0 & 0 & -1\\ 0 & 0 & 0\\ 1 & 0 & 0 \end{smallmatrix}\bigr)-\omega ^{3}\cdot \bigl(\begin{smallmatrix} 0 & 1 & 0\\ -1 & 0 & 0\\ 0 & 0 & 0 \end{smallmatrix}\bigr)=\bigl(\begin{smallmatrix} 0 & -\omega ^{3} & \omega ^{2}\\ \omega ^{3} & 0 & -\omega ^{1}\\ -\omega ^{2} & \omega ^{1} & 0 \end{smallmatrix}\bigr)[/math]

Por consiguiente se puede demostrar que el vector axial asociado a [math]\overrightarrow{\Omega}\times [/math] es [math]\dot{\Theta }(t)\cdot \overrightarrow{\omega } [/math] multiplicando la matriz de coordenadas de [math]\Omega [/math] por [math]\dot{\Theta }(t)[/math].

5 VISUALIZACIÓN DE LOS VECTORES VELOCIDAD

Visualización de los vectores velocidad

véase que los vectores son tangentes a la trayectoria

%pintamos los las partículas y calculamos
%los vectores de velocidad
plot3(coorX,coorY,coorZ, '-. b')
axis([-2,2,-2,2,0,2]);
hold on
%calculamos las coordenadas de los vectores velocidad
%de los puntos materiales
coorVx=zeros(1,20);
coorVy=zeros(1,20);
coorVz=zeros(1,20);
omega=[0,0,1];
radioaux=zeros(1,3);
veloaux=zeros(1,3);
for i=1:20
radioaux(1,1)=coorX(1,i);
radioaux(1,2)=coorY(1,i);
radioaux(1,3)=coorZ(1,i);
veloaux=cross(omega,radioaux);
coorVx(1,i)=veloaux(1);
coorVy(1,i)=veloaux(2);
coorVz(1,i)=veloaux(3);
clear('veloaux','radioaux');
radioaux=zeros(1,3);
veloaux=zeros(1,3);
end
for i=1:20
xaux=coorX(1,i);
yaux=coorY(1,i);
zaux=coorZ(1,i);
vxaux=coorVx(1,i);
vyaux=coorVy(1,i);
vzaux=coorVz(1,i);
%pintamos el vector velocidad con el comando quiver
quiver3(xaux,yaux,zaux,vxaux,vyaux,vzaux)
clear('xaux','yaux','zaux','vxaux','vyaux','vzaux');
end
hold off

6 TENSOR DE INERCIA

El momento angular del sistema se define por [math]L=\sum_{i=1}^{20}\overrightarrow {r_{i}}\times m_{i}\cdot\overrightarrow {v_{i}}[/math]. Sustituyendo el valor de la velocidad [math]\overrightarrow{v_{i}}=\overrightarrow{\omega}\times \overrightarrow{r_{i}} [/math], donde [math]\overrightarrow{\omega}[/math] es la velocidad angular, comprobamos analíticamente que podemos escribir [math]L=I\cdot \overrightarrow{\omega}[/math] donde I es un tensor de orden 2 conocido como tensor de inercia y sus componentes momentos de inercia.

[math]L=\sum_{i=1}^{20}\overrightarrow {r_{i}}\times m_{i}\cdot\overrightarrow {v_{i}}[/math] [math]=\sum_{i=1}^{20}\overrightarrow{r_{i}}\cdot (m_{i}\cdot \overrightarrow{\omega}\times \overrightarrow{r_{i}})[/math] [math]=\sum_{i=1}^{20}\begin{vmatrix} \overrightarrow{\omega}\cdot m_{i} & \overrightarrow{r_{i}}\\ \overrightarrow{r_{i}}\cdot \overrightarrow{\omega}\cdot m_{i} & \overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}} \end{vmatrix}[/math] [math]=\sum_{i=1}^{20}\left [ \overrightarrow{\omega}\cdot m_{i}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}})\cdot \overrightarrow{r_{i}}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{\omega}\cdot m_{i}) \right ][/math] [math]=\sum_{i=1}^{20}\left [m_{i} \cdot (\overrightarrow{\omega}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}}))-\overrightarrow{r_{i}}(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})\right ][/math] [math]=\sum_{i=1}^{20}[m_{i}[\overrightarrow{\omega}\cdot \left | \overrightarrow{r_{i}} \right |^{2}-\overrightarrow{r_{i}}(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})]][/math] [math]=\sum_{i=1}^{20}m_{i}\cdot [\left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-\overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}} ]\cdot \overrightarrow{\omega}[/math] [math]=\sum_{i=1}^{20}[m_{i}\cdot \left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-m_{i}\cdot \overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}]\overrightarrow{\omega}=L[/math] [math]=\frac{1}{2}m_{i}[(\overrightarrow{\omega}\cdot \overrightarrow{\omega})\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}})-(\overrightarrow{\omega}\cdot \overrightarrow{r_{i}})\cdot (\overrightarrow{\omega}\cdot \overrightarrow{r_{i}})]=\frac{1}{2}m_{i}\cdot \overrightarrow{\omega}[\overrightarrow{\omega}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}})-\overrightarrow{r_{i}}\cdot (\overrightarrow{\omega}\cdot \overrightarrow{r_{i}})]=\frac{1}{2}m_{i}\cdot \overrightarrow{\omega}[\left | \overrightarrow{r_{i}} \right |^{2}\cdot \overrightarrow{\omega}-(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})\cdot \overrightarrow{r_{i}}][/math] [math]=\frac{1}{2}\cdot \overrightarrow{\omega}[m_{i}\cdot (\left | \overrightarrow{r_{i}} \right |^{2}\cdot \overrightarrow{\omega})-m_{i}\cdot [(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})\cdot \overrightarrow{r_{i}}]]=\frac{1}{2}\overrightarrow{\omega}\cdot [m_{i}\cdot \left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-m_{i}\cdot \overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}]\cdot \overrightarrow{\omega}=\frac{1}{2}\cdot \overrightarrow{\omega}\cdot I\cdot \overrightarrow{\omega}[/math]

6.1 Cálculo de las componentes del tensor I de inercia en la base {e₁,e₂,e₃}

%calculamos el tensor de inercia
%primera matriz del tensor% %base i,j,k
%calculamos la primera matriz
clear ('vaux','i');
primeramatriz=zeros(3,3);
identidad=eye(3);
for i=1:20
primeramatrizaux=zeros(3,3);
vaux=zeros(1,3);
vaux(1,1)=coorX(1,i);
vaux(1,2)=coorY(1,i);
vaux(1,3)=coorZ(1,i);
modulocuadrado=vaux*vaux';
masa=masas(i);
momento=masa*modulocuadrado;
primeramatrizaux=momento*identidad;
primeramatriz=primeramatriz+primeramatrizaux;
clear('masa','primeramatrizaux','vaux','modulocuadrado','momento');
end
clear('i');
%calculamos la segunda matriz
segundamatriz=zeros(3,3);
clear ('vaux');
clear ('masa');
for i=1:20
segundamatrizaux=zeros(3,3);
segundamatrizauxaux=zeros(3,3);
vaux=zeros(1,3);
vaux(1,1)=coorX(1,i);
vaux(1,2)=coorY(1,i);
vaux(1,3)=coorZ(1,i);
masa=masas(i);
segundamatrizaux(1,1)=vaux(1,1)*vaux(1,1);
segundamatrizaux(1,2)=vaux(1,1)*vaux(1,2);
segundamatrizaux(1,3)=vaux(1,1)*vaux(1,3);
segundamatrizaux(2,1)=vaux(1,2)*vaux(1,1);
segundamatrizaux(2,2)=vaux(1,2)*vaux(1,2);
segundamatrizaux(2,3)=vaux(1,2)*vaux(1,3);
segundamatrizaux(3,1)=vaux(1,3)*vaux(1,1);
segundamatrizaux(3,2)=vaux(1,3)*vaux(1,2);
segundamatrizaux(3,3)=vaux(1,3)*vaux(1,3);
segundamatrizauxaux=masa*segundamatrizaux;
segundamatriz=segundamatriz+segundamatrizauxaux;
clear ('vaux','segundamatrizaux','segundamatrizauxaux');
end
clear('i');
tensordeinercia=zeros(3,3);
tensordeinercia=primeramatriz+segundamatriz

Ahora lo calculamos suponiendo que ω=e₃:

w=[0;0;1];
momentoejez=tensordeinercia*w;
momentoejez

7 RELACIÓN DEL TENSOR DE INERCIA Y LA ENERGÍA CINÉTICA

La energía cinética total del sistema se define como [math]E_{c}=\sum_{i=1}^{20}\frac{1}{2}\cdot m_{i}\cdot \left | \overrightarrow{v_{i}} \right |^{2}[/math]. Sustituyendo el valor de la velocidad [math]\overrightarrow{v_{i}}=\overrightarrow{\omega}\times \overrightarrow{r_{i}} [/math], comprobamos que podemos escribir [math]E_{c}=\frac{1}{2}\cdot \overrightarrow{\omega}\cdot I\cdot \overrightarrow{\omega}[/math], donde I es el tensor de inercia.

[math]L=\sum_{i=1}^{20}\overrightarrow {r_{i}}\times m_{i}\cdot\overrightarrow {v_{i}}[/math] [math]=\sum_{i=1}^{20}\overrightarrow{r_{i}}\cdot (m_{i}\cdot \overrightarrow{\omega}\times \overrightarrow{r_{i}})[/math] [math]=\sum_{i=1}^{20}\begin{vmatrix} \overrightarrow{\omega}\cdot m_{i} & \overrightarrow{r_{i}}\\ \overrightarrow{r_{i}}\cdot \overrightarrow{\omega}\cdot m_{i} & \overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}} \end{vmatrix}[/math] [math]=\sum_{i=1}^{20}\left [ \overrightarrow{\omega}\cdot m_{i}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}})\cdot \overrightarrow{r_{i}}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{\omega}\cdot m_{i}) \right ][/math] [math]=\sum_{i=1}^{20}\left [m_{i} \cdot (\overrightarrow{\omega}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}}))-\overrightarrow{r_{i}}(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})\right ][/math] [math]=\sum_{i=1}^{20}[m_{i}[\overrightarrow{\omega}\cdot \left | \overrightarrow{r_{i}} \right |^{2}-\overrightarrow{r_{i}}(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})]][/math] [math]=\sum_{i=1}^{20}m_{i}\cdot [\left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-\overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}} ]\cdot \overrightarrow{\omega}[/math] [math]=\sum_{i=1}^{20}[m_{i}\cdot \left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-m_{i}\cdot \overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}]\overrightarrow{\omega}=L[/math] [math]=\frac{1}{2}m_{i}[(\overrightarrow{\omega}\cdot \overrightarrow{\omega})\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}})-(\overrightarrow{\omega}\cdot \overrightarrow{r_{i}})\cdot (\overrightarrow{\omega}\cdot \overrightarrow{r_{i}})]=\frac{1}{2}m_{i}\cdot \overrightarrow{\omega}[\overrightarrow{\omega}\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{r_{i}})-\overrightarrow{r_{i}}\cdot (\overrightarrow{\omega}\cdot \overrightarrow{r_{i}})]=\frac{1}{2}m_{i}\cdot \overrightarrow{\omega}[\left | \overrightarrow{r_{i}} \right |^{2}\cdot \overrightarrow{\omega}-(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})\cdot \overrightarrow{r_{i}}][/math] [math]=\frac{1}{2}\cdot \overrightarrow{\omega}[m_{i}\cdot (\left | \overrightarrow{r_{i}} \right |^{2}\cdot \overrightarrow{\omega})-m_{i}\cdot [(\overrightarrow{r_{i}}\cdot \overrightarrow{\omega})\cdot \overrightarrow{r_{i}}]]=\frac{1}{2}\overrightarrow{\omega}\cdot [m_{i}\cdot \left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-m_{i}\cdot \overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}]\cdot \overrightarrow{\omega}=\frac{1}{2}\cdot \overrightarrow{\omega}\cdot I\cdot \overrightarrow{\omega}[/math]

%calculamos la energía cinética del sistema de partículas
energiaux=tensordeinercia*w;
W=(1/2)*w;
energiacinetica=(W')*energiaux;
energiacinetica

8 TEOREMA DE STEINER

El Teorema de Steiner permite relacionar el tensor de inercia respecto a un eje de giro que pasa por el centro de masas I_G con el tensor de inercia correspondiente a un eje de giro que pasa por cualquier otro punto I_P [math]\Rightarrow I_{P}=I_{G}+m\left \| \overrightarrow{a} \right \|^{2}\cdot 1-\overrightarrow{a}\otimes \overrightarrow{a}[/math] siendo [math]\overrightarrow{a}[/math] el vector que une el centro de masas G con el punto P.

Para empezar la demostración definimos I_P=I_G+T, donde T es un tensor de orden 2. El teorema quedará demostrado si [math]T=m(\left \| \overrightarrow{a} \right \|^{2}\cdot 1-\overrightarrow{a}\otimes \overrightarrow{a})[/math].

Definimos r'=(r-a) para demostrar lo siguiente: [math]I_{P}=\sum_{i=1}^{20}m_{i}(\left | \overrightarrow{{r_{i}}'} \right |^{2}\cdot 1-\overrightarrow{{r_{i}}'}\otimes \overrightarrow{{r_{i}}'})=\sum_{i=1}^{20}m_{i}(\left | \overrightarrow{r_{i}}- \overrightarrow{a}\right |^{2}\cdot 1-(\overrightarrow{r_{i}}-\overrightarrow{a})\otimes (\overrightarrow{r_{i}}-\overrightarrow{a}))=\sum_{i=1}^{20}m_{i}([\left | \overrightarrow{r_{i}} \right |^{2}+\left | \overrightarrow{a} \right |^{2}-2\cdot (\overrightarrow{r_{i}}-\overrightarrow{a})]\cdot 1-[\overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}-\overrightarrow{r_{i}}\otimes \overrightarrow{a}-\overrightarrow{a}\otimes \overrightarrow{r_{i}}+\overrightarrow{a}\otimes \overrightarrow{a}])[/math] [math]I_{P}=\sum_{i=1}^{20}(m_{i}\cdot \left | \overrightarrow{r_{i}} \right |^{2}\cdot 1+m_{i}\left | \overrightarrow{a} \right |^{2}\cdot 1-m_{i}\cdot 2\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{a})\cdot 1-m_{i}\cdot \overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}+m_{i}\overrightarrow{r_{i}}\otimes \overrightarrow{a}+m_{i}\cdot \overrightarrow{a}\otimes \overrightarrow{r_{i}}-m_{i}\cdot \overrightarrow{a}\otimes \overrightarrow{a})[/math] [math]I_{P}=\sum_{i=1}^{20}(m_{i}\cdot [\left | \overrightarrow{r_{i}} \right |^{2}\cdot 1-\overrightarrow{r_{i}}\otimes \overrightarrow{r_{i}}])+\sum_{i=1}^{20}(m_{i}\cdot [\left | \overrightarrow{a} \right |^{2}\cdot 1-2\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{a})\cdot 1+\overrightarrow{r_{i}}\otimes \overrightarrow{a}+\overrightarrow{a}\otimes \overrightarrow{r_{i}}-\overrightarrow{a}\otimes \overrightarrow{a}])[/math]

Steiner quedará demostrado si [math]\sum_{i=1}^{20}m_{i}\cdot (\left | \overrightarrow{a} \right |^{2}\cdot 1+\overrightarrow{r_{i}}\otimes \overrightarrow{a}+\overrightarrow{a}\otimes \overrightarrow{r_{i}}-\overrightarrow{a}\otimes \overrightarrow{a})=\sum_{i=1}^{20}m_{i}\cdot (\left | \overrightarrow{a} \right |^{2}\cdot 1+\overrightarrow{a}\otimes \overrightarrow{a})[/math] [math]\sum_{i=1}^{20}m_{i}\cdot (\left | \overrightarrow{a} \right |^{2}\cdot 1-2\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{a} )\cdot 1+\overrightarrow{r_{i}}\otimes\overrightarrow{a}+\overrightarrow{a}\otimes \overrightarrow{r_{i}}-\overrightarrow{a}\otimes \overrightarrow{a} - \left | \overrightarrow{a} \right |^{2}\cdot 1+\overrightarrow{a}\otimes \overrightarrow{a})=0[/math] [math]\sum_{i=1}^{20}m_{i}\cdot(-2\cdot (\overrightarrow{r_{i}}\cdot \overrightarrow{a} )\cdot 1+\overrightarrow{r_{i}}\otimes \overrightarrow{a}+\overrightarrow{a}\otimes \overrightarrow{r_{i}})=0[/math] [math](-2\cdot [(\sum_{i=1}^{20}m_{i}\cdot \overrightarrow{r_{i}})\cdot \overrightarrow{a})]\cdot 1+(\sum_{i=1}^{20}m_{i}\cdot \overrightarrow{r_{i}})\otimes \overrightarrow{a}+\overrightarrow{a}\otimes (\sum_{i=1}^{20}m_{i}\cdot\overrightarrow{r_{i}} )=0[/math]. Como por definición del centro de masas [math]\sum_{i=1}^{20}m_{i}\cdot\overrightarrow{r_{i}}=0[/math]. [math](-2\cdot (0\cdot \overrightarrow{a})\cdot 1+0\otimes\overrightarrow{a}+\overrightarrow{a}\otimes 0=0 [/math] tal como se quería demostrar.