1

fsolve を配列に適用しようとしています:

from __future__ import division
from math import fsum
from numpy import *
from scipy.optimize import fsolve
from scipy.constants import pi

nu = 0.05
cn = [0]
cn.extend([pi*n - pi/4 for n in range(1, 5 +1)])
b = linspace(200, 600, 400)
a = fsolve(lambda a: 1/b + fsum([1/(b+cn[i]) for i in range(1, 5 +1)]) - nu*(sqrt(a) - 1)**2, 3)

デフォルトでは許可されていません:

TypeError: only length-1 arrays can be converted to Python scalars

配列に fsolve を適用する方法はありますか?

編集：

#!/usr/bin/env python

from __future__ import division
import numpy as np
from scipy.optimize import fsolve

nu = np.float64(0.05)

cn = np.pi * np.arange(6) - np.pi / 4.
cn[0] = 0.

b = np.linspace(200, 600, 400)

cn.shape = (6,1)
cn_grid = np.repeat(cn, b.size, axis=1)
K = np.sum(1/(b + cn_grid), axis=0)

f = lambda a: K - nu*(np.sqrt(a) - 1)**2
a0 = 3. * np.ones(K.shape)
a = fsolve(f, a0)

print a

それを解決します。

4

2 に答える 2

3

fsumはPythonスカラー用であるため、ベクトル化についてはnumpyを参照する必要があります。5つの数値や単一のnumpy配列ではなく、5つのnumpy配列のリストを合計しようとしているため、メソッドが失敗している可能性があります。

まずcn、numpyを使用して再計算します。

import numpy as np

cn = np.pi * np.arange(6) - np.pi / 4.
cn[0] = 0.

次に、定数ベクトルであるため、前のfsumの結果を個別に計算します。これは1つの方法ですが、より効率的な方法がある場合があります。

cn.shape = (6,1)
cn_grid = np.repeat(cn, b.size, axis=1)
K = np.sum(1/(b + cn_grid), axis=0)

関数を再定義すると、次のように機能するKはずです。

f = lambda a: K - nu*(np.sqrt(a) - 1)**2

解を見つけるために使用fsolveするには、反復する適切な初期ベクトルを提供します。これはゼロベクトルを使用します：

a0 = np.zeros(K.shape)
a = fsolve(f, a0)

または使用できますa0 = 3：

a0 = 3. * np.ones(K.shape)
a = fsolve(f, a0)

この関数は可逆であるため、次f(a) = 0の2つの正確な解を確認できます。

a = (1 - np.sqrt(K/nu))**2

また

a = (1 + np.sqrt(K/nu))**2

fsolveから開始するときに最初の解決策を選択しているようで、。a0 = 0の2番目の解決策を選択しているようですa0 = 3。

于 2012-05-31T08:48:30.947 に答える
1

最小化する関数 (元の関数の 2 乗である必要があります) を定義してから、単純な最小化を使用できます (関数の導関数も定義することをお勧めします)。

funcOrig = lambda a: (K - nu*(np.sqrt(a) - 1)**2)
func2 = lambda a: funcOrig(a)**2
dotfunc2 = lambda a: 2*funcOrig(a)* (-nu * 2 * ( np.sqrt(a)-1) * 1./2./np.sqrt(a))
ret = scipy.optimize.fmin_l_bfgs_b(func2, np.ones(400)+1, fprime=dotfunc2, pgtol=1e-20)
于 2012-05-31T11:49:26.120 に答える