• Author / Uploaded
• owner salvador

• Categories
• Sistema de ecuaciones lineales
• Ecuaciones
• Ecuaciones diferenciales
• Calor
• Física y matemáticas

Citation preview

UNIVERSIDAD NACIONAL “HERMILIO VALDIZAN”- HUANUCO FACULTAD DE CIVIL Y ARQUITECTURA CURSO : INGENIERIA MÉTODOS NUMÉRICOS

PROFESORA

:

Jaimes Reátegui,

Sumaya INTEGRANTES

:

Domínguez

Muños, Jhon H. Mendoza Santiago, Brigner Salvador Salazar, Owner H

SOLUCIONES NUMÉRICAS PARA LAS ECUACIONES DIFERENCIALES PARCIALES CUERPO ISOTROPICO Si la conductividad térmica en cada uno de los puntos es independiente de la dirección del flujo de calor a través del punto. En un punto la Tº se obtiene resolviendo la siguiente ecuación

donde k, c y p son funciones de (x, y, z) y representan, respectiva­mente, la conductividad térmica, el calor específico y la densidad del cuerpo en el punto (x, y, z). Cuando k, c y p son constantes, a esta ecuación se le denomina ecuación simple tridimensional del calor, y se expresa como

Si la frontera del cuerpo es relativamente simple, la solución de esta ecuación se obtiene usando la serie de Fourier. En la general si k, c y no son constantes o cuando la frontera es irregular, la solución se obtiene mediante métodos de aproximación.



consideraremos la ecuación diferencial parcial elíptica denominada ecuación de Poisson:

 

0 

En esta ecuación suponemos que la función f describe los datos del problema en una región plana R cuya frontera denotamos con S. Este tipo de ecuaciones aparece de manera natural en el estudio de diversos problemas físicos dependientes del tiempo; por ejemplo la distribución de calor para estado estable en una región plana, la energía potencial de un punto en un plano sobre el que operan fuerzas gravitacionales y los problemas bidimencionales del estado estable que incluyen fluidos incompresibles.

Si la temperaturadentro de la region está determinada por su distribución en la frontera de la región, a las restricciones se les llama condiciones de frontera de Dirichlet. Éstas están dadas por  u(x,y) = g(x, y),  para toda (x, y) en S, o sea, la frontera de la región R. Véase la Fig.   

y

S R

(x,y):la temperatura se mantiene constante a g(x,y) grados

x



Este método permite el cálculo de la derivada de orden r de señales discretas, a partir de la determinación de las diferencias centrales existentes entre datos consecutivos, la cual viene dada de forma general mediante:



Para el cálculo de la derivada de cualquier orden se sustituye r en la ecuación anterior por el orden de la derivada de interes, obteniéndose para las derivadas de primer, segundo, tercer y cuarto orden las siguientes expresiones que corresponden a los valores de r=1; r=2; r=3 y r=4 respectivamente.



En la forma de diferencias, esto da como resultado el método de las diferencias centrales con un error local de truncamiento del orden ;

2   h 2   h 2    + 1 w i,j ­ (w i+1,j  w i­1,j )­   (w i,j+1  w i,j­1 ) = ­h2f (xi , y j )  k   k  

Los métodos de Gauss y Cholesky hacen parte de los métodos directos o finitos. Al cabo de un número finito de operaciones, en ausencia de errores de redondeo, se obtiene x* solución del sistema Ax = b.  El método de Gauss-Seidel hace parte de los métodos llamados indirectos o Iterativos. En nellos se comienza con x0 = 1 2 (x0 ; x0 ;…; x0 ), una aproximación inicial de la solución. A partir de x0 se construye una nueva aproximación de la solución, x1 1 2 n 1 = (x1 ; x1 ;…; x1 ). A partir de x se construye x2 (aquí el superíndice indica la iteración y no indica una potencia).  Asi sucesivamente se construye una sucesión de vectores {XK}, con el objetivo, no siempre garantizado, de que   

Generalmente los métodos indirectos son una buena opción cuando la matriz es muy grande y dispersa o rala (sparse), es decir, cuando el numero de elementos no nulos es pequeño comparado con n2, numero total de elementos de A. En estos casos se debe utilizar una estructura de datos adecuada que permita almacenar únicamente los elementos no nulos.  En cada iteración del método de Gauss-Seidel, hay n sub iteraciones. En la primera sub iteración se modifica únicamente x1. Las demás coordenadas x2, x3, ..., xn no se modifican. El calculo de x1 se hace de tal manera que se satisfaga la primera ecuación. 

 En

la segunda sub iteración se modifica únicamente x2. Las demás coordenadas x1, x3, ..., xn no se modifican. El calculo de x2 se hace de tal manera que se satisfaga la segunda ecuación.

 Así

sucesivamente, en la n-esima sub iteración se modifica únicamente xn. las demás coordenadas x1, x2, ..., xn-1 no se modifican. El calculo de xn se hace de tal manera que se satisfaga la n-esima ecuación.

 La

ecuación diferencial parcial elíptica que estudiaremos es la ecuación de

Poisson:2

2  u(x,y)  u(x,y) 2  u(x,y) = + = f(x,y) ; en (x,y)  R y 2 2 x y

 U(x,





y) = g(x, y); para (x,y) S, donde: (x,y) R= a < x < b, c < y < d

Y S denota la frontera de R. Para este análisis, suponemos que tanto f como g son continuas en sus dominios y que se garantiza una solución única.



El método usado es una adaptación de la técnica de diferencias para problemas con valor de frontera. El primer paso h = ( b – a ) / n y k = ( d – c ) / m. La división del intervalo [ a , b ] en n partes iguales de ancho h, y del intervalo [ c , d ] en m partes iguales i ,y j )como resultado una cuadricula en el de ancho k (xda rectángulo R al trazar líneas verticales y horizontales a través de los puntos con coordenadas , donde ym = d

.

.

.

.

y2 y1

. . . . . . . . .

xi = a + ih, para cada i = 0, 1, 2, 3, . . . , n,

. . .

yi = c + jk, para cada j = 0, 1, 2, 3, . . . , m.

. . .

y0 = c

x 0= a

x1

x2

x3

x4

x n= b



Las líneas son líneas de cuadrícula, y sus intersecciones son los puntos de red de la cuadrícula. En cada punto de red del interior de la cuadrícula con i = 1, 2, …, n-1 y con j = 1, 2, …, m-1, utilizando la serie de Taylor en la variable x alrededor de Xi para generar la fórmula de las diferencias centrales:  2u(x i , y j ) x 2

=

u(x i + 1, y j ) ­ 2*u(x i , y j ) + u(x i ­ 1, y j ) h2

h2  4 u ­ ( i , y j ), 12 x 4

Donde i  (x i ­ 1, x i +1 ) 

También usamos la serie de Taylor en la variable y alrededor de Yi para generar la fórmula de las diferencias centrales:  2u(xi , y j ) x

2

Donde

=

u(xi , y j+1 ) ­ 2*u(xi , y j ) + u(xi , y j ­ 1 ) k2

i  ( y i ­ 1, yi +1 )

h2  4 u ­ ( xi ,  j ), 4 12 y

El uso de estas fórmulas nos permite expresar la ecuación de Poisson en los puntos como: (x i ,y j ) u(xi+1, y j )­2*u(xi , y j )+u(xi­1, y j ) u(xi , y j+1 )­2*u(xi , y j )+u(xi , y j­1 ) h2  4 u h2  4 u + = f(xi ,y j )+ ( i , y j )+ ( xi ,  j ) h2 k2 12 x 4 12 y 4

Para toda i = 1, 2, … , n-1 y j = 1, 2, …, m-1 y las condiciones de frontera como: u(x o , y j ) = g(x o , y j ), para cada j = 0, 1, ... ,m, u(xn , y j ) = g(xn , y j ), para cada j = 0, 1, ... ,m, u(xi , y 0 ) = g(x i , y 0 ), para cada i = 0, 1, ... ,n ­ 1, u(xi , y m ) = g(xi , ym ), para cada i = 0, 1, ... ,n ­ 1.



En la forma de diferencias, esto da como resultado el método de las diferencias centrales con un error local de truncamiento del orden o(h2 + k 2 ) ;

2   h 2   h 2    + 1 w i,j ­ (w i+1,j  w i­1,j )­   (w i,j+1  w i,j­1 ) = ­h2f (xi , y j )  k   k  

Para toda i = 1, 2, … , n-1 y j = 1, 2, …, m-1, y

uo j = g(x o , y j ), para cada j = 0, 1, ... ,m, un j = g(xn , y j ), para cada j = 0, 1, ... ,m, ui0 = g(xi , y 0 ), para cada i = 0, 1, ... ,n ­ 1, uim = g(xi , y m ), para cada i = 0, 1, ... ,n ­ 1. 

Donde w i,j aproximau(xi ,y j )

.



u(x,y) La ecuación común contiene aproximaciones a en los puntos :

u(xi+1, y j ),u(xi , y j ),u(xi­1, y j ),u(xi , y j+1 ),u(xi , y j­1 ) Al reducir la parte de la cuadrícula donde estos puntos están situados, se observa que cada ecuación contiene aproximaciones en una región en forma de estrella(xalrededor de i, y j ) d yj + 1

X

yj

X

yj ­1

X

X

X

c

a

xi ­ 1

xi

xi + 1

b

Si utilizamos la información de las condiciones de frontera siempre que sea conveniente en el sistema dado, es (xipuntos , yj ) decir, en todos los adyacentes al punto de red de la frontera, tendremos un sistema lineal (n1)(m-1) por w i,j (n-1)(m-1) incógnitas son las aproximaciones a u(xi cuyas ,y j ) en el interior de los puntos de red.   El sistema lineal que contiene estas incógnitas se expresa más eficientemente en cálculos matriciales, si se introduce un remarcaje de los puntos interiores de la red. Un sistema de marcaje de estos puntos consiste en Pl = (xi , y j ) y Wl = w i,j , utilizar Donde l = i + (m - 1 – j)(n – 1), para toda i = 1, 2, … , n1 y j = 1, 2, …, m-1. Y así se marcan consecutivamente los puntos de red de izquierda a derecha y de arriba abajo. Por ejemplo, con n = 4 y m = 5 con el remarque se obtiene una cuadrícula cuyos puntos se muestran en el gráfico. Al marcar los puntos de este modo, se garantiza que el sistema necesario para determinar sea una matriz de banda con un ancho de banda máximo de 2n-1.

y5 P1

P2

P3

P4

P5

P6

y2

P7

P8

P9

y1

P10

P11

P12

y4 y3

y0

xo

x1

x2

x3

x4

a, b, c, d, m, n, f, g(condiciones de frontera), tol, N

PROCESO

X, Y, v,t, F, w

a: b: c: d: n: m: fun1: fun2: tol: N: w(i,j): h: k: l: norm: X: Y: F: G:

Limite inferior de las abscisas Limite superior de las abscisas Limite inferior de las ordenadas Limite superior de las ordenadas Numero de particiones del intervalo [a,b] Numero de particiones del intervalo [c,d] La función general, g(x,y) La función, f(x,y) Tolerancia para las iteraciones El numero de iteraciones Variable matricial que aproxima a g(x,y) Equidistancia entre los nudos horizontales Equidistancia entre los nudos verticales Contador que asigna el número de iteraciones Variable estático o valor constante Vector fila que representa a las abscisas Vector fila que representa a las ordenadas Matriz que representa a los puntos de f(x,y) Matriz que representa a los puntos de G(x,y)

Para el comprender mejor el programa se ha realizado un diagrama de flujo, mediante el cual este programa puede ser codificado en cualquier lenguaje de programación, facilitando así a la persona interesada en programar este tema.

INICIO LEER:a,b,c,d,n,m,fun,fun1,fun2,fun 3,fun4,tol,N

F=zeros(n+1,m+1); FOR i=1: n+1

h=(b-a)/n; X=zeros(1,n+1); v=zeros(1,n+1);

k=(d-c)/m; Y=zeros(1,m+1); t=zeros(1,m+1);

FOR j=1:m+1 y=Y(j); x=X(i); F(i,j)=eval(fun);

FOR i=1: n+1

NEXT (j)

X(i)=a+(i-1)*h; NEXT (i)

NEXT (i) w=zeros(n-1,m-1); la=(h^2)/(k^2); mu=2*(1+la);

FOR j=1: m+1 Y(i)=a+(i-1)*k; NEXT (j)

WHILE lnorm

V norm=abs(w(i-1,m-1)-z); w(i-1,m-1)=z;

NEXT (i) A

A

B

z=(-h^2*F(n,m)+t(m)+la*v(n)+w(n-2,m-1)+la*w(n-1,m-2))/mu

z=(-h^2*F(n,j)+t(j)+w(n-2,j-1)+la*w(n-1,j)+la*w(n-1,j-2))/mu;

abs(w(n-1,m-1)-z)>norm

F

abs(w(n-1,j-1)-z)>norm

F V

V

norm=abs(w(n-1,m-1)-z); w(n-1,m-1)=z;

norm=abs(w(n-1,j-1)-z); w(n-1,j-1)=z;

for j=m-1:-1:3

NEXT (j)

z=(-h^2*F(2,j)+Y(j)+la*w(1,j)+la*w(1,j-2)+w(2,j-1))/mu

z=(-h^2*F(2,2)+G(1,2)+la*G(2,1)+la*w(1,2)+w(2,1))/mu;

abs(w(1,j-1)-z)>norm

F

abs(w(1,1)-z)>norm

F V

V

norm=abs(w(1,j-1)-z); w(1,j-1)=z;

norm=abs(w(1,1)-z); w(1,1)=z;

FOR i=3:n-1 FOR i=3:n-1

z=(-h^2*F(i,j)+w(i-2,j-1)+la*w(i-1,j)+w(i,j-1)+la*w(i-1,j-2))/mu

z=(-h^2*F(i,2)+la*X(i)+w(i-2,1)+la*w(i-1,2)+w(i,1))/mu;

abs(w(i-1,j-1)-z)>norm

abs(w(i-1,1)-z)>norm

F F V

V

norm=abs(w(i-1,j-1)-z); w(i-1,j-1)=z;

NEXT (i)

norm=abs(w(i-1,1)-z); w(i-1,1)=z;

B

C

NEXT (i)

C

z=(-h^2*F(n,2)+t(2)+la*X(n)+w(n-2,1)+la*w(n-1,2))/mu

abs(w(n-1,j1)-z)>norm

F V norm=abs(w(n-1,1)-z); w(n-1,1)=z;

norm >