1. 要求
在平面XOZ上,椭圆的表达式为,椭圆曲线绕Z轴旋转一周得到椭球。点B在椭球上运动,已知A、C两点坐标。求B点坐标,满足B到A、C的距离之和最短。
2. 分析
参考
如何求椭圆绕x轴旋转一周后的曲面方程? - 知乎 (zhihu.com)https://www.zhihu.com/question/334080006?可以计算椭球方程如下
记A(ax, ay, az),B(bx, by, bz),C(cx, cy, cz), 则AB、BC的距离可以表示为:
目标函数为
其中满足
3. 代码
首先通过如下代码设置A、C两点的坐标
%% 设置A、C坐标
cx = 0.5; cy = 2; cz = 1;
ax = 0; ay = 0; az = 0;
A = [ax, ay, az];
C = [cx, cy, cz];
接着设置椭圆的参数a和c
%% 设置椭圆参数
a = 1; c = 3;
确定目标函数。用来求解两点之间的距离, fun是匿名函数,输入参数为B,输出为AB和AC的距离之和
%% 目标函数
fun = @(B) norm(A-B) + norm(C-B);
使用fmincon寻找约束非线性多变量函数的最小值,设置约束条件
%% 设置约束条件和初始点
B0 = [0 a 0];
AA = [];
b = [];
Aeq = [];
beq = [];
lb = [];
ub = [];
cc = @(B) [];
ceq = @(B) (B(1)*B(1)+B(2)*B(2))/a^2 + B(3)*B(3)/c^2 - 1;
nonlcon = @(B) deal(cc(B), ceq(B));
最后求解点B的最优位置
%% 求解
B = fmincon(fun,B0,AA,b,Aeq,beq,lb,ub,nonlcon)
完整代码如下
%% 设置A、C坐标
cx = 0.5; cy = 2; cz = 1;
ax = 0; ay = 0; az = 0;
A = [ax, ay, az];
C = [cx, cy, cz];
%% 设置椭圆参数
a = 1; c = 3;
%% 目标函数
fun = @(B) norm(A-B) + norm(C-B);
%% 设置约束条件和初始点
B0 = [0 a 0];
AA = [];
b = [];
Aeq = [];
beq = [];
lb = [];
ub = [];
cc = @(B) [];
ceq = @(B) (B(1)*B(1)+B(2)*B(2))/a^2 + B(3)*B(3)/c^2 - 1;
nonlcon = @(B) deal(cc(B), ceq(B));
%% 求解
B = fmincon(fun,B0,AA,b,Aeq,beq,lb,ub,nonlcon)