主程序main.m
close
clc
clear
syms x y
x=0:0.01:2*pi;
y=solve('(3.6*cos(x)+17*cos(y)+1)^2+(3.6*sin(x)+17*sin(y)-16.5)^2-36=0','y');
y=subs(y);
n=length(y)
for i=1:1:2
for j=1:1:n
if(~isreal(y(i,j)))
y(i,j)=0;
end
end
end
x1=x;
y1=y(1,:);
x1(y1==0)=[]
y1(y1==0)=[]
omg1= omega(x1,y1)
plot(x1,y1,'*')
grid on
figure
x2=x
y2=y(2,:)
x2(y2==0)=[]
y2(y2==0)=[]
omg2=omega(x2,y2)
plot(x2,y2,'*')
grid on
figure
plot(x1,omg1,'.',x2,omg2,'.')
grid on
omega函数:
function [omg]=omega(x,y)
%caculate omega
A=3.6*cos(x)+17*cos(y)+1
B=3.6*sin(x)+17*sin(y)-16.5
omg=3.6/17*(B.*cos(x)-A.*sin(x)) ./ (A.*sin(y)-B.*cos(y))
end