Задание:
1.
Рассмотреть стандартные типы окон, включая диалоговые.
Адаптировать окно
(quidetemplate1.fig в каталоге toolbox\matlab\uitools\guitemplates) под какую-либо простую программу.
2.
Создать GUI-окно, обладающее панелью инструментов, меню с тремя уровнями вложенности, которое возможно было бы масштабировать (включая пропорциональное увеличение осей и шрифтов). На панель инструментов вынести по крайней мере три пиктограммы. Предусмотреть небольшой логотип на окне
Указание: Возможности GUIDE здесь ограничены, в функции
OpeningFcn можно предусмотреть команды постановки дополнительных элементов в окно.
3.
Реализовать на базе Примера игру «О, счастливчик» с собственным GUI-интерфейсом.
84
Дополнительные примеры
Матричные комплекснозначные вычисления
Создаем случайную комплекснозначную матрицу и вычисляем несколько значений функции Бесселя.
>> x=rand(5,5)+(0+1i)*rand(5,5) x =
0.6649 + 0.6020i 0.6739 + 0.4222i
0.5485 + 0.8580i 0.7009 + 0.4983i 0.6343 +
0.8983i
0.3654 + 0.2536i 0.9994 + 0.9614i
0.2618 + 0.3358i 0.9623 + 0.4344i 0.8030 +
0.7546i
0.1400 + 0.8735i 0.9616 + 0.0721i
0.5973 + 0.6802i 0.7505 + 0.5625i 0.0839 +
0.7911i
0.5668 + 0.5134i 0.0589 + 0.5534i
0.0493 + 0.0534i 0.7400 + 0.6166i 0.9455 +
0.8150i
0.8230 + 0.7327i 0.3603 + 0.2920i
0.5711 + 0.3567i 0.4319 + 0.1133i 0.9159 +
0.6700i
>> bessel(0.6,x) ans =
0.6513 + 0.2129i 0.5987 + 0.1500i
0.7188 + 0.3544i 0.6288 + 0.1674i 0.7707 +
0.3259i
0.4243 + 0.1477i 0.8920 + 0.1439i
0.3871 + 0.2205i 0.6758 + 0.0800i 0.7521 +
0.2024i
0.5298 + 0.5571i 0.6231 + 0.0137i
0.6576 + 0.2658i 0.6641 + 0.1711i 0.4623 +
0.5360i
85 0.5860 + 0.2142i 0.3515 + 0.4157i
0.1348 + 0.0727i 0.6799 + 0.1904i 0.8091 +
0.1503i
0.7475 + 0.1884i 0.4303 + 0.1697i
0.5443 + 0.1510i 0.4387 + 0.0612i 0.7420 +
0.1374i
Немного меняем ее минор поэлементным умножением:
>> x(2:3,2:3)=[1 2;3 4];x(2:3,2:3)=x(2:3,2:3).*x(2:3,2:3) x =
0.6649 + 0.6020i 0.6739 + 0.4222i
0.5485 + 0.8580i 0.7009 + 0.4983i 0.6343 +
0.8983i
0.3654
+
0.2536i
1.0000 4.0000 0.9623 + 0.4344i 0.8030 +
0.7546i
0.1400
+
0.8735i
9.0000 16.0000 0.7505 + 0.5625i 0.0839 +
0.7911i
0.5668 + 0.5134i 0.0589 + 0.5534i
0.0493 + 0.0534i 0.7400 + 0.6166i 0.9455 +
0.8150i
0.8230 + 0.7327i 0.3603 + 0.2920i
0.5711 + 0.3567i 0.4319 + 0.1133i 0.9159 +
0.6700i
Решения
Лабораторная работа №2, п.1.
3.
% Книга расположена в текущем каталоге 21.xls
%Для корректного задания ссылок проверить в
Excel(2003) стиль ссылок
86
%Убрать флаг
Excel->Сервис->Параметры->Общие-
>Стиль ссылок R1C1
>> W=xlsread('21.xls', 1, 'A1:I5')
W =
1 1 2 NaN 1 3 -1 0 1
3 5 5 NaN 1 -5 2 3 4
-1 2 -2 NaN 2 5 -2 1 11 0 3 1 NaN NaN NaN NaN NaN
NaN
1 4 -11 NaN NaN NaN NaN NaN
NaN
>> B=W(:,1:3),A=W(1:3,5:9)
B =
1 1 2 3 5 5
-1 2 -2 0 3 1 1 4 -11
A =
1 3 -1 0 1 1 -5 2 3 4 2 5 -2 1 11 4.
>> A(2,3)=x*A(2,3)-1, B(3,1)=sin(y)*B(3,1)+1 5.
>> A.*B'
87 ans =
1.0000 9.0000 -0.3530 0 1.0000 1.0000 -25.0000 -5.4000 9.0000 16.0000 4.0000 25.0000 4.0000 1.0000 -
121.0000 6.
>> v=linspace(1,1,2);C=diag([v
1])+diag(x*v,1)+diag(y*v,-1),
C =
1.0000 -0.8500 0 9.7856 1.0000 -0.8500 0 9.7856 1.0000
>> A*B+C ans =
11.6470 17.1500 8.0000
-1.1676 -3.4000 -59.4500 27.2939 79.7856 -86.0000 7.
>> dlmwrite('myfile.txt', ans,
'delimiter',
'\t', 'precision', 6)
8.
>> trace(inv(ans)) ans =
0.1305
88
Лабораторная работа №2, п.2.
1.
% Ваши численные результаты в лабораторной работе могут отличаться
% ввиду случайности значений элементов матриц. Для самопроверки можно
% принудительно ввести именно данные матрицы А и В
>>
A=rand(3)+i*rand(3),
B=expm(A), z=B\ones(3,1)
A =
0.9501 + 0.4447i 0.4860 + 0.9218i
0.4565 + 0.4057i
0.2311 + 0.6154i 0.8913 + 0.7382i
0.0185 + 0.9355i
0.6068 + 0.7919i 0.7621 + 0.1763i
0.8214 + 0.9169i
B =
0.7557 + 1.3942i -1.0258 + 2.5609i -
1.2289 + 0.8812i
-1.7591 + 1.0176i 0.0147 + 1.6752i -
2.2026 + 0.9689i
-1.0854 + 2.3869i -0.8267 + 1.9390i -
0.1626 + 2.1086i z =
-0.0997 - 0.1304i
0.0039 - 0.3194i
-0.0813 + 0.0232i
89 2.
% Тех же результатов можно было бы достичь через команду diag или даже без нее
>>
C=cat(1,cat(2,A*B,B),cat(2,A,A.*B));D1=real(C);D2=
imag(C)
D1 =
-3.1588 -4.8146 -4.4527 0.7557 -
1.0258 -1.2289
-5.2554 -4.8658 -5.4803 -1.7591 0.0147 -2.2026
-5.2456 -5.3915 -5.3599 -1.0854 -
0.8267 -0.1626 0.9501 0.4860 0.4565 0.0980 -
2.8592 -0.9185 0.2311 0.8913 0.0185 -1.0329 -
1.2235 -0.9472 0.6068 0.7621 0.8214 -2.5489 -
0.9718 -2.0669
D2 =
1.1829 3.3544 -0.3721 1.3942 2.5609 0.8812
-0.5754 0.7272 -1.4280 1.0176 1.6752 0.9689 2.8754 2.8557 1.4947 2.3869 1.9390 2.1086 0.4447 0.9218 0.4057 1.6607 0.2990 -0.0963 0.6154 0.7382 0.9355 -0.8474 1.5040 -2.0425 0.7919 0.1763 0.9169 0.5889 1.3320 1.5829 3.
>>
Electr=@(x,y) sum(sum(((x-D1).^2+(y-
90
D2).^2).^(-0.5))); Electr(0,0) ans =
21.4684 4.
% при нехватке памяти/мощности подождать, уменьшить число точек до 100 или 500
>> [X,Y]=meshgrid(linspace(-2,2,1000));
>> Z=[];for k=1:length(X), for n=1:length(Y),
Z(k,n)=Electr(X(1,k),Y(n,1)); end ; end
>> contourf(X,Y,Z, logspace(0,1.8,10));colorbar;
Рисунок 19 – «Яичница-глазунья» электрического потенциала
5.
>> h=4/1000;[U,V] = gradient(Z,h,h);
%
Делаем матрицы разреженными, чтобы стрелки
91 не сливались
>> N=20;NN=1000/N;x=[];y=[];z=[];u=[];v=[];
>> for i=1:N, for j=1:N, m=NN*(i-1)+1;n=NN*(j-
1)+1;… x(i,j)=X(m,n);y(i,j)=Y(m,n);z(i,j)=Z(m,n);u(i,
j)=U(m,n);v(i,j)=V(m,n); end; end
% Строим лишний график ради осей полярных координат h = polar([0 2*pi], [0 2]); delete(h);hold on;
%
График скоростей и контурный, т.е. напряженностей и потенциала (рисунок 20). h=quiver(gca,x,y,z,v,8.0);set(h,'Color','r');c ontour(x,y,z);legend show;
Рисунок 20 - График скоростей, напряженностей и потенциала
6.
%
Строим логарифм потенциала, а цвет рациональнее изменять по обратной величине
(рисунок 21)
>>surf(X,Y,log10(Z),Z.^(-
1),'LineStyle','none');colormap hsv;colorbar
92
Рисунок 21 - Логарифм потенциала
7.
% Оформим в виде М-файла.
% Обратите внимание, что функция Electr вводится внутри цикла
% Потенциал нормируется на минимальное значение, и строится его логарифм function
M=MyCin(a) tt=0:0.05:2*pi;[x,y]=meshgrid(linspace(-
2,2,501));z=[];MM=[]; a=0.5;I=ones(6,1)*linspace(1,6,6); J=I';
D1=6*rand(6)-3;D2=6*rand(6)-3; D10=D1;D20=D2; flex1=@(t) D10.*(1+a*cos(I+t*J)); flex2=@(t) D20.*(1+a*sin(I+t*J));
Hfig=figure(1);Hax=axes();Hsur=surf(x,y,ones(lengt h(x)),
'LineStyle'
,
'none'
); axis manual
;axis([-2 2 -2 2 0 6]);caxis([0 6]);title(
''
); colorbar; colormap jet
; for t=tt,
D1=flex1(t);D2=flex2(t);
93
Electr=@(x,y) sum(sum(((x-D1).^2+(y-
D2).^2).^(-0.5))); for k=1:length(x), for n=1:length(y), z(k,n)=Electr(x(1,k),y(n,1)); end
; end
; z=z-min(min(z))+1;title(strcat(
'Время:
'
,num2str(t),
' сек'
)); set(Hsur,
'ZData'
,log(z),
'CData'
,log(z)); drawnow;
MM=[MM,getframe(Hfig)]; end
; movie2avi(MM,
'mym1.avi'
,
'fps'
,4);
M=1; end
Лабораторная работа №3, п.1.
1. function
Lab3=Lab31(A) rrr=@(t) cos(A*t).^2; zzz=@(t) 1+sqrt(A*rrr(t)).*sin(t); fi=linspace(0,2*pi,201);ro=rrr(fi);Z=zzz(fi);
[X,Y,Z]=pol2cart(fi,ro,Z);figure(1);plot3(X,Y,Z);g rid on
; figure(2);comet3(X,Y,Z);Lab3=0; end
Рисунок 22 – Пример графиков
94 2. function res=myfun1(a,s) if
((isvector(a))&(ischar(s))) switch s(1) case
'a'
w=sprintf(
'Среднее значение:
%f'
,mean(a)); case
'p'
w=sprintf(
'Полином в точке:
%f'
,polyval(a,1.1)); case
'g'
w=sprintf(
'График'
);plot(a); otherwise w=sprintf(
'Неверные данные!'
); end
; disp(w);res=0; else display(
'Неверный формат данных!'
);res=1; end
; end
>> t=10*rand(10,1)-6;myfun1(t,'pbx');
Полином в точке: -50.535823 3.
%Подготовьте xls-файл, столбцы
А,В,С
– предмет, имя, оценка
% Массив ячеек D имеет ту же структуру. Шапки не делать function
M=myfun2(flag) global
Dbase global
Dip function
R=reading()
%Заполните D сами
R=0; end
D={};
%Определение способа ввода - ручного или из Ексел- файла if
(isnumeric(flag))
95 reading(); else
[D1,D2]=xlsread(
'1234.xls'
,1,
'A:C'
);D=cat(2,
D2, num2cell(D1));clear
'D2'
'D1'
; end
%D-столбцы данных, их в массив записей Dbase
Dbase=cell2struct(D,{
'Subject'
,
'Name'
,
'Value'
},2);
%Данные из массива ячеек в вектора
Subject=char(D{:,1});Name=char(D{:,2});
Value=cell2mat(D(:,3)); clear(
'D'
);
[Name, Permut]=sortrows(Name);E=Dbase; for i=1:length(E), Dbase(i)=E(Permut(i)); end
; clear
'E'
;
%База данных пересортирована
%Теперь вычисляем среднее, делаем структуру Dip
[UnN, NN]=unique(Name,
'rows'
);
NN=[0;NN];dip=struct(
'Name'
,{
'xxx'
},
'Merit'
,{5});D
ip=[]; for i=1:size(UnN,1), dip.Name=UnN(i,:);db=0; for j=(NN(i)+1):NN(i+1), db=db+Dbase(j).Value; end
; dip.Merit=db/(NN(i+1)-NN(i));Dip=[Dip;dip]; end
;
%Усложняя:Структуру в вектора, их - на печать
DipCell=(struct2cell(Dip))',
DipCellName=char(DipCell{:,1});DipCellVal=cell2mat
(DipCell(:,2)); for i=1:size(Dip), disp(strcat(
'Студент: '
,DipCellName(i,:),
'
Балл: '
,num2str(DipCellVal(i)))), end
;
M=0; end
Лабораторная работа №3, п.2.
% Перед выполнением работы составить файл
MyMath.m по шаблону п.1 Примера
96 1.
% Внесено изменение в алгоритм. Итерации идут группами, примерно по сто (рисунок 23). function
M1=M1(x,p,flag)
MassOfPoint=@(x) 1; function
Region=Region(x) if
(sum(x.*x)>1) Region=0; else
Region=MassOfPoint(x); end end
N=1;NN=100*2^(p);eps=p*(1e-6); med1=0;med=1;media=[];summa=0; while
((abs(med1- med)>eps)&&(N*sqrt(eps)<1)) for i=1:NN, ksi=rand(1,p);summa=summa+Region(ksi); end
; media(N)=summa/NN;med=med1;med1=mean(media); summa=0;
%если захотим увидеть ход процесса, каждый десятый шаг if
(flag&&(mod(N,10))) [
'процессинг:
'
num2str(med1)
' '
num2str(med)], end
;
N=N+1; end
%вторым
аргументом идет погрешность расчета M1=[med1*((2*x).^p) 1/sqrt(N*NN)]; end
…. case
'1'
tic;Z=[];Z1=[]; for n=1:K, disp(n),res=M1(1,n,false);
Z=[Z;res(1),res(2)];
end
; for n=1:K, disp(n),res=M1(1,n,false);
Z1=[Z1;res(1),res(2)];
end
;
VolSp=@(x) pi^(x/2).*(1./gamma(0.5*x+1));
97
E=Z(:,2)'*ones(K);E1=Z1(:,2)'*ones(K);norm(Z(:,1)-
Z1(:,1)), x=1:K;figure(1);hold on
; plot(x,Z(:,1),
'-r'
);errorbar(x,Z(:,1),E,
'- r'
); plot(x,Z1(:,1),
'- g'
);errorbar(x,Z1(:,1),E1,
'-g'
); fplot(VolSp,[1 K],
'-m'
); legend(
'Монте-Карло1'
,
'dy'
,
'Монте-
Карло2'
,
'dy'
,
'Теория'
); grid on
;legend show
;
M=Z;toc,
>> MyMath('1',12)
… ans =
0.3075
Рисунок 23 – Группы итераций
Elapsed time is 222.764131 seconds. ans =
2.0000 0.0408 3.1900 0.0151
98 4.1597 0.0060 4.9331 0.0025 5.2497 0.0021 5.2142 0.0022 4.6321 0.0023 4.0785 0.0017 3.3158 0.0012 2.6167 0.0016 1.8250 0.0013 1.3400 0.0009 2. function
M2=M2(a,b) w=2*pi/(b-a); function
Fsin=Fsin(x)
Fsin=P5(x).*sin(w*n*x); end
; function
Fcos=Fcos(x)
Fcos=P5(x).*cos(w*n*x); end
; for n=0:10,
A(n+1)=quadl(@Fsin,a,b);B(n+1)=quadl(@Fcos,a,b); end
;
B(1)=B(1)/2;M2=abs(A+i*B)*2/(b-a); end
…
P5=@(x) polyval([5 0 4 3 2 1],x); switch
… case
'2'
roots([5 0 4 3 2 1]),x0=[-1 -0.5 0.5 1];[x0;P5(x0)], cf=abs(fft(P5(linspace(-
1,1,11))'));tf=M2(-1,1)';
[cf tf],
M=2;
>> MyMath('2',12)
99 ans =
0.4421 + 0.9076i
0.4421 - 0.9076i
-0.1798 + 0.5845i
-0.1798 - 0.5845i
-0.5247 ans =
-1.0000 -0.5000 0.5000 1.0000
-7.0000 0.0938 3.4063 15.0000 ans =
24.2000 2.0000 22.9582 3.1683 21.3636 2.6417 18.1306 2.0586 15.9210 1.6313 14.9136 1.3387 14.9136 1.1311 15.9210 0.9776 18.1306 0.8601 21.3636 0.7673 22.9582 0.6924
%
Расхождение левой и правой колонок обусловлено несовпадением частоты
% т.е. особенностью вычисления быстрого Фурье- преобразования
3.
)
x a
6
y
(
z
2
U
))
x
1
(
a
2
xy
(
xz
U
1 2
1
%%
%Тестовая функция
100 function
M3=M3(x,y,z)
A=K(1,1);h=1e-4; m3=((x.^2).*y).*z; function
Dz=Dz(x,y,z) part=@(x,z) A*x.*(z.^2); partZ=cat(4,part(x,z- h),part(x,z),part(x,z+h));
% Матрица по 4-му измерению, ведь x,y,z
- могут быть трехмерными
% Градиент по переменной берется по 2- му измерению partZ=permute(partZ,[1 4 3 2]);
Dz=gradient(partZ,h);Dz=permute(Dz,[1 4 3 2]);
Dz=Dz(:,:,:,2);Dz=reshape(Dz,size(m3)); end
M3=m3-(1+x.^2).*Dz(x,y,z); end
%%
%Рабочая функция function
M4=M4(x,y,z)
A=K(2,1);B=K(2,2);h=1e-4; m4=A*P5(x.*y)+erf(x-y+B.*z); function
Dz=Dz(x,y,z) function part0=part0(x,y,z) p=legendre(5,sin(x.*z+A*y));p=p(1,:,:,:); part0=reshape(p,size(m4)); end
; part=@(x,y,z) z.*part0(x,y,z); partZ=cat(4,part(x,y,z- h),part(x,y,z),part(x,y,z+h)); partZ=permute(partZ,[1 4 3 2]);
Dz=gradient(partZ,h);Dz=permute(Dz,[1 4 3 2]);
Dz=reshape(Dz(:,:,:,2),size(m4)); end
;
M4=m4-B*Dz(x,y,z); end
%%
101
%Строим графики по заданному указателю функции
ММ
function
M5=M5(MM,Lim,Name,flg)
%Ниже - процедура 3Д-графики, flg, если режим теста function
Mygraph3d=Mygraph3d(X,Y,Z,U,aver,s1,flag) if flag figure(
'Name'
,[
'3D-графики: '
Name]);axis(Lim);grid on
;box on
; xlabel(
'x'
);ylabel(
'y'
);zlabel(
'z'
); view(3);light(
'Position'
,[Lim(2)
Lim(4) Lim(6)]);camlight;lighting flat
; else hold on
; end
; p = patch(isosurface(X,Y,Z,U,aver));isonormals(X,Y,Z,U
,p); set(p,
'FaceColor'
,s1,
'EdgeColor'
,
'none'
);
Mygraph3d=gca; end
;
%Построение линий уровня на 2D
function
Mygraph2d=Mygraph2d(level)
%level - либо задаем число линий, либо вектор значений уровней
Hfig=figure(
'Name'
,[
'Семейство линий уровня для лаплассиана '
Name]); axis([Lim(1) Lim(2) Lim(5)
Lim(6)]);hold on
; xlabel(
'x'
);ylabel(
'z'
);LSpec={
'- r'
,
':b'
,
'-.g'
};
%Фиксируем y, строим линии уровня
LU(x,const,z)=const
[X,Z]=meshgrid(x,z); for yy=1:3, n=1+ceil((yy-
1)*N/2);Ly=squeeze(L(n,:,:))';
102
[C,Hcg]=contour(gca,X,Z,Ly,level,LSpec{yy}); set(Hcg,
'DisplayName'
,[
'at y='
num2str(y(n))],
'ShowText'
,
'off'
); if
(yy==1) set(Hcg,
'ShowText'
,
'on'
); end
; end
; grid on
;legend show
;
Mygraph2d=Hfig; end
;
%Расчеты
N=100;x=linspace(Lim(1),Lim(2),N+1);y=linspace(Lim
(3),Lim(4),N+1); z=linspace(Lim(5),Lim(6),N+1);[X,Y,Z]=meshgrid(x,y
,z);U=MM(X,Y,Z);
L=6*del2(U,(Lim(2)-Lim(1))/N,(Lim(4)-
Lim(3))/N,(Lim(6)-Lim(5))/N);
Um=mean(mean(mean(U)));Lm=mean(mean(mean(L))); s{1}=[
'U(r)='
,num2str(Um)];s{2}=[
'd2U(r)='
,num2str
(Lm)];
%Данные и лаплассиан сформированы. Далее графика if flg
Mygraph3d(X,Y,Z,U,Um,
'red'
,true);title(
'Тестовый пример - функция и ее лаплассиан'
);
Mygraph3d(X,Y,Z,L,Lm,
'blue'
,false); legend show
;legend(s,
'Location'
,
'NorthEast'
); hold off
; else
Mygraph3d(X,Y,Z,U,Um,
'green'
,true);title(
'Рабочая функция'
); legend show
;legend(s{1},
'Location'
,
'NorthEast'
);
103
Mygraph3d(X,Y,Z,L,Lm,
'magenta'
,true);title(
'Лаплас сиан рабочей функции'
); legend show
;legend(s{2},
'Location'
,
'NorthEast'
); end
;
%Построены две поверхности.
Теперь строим линии уровня лаплассиана Mygraph2d(3);
M5=0; end
%%
%Тестовая функция function
M3=M3(x,y,z)
A=K(1,1);h=1e-4; m3=((x.^2).*y).*z; function
Dz=Dz(x,y,z) part=@(x,z) A*x.*(z.^2); partZ=cat(4,part(x,z- h),part(x,z),part(x,z+h));
% Матрица по 4-му измерению, ведь x,y,z
- могут быть трехмерными
% Градиент по переменной берется по 2- му измерению partZ=permute(partZ,[1 4 3 2]);
Dz=gradient(partZ,h);Dz=permute(Dz,[1 4 3 2]);
Dz=Dz(:,:,:,2);Dz=reshape(Dz,size(m3)); end
M3=m3-(1+x.^2).*Dz(x,y,z);
%M3=A*sin(z)-x.^3-y.^4;
end
%%
%Рабочая функция function
M4=M4(x,y,z)
A=K(2,1);B=K(2,2);h=1e-4; m4=A*P5(x.*y)+erf(x-y+B.*z); function
Dz=Dz(x,y,z)
104 function part0=part0(x,y,z) p=legendre(5,sin(x.*z+A*y));p=p(1,:,:,:); part0=reshape(p,size(m4)); end
; part=@(x,y,z) z.*part0(x,y,z); partZ=cat(4,part(x,y,z- h),part(x,y,z),part(x,y,z+h)); partZ=permute(partZ,[1 4 3 2]);
Dz=gradient(partZ,h);Dz=permute(Dz,[1 4 3 2]);
Dz=reshape(Dz(:,:,:,2),size(m4)); end
;
M4=m4-B*Dz(x,y,z);
%M4=P5(x.*y);
end
%%
%Строим графики по заданному указателю функции
ММ
function
M5=M5(MM,Lim,Name,flg)
%Ниже - процедура 3Д-графики, flg, если режим теста function
Mygraph3d=Mygraph3d(X,Y,Z,U,aver,s1,flag) if flag figure(
'Name'
,[
'3D-графики: '
Name]);axis(Lim);grid on
;box on
; xlabel(
'x'
);ylabel(
'y'
);zlabel(
'z'
); view(3);light(
'Position'
,[Lim(2)
Lim(4) Lim(6)]);camlight;lighting flat
; else hold on
; end
; p = patch(isosurface(X,Y,Z,U,aver));isonormals(X,Y,Z,U
,p); set(p,
'FaceColor'
,s1,
'EdgeColor'
,
'none'
);
Mygraph3d=gca;
105 end
;
%Построение линий уровня на 2D
function
Mygraph2d=Mygraph2d(level)
%level - либо задаем число линий, либо вектор значений уровней
Hfig=figure(
'Name'
,[
'Семейство линий уровня для лаплассиана '
Name]); axis([Lim(1) Lim(2) Lim(5)
Lim(6)]);hold on
; xlabel(
'x'
);ylabel(
'z'
);LSpec={
'- r'
,
':b'
,
'-.g'
};
%Фиксируем y, строим линии уровня
LU(x,const,z)=const
[X,Z]=meshgrid(x,z); for yy=1:3, n=1+ceil((yy-
1)*N/2);Ly=squeeze(L(n,:,:))';
[C,Hcg]=contour(gca,X,Z,Ly,level,LSpec{yy}); set(Hcg,
'DisplayName'
,[
'at y='
num2str(y(n))],
'ShowText'
,
'off'
); if
(yy==1) set(Hcg,
'ShowText'
,
'on'
); end
; end
; grid on
;legend show
;
Mygraph2d=Hfig; end
;
%Расчеты
N=100;x=linspace(Lim(1),Lim(2),N+1);y=linspace(Lim
(3),Lim(4),N+1); z=linspace(Lim(5),Lim(6),N+1);[X,Y,Z]=meshgrid(x,y
,z);U=MM(X,Y,Z);
L=6*del2(U,(Lim(2)-Lim(1))/N,(Lim(4)-
Lim(3))/N,(Lim(6)-Lim(5))/N);
Um=mean(mean(mean(U)));Lm=mean(mean(mean(L)));
106 s{1}=[
'U(r)='
,num2str(Um)];s{2}=[
'd2U(r)='
,num2str
(Lm)];
%Данные и лаплассиан сформированы. Далее графика if flg
Mygraph3d(X,Y,Z,U,Um,
'red'
,true);title(
'Тестовый пример - функция и ее лаплассиан'
);
Mygraph3d(X,Y,Z,L,Lm,
'blue'
,false); legend show
;legend(s,
'Location'
,
'NorthEast'
); hold off
; else
Mygraph3d(X,Y,Z,U,Um,
'green'
,true);title(
'Рабочая функция'
); legend show
;legend(s{1},
'Location'
,
'NorthEast'
);
Mygraph3d(X,Y,Z,L,Lm,
'magenta'
,true);title(
'Лаплас сиан рабочей функции'
); legend show
;legend(s{2},
'Location'
,
'NorthEast'
); end
;
%Построены две поверхности. Теперь строим линии уровня лаплассиана
Mygraph2d(3);
M5=0; end
… case
'3'
tic;v=[1 2 5 10 0 2];M5(@M3,v,
'Test'
,true); v=[-1 1 0 1 -1 0];M5(@M4,v,
'Work'
,false);
M=3;toc,
>> K=[1 -1;0.1 -1];
>> MyMath('3',K);
Elapsed time is 30.782295 seconds.
107
Рисунок 24 – Тестовая функция
Рисунок 25 – Линии уровня лапласиана тестовой функции
108
Рисунок 26 – Рабочая фукция
Рисунок 27 – Лапласиан рабочей функции
109
Рисунок 28 - Линии уровня лапласиана рабочей функции
4.
%%
function
M6=M6(MM)
N=101;scrsz = [1 1 0.8 0.8].*get(0,
'ScreenSize'
);s=
''
; figure(
'Name'
,
'Simulation Plot
Window'
,
'NumberTitle'
,
'off'
,
'Position'
,scrsz);
%1
subplot(2,2,1);x=linspace(K(3,1)/2,2*K(3,1),N);
[s,error]=sprintf(
'U(x,y0,z0) at y0=%3.2f, z0=%3.2f'
,K(3,2),K(3,3)); if
(func2str(MM)==
'MyMath/M3'
)
[leg,error]=sprintf(
'a=%3.2f'
,K(1,1)); else
[leg,error]=sprintf(
'(a,b)=(%3.2f,%3.2f)'
,K(2,1),K
(2,2)); end
110
U=MM(x,K(3,2),K(3,3));plot(x,U,
'r'
,
'DisplayName'
,l eg); title(s);grid on
;legend(
'Location'
,
'NE'
);legend show
;xlabel(
'x'
);ylabel(
'U(x)'
);
%2
subplot(2,2,2);y=linspace(K(3,2)/2,2*K(3,2),N);[X,
Y]=meshgrid(x,y);U=MM(X,Y,K(3,3));
[s,error]=sprintf(
'U(x,y,z0) at z0=%3.2f'
,K(3,3)); contourf(x,y,U,
'DisplayName'
,leg); title(s);xlabel(
'x'
);ylabel(
'y'
);grid on
;legend(
'Location'
,
'NE'
);legend show
;colorbar;
%3
subplot(2,2,3);z=linspace(K(3,3)/2,2*K(3,3),N);[X,
Y,Z]=meshgrid(x,y,z);U=MM(X,Y,Z); d=linspace(0.5,2,4);xd=K(3,1)*d;yd=K(3,2)*d;zd=K(3
,3)*d(2:3); h = slice(x,y,z,U,xd,yd,zd);set(h,
'FaceColor'
,
'interp'
,
'EdgeColor'
,
'none'
); axis tight
;title(
'Объемный график сечений для U(x,y,z)'
); xlabel(
'x'
);ylabel(
'y'
);zlabel(
'z'
);grid on
;colorbar;
%4
PU=[]; if
(func2str(MM)==
'MyMath/M3'
)
[s,error]=sprintf(
'f(x,a)=U(x,y0,z0) at y0=%3.2f, z0=%3.2f'
,K(3,2),K(3,3)); a=K(1,1);y=K(3,2);z=K(3,3);aa=linspace(a/2,2*a,N); for k=aa,
K(1,1)=k;d=MM(x,y,z);PU=[PU;d]; end
;K(1,1)=a;
[X,Y]=meshgrid(x,aa);subplot(2,2,4);surfc(X,Y,PU,
'
LineStyle'
,
'None'
);
111 xlabel(
'x'
);ylabel(
'a'
);zlabel(
'U(x,a)'
);colorbar;
axis tight
;axis on
;box on
; title(s); else
[s,error]=sprintf(
'f(a,b)=U(x0,y0,z0) at x0=%3.2f, y0=%3.2f, z0=%3.2f'
,K(3,1),K(3,2),K(3,3)); a=K(2,1);b=K(2,2);x=K(3,1);y=K(3,2);z=K(3,3); aa=linspace(a/2,2*a,N);bb=linspace(b/2,2*b,N); for m=aa, d=[];K(2,1)=m;
for n=bb,
K(2,2)=n;d=[d;MM(x,y,z)]; end
; PU=[PU d]; end
;
K(2,1)=a;K(2,2)=b;
[X,Y]=meshgrid(aa,bb);subplot(2,2,4);surfc(X,Y,PU,
'LineStyle'
,
'None'
); xlabel(
'a'
);ylabel(
'b'
);zlabel(
'U(a,b)'
);colorbar;
axis tight
;axis on
;box on
; title(s); end
M6=0; end
… case
'4'
tic;M6(@M3);M6(@M4);M=4;toc,
>> K=[2 0 0;0.2 1 0;-0.5 0.5 1];
>> MyMath('4',K);
Elapsed time is 36.239721 seconds.
112
Рисунок 29 – 3-D график функции
Рисунок 30 -3-D график функции
113 5. case
'5'
N=101;x=linspace(K(2,1),K(2,2),N);U=M3(x,K(3,2),K(
3,3))-0.1*P5(x); plot(x,U);grid on
;
%повторить несколько раз предыдущие две строки
M=5;
>> K=[2 0 0;-1 2 0;0 0.5 1]
K =
2.0000 0 0
-1.0000 2.0000 0 0 0.5000 1.0000
>> MyMath('5',K);
>> K(2,1)=-0.1;K(2,2)=0.1;
>> MyMath('5',K);
>> K(2,1)=-0.03;K(2,2)=-0.02;
>> MyMath('5',K);
>> K(2,1)=-0.024;K(2,2)=-0.023;
>> K(2,1)=-0.024;K(2,2)=-0.023;
>>
>> MyMath('5',K);
>> K(2,1)=-0.0238;K(2,2)=-0.0236;
>> MyMath('5',K);
>> x=-0.0238;
…
U=M3(x,K(3,2),K(3,3));q=polyfit(x,U,5),q=q-0.1*[5 0 4 3 2 1];roots(q),
M=5;
>> K=[2 0 0;-1 2 0;0 0.5 1]
K =
114 2.0000 0 0
-1.0000 2.0000 0 0 0.5000 1.0000
>> MyMath('5',K); q =
-0.0000 0.0000 -4.0000 0.5000 -
4.0000 0.0000 ans =
-0.0321 + 2.7775i
-0.0321 - 2.7775i
0.0440 + 1.0434i
0.0440 - 1.0434i
-0.0238
% Нетрудно видеть, что корни и коэффициенты полинома найдены верно
Скачать 1.64 Mb.