EDP El\u00edpticas<\/p>\n\n\n\n
ejercicio<\/a><\/p>\n\n\n\n cont\u00ednua a discreta<\/a><\/p>\n\n\n\n algoritmo discreta<\/a><\/p>\n\n\n\n algoritmo iterativo<\/a><\/p>\n\n\n\n gr\u00e1fica<\/a><\/p>\n<\/div>\n\n\n\n Referencia<\/strong><\/em>: Chapra 29.1 p866, Rodr\u00edguez 10.3 p425, Burden 12.1 p694<\/p>\n\n\n \\frac{\\partial ^2 u}{\\partial x^2} + \\frac{\\partial ^2 u}{ \\partial y^2} = 0<\/span>\n\n\n\n Valores de frontera: T<\/strong>a = 60, T<\/strong>b = 25, T<\/strong>c = 50, T<\/strong>d = 70 (ecuaci\u00f3n de Laplace, Ecuaci\u00f3n de Poisson con f(x,y)=0)<\/p>\n\n\n\n EDP El\u00edpticas<\/p>\n\n\n\n ejercicio<\/a><\/p>\n\n\n\n cont\u00ednua a discreta<\/a><\/p>\n\n\n\n algoritmo discreta<\/a><\/p>\n\n\n\n algoritmo iterativo<\/a><\/p>\n\n\n\n gr\u00e1fica<\/a><\/p>\n<\/div>\n\n\n\n El resultado de las operaciones en la expresi\u00f3n usando el algoritmo es:<\/p>\n\n\n\n El desarrollo anal\u00edtico inicia convirtiendo la ecuaci\u00f3n diferencial parcial el\u00edptica<\/strong> de la forma cont\u00ednua<\/strong> a su representaci\u00f3n discreta<\/strong> usando las expresiones de derivadas en diferencias divididas<\/strong>. El cambio de la expresi\u00f3n se realiza usando Sympy<\/strong>.<\/p>\n\n\n\n La EDP se escribe en formato Sympy en el bloque de ingreso. La ecuaci\u00f3n EDP se compone del lado izquierdo (LHS) y lado derecho (RHS), indicando u<\/strong> como una funci\u00f3n de las variables x,y.<\/p>\n\n\n las diferencias divididas se ingresan como un diccionario:<\/p>\n\n\n Las condiciones iniciales en los bordes, pueden ser funciones matem\u00e1ticas, por lo que se usa el formato lambda. <\/p>\n\n\n\n En el ejercicio b\u00e1sico presentado en la parte te\u00f3rica, los bordes tienen temperaturas constantes, por lo que en \u00e9sta secci\u00f3n se escriben las condiciones como la constante mas cero por la variable independiente, facilitando la evaluaci\u00f3n de un vector xi por ejemplo.<\/p>\n\n\n Los dem\u00e1s par\u00e1metros del ejercicio se ingresan m\u00e1s adelante de forma semejante al ejercicio con Numpy.<\/p>\n\n\n\n Todas las expresiones se escriben al lado izquierdo (LHS) de la igualdad, la parte del lado derecho(RHS) se la prefiere mantener en cero. La expresi\u00f3n se organiza manteniendo el coeficiente en 1 para el termino de Dx de mayor orden. Se sustituyen las derivadas por diferencias divididas, para luego simplificar la expresi\u00f3n usando \u03bb como referencia.<\/p>\n\n\n\n EDP El\u00edpticas<\/p>\n\n\n\n ejercicio<\/a><\/p>\n\n\n\n cont\u00ednua a discreta<\/a><\/p>\n\n\n\n algoritmo discreta<\/a><\/p>\n\n\n\n algoritmo iterativo<\/a><\/p>\n\n\n\n gr\u00e1fica<\/a><\/p>\n<\/div>\n\n\n\n Los resultados al aplicar una operaci\u00f3n a la expresi\u00f3n se guardan en un diccionario, para revisar cada paso intermedio al final<\/p>\n\n\n EDP El\u00edpticas<\/p>\n\n\n\n ejercicio<\/a><\/p>\n\n\n\n cont\u00ednua a discreta<\/a><\/p>\n\n\n\n algoritmo discreta<\/a><\/p>\n\n\n\n algoritmo iterativo<\/a><\/p>\n\n\n\n gr\u00e1fica<\/a><\/p>\n<\/div>\n\n\n\n El desarrollo del m\u00e9todo iterativo <\/strong>para una ecuaci\u00f3n diferencial parcial el\u00edptica<\/strong> usando Sympy<\/strong>, sigue a continuaci\u00f3n del anterior.<\/p>\n\n\n\n Se a\u00f1aden las instrucciones para iterar la expresi\u00f3n discreta para cada i,j dentro de la matriz U creada para el ejercicio. Se aplican las condiciones iniciales a los bordes, para luego proceder a iterar.<\/p>\n\n\n\n Como en las actividades del curso se requiere realizar al menos 3 iteraciones para las expresiones del algoritmo con papel y l\u00e1piz, solo se presentar\u00e1n los resultados para la primera fila en j=0<\/p>\n\n\n\n Los resultados del algoritmo luego del resultado del algoritmo anterior se presentan como:<\/p>\n\n\n\n Se pueden sustituir los valores de las expresiones, evaluando cada u[i,j] y completando la matriz y para la gr\u00e1fica. Esta parte queda como tarea<\/strong>.<\/p>\n\n\n\n Las instrucciones adicionales al algoritmo de EDP el\u00edptica cont\u00ednua a discreta empiezan con la creaci\u00f3n de la matriz u_xy para los valores, los vectores xi,yj y una matriz u_mask que indica si se dispone de un valor calculado o conocido en ese punto (x,y) . Para el m\u00e9todo iterativo, la matriz se rellena con el promedio de los valores m\u00e1ximos de las funciones dadas para los bordes de la placa.<\/p>\n\n\n\n Tarea<\/strong>: El algoritmo desarrolla el c\u00e1lculo para j_k=1, solo la primera fila. Se deja como tarea desarrollar para toda la matriz.<\/p>\n\n\n EDP El\u00edpticas<\/p>\n\n\n\n ejercicio<\/a><\/p>\n\n\n\n cont\u00ednua a discreta<\/a><\/p>\n\n\n\n algoritmo discreta<\/a><\/p>\n\n\n\n algoritmo iterativo<\/a><\/p>\n\n\n\n
\n\n\n\n1. Ejercicio<\/h2>\n\n\n\n
Longitud en x<\/strong>0 = 0, x<\/strong>n = 2, y<\/strong>0 = 0, y<\/strong>n = 1.5
Tama\u00f1o de paso dx = 0.25, dy = dx
iteramax=100, tolera = 0.0001<\/p>\n\n\n\n
\n\n\n\n
\n\n\n\n2. EDP El\u00edpticas cont\u00ednua a discreta<\/h2>\n\n\n\n
EDP El\u00edptica cont\u00ednua a discreta\n ordenDx : 2\n ordenDy : 2\n\n edp=0 :\n 2 2 \n \u2202 \u2202 \n\u2500\u2500\u2500(u(x, y)) + \u2500\u2500\u2500(u(x, y)) = 0\n 2 2 \n\u2202x \u2202y \n K_ : 1\n\n discreta=0 :\n-2\u22c5U(i, j) + U(i, j - 1) + U(i, j + 1) -2\u22c5U(i, j) + U(i - 1, j) + U(i + 1, j)\n\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 2 2 \n Dy Dx \n Lambda_L :\n 2\nDy \n\u2500\u2500\u2500\n 2\nDx \n Lambda L_k : 1\n\n (discreta=0)*Dy**ordeny :\n 2 2 2\n 2\u22c5Dy \u22c5U(i, j) Dy \u22c5U(i - 1, j) Dy \u22c5 U(i + 1, j)\n-2\u22c5U(i, j) + U(i, j - 1) + U(i, j + 1) - \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 2 2 2 \n Dx Dx Dx \n\n discreta_L = 0 :\nL\u22c5U(i - 1, j) + L\u22c5U(i + 1, j) + (-2\u22c5L - 2)\u22c5U(i, j) + U(i, j - 1) + U(i, j + 1) = 0\n\n discreta :\n-4\u22c5U(i, j) + U(i, j - 1) + U(i, j + 1) + U(i - 1, j) + U(i + 1, j) = 0\n\ndiscreta_iterativa\n U(i, j - 1) + U(i, j + 1) + U(i - 1, j) + U(i + 1, j)\nU(i, j) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 \n<\/code><\/pre>\n\n\n\n\n# INGRESO\nfxy = 0*x+0*y # f(x,y) = 0 , ecuacion de Poisson\n# ecuacion: LHS=RHS\nLHS = sym.diff(u,x,2) + sym.diff(u,y,2)\nRHS = fxy\nEDP = sym.Eq(LHS-RHS,0)\n\n<\/pre><\/div>\n\n\n
\n# centrada, centrada, atras\ndif_dividida ={sym.diff(u,x,2): (U.subs(i,i-1)-2*U+U.subs(i,i+1))\/(Dx**2),\n sym.diff(u,y,2): (U.subs(j,j-1)-2*U+U.subs(j,j+1))\/(Dy**2),\n sym.diff(u,y,1): (U - U.subs(j,j-1))\/Dy}\n<\/pre><\/div>\n\n\n\n# Condiciones iniciales en los bordes\nfya = lambda y: 60 +0*y # izquierda\nfyb = lambda y: 25 +0*y # derecha\nfxc = lambda x: 50 +0*x # inferior, funci\u00f3n inicial \nfxd = lambda x: 70 +0*x # superior, funci\u00f3n inicial \n<\/pre><\/div>\n\n\n
\n\n\n\n
\n\n\n\n3. Algoritmo cont\u00ednua a discreta en Python<\/h2>\n\n\n\n
\n# Ecuaciones Diferenciales Parciales Elipticas\n# EDP El\u00edpticas cont\u00ednua a discreta con Sympy\nimport numpy as np\nimport sympy as sym\n\n# u(x,y) funciones continuas y variables simb\u00f3licas usadas\nx = sym.Symbol('x',real=True)\ny = sym.Symbol('y',real=True)\nu = sym.Function('u')(x,y) # funcion\nf = sym.Function('f')(x,y) # funcion complemento\n# U[i,j] funciones discretas y variables simb\u00f3licas usadas\ni = sym.Symbol('i',integer=True,positive=True) # indices\nj = sym.Symbol('j',integer=True,positive=True)\nDx = sym.Symbol('Dx',real=True,positive=True)\nDy = sym.Symbol('Dy',real=True,positive=True)\nL = sym.Symbol('L',real=True)\nU = sym.Function('U')(i,j)\n\n# INGRESO\nfxy = 0*x+0*y # f(x,y) = 0 , ecuacion de Poisson\n# ecuacion edp : LHS=RHS\nLHS = sym.diff(u,x,2) + sym.diff(u,y,2)\nRHS = fxy\nedp = sym.Eq(LHS-RHS,0)\n\n# centrada, centrada, atras\ndif_dividida ={sym.diff(u,x,2): (U.subs(i,i-1)-2*U+U.subs(i,i+1))\/(Dx**2),\n sym.diff(u,y,2): (U.subs(j,j-1)-2*U+U.subs(j,j+1))\/(Dy**2),\n sym.diff(u,y,1): (U - U.subs(j,j-1))\/Dy}\n\n# Condiciones iniciales en los bordes\nfya = lambda y: 60 +0*y # izquierda\nfyb = lambda y: 25 +0*y # derecha\nfxc = lambda x: 50 +0*x # inferior \nfxd = lambda x: 70 +0*x # superior \n\n# dimensiones de la placa\nx0 = 0 # longitud en x\nxn = 2\ny0 = 0 # longitud en y\nyn = 1.5\n# muestreo en ejes, discreto, supone dx=dy\ndx = 0.25 # Tama\u00f1o de paso\ndy = dx # supone dx=dy\niteramax = 100 # revisa convergencia\ntolera = 0.0001\n\nverdigitos = 2 # para mostrar en tabla de resultados\ncasicero = 1e-15 # para redondeo de t\u00e9rminos en ecuacion\n\n# PROCEDIMIENTO\ndef edp_discreta(edp,dif_dividida,x,y,u):\n ''' EDP cont\u00ednua a discreta, usa diferencias divididas\n proporcionadas en par\u00e1metros, indica las variables x,y\n con funci\u00f3n u de (x,y)\n '''\n resultado={}\n # expresi\u00f3n todo a la izquierda LHS (Left Hand side)\n LHS = edp.lhs\n RHS = edp.rhs\n if not(edp.rhs==0):\n LHS = LHS-RHS\n RHS = 0\n edp = sym.Eq(LHS,RHS)\n # orden de derivada por x, y\n edp_x = edp.subs(x,0)\n edp_y = edp.subs(y,0)\n ordenDx = sym.ode_order(edp_x,u)\n ordenDy = sym.ode_order(edp_y,u)\n resultado['ordenDx'] = ordenDx # guarda en resultados\n resultado['ordenDy'] = ordenDy\n # coeficiente derivada orden mayor a 1 (d2u\/dx2)\n coeff_x = edp_coef_Dx(edp,x,ordenDx)\n if not(coeff_x==1):\n LHS = LHS\/coeff_x\n RHS = RHS\/coeff_x\n edp = sym.Eq(LHS,RHS)\n K_ = edp_coef_Dx(edp,y,ordenDy)\n if abs(K_)%1<casicero: # si es entero\n K_ = int(K_)\n resultado['edp=0'] = edp\n resultado['K_'] = K_\n \n discreta = edp.lhs # EDP discreta\n for derivada in dif_dividida: # reemplaza diferencia dividida\n discreta = discreta.subs(derivada,dif_dividida[derivada])\n resultado['discreta=0'] = discreta\n return (resultado)\n\ndef edp_coef_Dx(edp,x,ordenx):\n ''' Extrae el coeficiente de la derivada Dx de ordenx\n edp es la ecuaci\u00f3n como lhs=rhs\n '''\n coeff_x = 1.0 # valor inicial\n # separa cada t\u00e9rmino de suma\n term_suma = sym.Add.make_args(edp.lhs)\n for term_k in term_suma:\n if term_k.is_Mul: # mas de un factor\n factor_Mul = sym.Mul.make_args(term_k)\n coeff_temp = 1; coeff_usar=False\n # separa cada factor de t\u00e9rmino \n for factor_k in factor_Mul:\n if not(factor_k.is_Derivative):\n coeff_temp = coeff_temp*factor_k\n else: # factor con derivada de ordenx\n partes = factor_k.args\n if partes[1]==(x,ordenx):\n coeff_usar = 1\n if coeff_usar==True:\n coeff_x = coeff_x*coeff_temp\n return(coeff_x)\n\ndef redondea_coef(ecuacion, precision=6,casicero = 1e-15):\n ''' redondea coeficientes de t\u00e9rminos suma de una ecuacion\n ecuaci\u00f3n como lhs=rhs\n '''\n tipo = type(ecuacion)\n tipo_eq = False\n if tipo == sym.core.relational.Equality:\n RHS = ecuacion.rhs\n ecuacion = ecuacion.lhs\n tipo = type(ecuacion)\n tipo_eq = True\n\n if tipo == sym.core.add.Add: # t\u00e9rminos suma de ecuacion\n term_sum = sym.Add.make_args(ecuacion)\n ecuacion = sym.S.Zero # vacia\n for term_k in term_sum:\n # factor multiplicativo de termino suma\n term_mul = sym.Mul.make_args(term_k)\n producto = sym.S.One\n for factor in term_mul:\n if not(factor.has(sym.Symbol)): # es numerico\n factor = np.around(float(factor),precision)\n if (abs(factor)%1)<casicero: # si es entero\n factor = int(factor)\n producto = producto*factor\n ecuacion = ecuacion + producto\n \n if tipo == sym.core.mul.Mul: # termino \u00fanico, busca factores\n term_mul = sym.Mul.make_args(ecuacion)\n producto = sym.S.One\n for factor in term_mul:\n if not(factor.has(sym.Symbol)): # es numerico\n factor = np.around(float(factor),precision)\n if (abs(factor)%1)<casicero: # si es entero\n factor = int(factor)\n producto = producto*factor\n ecuacion = producto\n \n if tipo == float: # # solo un numero\n if (abs(ecuacion)%1)<casicero: \n ecuacion = int(ecuacion)\n if tipo_eq==True: # era igualdad, integra lhs=rhs\n ecuacion = sym.Eq(ecuacion,RHS)\n \n return(ecuacion)\n\ndef edp_simplificaLamba(resultado,dx,dy):\n '''simplifica ecuacion con valores de lambda, dx y dy\n entregando la edp discreta simplificada\n '''\n discreta = resultado['discreta=0']\n ordenDy = resultado['ordenDy']\n ordenDx = resultado['ordenDx']\n K_ = resultado['K_']\n lamb = (Dy**ordenDy)\/(Dx**ordenDx)\n if ordenDy==1 and ordenDx==2:\n lamb = lamb\/K_\n resultado['Lambda_L'] = lamb\n # valor de Lambda en ecuacion edp\n L_k = lamb.subs([(Dx,dx),(Dy,dy)])\n if abs(L_k)%1<casicero: # si es entero\n L_k = int(L_k)\n resultado['Lambda L_k'] = L_k\n # simplifica con lambda L\n discreta_L = sym.expand(discreta*(Dy**ordenDy),mul=True)\n resultado['(discreta=0)*Dy**ordeny'] = discreta_L\n discreta_L = edp_sustituye_L(resultado)\n discreta_L = discreta_L.subs(lamb,L)\n discreta_L = sym.collect(discreta_L,U)\n discreta_L = sym.Eq(discreta_L,0)\n resultado['discreta_L = 0'] = discreta_L\n # sustituye constantes en ecuaci\u00f3n a iterar\n discreta_L = discreta_L.subs([(Dx,dx),(Dy,dy),(L,L_k)])\n discreta_L = redondea_coef(discreta_L)\n resultado['discreta'] = discreta_L\n return(resultado)\n\ndef edp_sustituye_L(resultado):\n ''' sustituye lambda con Dy**ordeny\/Dx**x\/K_\n por L, al simplificar Lambda\n '''\n discreta = resultado['(discreta=0)*Dy**ordeny']\n ordenDy = resultado['ordenDy']\n ordenDx = resultado['ordenDx']\n discreta_L = 0\n # separa cada t\u00e9rmino de suma\n term_suma = sym.Add.make_args(discreta)\n for term_k in term_suma:\n # busca partes de L y cambia por valor L\n cambiar = False # por orden de derivada\n if term_k.has(Dx) and term_k.has(Dy):\n partes = term_k.args\n ordeny=1\n ordenx=1\n for unaparte in partes: \n if unaparte.has(Dy):\n if unaparte.is_Pow:\n partey = unaparte.args\n ordeny = partey[1]\n if unaparte.has(Dx):\n if unaparte.is_Pow:\n partey = unaparte.args\n ordenx = partey[1]\n if (ordeny<=ordenDy and ordenx<=-ordenDx):\n cambiar=True\n if cambiar==True:\n term_k = term_k*L\/resultado['Lambda_L']\n discreta_L = discreta_L + term_k\n # simplifica unos con decimal a entero \n discreta_L = redondea_coef(discreta_L)\n return(discreta_L)\n\n# PROCEDIMIENTO -------------------------------\n# transforma edp continua a discreta\nresultado = edp_discreta(edp,dif_dividida,x,y,u)\nresultado = edp_simplificaLamba(resultado,dx,dy)\ndiscreta = resultado['discreta']\n\n# SALIDA\nnp.set_printoptions(precision=verdigitos)\nalgun_numero = [int,float,str,'Lambda L_k']\nprint('EDP El\u00edptica cont\u00ednua a discreta')\nfor entrada in resultado:\n tipo = type(resultado[entrada])\n if tipo in algun_numero or entrada in algun_numero:\n print('',entrada,':',resultado[entrada])\n else:\n print('\\n',entrada,':')\n sym.pprint(resultado[entrada])\n<\/pre><\/div>\n\n\n
\n\n\n\n
\n\n\n\n4. EDP El\u00edpticas - M\u00e9todo iterativo con Sympy-Python<\/h2>\n\n\n\n
discreta_iterativa\n U(i, j - 1) + U(i, j + 1) + U(i - 1, j) + U(i + 1, j)\nU(i, j) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 \n\n iterar i,j:\nj: 1 ; i: 1\n U(0, 1) U(1, 0) U(1, 2) U(2, 1)\nU(1, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(1, 1) = 53.125\nj: 1 ; i: 2\n U(1, 1) U(2, 0) U(2, 2) U(3, 1)\nU(2, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(2, 1) = 51.40625\nj: 1 ; i: 3\n U(2, 1) U(3, 0) U(3, 2) U(4, 1)\nU(3, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(3, 1) = 50.9765625\nj: 1 ; i: 4\n U(3, 1) U(4, 0) U(4, 2) U(5, 1)\nU(4, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(4, 1) = 50.869140625\nj: 1 ; i: 5\n U(4, 1) U(5, 0) U(5, 2) U(6, 1)\nU(5, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(5, 1) = 50.84228515625\nj: 1 ; i: 6\n U(5, 1) U(6, 0) U(6, 2) U(7, 1)\nU(6, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(6, 1) = 50.83<\/strong>55712890625\nj: 1 ; i: 7\n U(6, 1) U(7, 0) U(7, 2) U(8, 1)\nU(7, 1) = \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u2500\u2500\n 4 4 4 4 \nU(7, 1) = 44.27<\/strong>1392822265625\n\nM\u00e9todo iterativo EDP El\u00edptica\ncontinuar el desarrollo con: \nxi: [0. 0.25 0.5 0.75 1. 1.25 1.5 1.75 2. ]\nyj: [0. 0.25 0.5 0.75 1. 1.25 1.5 ]\n j, U[i,j]\n2 [60. 51.25 51.25 51.25 51.25 51.25 51.25 51.25 25. ]\n1 [60. 53.12 51.41 50.98 50.87 50.84 50.84<\/strong> 44.27<\/strong> 25. ]\n0 [50. 50. 50. 50. 50. 50. 50. 50. 50.]\n<\/code><\/pre>\n\n\n\nInstrucciones en Python<\/h3>\n\n\n\n
\n# ITERAR para cada i,j dentro de U ------------\n# x[i] , y[j] valor en posici\u00f3n en cada eje\nxi = np.arange(x0,xn+dx\/2,dx)\nyj = np.arange(y0,yn+dy\/2,dy)\nn = len(xi)\nm = len(yj)\n\n# Matriz U \nu_xy = np.zeros(shape=(n,m),dtype = float)\nu_xy = u_xy*np.nan # valor inicial dentro de u\n# llena u con valores en fronteras\nu_xy[0,:] = fya(yj) # izquierda Ta\nu_xy[n-1,:] = fyb(yj) # derecha Tb\nu_xy[:,0] = fxc(xi) # inferior Tc\nu_xy[:,m-1] = fxd(xi) # superior Td\nu0 = np.copy(u_xy) # matriz u inicial\n\n# u_mask[i,j] con valores iniciales o calculados:True\nu_mask = np.zeros(shape=(n,m),dtype=bool)\nu_mask[0,:] = True # izquierda\nu_mask[-1,:] = True # derecha\nu_mask[:,0] = True # inferior\nu_mask[:,-1] = True # superior\n\n# valor inicial de iteraci\u00f3n dentro de u\n# promedio = (Ta+Tb+Tc+Td)\/4\npromedio = (np.max(fya(yj))+np.max(fyb(yj))+\\\n np.max(fxc(xi))+np.max(fxd(xi)))\/4\nu_xy[1:n-1,1:m-1] = promedio\nu_mask[1:n-1,1:m-1] = True\nu0 = np.copy(u_xy) # copia para revisi\u00f3n\n\ndef edp_sustituyeValorU(discreta,xi,yj,u_xy,u_mask):\n '''Sustituye en edp discreta los valores conocidos de U[i,j]\n tomados desde u_xy, marcados con u_mask\n u_mask indica si el valor se ha calculado con edp.\n '''\n LHS = discreta.lhs # lado izquierdo de ecuacion\n RHS = discreta.rhs # lado derecho\n # sustituye U[i,j] con valores conocidos\n A_diagonal = [] # lista de i,j para matriz de coeficientes A\n # Separa t\u00e9rminos suma\n term_suma = sym.Add.make_args(LHS)\n for term_k in term_suma:\n # busca U[i,j] y cambia por valor uxt[i,j]\n cambiar = False ; cambiar_valor = 0 ; cambiar_factor = 0\n # separa cada factor de t\u00e9rmino\n factor_Mul = sym.Mul.make_args(term_k)\n for factor_k in factor_Mul:\n # busca U[i,j] en matriz uxt[i,j]\n if factor_k.is_Function:\n [i_k,j_k] = factor_k.args\n if not(np.isnan(u_xy[i_k,j_k])):\n cambiar = u_mask[i_k,j_k]\n cambiar_factor = factor_k\n cambiar_valor = u_xy[i_k,j_k]\n else:\n A_diagonal.append([i_k,j_k,term_k\/factor_k])\n # encontr\u00f3 valor U[i,j],term_k va a la derecha de ecuaci\u00f3n\n if cambiar==True:\n LHS = LHS - term_k\n term_ki = term_k.subs(cambiar_factor,cambiar_valor)\n RHS = RHS - term_ki\n discreta = sym.Eq(LHS,RHS)\n B_diagonal = RHS\n resultado = [discreta,A_diagonal,B_diagonal]\n return (resultado)\n\n# M\u00e9todo iterativo para EDP El\u00edptica\n# separa t\u00e9rmino U[j,j] del centro,con valor no conocido\nbuscar = U.subs(j,j)\ndiscreta = sym.solve(discreta,buscar)\ndiscreta = sym.factor_terms(discreta[0])\ndiscreta = sym.Eq(buscar,discreta)\nresultado['discreta_iterativa'] = discreta\n\nprint('\\ndiscreta_iterativa')\nsym.pprint(discreta)\n\n# Desarrollo de iteraciones en cada nodo i,j, para la fila j_k\n# genera el valor en cada punto\n# Nota: Solo la primera fila, Tarea desarrollar para la matriz\nprint('\\n iterar i,j:')\nj_k = 1\nfor i_k in range(1,n-1,1):\n print('j: ',j_k,' ; ','i: ',i_k)\n discreta_ij = discreta.subs({i:i_k,j:j_k,\n x:xi[i_k],y:yj[j_k]})\n sym.pprint(discreta_ij)\n RHS = discreta_ij.rhs\n discreta_ij = sym.Eq(RHS,0)\n # usa valores de frontera segun u_mask con True\n discreta_k = edp_sustituyeValorU(discreta_ij,\n xi,yj,u_xy,u_mask)\n u_xy[i_k,j_k] = sym.Float(-discreta_k[2])\n print(buscar.subs({i:i_k,j:j_k}),\n '=',u_xy[i_k,j_k])\n \n# SALIDA\nnp.set_printoptions(precision=verdigitos)\nprint('\\nM\u00e9todo iterativo EDP El\u00edptica')\nprint('continuar el desarrollo con: ')\nprint('xi:',xi)\nprint('yj:',yj)\nprint(' j, U[i,j]')\nfor j_k in range(2,-1,-1):\n print(j_k, (u_xy[:,j_k]))\n<\/pre><\/div>\n\n\n
\n\n\n\n