Elementos Finitos.docx

Universidad de Cuenca Elementos finitos

FACULTAD DE INGENIERÍA ESCUELA DE INGENIERÍA CIVIL ELEMENTOS FINITOS

“MÉTODO DE LOS ELEMENTOS FINITOS: Ecuación del calor” INFORME N° 2 INTEGRANTES:

Acero Mainato Paiwa Paccha

0302917521

Uyaguari Perea Johanna Alexandra

0105701346

Valla Poma Shirley Johanna

0106698053

Villa Clavijo Damián Mateo

0106488836

DOCENTE:

Ing. Esteban Samaniego

FECHA DE ENTREGA DE REPORTE: Viernes, 6 de julio de 2018

Universidad de Cuenca Elementos finitos

 Metodología El método de los elementos finitos es aplicado en 2 dimensiones para el análisis de la conducción de calor a través de superficies con distintas condiciones de borde, que definen el comportamiento y parámetros iniciales que permiten calcular la distribución interna del calor. Para la modificación e implementación de nuevos parámetros, tomamos como base el código existente del documento y se le fue reescribiendo, con algunas variaciones de forma que no alteran el funcionamiento. Para los tres ejemplos presentados en este informe, la malla se generó mediante el paquete PdeTool de Matlab. Para ingresar las coordenadas del mallado generado por PdeTool a la ecuación de calor, utilizamos el código “conversión.mat” (véase anexo 2). En dicho código se deben realizar ciertas modificaciones, en las condiciones de Dirichlet dependiendo del problema que se plantee.

 Resultados y Discusión Para la validación de las modificaciones realizadas al código, se utilizó la herramienta Pdetool de Matlab. Teniendo una superficie con sus condiciones de contorno, su flujo y su conductividad definidas, se puede obtener mediante la opción Plot, la distribución de calor en la misma; de esta manera para la misma superficie y bajo las mismas condiciones se emplea el código modificado y se obtiene la misma distribución de calor, como se observa en los siguientes ejemplos. Así, al comparar estas superficies se puede comprobar la funcionalidad del código.

Ejemplo 1 Datos Para este ejemplo se utilizo una constante de conductividad térmica de 0.5 este valor se utilizó en el código principal (ver anexo 1), en el cual K multiplica a la matriz de rigidez “A”. Se considero que las condiciones de Neumann en el contorno igual a cero, este valor se lo utilizo en el código “g.m”, el cual devuelve un vector de ceros. En cuanto a las condiciones de Dirichlet, se utilizo el valor de 25, el cual se introdujo en el código “u_d.,m”, este valor multiplica a un vector de unos. En este mismo código se cambio el primer argumento de “size” ya que trabajamos con una malla triangular. A continuación, se presenta el ingreso de datos en la herramienta PdeTool, para generar la distribución térmica en la figura realizada.

Universidad de Cuenca Elementos finitos

Ilustración 1: coeficiente de conductividad térmica

Ilustración 2: condiciones de Dirichlet y Neumann

Geometría Se compara las distribuciones de temperaturas obtenidas con la herramienta Pdetool y el código modificado, bajo las siguientes condiciones:  

Número de nodos = 136 (Las coordenadas no se colocan por ser demasiado extensas) Número de triángulos = 206 (Los nodos no se colocan por ser demasiado extensos)



Condiciones de Contorno

Ilustración 3: mallado de la figura propuesta

Universidad de Cuenca Elementos finitos

Ilustración 4: Regiones y condiciones Dirichlet (rojo)

Resultados Mediante PdeTool:

Ilustración 5: resultado usando Pdetool

Mediante el código:

Universidad de Cuenca Elementos finitos

Ilustración 6: resultado mediante el código de Ecuación de Calor

En el grafico se observa que el calor se concentra en todo el borde de la figura (color amarillo) debido a que las condiciones de Dirichlet fueron aplicadas en todo el contorno. En el interior de la figura se observa un descenso de temperatura. Vista tridimensional

Ilustración 7: grafica en 3D

Comentarios Al aumentar el coeficiente de conductividad térmica (K), el material será mejor conductor del calor. Se observa que la distribución térmica de la gráfica obtenida de PdeTool coincide con la gráfica obtenida del código principal.

Universidad de Cuenca Elementos finitos

Ejemplo 2 Se hace un análisis de convergencia, utilizando diferentes mallas con mayor aproximación en una misma figura, es decir, en diferentes pruebas del código para una misma superficie, se va aumentando la cantidad de triángulos con la finalidad de mejorar la aproximación. De esta manera es posible incrementar la precisión de la distribución de temperaturas. Datos Para este ejemplo se utilizó una constante de conductividad térmica de 1, este valor se utilizó en el código principal (ver anexo 1), en el cual K multiplica a la matriz de rigidez “A”. Se considero que las condiciones de Neumann en el contorno son igual a cero, este valor se utilizó en el código “g.m”, el cual devuelve un vector de ceros. En cuanto a las condiciones de Dirichlet, se utilizó el valor de 10, el cual se introdujo en el código “u_d.,m”, este valor multiplica a un vector de unos. En este mismo código se cambió el primer argumento de “size” ya que trabajamos con una malla triangular. A continuación, se presenta el ingreso de datos en la herramienta PdeTool, para generar la distribución térmica en la figura realizada.

Ilustración 8: coeficiente de conductividad térmica

Ilustración 9: condiciones de Dirichlet y Neumann

Geometría

Universidad de Cuenca Elementos finitos

Ilustración 10: regiones del elemento por analizar

Las condiciones de borde de Dirichlet, no se colocaron en todo el contorno sino en algunos nodos de la región R1 y R3 en los siguientes nodos: Nodos

Condiciones de Dirichlet 10, 9, 8, 10 19, 22 Tabla 1: nodos y Dirichlet

Resultados Primera triangulación  Número de nodos = 514 (Las coordenadas no se colocan por ser demasiado extensas)  Número de triángulos = 896 (Los nodos no se colocan por ser demasiado extensos)

Ilustración 11: resultado usando el código de la ecuación del calor

Segunda triangulación

Universidad de Cuenca Elementos finitos

 

Número de nodos = 1923(Las coordenadas no se colocan por ser demasiado extensas) Número de triángulos = 3584 (Los nodos no se colocan por ser demasiado extensos)

Ilustración 12: Resultados utilizando PdeTool

Comentarios Con respecto a la convergencia, se puede notar que a medida que aumenta el número de triángulos en la malla, la diferencia entre temperaturas máximas va disminuyendo; esto quiere decir que los valores convergen hacia un valor determinado, lo cual depende de la aproximación que se le dé a la figura estudiada.

Ejemplo 3 Datos Se presenta el ingreso de datos en la herramienta PdeTool, para generar la distribución térmica en la figura realizada.

Ilustración 13: condiciones de Dirichlet y Neumann

Universidad de Cuenca Elementos finitos

Ilustración 14: coeficiente de conductividad térmica

Geometría

Ilustración 15: Nodos y regiones

Ilustración 16: triangulación de la figura - PdeTool

Las condiciones de borde de Dirichlet, no se colocaron en todo el contorno sino en los siguientes nodos: Nodos

Condiciones de Dirichlet

2-44-101-123-43-121122-41-112-116-39-48100 46-37-97-99-35-106108-34-4-3 Tabla 2: nodos y Dirichlet

Resultados Mediante código de conducción del calor

Universidad de Cuenca Elementos finitos

Ilustración 17: Resultados utilizando la ecuación del calor

Ilustración 17: resultados mediante código – ecuación del calor

Resultados de PdeTool

Ilustración 1918: resultados con PdeTool

Universidad de Cuenca Elementos finitos

En los gráficos se observa que el calor se distribuye por la base y las paredes de la figura. Al comparar los gráficos obtenidos del programa Pdetool y del código de la ecuación del calor, se observa que la distribución de calor es similar en ambos casos. El coeficiente de conductividad térmica es de 0.5, lo que significa que el material no es un buen conductor del calor, es porque se observa amplias zonas con temperaturas bajas rodeando a los 0oC. (zonas de color azul).

 Referencias 

Alberty, J., Carstensen, C., & A. Funken, S. (1999). Remarks around 50 lines of Matlab: short finite element. Numerical Algorithms, 117-137.



The MathWorks, Inc. (2016). MATLAB.



Notas de Clase

 ANEXOS Anexo 1 %FEM2D_HEAT finite element method for two-dimensional heat equation. %Initialisation load coordinates.dat; coordinates(:,1)=[]; load elements3.dat; elements3(:,1)=[]; eval('load neumann.dat; neumann(:,1) = [];','neumann=[];');3 load dirichlet.dat; dirichlet(:,1) = []; FreeNodes=setdiff(1:size(coordinates,1),unique(dirichlet)); % valores de A que no estan en B A = sparse(size(coordinates,1),size(coordinates,1)); % Guardo memoria con matriz de tipo sparse B = sparse(size(coordinates,1),size(coordinates,1)); T = 1; dt = 0.01; N = T/dt; U = zeros(size(coordinates,1),N+1); % Assembly: Montaje de la matriz de rigidez A en un bucle sobre todo triang for j = 1:size(elements3,1) A(elements3(j,:),elements3(j,:)) = 10*(A(elements3(j,:), ... elements3(j,:)) + (stima3(coordinates(elements3(j,:),:)))); end for j = 1:size(elements3,1) B(elements3(j,:),elements3(j,:)) = B(elements3(j,:),elements3(j,:)) + ... det([1,1,1;coordinates(elements3(j,:),:)'])*[2,1,1;1,2,1;1,1,2]/24; end % Initial Condition U(:,1) = zeros(size(coordinates,1),1); % time steps for n = 2:N+1 b = sparse(size(coordinates,1),1); % Volume Forces for j = 1:size(elements3,1) b(elements3(j,:)) = b(elements3(j,:)) + ... det([1,1,1; coordinates(elements3(j,:),:)']) * ... dt*f(sum(coordinates(elements3(j,:),:))/3,n*dt)/6; end % Neumann conditions for j = 1 : size(neumann,1) b(neumann(j,:)) = b(neumann(j,:)) + ... norm(coordinates(neumann(j,1),:)-coordinates(neumann(j,2),:)) * ... dt*g(sum(coordinates(neumann(j,:),:))/2,n*dt)/2; end % previous timestep b = b + B * U(:,n-1); % Dirichlet conditions u = sparse(size(coordinates,1),1); u(unique(dirichlet)) = u_d(coordinates(unique(dirichlet),:),n*dt);

Universidad de Cuenca Elementos finitos b = b - (dt * A + B) * u; % Computation of the solution u(FreeNodes) = (dt*A(FreeNodes,FreeNodes)+ ... B(FreeNodes,FreeNodes))\b(FreeNodes); U(:,n) = u; end % graphic representation show(elements3,[],coordinates,full(U(:,N+1)));

Anexo 2 function conversion(p,e,t) % asignamos los nodos para cada uno de los triangulos generados por la % malla, la esta definida con un sentido horario de la forma % |elemento|nodo1|nodo2|nodo3| for i=1:length(t) elements3(i,1)=i; elements3(i,2)=t(1,i); elements3(i,3)=t(2,i); elements3(i,4)=t(3,i); end %definimos las condiciones de neumann y dirichlet en funci?n al elemento contd=0; contn=0; for i=1:length(e) if(e(5,i)==10 || e(5,i)==9 || e(5,i)==8 || e(5,i)==19 ) contd=contd+1; dirichlet(contd,1)=contd; dirichlet(contd,2)=e(1,i); dirichlet(contd,3)=e(2,i); end if(e(5,i)==22) contd=contd+1; dirichlet(contd,1)=contd; dirichlet(contd,2)=e(1,i); dirichlet(contd,3)=e(2,i); end if(e(5,i)==3 || e(5,i)==4) contn=contn+1; neumann(contn,1)=contn; neumann(contn,2)=e(1,i); neumann(contn,3)=e(2,i); end end %asignamos las coordenadas para cada elemento a trav?s de la variable p (puntos) %entregada por toolbox de la forma |elemento|x|y| for i=1:length(p) coordinates(i,1)=i; coordinates(i,2)=p(1,i); coordinates(i,3)=p(2,i); end %reescribimos nuestras variables .dat para trabajar con el programa %entregado dlmwrite('coordinates.dat',coordinates) dlmwrite('elements3.dat',elements3) dlmwrite('dirichlet.dat',dirichlet) dlmwrite('neumann.dat',neumann) end