Skip to main content
se añadieron 4697 caracteres en el cuerpo
Origen Enlace
FJSevilla
  • 62.1k
  • 7
  • 43
  • 68

También es posible hacerlo de forma recursiva. El código en general

EDICIÓN:

Como te comenté hay algunas cosas que se puede hacer más eficiente, por ejemplo las funciones (p, UI, tpueden mejorar, etc) deberian estar declaradas fuera de calcular() para no volveraquí dejo el código con algunos cambios que comentaré después:

import itertools
import numpy as np
from scipy import optimize
from scipy import integrate

def p(x,b,n,y,s): return 1.0/(x*(b*n+y*np.log(x)-s*b*x))
def Ui(x,x1,b,n,y,s): return (integrate.fixed_quad(p, x, x1, args=(b,n,y,s), n=5)[0]-1.0)
def t(x,x1,b,n,y,s): return optimize.newton(Ui, x, args=(x1,b,n,y,s),  maxiter=5000)
def ii(x,b,n,y,s): return ((y/b)*np.log(x)-s*x+n)
def di(x,b,n,y,s): return (s*x-(y/b))*(y*np.log(x)+n*b-s*b*x)
def p2(x,b,n,y,s): return (s*x-n)/(x*((y*np.log(x)+n*b-s*b*x)**2))
def pe(x,b,n,y,s): return (-np.log(x))/(x*((y*np.log(x)+n*b-s*b*x)**2))

def calcular(B):
    LA= np.array(np.array(B)[:,0])
    y= LA[0]
    LB= np.array(B[1])[0]
    b= LB[0]

    x=0.99;  x1=1.0
    X=[1]                
    for _ in range(1,15): x1=t(x,x1, b, n, y, s); x=x1; X.append(x)
    Xu = X[1:]

    Y=X[0:2];C=X[1:3];D=X[2:4];F=X[3:5];G=X[4:6];H=X[5:7];J=X[6:8]
    K=X[7:9];L=X[8:10];Q=X[9:11];W=X[10:12];E=X[11:13];R=X[12:14];T=X[13:15]

    fi=[ii(x, b, n, y, s) for x in Xu]
    iteorica=np.matrix([[fi[0]],[fi[2]],[fi[3]],[fi[4]],[fi[5]],[fi[6]],[fi[6]],[fi[7]],[fi[8]],[fi[9]],[fi[10]],[fi[11]],[fi[12]],[fi[13]]])
    vdi=[di(x, b, n, y, s)  for x in Xu]
    
    def taobeta(N,M): return integrate.fixed_quad(p2, M, N, args=(b,n,y,s), n=10)[0]
    valtaobeta= list(itertools.starmap(taobeta, [(Y),(C),(D),(F),(G),(H),(J),(K),(L),(Q),(W),(E),(R),(T)]))
    sumtaobeta = np.cumsum(valtaobeta)
    multib=[vdi[i]*sumtaobeta[i] for i in range(len(vdi))]
    
    def taogamma(N,M): return integrate.fixed_quad(pe, M, N, args=(b,n,y,s), n=10)[0]
    valtaogamma=list(itertools.starmap(taogamma, [(Y),(C),(D),(F),(G),(H),(J),(K),(L),(Q),(W),(E),(R),(T)]))
    sumtaogamma=np.cumsum(valtaogamma)
    g=[vdi[i]*sumtaogamma[i] for i in range(len(vdi))]
    
    jacobiana = np.matrix(np.column_stack((multib,g)))
    return(iteorica, jacobiana)

n=763.0
s=762.0     
itecho=np.matrix([[3],[6],[25],[73],[222],[294],[258],[237],[191],[125],[69],[27],[11],[4]])
B = np.matrix([[0.5673],[0.0026]])
for i in range(3):
    print('Iteración {}, valor de B:\n {}'.format(i, B).decode('UTF-8'))
    iteorica, jacobiana = calcular(B)
    B = (np.dot(np.linalg.inv(jacobiana.T*jacobiana),jacobiana.T*(iteorica-itecho)))+B

Los cambios más significativo son:

  • Limpiar un poco los import como te comenté, es especialmete importante evitar from module import *, sobretodo con módulos extensos y complejos como NumPy o SciPy. Puedes terminar reescribiendo funciones de estos módulos sin querer y con resultados desagradables en general.

  • En muchas de tus operaciones es plausible usar list-comprehensions evitando variables no necesarias y mejorando la eficiencia, por ejemplo:

      fi=[]
      for x in Xu:
          qq=ii(x, b, n, y, s)
          fi.append(qq)
    

Puede quedar simplemente como:

    fi=[ii(x, b, n, y, s) for x in Xu]
  • Haces varias sumas acumuladas usando ciclos for, puedes hacer esto eficientemente usando numpy.cumsum. Si usaras una versión de Python >= 3.2 podrias usar también itertools.acumulate. De esta forma puedes simplificar eficientemente cosas como:

      suma=0
      sumtaobeta=[]
      for i in valtaobeta:
          suma=suma+i
          sumtaobeta.append(suma)
    

Por:

    sumtaobeta = np.cumsum(valtaobeta)
  • Para calcular jacobiana puedes usar np.column_stack, reduciendo:

      jacobiana = matrix([[multib[0],g[0]],[multib[1],g[1]],[multib[2],g[2]],[multib[3],g[3]],
      [multib[4],g[4]],[multib[5],g[5]],[multib[6],g[6]],[multib[7],g[7]],   [multib[8],g[8]],[multib[9],g[9]],
      [multib[10],g[10]],[multib[11],g[11]],[multib[12],g[12]],[multib[13],g[13]]])
    

a declararlas en cada iteraciónsimplemente:

    jacobiana = np.matrix(np.column_stack((multib,g)))
  • Todas las funciones que sean fáciles de sacar de calcular() están ahora fuera, solo quedan dentro taobeta y taogamma porque usas itertools.starmap y no es posible pasar los argumentos b,n,y y s a través de ella. Se podrian sacar también on algunos cambios pero no creo que merezca la pena. Esto evita que cada vez que llamemos a calcular de redefinan todas estas funciones, ahora son globales y no están dentro de calcular. Nota que la segunda definición de p que haces en tu código (la que usa taobeta) ahora se llama p2.

Esto son solo ideas por si te interesa aplicarlas. No quiero modificar demasiadoSeguramente se pueden hacer otras mejoras pero hay conceptos del procedimiento matemático de tu código sin saber que hace cada cosa pero hay algunas cosaslógicamente se me escapan por lo que no termino de entender, por ejemplome atrevo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella. Tampoco deberias importar NumPytocar más sin saber el significado de 4 formas distintas, esto termina siendo confuso y peligroso, sobretodo si usas la forma from module import *cada variable/función.

También es posible hacerlo de forma recursiva. El código en general se puede hacer más eficiente, por ejemplo las funciones (p, UI, t, etc) deberian estar declaradas fuera de calcular() para no volver a declararlas en cada iteración. No quiero modificar demasiado tu código sin saber que hace cada cosa pero hay algunas cosas que no termino de entender, por ejemplo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella. Tampoco deberias importar NumPy de 4 formas distintas, esto termina siendo confuso y peligroso, sobretodo si usas la forma from module import *.

También es posible hacerlo de forma recursiva.

EDICIÓN:

Como te comenté hay algunas cosas que se pueden mejorar, aquí dejo el código con algunos cambios que comentaré después:

import itertools
import numpy as np
from scipy import optimize
from scipy import integrate

def p(x,b,n,y,s): return 1.0/(x*(b*n+y*np.log(x)-s*b*x))
def Ui(x,x1,b,n,y,s): return (integrate.fixed_quad(p, x, x1, args=(b,n,y,s), n=5)[0]-1.0)
def t(x,x1,b,n,y,s): return optimize.newton(Ui, x, args=(x1,b,n,y,s),  maxiter=5000)
def ii(x,b,n,y,s): return ((y/b)*np.log(x)-s*x+n)
def di(x,b,n,y,s): return (s*x-(y/b))*(y*np.log(x)+n*b-s*b*x)
def p2(x,b,n,y,s): return (s*x-n)/(x*((y*np.log(x)+n*b-s*b*x)**2))
def pe(x,b,n,y,s): return (-np.log(x))/(x*((y*np.log(x)+n*b-s*b*x)**2))

def calcular(B):
    LA= np.array(np.array(B)[:,0])
    y= LA[0]
    LB= np.array(B[1])[0]
    b= LB[0]

    x=0.99;  x1=1.0
    X=[1]                
    for _ in range(1,15): x1=t(x,x1, b, n, y, s); x=x1; X.append(x)
    Xu = X[1:]

    Y=X[0:2];C=X[1:3];D=X[2:4];F=X[3:5];G=X[4:6];H=X[5:7];J=X[6:8]
    K=X[7:9];L=X[8:10];Q=X[9:11];W=X[10:12];E=X[11:13];R=X[12:14];T=X[13:15]

    fi=[ii(x, b, n, y, s) for x in Xu]
    iteorica=np.matrix([[fi[0]],[fi[2]],[fi[3]],[fi[4]],[fi[5]],[fi[6]],[fi[6]],[fi[7]],[fi[8]],[fi[9]],[fi[10]],[fi[11]],[fi[12]],[fi[13]]])
    vdi=[di(x, b, n, y, s)  for x in Xu]
    
    def taobeta(N,M): return integrate.fixed_quad(p2, M, N, args=(b,n,y,s), n=10)[0]
    valtaobeta= list(itertools.starmap(taobeta, [(Y),(C),(D),(F),(G),(H),(J),(K),(L),(Q),(W),(E),(R),(T)]))
    sumtaobeta = np.cumsum(valtaobeta)
    multib=[vdi[i]*sumtaobeta[i] for i in range(len(vdi))]
    
    def taogamma(N,M): return integrate.fixed_quad(pe, M, N, args=(b,n,y,s), n=10)[0]
    valtaogamma=list(itertools.starmap(taogamma, [(Y),(C),(D),(F),(G),(H),(J),(K),(L),(Q),(W),(E),(R),(T)]))
    sumtaogamma=np.cumsum(valtaogamma)
    g=[vdi[i]*sumtaogamma[i] for i in range(len(vdi))]
    
    jacobiana = np.matrix(np.column_stack((multib,g)))
    return(iteorica, jacobiana)

n=763.0
s=762.0     
itecho=np.matrix([[3],[6],[25],[73],[222],[294],[258],[237],[191],[125],[69],[27],[11],[4]])
B = np.matrix([[0.5673],[0.0026]])
for i in range(3):
    print('Iteración {}, valor de B:\n {}'.format(i, B).decode('UTF-8'))
    iteorica, jacobiana = calcular(B)
    B = (np.dot(np.linalg.inv(jacobiana.T*jacobiana),jacobiana.T*(iteorica-itecho)))+B

Los cambios más significativo son:

  • Limpiar un poco los import como te comenté, es especialmete importante evitar from module import *, sobretodo con módulos extensos y complejos como NumPy o SciPy. Puedes terminar reescribiendo funciones de estos módulos sin querer y con resultados desagradables en general.

  • En muchas de tus operaciones es plausible usar list-comprehensions evitando variables no necesarias y mejorando la eficiencia, por ejemplo:

      fi=[]
      for x in Xu:
          qq=ii(x, b, n, y, s)
          fi.append(qq)
    

Puede quedar simplemente como:

    fi=[ii(x, b, n, y, s) for x in Xu]
  • Haces varias sumas acumuladas usando ciclos for, puedes hacer esto eficientemente usando numpy.cumsum. Si usaras una versión de Python >= 3.2 podrias usar también itertools.acumulate. De esta forma puedes simplificar eficientemente cosas como:

      suma=0
      sumtaobeta=[]
      for i in valtaobeta:
          suma=suma+i
          sumtaobeta.append(suma)
    

Por:

    sumtaobeta = np.cumsum(valtaobeta)
  • Para calcular jacobiana puedes usar np.column_stack, reduciendo:

      jacobiana = matrix([[multib[0],g[0]],[multib[1],g[1]],[multib[2],g[2]],[multib[3],g[3]],
      [multib[4],g[4]],[multib[5],g[5]],[multib[6],g[6]],[multib[7],g[7]],   [multib[8],g[8]],[multib[9],g[9]],
      [multib[10],g[10]],[multib[11],g[11]],[multib[12],g[12]],[multib[13],g[13]]])
    

a simplemente:

    jacobiana = np.matrix(np.column_stack((multib,g)))
  • Todas las funciones que sean fáciles de sacar de calcular() están ahora fuera, solo quedan dentro taobeta y taogamma porque usas itertools.starmap y no es posible pasar los argumentos b,n,y y s a través de ella. Se podrian sacar también on algunos cambios pero no creo que merezca la pena. Esto evita que cada vez que llamemos a calcular de redefinan todas estas funciones, ahora son globales y no están dentro de calcular. Nota que la segunda definición de p que haces en tu código (la que usa taobeta) ahora se llama p2.

Esto son solo ideas por si te interesa aplicarlas. Seguramente se pueden hacer otras mejoras pero hay conceptos del procedimiento matemático de tu código que lógicamente se me escapan por lo que no me atrevo a tocar más sin saber el significado de cada variable/función.

se añadieron 413 caracteres en el cuerpo
Origen Enlace
FJSevilla
  • 62.1k
  • 7
  • 43
  • 68

También es posible hacerlo de forma recursiva. El código en general se puede hacer más eficiente, por ejemplo las funciones (p, UI, t, etc) deberian estar declaradas fuera de calcular() para no volver a declararlas en cada iteración. No quiero modificar demasiado tu código sin saber que hace cada cosa pero hay algunas cosas que no termino de entender, por ejemplo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella. Tampoco deberias importar NumPy de 4 formas distintas, esto termina siendo confuso y peligroso, sobretodo si usas la forma from module import *.

También es posible hacerlo de forma recursiva. El código en general se puede hacer más eficiente, por ejemplo las funciones (p, UI, t, etc) deberian estar declaradas fuera de calcular() para no volver a declararlas en cada iteración. No quiero modificar demasiado tu código sin saber que hace cada cosa pero hay algunas cosas que no termino de entender, por ejemplo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella.

También es posible hacerlo de forma recursiva. El código en general se puede hacer más eficiente, por ejemplo las funciones (p, UI, t, etc) deberian estar declaradas fuera de calcular() para no volver a declararlas en cada iteración. No quiero modificar demasiado tu código sin saber que hace cada cosa pero hay algunas cosas que no termino de entender, por ejemplo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella. Tampoco deberias importar NumPy de 4 formas distintas, esto termina siendo confuso y peligroso, sobretodo si usas la forma from module import *.

se añadieron 413 caracteres en el cuerpo
Origen Enlace
FJSevilla
  • 62.1k
  • 7
  • 43
  • 68

También es posible hacerlo de forma recursiva. El código en general se puede hacer más eficiente, por ejemplo las funciones (p, UI, t, etc) deberian estar declaradas fuera de calcular() para no volver a declararlas en cada iteración. No quiero modificar demasiado tu código sin saber que hace cada cosa pero hay algunas cosas que no termino de entender, por ejemplo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella.

También es posible hacerlo de forma recursiva.

También es posible hacerlo de forma recursiva. El código en general se puede hacer más eficiente, por ejemplo las funciones (p, UI, t, etc) deberian estar declaradas fuera de calcular() para no volver a declararlas en cada iteración. No quiero modificar demasiado tu código sin saber que hace cada cosa pero hay algunas cosas que no termino de entender, por ejemplo a la función UI necesita 5 parámetros de entrada pero solo usas uno dentro de ella.

Origen Enlace
FJSevilla
  • 62.1k
  • 7
  • 43
  • 68
Loading