Python fsolve() жалуется на форму. Зачем?

Question

Имея функцию f (x, y, z), мне нужно решить ограничение f (x, y, z) = 0, а затем построить его. Я пытался найти для каждой пары (у, г) значения х, для которых F (х, у, г) = 0:

from numpy import * 
from scipy.optimize import fsolve 

def func(x,y,z): 
    return x+y+z 

y = linspace(0,1,100) 
z = linspace(0,1,100) 
x0 = zeros((y.size,z.size)) + 0.5 # the initial guess 
yz = (y[:,newaxis],z[newaxis,:]) # the other parameters 
x, info, iterations, message = fsolve(func,x0,yz) 
contour(y,z,x) 

Python (2.7.5) говорит «TypeError: fsolve: есть несоответствие между формой ввода и вывода аргумента «func» «func».

Но если я проверить это сам, это дает ту же форму:

func(x0,y[:,newaxis],z[:,newaxis]).shape == x0.shape 

возвращает True.

Почему fsolve() жалуется?

Answer 1

fsolve ожидает аргумент x, а возвращаемое значение func является скалярным или одномерным массивом. Вам нужно будет изменить свой код для работы с плоскими значениями x. Например.

def func(x, y, z): 
    x = x.reshape(y.size, z.size) 
    return (x + y + z).ravel() 

и что-то подобное для вызова fsolve:

sol, info, ier, mesg = fsolve(func, x0.ravel(), args=yz, full_output=True) 
x = sol.reshape(y.size, z.size) 

Answer 2

Вот сравнение с krylov метода рекламируемый в scipy.optimize tutorial:

from numpy import linspace, zeros, newaxis 
import time 
from scipy.optimize import root 

def func(x,y,z): 
    x = x.reshape(y.size, z.size) 
    f = x + y + z 
    f = f.ravel() 
    return f 

n = 50 
y = linspace(0,1,n) 
z = linspace(0,1,n) 
x0 = zeros((y.size,z.size)) + 0.5 # the initial guess 
yz = (y[:,newaxis],z[newaxis,:]) # the other parameters 

start = time.time() 
sol1 = root(func, x0.ravel(), args=yz, method='hybr', tol=1e-7) # same as fsolve 
x1 = sol1.x.reshape(y.size, z.size) 
print("(fsolve) time taken (sec): %g" % (time.time() - start,)) 
print("(fsolve) successful: %r (%s)" % (sol1.success, sol1.message)) 
print("(fsolve) max error: %g" % (abs(func(x1, *yz)).max(),)) 

start = time.time() 
sol2 = root(func, x0.ravel(), args=yz, method='krylov', tol=1e-9) 
x2 = sol2.x.reshape(y.size, z.size) 
print("(krylov) time taken (sec): %g" % (time.time() - start,)) 
print("(krylov) successful: %r (%s)" % (sol2.success, sol2.message)) 
print("(krylov) max error: %g" % (abs(func(x2, *yz)).max(),)) 

Печать

 
(fsolve) time taken (sec): 26.9296 
(fsolve) successful: False (The iteration is not making good progress, as measured by the 
    improvement from the last ten iterations.) 
(fsolve) max error: 1.52656e-16 
(krylov) time taken (sec): 0.0173709 
(krylov) successful: True (A solution was found at the specified tolerance.) 
(krylov) max error: 1.11022e-16