Matlab ejemplos

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 y aproximamos lassegundas 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 lograobtener 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 u15 0 0 0 0 0 1 −4 1 0 0 0 0 0 0 0 0 1 00 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 0 0..... 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 00 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 puntos z va a valer cero, por condición de bordeantedicha, 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; if posicion>0 A(i,posicion)=1; end posicion=i+9; if posicion0A(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 ser nuestro vector de solución, obtenemos los valores solución de z expresados en elvector 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 la distancia a la solución real (X) está acotada por la norma de Q y la tolerancia....