DIFERENCIAS FINITAS

Consideremos el problema de hallar la deflexión, en 100 puntos equi−espaciados, de una placa cuadrada, de
2m x 2m y con un espesor de z = 10−2, que se encuentra simplemente apoyada en sus orillas.
La misma esta sometida a una carga q = 33,60 [KN/m2], tiene un módulo de elasticidad de = 2 * 1011 y la
relación de Poisson es = 0,3.
La deflexión en la dimensión Z es regida por la E.D.P. elíptica:z + 2z + z = q
x4 x2y2 y4 D
Donde D = z3 es la rigidez flexional.
12(1−)
El problema puede descomponerse en
1) U = z + z con z = 0 sobre el borde
x2 y2
2) u + u = q con u = 0 sobre el borde
x2 y2 D
Utilizando el método de diferencias finitas, aproximamos primero u = u(x,y) a partir de 2, y luego con el
resultado, obtendremos una aproximación z para la deflexión z = z(x,y).
Tomamos 2 yaproximamos las segundas derivadas parciales en el punto de la malla (x, y); de esta forma
obtenemos:
Uxx (x,y) " U(xi+h,yj) − 2U(xi,yj) + U(xi−h,yj)
h2
Uyy (x,y) " U(xi,yj+h) − 2U(xi,yj) + U(xi,yj−h)
h2
con xi = ih, i = 1, .........., 11
yj = jh, j = 1, .........., 11
h = 0,2
Por lo que a partir de la ecuación de Poisson: U(x,y) = f(x,y), se obtiene
Uxx (xi,yj) + Uyy (xi,yj) = f(x,y) "i,j, de donde se logra obtener el algoritmo de diferencias finitas:

1

Ui+1,j − 2Ui,j + Ui−1,j + Ui,j+1 − 2Uij + Ui,j−1 = q lo que es igual a:
h2 h2 D
Ui−1,j + Ui,j−1 − 4Ui,j +Ui,j+1 + Ui+1,j = qh2
D
De aquí se puede obtener un sistema lineal del tipo AX = b, donde el vector de las incógnitas es de dimensión
121. A es una matriz 121 x 121 compuesto por los factores de U, y b = [qh2/D],también es una matriz 1 x
121.
Por lo que A será la siguiente Matriz:
−4 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 ................0 u11 0
1 −4 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 ................0 u12 0
0 1 −4 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 ................0 u13 0
0 0 1 −4 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 ................0 u14 0
0 0 0 0 −4 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 ................0 u150
0 0 0 0 1 −4 1 0 0 0 0 0 0 0 0 1 0 0 0 0 ................0 u16 0
......... ..... .....
1 0 0 0 0 0 0 −4 1 0 0 0 0 0 0 0 0 0 1 0 ............... 0 uij qh2/D
0 1 0 0 0 0 0 1 −4 1 0 0 0 0 0 0 0 0 0 1 ............... 0 uij+1 qh2/D
............ ....... .......
............... ...... .........
................... ....... ......
0..... 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 −4 1 unn−1 00..... 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 −4 unn 0
Como se puede apreciar la condición de borde u = 0, hace que existan puntos donde el vector incógnita toma
valores iguales a 0, los cuales pueden ser eliminados del sistema, para hacerla más pequeño y por tanto se
necesitarán menos operaciones para resolverlo. Estos puntos son eliminados de la matriz A, ya que en estos
mismos puntosz va a valer cero, por condición de borde antedicha, es decir la deflexión de la placa en sus
bordes es igual a 0.
Así que resumiendo, nuestro problema se reduce a encontrar la solución del sistema:
A81x81 U81x1 = b81x1.

2

Cargando la matriz A y el vector B en Matlab, de la siguiente manera:
A=zeros(81);
for i=1:81
A(i,i)=−4;
resto = mod(i,9);
if resto == 0
posicion=i−1;
ifposicion>0
A(i,posicion)=1;
end
posicion=i+9;
if posicion0
A(i,posicion)=1;
end
end
if resto == 1
posicion=i−9;
if posicion>0
A(i,posicion)=1;
end
posicion=i+1;
if posiciontol)
yanterior = U;
for i=1:81;
y(i)=(L(i,1:81)*U+B(i))./d(i);
U(i)=y(i); %el nuevo valor calculado ya lo sustituyo en el anterior
end
Donde realizando el cambio de y(i) por z(i), y pasando el vector U a sernuestro vector de solución,
obtenemos los valores solución de z expresados en el vector S.
Si graficamos las soluciones obtenidas en el Vector S, obtenemos la siguiente gráfica:

6

Sabemos que Norm(XK+1 − X) " (Norm (Q)/1−Norm(Q)) (Norm(Xk+1 − XK)
Como tenemos Condición de Parada (Norm(Xk+1 − XK) < tol, entonces podemos decir que la Norm(Xk+1 −
X) " (Norm(Q)/1−Norm(Q)) * tol. Es decir que...