%
%%
function varargout = Threebody(varargin)
% THREEBODY MATLAB code for Threebody.fig
% THREEBODY, by itself, creates a new THREEBODY or raises the existing
% singleton*.
%
% H = THREEBODY returns the handle to a new THREEBODY or the handle to
% the existing singleton*.
%
% THREEBODY('CALLBACK',hObject,eventData,handles,...) calls the local
% function named CALLBACK in THREEBODY.M with the given input arguments.
%
% THREEBODY('Property','Value',...) creates a new THREEBODY or raises the
% existing singleton*. Starting from the left, property value pairs are
% applied to the GUI before Threebody_OpeningFcn gets called. An
% unrecognized property name or invalid value makes property application
% pause. All inputs are passed to Threebody_OpeningFcn via varargin.
%
% *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one
% instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES
% Edit the above text to modify the response to menu_help Threebody
% Last Modified by GUIDE v2.5 28-Oct-2017 15:50:30
% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name', mfilename, ...
'gui_Singleton', gui_Singleton, ...
'gui_OpeningFcn', @Threebody_OpeningFcn, ...
'gui_OutputFcn', @Threebody_OutputFcn, ...
'gui_LayoutFcn', [] , ...
'gui_Callback', []);
if nargin && ischar(varargin{1})
gui_State.gui_Callback = str2func(varargin{1});
end
if nargout
[varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT
%%
% --- Executes just before Threebody is made visible.
function Threebody_OpeningFcn(hObject, eventdata, handles, varargin)
%%
% This function has no output args, see OutputFcn.
% hObject handle to figure
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% varargin command line arguments to Threebody (see VARARGIN)
%初始化天体的质量坐标和速度
handles.Star=Model_Init_1();
handles.type=1;
handles.pretype=handles.type;
%运行标志位，判断是否在正在运行
global Is_Running;
%初始化标志位，判断是否在已初始化
global Is_Init;
%默认状态为已初始化且正在运行
Is_Running=true;
Is_Init=true;
handles.output = hObject;
% Update handles structure
guidata(hObject, handles);
% UIWAIT makes Threebody wait for user response (see UIRESUME)
% uiwait(handles.figure1);
%%
% --- Outputs from this function are returned to the command line.
function varargout = Threebody_OutputFcn(hObject, eventdata, handles)
%%
% varargout cell array for returning output args (see VARARGOUT);
% hObject handle to figure
% eventdata reserved - to be defined in a future version of MATLAB
% handles structure with handles and user data (see GUIDATA)
% Get default command line output from handles structure
varargout{1} = handles.output;
%持续运行run()函数，等待中断
run();
%%
%星体运动函数
function run()
%%
global Is_Running;
global Is_Init;
%万有引力常数
G=6.67*10^(-11);
%以0.00001的时间步长迭代
T=0.00001;
%注意死循环
while (1)
handles=guidata(gcf);
%退出
if(isempty(handles))
close all;
return;
end
%已初始化意为一个阶段内的运算结束，将结果显示
if(Is_Init==true)
%获取参数
star=handles.Star;
%--------------------一号星体--------------------
M1=star(1).M;
X1=star(1).X;
Y1=star(1).Y;
Z1=star(1).Z;
U1=star(1).U;
V1=star(1).V;
W1=star(1).W;
%--------------------二号星体--------------------
M2=star(2).M;
X2=star(2).X;
Y2=star(2).Y;
Z2=star(2).Z;
U2=star(2).U;
V2=star(2).V;
W2=star(2).W;
%--------------------三号星体--------------------
M3=star(3).M;
X3=star(3).X;
Y3=star(3).Y;
Z3=star(3).Z;
U3=star(3).U;
V3=star(3).V;
W3=star(3).W;
%清空
cla;
axis([-1000 1010 -1000 1000 -1000 1000]);
%不采用默认坐标样式，后面自己画线
axis off;
%默认视角下的右手系
set(gca,'XDir','reverse');
set(gca,'YDir','reverse');
%画网格
x=-1000:125:1000;
y=-1000:125:1000;
for i=1:length(x)
line([x(i) x(i)],[1000 -1000],[0 0],'Color',[0.1 0.1 0.1],'LineWidth',0.5);
line([1000 -1000],[y(i) y(i)],[0 0],'Color',[0.1 0.1 0.1],'LineWidth',0.5);
end
%画箭头
line([1000 0],[0 0],[0 0],'Color',[1 0 1],'LineWidth',2);
line([0 0],[1000 0],[0 0],'Color',[1 0 1],'LineWidth',2);
line([0 0],[0 0],[1000 0],'Color',[1 0 1],'LineWidth',2);
%标记坐标轴
text(1050,0,0,'X','Color',[1 0 1],'FontSize',12);
text(0,1050,0,'Y','Color',[1 0 1],'FontSize',12);
text(0,0,1050,'Z','Color',[1 0 1],'FontSize',12);
%标记单位长度
text(125,0,0,'125','Color',[1 0 1],'FontSize',12);
text(0,125,0,'125','Color',[1 0 1],'FontSize',12);
%初始设置天体颜色、点型、大小等参数，每次循环都重新设置整个画面
if(M1~=0)
h=line('Color',[1 0 0],'Marker','.','MarkerSize',25,'erasemode','normal');
line(X1,Y1,Z1,'Color',[1 0 0],'Marker','.','MarkerSize',10);
end
if(M2~=0)
i=line('Color',[0 1 0],'Marker','.','MarkerSize',25,'erasemode','normal');
line(X2,Y2,Z2,'Color',[0 1 0],'Marker','.','MarkerSize',10);
end
if(M3~=0)
l=line('Color',[0 0 1],'Marker','.','MarkerSize',25,'erasemode','normal');
line(X3,Y3,Z3,'Color',[0 0 1],'Marker','.','MarkerSize',10);
end
%刷新画面
drawnow;
Is_Init=false;
else
[y,fs] = audioread('Windows Background.wav');
%每计算200次重绘一次
for k=1:200;
%计算距离
R12=sqrt((X1-X2)^2+(Y1-Y2)^2+(Z1-Z2)^2);
R13=sqrt((X1-X3)^2+(Y1-Y3)^2+(Z1-Z3)^2);
R23=sqrt((X2-X3)^2+(Y2-Y3)^2+(Z2-Z3)^2);
%先判断星球会不会碰撞爆炸，相对距离以8位限
if ((R12<=8)&&(M1~=0)&&(M2~=0))
M1=0;
M2=0;
delete(h);
delete(i);
%碰撞音效
sound(y,fs);
end
if ((R23<=8)&&(M3~=0)&&(M3~=0))
M2=0;
M3=0;
delete(i);
delete(l);
sound(y,fs);
end
if ((R13<=8)&&(M1~=0)&&(M3~=0))
M1=0;
M3=0;
delete(h);
delete(l);
sound(y,fs);
end
%万有引力定律
f12=G*M1*M2/R12^3*[X1-X2,Y1-Y2,Z1-Z2];
f23=G*M2*M3/R23^3*[X2-X3,Y2-Y3,Z2-Z3];
f31=G*M3*M1/R13^3*[X3-X1,Y3-Y1,Z3-Z1];
if(Is_Running==true)
t=T;
%停止运动
else
t=0;
end
%星体1
if (M1~=0)
Ax1=(f31(1)-f12(1))/M1;
Ay1=(f31(2)-f12(2))/M1;
Az1=(f31(3)-f12(3))/M1;
%注意，需先计算坐标，后计算速度
X1=X1+U1*t+1/2*Ax1*t^2;
Y1=Y1+V1*t+1/2*Ay1*t^2;
Z1=Z1+W1*t+1/2*Az1*t^2;
U1=U1+Ax1*t;
V1=V1+Ay1*t;
W1=W1+Az1*t;
end
%星体2
if (M2~=0)
Ax2=(f12(1)-f23(1))/M2;
Ay2=(f12(2)-f23(2))/M2;
Az2=(f12(3)-f23(3))/M2;
X2=X2+U2*t+1/2*Ax2*t^2;
Y2=Y2+V2*t+1/2*Ay2*t^2;
Z2=Z2+W2*t+1/2*Az2*t^2;
U2=U2+Ax2*t;
