Scipy fmin_cg не повторяется - PullRequest
0 голосов
/ 11 июля 2020

Я читал много сообщений о том, как использовать эту процедуру, но не могу понять, что делаю неправильно, я проверил размеры и входные данные, но я все еще в той же ситуации.

Мне нужно оптимизировать целевую функцию с использованием и без использования ее градиентной (якобианской) матрицы, когда я запускаю код ниже, процедура не повторяется, он говорит:

SciPy
Optimization terminated successfully.
         Current function value: 0.000000
         Iterations: 0
         Function evaluations: 1
         Gradient evaluations: 1

m0  : [  100   500 28000   500]
m: [  100   500 28000   500]

Это код, который у меня есть прямо сейчас:

Первая попытка: это без градиента функции

*
#Data
data = np.loadtxt('http://carina.fcaglp.unlp.edu.ar/~jgomez/academic/mpp/data/grav-cilindro2D-invertir.txt') 

x,d0= data[:,0],data[:,1]  # columnas del dato en cada vector (coordinate,anomaly)
d0=d0*1E-5 # m/s2

def gz_cilindro(m,x):
  #Constant
  G=6.67408e-11  #m^3/kg.s^2
  #
  ctte= 2 * np.pi * G 
  #
  gz=np.zeros(len(x))
  for i in range (len(x)):
    gz[i]= (ctte* m[0] * (m[1]**2) * m[3] )/ ((x[i] - m[2])**2 + m[3]**2) #m/s2
  #
  return gz

def cost(m,d0,x):
  #
  dc=gz_cilindro(m,x)
  #
  #Error vector
  e=np.zeros(len(dc))
  for i in range(len(dc)):
    e[i]=d0[i]-dc[i]
  #
  q = e.T@e
  #
  return q/2  #scalar

#
m0 = np.array([100, 500 , 28000, 500])   #Initial values
print('SciPy')
m = optimize.fmin_cg(cost,m0,fprime= None,args=(d0,x))
print(m)

Вторая попытка: эта вторая версия использует функцию градиента, которая использует те же функции, что и выше :

def jacob(m,x):
    #Constant
    G=6.67408e-11  #m^3/kg.s^2
    fac=2 * np.pi *G         #Factor constante
    A=np.array([fac * m[1]**2 * m[3] / ( (x-m[2])**2 + m[3]**2 ),
              2* fac * m[0]* m[1] * m[3] / ( (x-m[2])**2 + m[3]**2 ), 
              2* fac * m[0]* m[1]**2 * m[3]*(x-m[2]) / (( (x-m[2])**2 + m[3]**2 )**2),
             fac * m[0]* m[1]**2 * ( (x-m[2])**2 - m[3]**2  ) / (( (x-m[2])**2 + m[3]**2 )**2)]).T
    return A


def grad_cost (m,d0,x):
  #
  dc=gz_cilindro(m,x)
  #Data
  N=len(d0)
  #
  #Error
  e=np.zeros(N)
  for i in range (N):
      e[i] = d0[i]-dc[i]
  #
  A=jacob(m,x)
  #Defino grad_q
  grad_q= -A.T@e  
  #
  return grad_q  #1D vector
#
m0 = np.array([100, 500 , 28000, 500])

print('SciPy')
m = optimize.fmin_cg(cost,m0, fprime= grad_cost, args=(d0,x))

Я сделал свой собственный код о самом крутом градиенте, и я хочу запустить этот цикл, чтобы сравнить результаты.

Это мой первый пост здесь, я не знайте, если вам нужна дополнительная информация, если да, дайте мне знать, и я обновлю вопрос.

Снова извините за вопрос, я пытаюсь решить свои собственные проблемы, но этот ...

Большое спасибо, и я прыгаю e у вас хороший день.

...