如何求解一个具有两个变量的方程

Question

我想找到满足方程的点集(x，y(

A((x-x0(^2(+2B(x-x0((y-y0(+C((y-y0(^2(=deltachi2

与蟒蛇。这是我的代码

import math
deltachi2=2.3 #2.3 for 1-sigma and 6.18 for 2-sigma
A=3.
B=2.
C=1.1
x0=1.5
y0=3.
x_ini=1.
x_fin=2.
y_ini=2.
y_fin=4.
number_of_x_devisions=100
delta_x= float(x_fin-x_ini)/float(number_of_x_devisions)
number_of_y_devisions=100
delta_y= float(math.log10(y_fin)-math.log10(y_ini))/float(number_of_y_devisions)
filename="test.txt"
outfile=open(filename,"w") 
for i in range(number_of_x_devisions+1):
x_i=x_ini+ i*delta_x
for j in range(number_of_y_devisions+1): 
y_j=y_ini * pow(10,(delta_y*j))
XX=A*pow((x_i-x0),2)+2.*B*(x_i-x0)*(y_j-y0)+C*pow((y_j-y0),2)
if XX==deltachi2:
outfile.write("%1.12e       %1.12e       %1.12en"%(x_i,y_j,XX))

outfile.close()

正如我所说，我需要找到使等高线图为 1 西格玛和 2 西格玛的点集 (x，y((deltachi2 对于 1-sigam 为 2.3，对于 2-sigma 为 6.18(。但是，上面的代码给出了空白输出。如何正确找到轮廓？

Answer 1

上面的代码存在缩进问题;for j循环应嵌套在for i循环内，以对两个"维度"进行扫描。

此外，由于您使用的是浮点数学(并且使用朴素网格搜索而不是更好的算法(，因此可能无法找到确切的解决方案，因此需要使用 epsilon 或精度阈值来查找近似解决方案。

此外，我稍微重构了您的代码，以将x_i和y_j值的生成移动到生成器函数中，因此问题本身更容易阅读。

import math

def x_gen(x_ini, x_fin, num_div):
delta_x = float(x_fin - x_ini) / float(num_div)
for i in range(num_div):
yield x_ini + i * delta_x

def y_gen(y_ini, y_fin, num_div):
delta_y = float(math.log10(y_fin) - math.log10(y_ini)) / float(num_div)
for j in range(num_div + 1):
yield y_ini * pow(10, (delta_y * j))

EPSILON = .01
deltachi2 = 2.3  # 2.3 for 1-sigma and 6.18 for 2-sigma
A = 3.
B = 2.
C = 1.1
x0 = 1.5
y0 = 3.
for x_i in x_gen(1, 2, 100):
for y_j in y_gen(2, 4, 100):
XX = (
A * pow((x_i - x0), 2)
+ 2. * B * (x_i - x0) * (y_j - y0)
+ C * pow((y_j - y0), 2)
)
err = abs(XX - deltachi2)
if err < EPSILON:
print((x_i, y_j, XX, err))

这段代码打印出几个可能的解决方案(我也添加了打印错误(：

(1.0, 2.411615655381521, 2.307584439928287, 0.007584439928287168)
(1.01, 2.394957409237857, 2.3088676681935825, 0.008867668193582645)
(1.02, 2.378414230005442, 2.3096504347952513, 0.009650434795251517)
(1.03, 2.3619853228590606, 2.309936594096949, 0.009936594096949136)
(1.04, 2.3456698984637576, 2.309730056780758, 0.009730056780758112)
(1.05, 2.3294671729369116, 2.309034788099698, 0.009034788099698154)
(1.06, 2.3133763678105748, 2.3078548061624873, 0.00785480616248746)
(1.07, 2.29739670999407, 2.3061941802500714, 0.006194180250071568)
(1.08, 2.2815274317368472, 2.3040570291634124, 0.004057029163412551)
(1.09, 2.265767770591597, 2.3014475196020174, 0.0014475196020176107)
(1.1, 2.250116969377619, 2.298369864572757, 0.0016301354272427737)
(1.11, 2.23457427614444, 2.2948283218284615, 0.0051716781715382965)
(1.12, 2.21913894413569, 2.290827192335828, 0.009172807664171945)
(1.18, 2.1140360811227605, 2.304659088270414, 0.00465908827041428)
(1.19, 2.0994333672461347, 2.2971249106473923, 0.002875089352607496)
(1.23, 2.027918959580058, 2.307983227711846, 0.00798322771184612)
(1.24, 2.0139111001134378, 2.2979409062094645, 0.0020590937905353712)
(1.77, 3.9723699817481437, 2.3089132998333692, 0.008913299833369415)
(1.8, 3.9176811903477073, 2.2975700722470327, 0.0024299277529671315)
(1.88, 3.784230586902383, 2.301749866868204, 0.0017498668682041085)
(1.9100000000000001, 3.7321319661472296, 2.2946153619215264, 0.005384638078473447)
(1.99, 3.6050018504433208, 2.3087335898127352, 0.00873358981273542)