又遇上了mathematica的bug
本帖最后由 mathematica 于 2020-6-2 13:34 编辑一直都很相信mathematica软件，但是最近还是遇上了这个软件的bug
Clear["Global`*"];(*Clear all variables*)
{xa,ya}={1,3}
{xb,yb}={0,0}
{xc,yc}={4,0}
{xd,yd}={4,2}
k1=(yd-ye)/(xd-xe)
k2=(ya-ye)/(xa-xe)
cons=(k2-k1)/(1+k1*k2)+1//Together//Numerator
ans=NMinimize[{MB+MC+ME,
MB^2==(x-xb)^2+(y-yb)^2&&
MC^2==(x-xc)^2+(y-yc)^2&&
ME^2==(x-xe)^2+(y-ye)^2&&
cons==0&&
MB>=0&&MC>=0&&ME>=0
},{xe,ye,x,y,MB,MC,ME},AccuracyGoal->20,PrecisionGoal->20,WorkingPrecision->30]
这个求解结果不对！
问题背景如下：
{6.76386731323477275775716172609, {xe ->
1.65895224582559307947875038393,
ye -> 2.21070099731101484070201512004,
x -> 0.0000379794082457822341383996359544,
y -> 5.64635171044463173644820455826*10^-6,
MB -> 3.59744981889923378613703373725*10^-7,
MC -> 3.99996202059573747723047807369,
ME -> 2.76390493289405339060330503870}}
这个求解结果根本不对！ mathematica 发表于 2020-6-2 13:13
{6.76386731323477275775716172609, {xe ->
1.65895224582559307947875038393,
ye -> 2.210700997 ...
我用lingo17 64bit win7操作系统，计算
!初始赋值;
xa=1;ya=3;
xb=0;yb=0;
xc=4;yc=0;
xd=4;yd=2;
!E点在圆上;
20 - 6*xe + xe^2 - 8*ye + ye^2=0;
!求三条线段的长度;
MB^2=(x-xb)^2+(y-yb)^2;
MC^2=(x-xc)^2+(y-yc)^2;
ME^2=(x-xe)^2+(y-ye)^2;
!三条线段和的最小值;
min=MB+MC+ME;
求解结果如下：
Global optimal solution found.
Objective value: 5.2947230131116
Objective bound: 5.2947230131116
Infeasibilities: 0.12612133559742E-12
Extended solver steps: 0
Total solver iterations: 405
Elapsed runtime seconds: 0.63
Model Class: QP
Total variables: 7
Nonlinear variables: 7
Integer variables: 0
Total constraints: 5
Nonlinear constraints: 4
Total nonzeros: 16
Nonlinear nonzeros: 15
Variable Value Reduced Cost
XA 1.0000000000000 0.0000000000000
YA 3.0000000000000 0.0000000000000
XB 0.0000000000000 0.0000000000000
YB 0.0000000000000 0.0000000000000
XC 4.0000000000000 0.0000000000000
YC 0.0000000000000 0.0000000000000
XD 4.0000000000000 0.0000000000000
YD 2.0000000000000 0.0000000000000
XE 2.7031122435702 -0.22427084678506E-04
YE 1.7837288838948 0.0000000000000
MB 2.8201390812269 0.0000000000000
X 2.6079489165186 0.0000000000000
Y 1.0732133461211 0.0000000000000
MC 1.7577238421647 -0.57155544394449E-04
ME 0.71686008971999 -0.87320651063205E-05
Row Slack or Surplus Dual Price
1 0.0000000000000 0.0000000000000
2 0.0000000000000 0.0000000000000
3 0.0000000000000 0.92471568638404
4 0.0000000000000 0.38051005460772
5 0.0000000000000 -0.79207705459366
6 0.0000000000000 0.61053538883563
7 0.0000000000000 0.0000000000000
8 0.0000000000000 0.0000000000000
9 -0.10658141036402E-13 0.22360941607654
10 -0.26645352591004E-14 -0.17729622036317
11 -0.12612133559742E-12 -0.28447504993525
12 -0.12323475573339E-13 -0.69749226272013
13 5.2947230131116 -1.0000000000000
本帖最后由 lsr314 于 2020-6-2 14:33 编辑
这道题和https://bbs.emath.ac.cn/thread-17297-1-1.html完全一样，就是求一点使之到两点和一圆的圆心距离之和减去圆的半径取最小值，直线是圆的特殊情况。
作图也完全一致，先做出过ABE（字母以下图为准，下同）的圆D，再做以GF为边的等边三角形GFK，DK与圆D的交点E就是所求点。
注：过ABE的圆D是指满足∠BDA=90°的圆。 lsr314 发表于 2020-6-2 14:23
这道题和https://bbs.emath.ac.cn/thread-17297-1-1.html完全一样，就是求一点使之到两点和一圆的圆心距离 ...
能看出来mathematica为什么bug吗？ lsr314 发表于 2020-6-2 14:23
这道题和https://bbs.emath.ac.cn/thread-17297-1-1.html完全一样，就是求一点使之到两点和一圆的圆心距离 ...
还是老牌优化软件lingo好，mathematica有bug 本帖最后由 chyanog 于 2020-6-2 19:46 编辑
稍作修改即可解决问题。
NMinimize[{MB+MC+ME,MB==Sqrt[(x-xb)^2+(y-yb)^2]&&MC==Sqrt[(x-xc)^2+(y-yc)^2]&&ME==Sqrt[(x-xe)^2+(y-ye)^2]&&cons==0},{xe,ye,x,y,MB,MC,ME},WorkingPrecision->40]
不管是解方程组或求最值，没有必要最好不要用太多的变量，不然算起来可能比较慢
NMinimize[{Sqrt[(x-xb)^2+(y-yb)^2]+Sqrt[(x-xc)^2+(y-yc)^2]+Sqrt[(x-xe)^2+(y-ye)^2],cons==0},{xe,ye,x,y},WorkingPrecision->40]
这样得到的精度更高
B1={0,0};
C1={4,0};
E1={3+√5Cos[θ],4+√5Sin[θ]};
NMinimize&[{x,y}-p], {p,{B1,C1,E1}}], {x,y,θ}, WorkingPrecision->40]
Output: {
5.294723011573140131084264863104979933979,
{
x->2.607891000963563319472106916168673875392,
y->1.073258586982084715783742532559337791977,
θ->-1.703977867866132654704426139257573526657
}
}
补充内容 (2020-6-3 10:52):
FindMinimum[{MB+MC+ME,MB^2==(x-xb)^2+(y-yb)^2&&MC^2==(x-xc)^2+(y-yc)^2&&ME^2==(x-xe)^2+(y-ye)^2&&cons==0&&MB>=0&&MC>=0&&ME>=0},{xe,ye,x,y,MB,MC,ME},MaxIterations->10^4]
这样也行 本帖最后由 王守恩 于 2020-6-3 07:57 编辑
lsr314 发表于 2020-6-2 14:23
这道题和https://bbs.emath.ac.cn/thread-17297-1-1.html完全一样，就是求一点使之到两点和一圆的圆心距离 ...
说得好：求1点使之到3点(2点+1圆弧)距离之和最小
HE+HF+HG=EK=DK-DE
=\(\D\sqrt{(3-2)^2+(4+2\sqrt{3})^2}-\sqrt{5}\)
比较这题与那题：
这题的最小值 EK 是垂直圆弧 BA
那题的最小值 EP 是垂直直线 BC
谢谢 mathematica！找了道好题。 chyanog 发表于 2020-6-2 19:16
稍作修改即可解决问题。
不管是解方程组或求最值，没有必要最好不要用太多的变量，不然算起来可能比较慢 ...
个人认为是mathematica博而不专导致的，
lingo专门搞最优化这一方面，所以可能算法经受的考验比较多一些。
经历的变量考验比较多一些，但是mathematica博而不专，
经历的优化的考验少一些，经历的大的优化考验比较少。
可能mathematica做做数学题还行，真的做专门的实际工程应用，
估计还是lingo更可靠！ chyanog 发表于 2020-6-2 19:16
稍作修改即可解决问题。
不管是解方程组或求最值，没有必要最好不要用太多的变量，不然算起来可能比较慢 ...
chyanog一口气改写了3个正确版本代码,好教材为了研究方便,配了一张附图
计算结果
min=5.294723011573140131041587860996532284563=5.29472,
xe=2.703076611599442267363134882782538481398,
ye=1.783733657382143174433811931411245739003,
xm=2.607891000922405058296016809677233768157,
ym=1.073258587056624216692495174703803557058,
MBm=2.820102740572202118217710274671448817505,
MCm=1.757797331891138312630239889146542099598,
MEm=0.7168229391097997001936376971785413674598