我试图在Python中找到以下方程组的最优解:
(x-x1)^2 + (y-y1)^2 - r1^2 = 0 (x-x2)^2 + (y-y2)^2 - r2^2 = 0 (x-x3)^2 + (y-y3)^2 - r3^2 = 0
给定点(x,y)和半径(r)的值:
x1, y1, r1 = (0, 0, 0.88) x2, y2, r2 = (2, 0, 1) x3, y3, r3 = (0, 2, 0.75)
求点(x,y)的最优解的最佳方法是什么?使用上面的例子,它将是:~(1,1)
lyfkaqu11#
如果我没理解错你的问题,我想这就是你想要的:
from scipy.optimize import minimize import numpy as np def f(coord, x, y, r): return np.sum(((coord[0] - x)**2) + ((coord[1] - y)**2) - (r**2)) x = np.array([0, 2, 0]) y = np.array([0, 0, 2]) r = np.array([.88, 1, .75]) # initial (bad) guess at (x,y) values initial_guess = np.array([100, 100]) res = minimize(f, initial_guess, args=(x, y, r))
得到:
>>> print res.x [0.66666665 0.66666665]
您也可以尝试最小二乘法,该方法期望目标函数返回一个向量。它希望最小化该向量的平方和。使用最小二乘法,您的目标函数将如下所示:
def f2(coord, x, y, r): # notice that we're returning a vector of dimension 3 return ((coord[0] - x)**2) + ((coord[1] - y)**2) - (r**2)
你可以这样最小化它:
from scipy.optimize import leastsq res = leastsq(f2, initial_guess, args=(x, y, r))
>>> print res[0] >>> [0.77958134 0.8580946 ]
这基本上与使用minimize并将原始目标函数重写为相同:
minimize
def f(coord, x, y, r): vec = ((coord[0] - x)**2) + ((coord[1] - y)**2) - (r**2) # return the sum of the squares of the vector return np.sum(vec**2)
这将产生:
>>> print res.x >>> [0.77958326 0.85809648]
请注意,args与leastsq的处理方式稍有不同,并且这两个函数返回的数据结构也不同。有关详细信息,请参阅scipy.optimize.minimize和scipy.optimize.leastsq的文档。请参阅scipy.optimize文件以取得更多最佳化选项。
args
leastsq
scipy.optimize.minimize
scipy.optimize.leastsq
scipy.optimize
zy1mlcev2#
我注意到,在被接受的解决方案中的代码不再工作了......我想也许scipy.optimize在答案发布后改变了它的接口。我可能是错的。无论如何,我支持使用scipy.optimize中的算法的建议,并且被接受的答案确实演示了如何(或者曾经这样做过,如果接口已经改变的话)。我在这里添加了一个额外的答案,纯粹是为了建议一个替代的包,它在核心使用scipy.optimize算法,但对约束优化来说更健壮。这个包是mystic。其中一个重大的改进是mystic提供了约束全局优化。首先,下面是一个示例,其执行方式与scipy.optimize.minimize非常相似,但使用了全局优化器。
mystic
from mystic import reduced @reduced(lambda x,y: abs(x)+abs(y)) #choice changes answer def objective(x, a, b, c): x,y = x eqns = (\ (x - a[0])**2 + (y - b[0])**2 - c[0]**2, (x - a[1])**2 + (y - b[1])**2 - c[1]**2, (x - a[2])**2 + (y - b[2])**2 - c[2]**2) return eqns bounds = [(None,None),(None,None)] #unnecessary a = (0,2,0) b = (0,0,2) c = (.88,1,.75) args = a,b,c from mystic.solvers import diffev2 from mystic.monitors import VerboseMonitor mon = VerboseMonitor(10) result = diffev2(objective, args=args, x0=bounds, bounds=bounds, npop=40, \ ftol=1e-8, disp=False, full_output=True, itermon=mon) print result[0] print result[1]
结果如下所示:
Generation 0 has Chi-Squared: 38868.949133 Generation 10 has Chi-Squared: 2777.470642 Generation 20 has Chi-Squared: 12.808055 Generation 30 has Chi-Squared: 3.764840 Generation 40 has Chi-Squared: 2.996441 Generation 50 has Chi-Squared: 2.996441 Generation 60 has Chi-Squared: 2.996440 Generation 70 has Chi-Squared: 2.996433 Generation 80 has Chi-Squared: 2.996433 Generation 90 has Chi-Squared: 2.996433 STOP("VTRChangeOverGeneration with {'gtol': 1e-06, 'target': 0.0, 'generations': 30, 'ftol': 1e-08}") [ 0.66667151 0.66666422] 2.99643333334
如前所述,reduced中lambda的选择会影响优化器找到的点,因为方程没有实际的解。mystic还提供了将符号方程转换为函数的功能,其中得到的函数可以用作目标函数或罚函数。下面是相同的问题,但使用方程作为罚函数而不是目标函数。
reduced
lambda
def objective(x): return 0.0 equations = """ (x0 - 0)**2 + (x1 - 0)**2 - .88**2 == 0 (x0 - 2)**2 + (x1 - 0)**2 - 1**2 == 0 (x0 - 0)**2 + (x1 - 2)**2 - .75**2 == 0 """ bounds = [(None,None),(None,None)] #unnecessary from mystic.symbolic import generate_penalty, generate_conditions from mystic.solvers import diffev2 pf = generate_penalty(generate_conditions(equations), k=1e12) result = diffev2(objective, x0=bounds, bounds=bounds, penalty=pf, \ npop=40, gtol=50, disp=False, full_output=True) print result[0] print result[1]
结果:
[ 0.77958328 0.8580965 ] 3.6473132399e+12
结果与以前不同,因为应用的惩罚与我们之前在reduced中应用的惩罚不同。在mystic中,您可以选择要应用的惩罚。我们指出这个方程没有解。从上面的结果可以看出,这个结果是严重不利的,所以这是一个很好的迹象,表明没有解。然而,mystic有另一种方法,你可以看到那里没有解。而不是应用一个更传统的penalty,其惩罚违反约束的解... mystic提供constraint,其本质上是核变换,其移除不满足常数的所有潜在解。
penalty
constraint
def objective(x): return 0.0 equations = """ (x0 - 0)**2 + (x1 - 0)**2 - .88**2 == 0 (x0 - 2)**2 + (x1 - 0)**2 - 1**2 == 0 (x0 - 0)**2 + (x1 - 2)**2 - .75**2 == 0 """ bounds = [(None,None),(None,None)] #unnecessary from mystic.symbolic import generate_constraint, generate_solvers, simplify from mystic.symbolic import generate_penalty, generate_conditions from mystic.solvers import diffev2 cf = generate_constraint(generate_solvers(simplify(equations))) result = diffev2(objective, x0=bounds, bounds=bounds, \ constraints=cf, \ npop=40, gtol=50, disp=False, full_output=True) print result[0] print result[1]
[ nan 657.17740835] 0.0
其中,nan实际上表示没有有效的解。顺便说一句,我是作者,所以我有一些偏见。但是,mystic已经存在了几乎和scipy.optimize一样长的时间,是成熟的,并且在这段时间里有一个更稳定的接口。关键是,如果你需要一个更灵活和强大的约束非线性优化器,我建议mystic。
nan
hgtggwj03#
这些方程可以看作是描述了二维空间中三个圆的圆周上的所有点,解就是圆的交点。两个圆的半径之和小于圆心之间的距离,因此两个圆不会重叠。我按比例绘制了两个圆:
没有点满足这个方程组。
bvk5enib4#
我做了一个例子脚本如下.注意最后一行将找到一个最优解(a,b):
import numpy as np import scipy as scp import sympy as smp from scipy.optimize import minimize a,b = smp.symbols('a b') x_ar, y_ar = np.random.random(3), np.random.random(3) x = np.array(smp.symbols('x0:%d'%np.shape(x_ar)[0])) y = np.array(smp.symbols('y0:%d'%np.shape(x_ar)[0])) func = np.sum(a**2+b**2-x*(a+b)+2*y) print func my_func = smp.lambdify((x,y), func) print 1.0/3*my_func(x_ar,y_ar) ab = smp.lambdify((a,b),my_func(x_ar,x_ar)) print ab(1,2) def ab_v(x): return ab(*tuple(x)) print ab_v((1,2)) minimize(ab_v,(0.1,0.1))
输出为:
3*a**2 + 3*b**2 - x0*(a + b) - x1*(a + b) - x2*(a + b) + 2*y0 + 2*y1 + 2*y2 1.0*a**2 - 0.739792011558683*a + 1.0*b**2 - 0.739792011558683*b +0.67394435712335 12.7806239653 12.7806239653 Out[33]: status: 0 success: True njev: 3 nfev: 12 hess_inv: array([[1, 0], [0, 1]]) fun: 3.6178137388030356 x: array([ 0.36989601, 0.36989601]) message: 'Optimization terminated successfully.' jac: array([ 5.96046448e-08, 5.96046448e-08])
4条答案
按热度按时间lyfkaqu11#
如果我没理解错你的问题,我想这就是你想要的:
得到:
您也可以尝试最小二乘法,该方法期望目标函数返回一个向量。它希望最小化该向量的平方和。使用最小二乘法,您的目标函数将如下所示:
你可以这样最小化它:
得到:
这基本上与使用
minimize
并将原始目标函数重写为相同:这将产生:
请注意,
args
与leastsq
的处理方式稍有不同,并且这两个函数返回的数据结构也不同。有关详细信息,请参阅scipy.optimize.minimize
和scipy.optimize.leastsq
的文档。请参阅
scipy.optimize
文件以取得更多最佳化选项。zy1mlcev2#
我注意到,在被接受的解决方案中的代码不再工作了......我想也许
scipy.optimize
在答案发布后改变了它的接口。我可能是错的。无论如何,我支持使用scipy.optimize
中的算法的建议,并且被接受的答案确实演示了如何(或者曾经这样做过,如果接口已经改变的话)。我在这里添加了一个额外的答案,纯粹是为了建议一个替代的包,它在核心使用
scipy.optimize
算法,但对约束优化来说更健壮。这个包是mystic
。其中一个重大的改进是mystic
提供了约束全局优化。首先,下面是一个示例,其执行方式与
scipy.optimize.minimize
非常相似,但使用了全局优化器。结果如下所示:
如前所述,
reduced
中lambda
的选择会影响优化器找到的点,因为方程没有实际的解。mystic
还提供了将符号方程转换为函数的功能,其中得到的函数可以用作目标函数或罚函数。下面是相同的问题,但使用方程作为罚函数而不是目标函数。结果:
结果与以前不同,因为应用的惩罚与我们之前在
reduced
中应用的惩罚不同。在mystic
中,您可以选择要应用的惩罚。我们指出这个方程没有解。从上面的结果可以看出,这个结果是严重不利的,所以这是一个很好的迹象,表明没有解。然而,
mystic
有另一种方法,你可以看到那里没有解。而不是应用一个更传统的penalty
,其惩罚违反约束的解...mystic
提供constraint
,其本质上是核变换,其移除不满足常数的所有潜在解。结果:
其中,
nan
实际上表示没有有效的解。顺便说一句,我是作者,所以我有一些偏见。但是,
mystic
已经存在了几乎和scipy.optimize
一样长的时间,是成熟的,并且在这段时间里有一个更稳定的接口。关键是,如果你需要一个更灵活和强大的约束非线性优化器,我建议mystic
。hgtggwj03#
这些方程可以看作是描述了二维空间中三个圆的圆周上的所有点,解就是圆的交点。
两个圆的半径之和小于圆心之间的距离,因此两个圆不会重叠。我按比例绘制了两个圆:
没有点满足这个方程组。
bvk5enib4#
我做了一个例子脚本如下.注意最后一行将找到一个最优解(a,b):
输出为: