⭐ 欢迎来到虫虫下载站! | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站

📄 mmos.m

📁 环境数学模型
💻 M
📖 第 1 页 / 共 2 页
字号:
function edit14_Callback(hObject, eventdata, handles)% hObject    handle to edit14 (see GCBO)% eventdata  reserved - to be defined in a future version of MATLAB% handles    structure with handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of edit14 as text%        str2double(get(hObject,'String')) returns contents of edit14 as a double

% --- Executes during object creation, after setting all properties.
function edit15_CreateFcn(hObject, eventdata, handles)
% hObject    handle to edit15 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end



function edit15_Callback(hObject, eventdata, handles)
% hObject    handle to edit15 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of edit15 as text
%        str2double(get(hObject,'String')) returns contents of edit15 as a double


% --- Executes on button press in pushbutton1.function pushbutton1_Callback(hObject, eventdata, handles)% hObject    handle to pushbutton1 (see GCBO)% eventdata  reserved - to be defined in a future version of MATLAB% handles    structure with handles and user data (see GUIDATA)
theta=str2num(get(handles.edit2,'String'));
T=str2num(get(handles.edit3,'String'));
K1=str2num(get(handles.edit4,'String'));
K2=str2num(get(handles.edit5,'String'));
x10=str2num(get(handles.edit13,'String'));
x20=str2num(get(handles.edit14,'String'));
sigma=str2num(get(handles.edit8,'String'));
Q_s=get(handles.edit9,'String');
Q_plus=str2num(get(handles.edit10,'String'));
Q_minus=str2num(get(handles.edit11,'String'));
Q_0=str2num(get(handles.edit12,'String'));
razb=10;

alpha=[2.4048
5.5201
8.6537
11.7915 
14.9309
18.0711
21.2116
24.3525
27.4935
30.6346
33.7758
36.9171
40.0584
43.1998
46.3412
49.4826
52.6241
55.7655
58.9070
62.0485
65.1900
68.3315
71.4730
74.6145
77.7560
80.8976
84.0391
87.1806
90.3222
93.4637
96.6053
99.7468
102.8884
106.0299
109.1715
112.3131
115.4546
118.5962
121.7377
124.8793
128.0209
131.1624
134.3040
137.4456
140.5872
143.7287
146.8703
150.0119
153.1535
156.2950
159.4366
162.5782
165.7198
168.8613
172.0029
175.1445
178.2861
181.4277
184.5692
187.7108];

ii=find((alpha.*alpha+sigma)*T<700);
M=ii(end);
set(handles.edit15,'String',num2str(M));
alpha=alpha(1:M);


eQ_s=sym('exp((sigma+alpha*alpha)*(t-T))')*sym(Q_s);
ieQ_s=int(eQ_s,'t','0','T');
f=inline(ieQ_s);

r0=sqrt(x10.*x10/K1+x20.*x20/K2);

[r,phi,Alpha]=meshgrid(0:1/50:1,0:2*pi/20:2*pi,alpha);
C=sum(2*pi*r0*besselj(0,Alpha*r0).*besselj(0,Alpha.*r).*f(T,Alpha,sigma),3);

[r,phi]=meshgrid(0:1/50:1,0:2*pi/20:2*pi);
surf(sqrt(K1)*r.*cos(phi),sqrt(K2)*r.*sin(phi),C);

XLabel('x_1');
YLabel('x_2');
ZLabel('C');

rotate3D on;
% --- Executes on button press in pushbutton2.function pushbutton2_Callback(hObject, eventdata, handles)% hObject    handle to pushbutton2 (see GCBO)% eventdata  reserved - to be defined in a future version of MATLAB% handles    structure with handles and user data (see GUIDATA)

theta=str2num(get(handles.edit2,'String'));
T=str2num(get(handles.edit3,'String'));
K1=str2num(get(handles.edit4,'String'));
K2=str2num(get(handles.edit5,'String'));
x10=str2num(get(handles.edit13,'String'));
x20=str2num(get(handles.edit14,'String'));
sigma=str2num(get(handles.edit8,'String'));
Q_s=get(handles.edit9,'String');
Q_plus=str2num(get(handles.edit10,'String'));
Q_minus=str2num(get(handles.edit11,'String'));
Q_0=str2num(get(handles.edit12,'String'));

alpha=[
2.4048
5.5201
8.6537
11.7915
14.9309
18.0711
21.2116
24.3525
27.4935
30.6346
33.7758
36.9171
40.0584
43.1998
46.3412
49.4826
52.6241
55.7655
58.9070
62.0485
65.1900
68.3315
71.4730
74.6145
77.7560
80.8976
84.0391
87.1806
90.3222
93.4637
96.6053
99.7468
102.8884
106.0299
109.1715
112.3131
115.4546
118.5962
121.7377
124.8793
128.0209
131.1624
134.3040
137.4456
140.5872
143.7287
146.8703
150.0119
153.1535
156.2950
159.4366
162.5782
165.7198
168.8613
172.0029
175.1445
178.2861
181.4277
184.5692
187.7108];

ii=find((alpha.*alpha+sigma)*T<700);
M=ii(end);
set(handles.edit15,'String',num2str(M));
alpha=alpha(1:M);

eQ_s=sym('exp((sigma+alpha*alpha)*(t-T))')*sym(Q_s);
ieQ_s=int(eQ_s,'t','0','T');
f=inline(ieQ_s);



r0=sqrt(x10.*x10/K1+x20.*x20/K2);

max00=0;

for T1=0:T/15:T

[r,Alpha]=meshgrid(0:1/10:1,alpha);
C=sum(2*pi*r0*bessel(0,Alpha*r0).*bessel(0,Alpha.*r).*f(T1,Alpha,sigma),2);

max0=max(C);
if max0>max00
    max00=max0;
end;
    
end;

max00=max00/10;

for T1=0:T/15:T

[r,phi,Alpha]=meshgrid(0:1/50:1,0:2*pi/20:2*pi,alpha);
C=sum(2*pi*r0*bessel(0,Alpha*r0).*bessel(0,Alpha.*r).*f(T1,Alpha,sigma),3);

[r,phi]=meshgrid(0:1/50:1,0:2*pi/20:2*pi);
plot3(0,0,max00)
hold on;
plot3(0,0,0)
surf(sqrt(K1)*r.*cos(phi),sqrt(K2)*r.*sin(phi),C);
hold off;
F(round(T1*15)+1)=getframe;
end;

movie(F);
grid on;

XLabel('x_1');
YLabel('x_2');
ZLabel('C');

rotate3D on;
% --- Executes on button press in pushbutton3.function pushbutton3_Callback(hObject, eventdata, handles)% hObject    handle to pushbutton3 (see GCBO)% eventdata  reserved - to be defined in a future version of MATLAB% handles    structure with handles and user data (see GUIDATA)theta=str2num(get(handles.edit2,'String'));
T=str2num(get(handles.edit3,'String'));
K1=str2num(get(handles.edit4,'String'));
K2=str2num(get(handles.edit5,'String'));
x10=str2num(get(handles.edit13,'String'));
x20=str2num(get(handles.edit14,'String'));
sigma=str2num(get(handles.edit8,'String'));
Q_s=get(handles.edit9,'String');
Q_plus=str2num(get(handles.edit10,'String'));
Q_minus=str2num(get(handles.edit11,'String'));
Q_0=str2num(get(handles.edit12,'String'));

alpha=[
2.4048
5.5201
8.6537
11.7915
14.9309
18.0711
21.2116
24.3525
27.4935
30.6346
33.7758
36.9171
40.0584
43.1998
46.3412
49.4826
52.6241
55.7655
58.9070
62.0485
65.1900
68.3315
71.4730
74.6145
77.7560
80.8976
84.0391
87.1806
90.3222
93.4637
96.6053
99.7468
102.8884
106.0299
109.1715
112.3131
115.4546
118.5962
121.7377
124.8793
128.0209
131.1624
134.3040
137.4456
140.5872
143.7287
146.8703
150.0119
153.1535
156.2950
159.4366
162.5782
165.7198
168.8613
172.0029
175.1445
178.2861
181.4277
184.5692
187.7108];

ii=find((alpha.*alpha+sigma)*theta<700);
M=ii(end);
set(handles.edit15,'String',num2str(M));
alpha=alpha(1:M);

r0=sqrt(x10.*x10/K1+x20.*x20/K2);

[r,tau,Alpha]=ndgrid(0:1/30:1,0:theta/30:theta,alpha);

ee=(1-exp((sigma+Alpha.*Alpha).*(tau-theta)))./(sigma+Alpha.*Alpha);
F=sum(2*pi*r0*bessel(0,Alpha*r0).*bessel(0,Alpha.*r).*ee,3);

beta=(theta*Q_plus*Q_plus-Q_0)/(Q_plus*Q_plus-Q_minus*Q_minus);

[r,tau]=ndgrid(0:1/30:1,0:theta/30:theta);
i_minus=find(tau>beta);

F1=F;
F1(i_minus)=0;
F2=F-F1;

Vm=(Q_plus*sum(F2,2)+Q_minus*sum(F1,2))*theta/10;

[r,phi]=ndgrid(0:1/30:1,0:2*pi/30:2*pi);
V=r*0;

for i=1:30+1
V(:,i)=Vm;    
end;    
figure(1);
surf(sqrt(K1)*r.*cos(phi),sqrt(K2)*r.*sin(phi),V);

XLabel('x_1');
YLabel('x_2');
ZLabel('V');

rotate3D on;

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -