Bài tập lớn Matlab Giải tích 2 - Phần code

Nội dung text: Bài tập lớn Matlab Giải tích 2 - Phần code

1. else, tmp1=[tmp1;tmp1;tmp1;tmp1]; tmp1=T(:,[ii ii ii ii]).*tmp1; tmp2=zeros(4,16*numii); tmp2(:)=tmp1(:); X0=zeros(4,4*numii); X0(:)=sum(tmp2)'; end; X0=X0./(ones(4,1)*X0(4,:)); % transform stop points tmp1 = [(2*stop(:,ii)-start(:,ii)-axm(:,ii))./axr(:,ii);ones(1,numii)]; tmp1 = [tmp1 tmp1 tmp1 tmp1]; if (oneax), Xf=T*tmp1; else, tmp1=[tmp1;tmp1;tmp1;tmp1]; tmp1=T(:,[ii ii ii ii]).*tmp1; tmp2=zeros(4,16*numii); tmp2(:)=tmp1(:); Xf=zeros(4,4*numii); Xf(:)=sum(tmp2)'; end; Xf=Xf./(ones(4,1)*Xf(4,:)); % compute perpendicular directions pixfact = ((limrange(1,ii)./limrange(2,ii)).*(ap(2,ii)./ap(1,ii))).^2; pixfact = [pixfact pixfact pixfact pixfact]; pixfact = [pixfact;1./pixfact]; [dummyval,jj] = max(abs(Xf(1:2,:)-X0(1:2,:))); jj1 = ((1:4)'*ones(1,length(jj))==ones(4,1)*jj); jj2 = ((1:4)'*ones(1,length(jj))==ones(4,1)*(3-jj)); jj3 = jj1(1:2,:); Xf(jj1)=Xf(jj1)+(Xf(jj1)-X0(jj1)==0); %eaj new 2/24/98 Xp = X0; Xp(jj2) = X0(jj2) + ones(sum(jj2(:)),1); Xp(jj1) = X0(jj1) - (Xf(jj2)-X0(jj2))./(Xf(jj1)-X0(jj1)) .* pixfact(jj3); % inverse transform the cross points if (oneax), Xp=invT*Xp; else, tmp1=[Xp;Xp;Xp;Xp]; tmp1=invT(:,[ii ii ii ii]).*tmp1; tmp2=zeros(4,16*numii); tmp2(:)=tmp1(:); Xp=zeros(4,4*numii); Xp(:)=sum(tmp2)'; end; Xp=(Xp(1:3,:)./(ones(3,1)*Xp(4,:))).*axr(:,[ii ii ii ii])+axm(:,[ii ii ii ii]); basecross(:,ii) = Xp(:,0*numii+(1:numii)); tipcross(:,ii) = Xp(:,1*numii+(1:numii)); sbasecross(:,ii) = Xp(:,2*numii+(1:numii)); stipcross(:,ii) = Xp(:,3*numii+(1:numii)); end; % compute all points % compute start points axm11 = [axm axm axm axm axm axm axm axm axm axm axm]; axr11 = [axr axr axr axr axr axr axr axr axr axr axr]; st = [stoppoint tippoint basepoint sbasepoint stippoint startpoint stippoint sbasepoint basepoint tippoint stoppoint]; tmp1 = (st - axm11) ./ axr11; tmp1 = [tmp1; ones(1,size(tmp1,2))]; if (oneax), X0=T*tmp1; else, tmp1=[tmp1;tmp1;tmp1;tmp1]; tmp1=[T T T T T T T T T T T].*tmp1; tmp2=zeros(4,44*narrows); tmp2(:)=tmp1(:); X0=zeros(4,11*narrows); X0(:)=sum(tmp2)'; end; X0=X0./(ones(4,1)*X0(4,:)); % compute stop points tmp1 = ([start tipcross basecross sbasecross stipcross stop stipcross sbasecross basecross tipcross start] - axm11) ./ axr11; tmp1 = [tmp1; ones(1,size(tmp1,2))]; if (oneax), Xf=T*tmp1; else, tmp1=[tmp1;tmp1;tmp1;tmp1]; tmp1=[T T T T T T T T T T T].*tmp1; tmp2=zeros(4,44*narrows); tmp2(:)=tmp1(:); Xf=zeros(4,11*narrows); Xf(:)=sum(tmp2)'; end; Xf=Xf./(ones(4,1)*Xf(4,:)); % compute lengths len0 = len.*((ends==1)|(ends==3)).*tan(tipangle/180*pi); slen0 = len.*((ends==2)|(ends==3)).*tan(tipangle/180*pi); le = [zeros(1,narrows) len0 wid/2 wid/2 slen0 zeros(1,narrows) -slen0 -wid/2 -wid/2 - len0 zeros(1,narrows)];
2. mask=any([ones(1,2+size(zz,2));diff([xx yy zz],[],1)],2); xx=xx(mask); yy=yy(mask); if ~isempty(zz), zz=zz(mask); end; end; % plot the patch or line xyz = {'XData',xx,'YData',yy,'ZData',zz,'Tag',ArrowTag}; if newpatch(k)|newline(k), if newpatch(k), H(k) = patch(xyz{:}); else, H(k) = line(xyz{:}); end; if ~isempty(oldh), arrow_copyprops(oldh(k),H(k)); end; else, if ispatch(k), xyz={xyz{:},'CData',[]}; end; set(H(k),xyz{:}); end; end; if ~isempty(oldh), delete(oldh(oldh~=H)); end; % % additional properties set(H,'Clipping','off'); set(H,{'UserData'},num2cell(ud,2)); if (length(extraprops)>0), set(H,extraprops{:}); end; % handle choosing arrow Start and/or Stop locations if unspecified [H,oldaxlims,errstr] = arrow_clicks(H,ud,x,y,z,ax,oldaxlims); if ~isempty(errstr), error([upper(mfilename) ' got ' errstr]); end; % set the output if (nargout>0), h=H; end; % make sure the axis limits did not change if isempty(oldaxlims), ARROW_AXLIMITS = []; else, lims = get(oldaxlims(:,1),{'XLim','YLim','ZLim'})'; lims = reshape(cat(2,lims{:}),6,size(lims,2)); mask = arrow_is2DXY(oldaxlims(:,1)); oldaxlims(mask,6:7) = lims(5:6,mask)'; ARROW_AXLIMITS = oldaxlims(find(any(oldaxlims(:,2:7)'~=lims)),:); if ~isempty(ARROW_AXLIMITS), warning(arrow_warnlimits(ARROW_AXLIMITS,narrows)); end; end; else, % don't create the patch, just return the data h=x; yy=y; zz=z; end; function out = arrow_defcheck(in,def,prop) % check if we got 'default' values out = in; if ~isstr(in), return; end; if size(in,1)==1 & strncmp(lower(in),'def',3), out = def; elseif ~isempty(prop), error([upper(mfilename) ' does not recognize ''' in(:)' ''' as a valid ''' prop ''' string.']); end; function [H,oldaxlims,errstr] = arrow_clicks(H,ud,x,y,z,ax,oldaxlims)
3. set(fig,'WindowButtonMotionFcn',[mfilename '(''callback'',''motion'');']); % wait for the mouse button to be released try waitfor(fig,'WindowButtonUpFcn',''); catch errstr = lasterr; wasInterrupted = 1; end; end; if ~wasInterrupted, feval(mfilename,'callback','motion'); end; % restore some things set(gcf,oldFigProps,oldFigValue); set(H,oldArrowProps,oldArrowValue); function arrow_callback(varargin) % handle redrawing callbacks if nargin==0, return; end; str = varargin{1}; if ~isstr(str), error([upper(mfilename) ' got an invalid Callback command.']); end; s = lower(str); if strcmp(s,'motion'), % motion callback global ARROW_CLICK_H ARROW_CLICK_PROP ARROW_CLICK_AX ARROW_CLICK_USE_Z feval(mfilename,ARROW_CLICK_H,ARROW_CLICK_PROP,arrow_point(ARROW_CLICK_AX,ARROW_CLICK_USE_Z) ); drawnow; else, error([upper(mfilename) ' does not recognize ''' str(:).' ''' as a valid Callback option.']); end; function out = arrow_point(ax,use_z) % return the point on the given axes if nargin==0, ax=gca; end; if nargin<2, use_z=~arrow_is2DXY(ax)|~arrow_planarkids(ax); end; out = get(ax,'CurrentPoint'); out = out(1,:); if ~use_z, out=out(1:2); end; function [wasKeyPress,wasInterrupted,errstr] = arrow_wfbdown(fig) % wait for button down ignoring object ButtonDownFcn's if nargin==0, fig=gcf; end; errstr = ''; % save ButtonDownFcn values objs = findobj(fig); buttonDownFcns = get(objs,'ButtonDownFcn'); mask=~strcmp(buttonDownFcns,''); objs=objs(mask); buttonDownFcns=buttonDownFcns(mask); set(objs,'ButtonDownFcn',''); % save other figure values figProps = {'KeyPressFcn','WindowButtonDownFcn'}; figValue = get(fig,figProps); % do the real work set(fig,'KeyPressFcn','set(gcf,''KeyPressFcn'','''',''WindowButtonDownFcn'','''');', 'WindowButtonDownFcn','set(gcf,''WindowButtonDownFcn'','''')'); lasterr(''); try waitfor(fig,'WindowButtonDownFcn',''); wasInterrupted = 0; catch wasInterrupted = 1; end wasKeyPress = ~wasInterrupted & strcmp(get(fig,'KeyPressFcn'),''); if wasInterrupted, errstr=lasterr; end;
4. msg = [upper(mfilename) ' changed the axis limits ' msg 'when adding the arrow']; if (narrows>1), msg=[msg 's']; end; out = [msg '.' sprintf('\n') ' Call ' upper(mfilename) ' FIXLIMITS to reset them now.']; function arrow_copyprops(fm,to) % copy line properties to patches props = {'EraseMode','LineStyle','LineWidth','Marker','MarkerSize', 'MarkerEdgeColor','MarkerFaceColor','ButtonDownFcn', 'Clipping','DeleteFcn','BusyAction','HandleVisibility', 'Selected','SelectionHighlight','Visible'}; lineprops = {'Color', props{:}}; patchprops = {'EdgeColor',props{:}}; patch2props = {'FaceColor',patchprops{:}}; fmpatch = strcmp(get(fm,'Type'),'patch'); topatch = strcmp(get(to,'Type'),'patch'); set(to( fmpatch& topatch),patch2props,get(fm( fmpatch& topatch),patch2props)); %p->p set(to(~fmpatch&~topatch),lineprops, get(fm(~fmpatch&~topatch),lineprops )); %l->l set(to( fmpatch&~topatch),lineprops, get(fm( fmpatch&~topatch),patchprops )); %p->l set(to(~fmpatch& topatch),patchprops, get(fm(~fmpatch& topatch),lineprops) ,'FaceColor','none'); %l->p function arrow_props % display further help info about ARROW properties c = sprintf('\n'); disp([c 'ARROW Properties: Default values are given in [square brackets], and other' c ' acceptable equivalent property names are in (parenthesis).' c c ' Start The starting points. For N arrows, B' c ' this should be a Nx2 or Nx3 matrix. /|\ ^' c ' Stop The end points. For N arrows, this /|||\ |' c ' should be a Nx2 or Nx3 matrix. //|||\\ L|' c ' Length Length of the arrowhead (in pixels on ///|||\\\ e|' c ' screen, points on a page). [16] (Len) ////|||\\\\ n|' c ' BaseAngle Angle (degrees) of the base angle /////|D|\\\\\ g|' c ' ADE. For a simple stick arrow, use //// ||| \\\\ t|' c ' BaseAngle=TipAngle. [90] (Base) /// ||| \\\ h|' c ' TipAngle Angle (degrees) of tip angle ABC. // || \\ |' c ' [16] (Tip) / base ||| \ V' c ' Width Width of the base in pixels. Not E angle || C' c ' the ''LineWidth'' prop. [0] (Wid) |||tipangle' c ' Page If provided, non-empty, and not NaN, |||' c ' this causes ARROW to use hardcopy |||' c ' rather than onscreen proportions. A' c ' This is important if screen aspect > ' c ' vastly different. []' c ' CrossDir A vector giving the direction towards which the fletches' c ' on the arrow should go. [computed such that it is perpen-' c ' dicular to both the arrow direction and the view direction' c ' (i.e., as if it was pasted on a normal 2-D graph)] (Note' c ' that CrossDir is a vector. Also note that if an axis is' c ' plotted on a log scale, then the corresponding component' c ' of CrossDir must also be set appropriately, i.e., to 1 for' c ' no change in that direction, >1 for a positive change, >0' c ' and <1 for negative change.)' c ' NormalDir A vector normal to the fletch direction (CrossDir is then' c ' computed by the vector cross product [Line]x[NormalDir]). []' c
5. 'EdgeColor','r','FaceColor','none','LineWidth',2); % Stick arrow h7 = feval(mfilename,[-1.6 -1.65 -6.5],[0 -1.65 -6.5],[],16,16); t4=text(-1.5,-1.65,-7.25,'global mininum'); set(t4,'HorizontalAlignment','center'); % Normal, black fill h8 = feval(mfilename,[-1.4 0 -7.2],[-1.4 0 -3],'FaceColor','k'); t5=text(-1.5,0,-7.75,'local minimum'); set(t5,'HorizontalAlignment','center'); % Gray fill, crossdir specified, 'LineStyle' h9 = feval(mfilename,[-3 2.2 -6],[-3 2.2 -.05],36,[],27,6,[],[0 -1 0], 'EdgeColor','k','FaceColor',.75*[1 1 1],'LineStyle',' '); % a series of normal arrows, linearly spaced, crossdir specified h10y=(0:4)'/3; h10 = feval(mfilename,[-3*ones(size(h10y)) h10y -6.5*ones(size(h10y))], [-3*ones(size(h10y)) h10y -.05*ones(size(h10y))], 12,[],[],[],[],[0 -1 0]); % a series of normal arrows, linearly spaced h11x=(1:.33:2.8)'; h11 = feval(mfilename,[h11x -3*ones(size(h11x)) 6.5*ones(size(h11x))], [h11x -3*ones(size(h11x)) -.05*ones(size(h11x))]); % series of magenta arrows, radially oriented, crossdir specified h12x=2; h12y=-3; h12z=axlim(5)/2; h12xr=1; h12zr=h12z; ir=.15;or=.81; h12t=(0:11)'/6*pi; h12 = feval(mfilename, [h12x+h12xr*cos(h12t)*ir h12y*ones(size(h12t)) h12z+h12zr*sin(h12t)*ir],[h12x+h12xr*cos(h12t)*or h12y*ones(size(h12t)) h12z+h12zr*sin(h12t)*or], 10,[],[],[],[], [-h12xr*sin(h12t) zeros(size(h12t)) h12zr*cos(h12t)], 'FaceColor','none','EdgeColor','m'); % series of normal arrows, tangentially oriented, crossdir specified or13=.91; h13t=(0:.5:12)'/6*pi; locs = [h12x+h12xr*cos(h13t)*or13 h12y*ones(size(h13t)) h12z+h12zr*sin(h13t)*or13]; h13 = feval(mfilename,locs(1:end-1,:),locs(2:end,:),6); % arrow with no line ==> oriented downwards h14 = feval(mfilename,[3 3 .100001],[3 3 .1],30); t6=text(3,3,3.6,'no line'); set(t6,'HorizontalAlignment','center'); % arrow with arrowheads at both ends h15 = feval(mfilename,[-.5 -3 -3],[1 -3 -3],'Ends','both','FaceColor','g', 'Length',20,'Width',3,'CrossDir',[0 0 1],'TipAngle',25); h=[h1;h2;h3;h4;h5;h6;h7;h8;h9;h10;h11;h12;h13;h14;h15]; function h = arrow_demo2(in) axlim = in.axlim; dolog = 1; if (dolog), set(in.hs,'YData',10.^get(in.hs,'YData')); end; shading('interp'); view(2); title(['Demo of the capabilities of the ARROW function in 2-D']); hold on; [C,H]=contour(in.x,in.y,in.z,20,'-'); hold off; for k=H',
6. z1=sqrt((r.^2-1)*c^2); z2=-sqrt((r.^2-1)*c^2); surf(x,y,z1) hold on surf(x,y,z2) xlabel('Truc Ox') ylabel('Truc Oy') zlabel('Truc Oz') title('Do thi hyperboloid 1 tang: x^2/a^2+y^2/b^2-z^2/c^2=1') rotate3d on hold off HYPERBOLOID 2 TẦNG clear all clc syms x y z real a=3; b=4; c=2; r=3; r=0:0.1:r; phi=0:0.1:2*pi+0.1; [r phi]=meshgrid(r,phi);%chia luoi x=a*r.*cos(phi); y=b*r.*sin(phi); z1=sqrt((r.^2+1))*c; z2=-sqrt((r.^2+1))*c; surf(x,y,z1) hold on surf(x,y,z2) xlabel('Truc Ox') ylabel('Truc Oy') zlabel('Truc Oz') title('Do thi hyperboloid 2 tang: x^2/a^2+y^2/b^2-z^2/c^2=-1') rotate3d on hold off NÓN 2 PHÍA clear all clc syms x y z a=3; b=4; c=4; r=0:0.1:1; phi=0:0.1:2*pi+0.1; [r phi]=meshgrid(r,phi);%chia luoi x=a*r.*cos(phi); y=b*r.*sin(phi); z1=r*c; z2=-r*c; surfc(x,y,z1) hold on surfc(x,y,z2) xlabel('Truc Ox') ylabel('Truc Oy') zlabel('Truc Oz') title('Do thi Paraboloid eliptic: y=x^2/a^2+z^2/b^2') rotate3d on hold off PARABOLOID ELIPTIC clear all
7. [x z]=meshgrid(x,z);%chia luoi y=z.^2/(a^2)-x.^2/(b^2); surfc(x,y,z) xlabel('Truc Ox') ylabel('Truc Oy') zlabel('Truc Oz') title('Do thi Paraboloid eliptic: y=z^2-y^2') rotate3d on TRỤ ELIPSE clear all clc syms x y z a=3; b=4; z=norm([a b]); r=-1:0.1:1; phi=0:0.1:2*pi+0.1; [z phi]=meshgrid(z*r,phi);%chia luoi x=a*cos(phi); y=b*sin(phi); surfc(x,y,z) xlabel('Truc Ox') ylabel('Truc Oy') zlabel('Truc Oz') title('Do thi Tru elipsoid: x^2/a^2+z^2/b^2=1') rotate3d on TRỤ PARABOLIC clear all clc clf syms x y z a=2; x=-a:0.1:2; z=-a:0.1:a; [z x]=meshgrid(z,x);%chia luoi y=x.^2; surfc(x,y,z) xlabel('Truc Ox') ylabel('Truc Oy') zlabel('Truc Oz') title('Do thi Tru parabol: y=x^2') rotate3d on MỘT SỐ HÀM VẼ CONE: function cone(a,b,c); % cone(a,b) ve non x^2/a^2+y^2/b^2=z^2/c^2 (non chinh tac)trong toa do tru clc clf phi=linspace(0,2*pi,50); r=linspace(0,1,50); [r phi]=meshgrid(r,phi);%tao luoi theo goc va ban kinh thiet dien. x=a*r.*cos(phi);%x tinh theo toa do tru (cuc) y=b*r.*sin(phi);%y tinh theo toa do tru (cuc) z=c*r;% pt mat non tren, -z la mat non duoi %lenh ve mat cong, mau luoi mat dinh la black,trong hinh mau luoi la white %('none' la ve khong luoi) %do trong suot la 0.5
8. CHÉO HÓA: function cheohoa A=input('nhap ma tran A='); [a,b]=size(A); while a~=b disp('Ma tran A khong vuong') A=input('nhap lai ma tran A='); [a,b]=size(A); end syms x lamda=solve(det(A-x*eye(size(A)))); lamda=double(lamda); i=1; while i 0.000001 lamda(i)=[]; i=i-1; else lamda(i)=lamda(i)-imag(lamda(i))*1i; end i=i+1; end lamda=unique(lamda); [n,~]=size(lamda); P=[]; for i=1:n P=[null(A-lamda(i)*eye(a)) P]; end [~,n]=size(P); if a~=n disp('khong cheo hoa ma tran A duoc') return end D=P\A*P; D(abs(D)<0.00001)=0; disp('ma tran cheo hoa cua A la:') disp(D) end CHUYỂN CƠ SỞ: function chuyencs E=input('nhap ma tran E: ');%nhap cs E [m,n]=size(E); %tim kich thuoc ma tran E while m~=n || rank(E)~=n %xem E co phai cs ko disp('E khong phai la co so') E=input('nhap lai ma tran E: '); [m,n]=size(E); end F=input('nhap ma tran F: ');%nhap ma tran F [a,b]=size(F); while b~=rank(F) || a~=b || a~=m disp('F khong phai la co so') F=input('nhap lai ma tran F: '); [a,b]=size(F); end X=F/E; %tim ma tran toa do cua F trong E X=X';% chuyen vi ma tran X disp('ma tran chuyen co so E sang F la:') disp(X) end ĐA THỨC: function dathuc
9. j=i+1; while j m k=n-m; end for i=1:m if A(i,:)==zeros(1,n) k=k+1;
10. end %ham sap xep cot cua ma tran tang dan function x=sapxep(a,i) [m,~]=size(a); j=i+1; while j<=m if a(i,i)==0 n=a(j,: ); a(i,: )=a(j,: ); a(i,: )=n; j=j+1; else break end end x=a; end HÌNH CHIẾU: function hinhchieu f=input('nhap khong gian con F: '); v=input('nhap vector x: '); [m,n]=size(f);[~,nv]=size(v); if n~=nv disp('ban nhap sai vector, xin kiem tra lai') else matrix=zeros(m); %tao ma tran 0 co m hang, m cot for i=1:m for j=1:m matrix(i,j)=tichvohuong(f(i,:),f(j,:)); end end exp=zeros(m,1); for i=1:m exp(i,1)=tichvohuong(v,f(i,:));%phan ma tran mo rong end x=matrix\exp; prfv=zeros(1,n); for i=1:m prfv=prfv+x(i,1)*f(i,:); end disp('hinh chieu cua x xuong F la vetor: ') disp(prfv) end end %ham tinh tich vo huong function a=tichvohuong(f1,f2) [~,n]=size(f1); a=0; for i=1:n a=a+f1(1,i)*f2(1,i); end end KHOẢNG CÁCH: function khoangcach syms x F=input('nhap khong gian vector con F: '); %nhap luon ca x vd:[x^4+7*x^3-3*x^2-x+7;x^3- 6*x^2+5*x-9] [m,~]=size(F); if m==1 F=F'; end
11. for i=1:m k(i)=poly2sym(f(i,:)); end end MA TRẬN BẬC THANG: function matranbt a=input('nhap ma tran vao '); [m,n]=size(a ); for i=1:min(m,n) [a]=sapxep(a,i); if a(i,i)~=0 a(i,: )=a(i,: )/a(i,i); for j=m:-1:i+1 a(j,: )=a(j,: )-a(i,: )*a(j,i); end end end disp('ma tran bac thang la: ') disp(a) end %ham sap xep a(i,i) ~=0 function x=sapxep(a,i) [m,~]=size(a); j=i+1; while j<=m if a(i,i)==0 n=a(i,: ); a(i,: )=a(j,: ); a(j,: )=n; j=j+1; else break end end x=a; end MA TRẬN TUYẾN TÍNH: function matrantt E=input('nhap ma tran cs E= '); [m,n]=size(E); while rank(E)~=n || m~=n disp('E khong phai la cs, yeu cau nhap lai') E=input('nhap lai ma tran cs E= '); [m,n]=size(E); end A=input('nhap ma tran tuyen tinh A trong cs E' ); while size(E)~=size(A) disp('A khong phai la ma tran tuyen tinh trong cs E, yeu cau nhap lai') A=input('nhap lai ma tran tuyen tinh A='); end B=input('nhap ma tran B= '); [m,n]=size(B); while rank(B)~=n || m~=n disp('B khong phai la cs, yeu cau nhap lai') B=input('nhap lai ma tran cs B= '); [m,n]=size(B); end while size(B)~=size(E) disp('ma tran B khong tuong thich voi ma tran A, yeu cau nhap lai') B=input('nhap ma tran cs B='); end
12. end disp('ma tran nghich dao la: ') disp(I) end end %ham sap xep cot cua ma tran tang dan function [x,k]=sapxep(a,i,I) [m,~]=size(a); j=i; while j<=m if a(i,i)==0 n=a(i,: ); N=I(i,: ); a(i,: )=a(j,: ); I(i,: )=I(j,: ); a(j,: )=n; I(j,: )=N; j=j+1; else break end end x=a; k=I; end PHỤ HỢP: function phuhop A=input('nhap ma tran A= '); C=zeros(size(A)); [m,n]=size(A); if m~=n disp('ma tran A khong vuong') return end for i=1:m for j=1:n B=A; B(:,j)=[]; B(i,:)=[]; C(i,j)=((-1)^(i+j))*det(B); end end C=C'; disp('ma tran phu hop cua A la:') disp(C) end TOÀN PHƯƠNG: function toanphuong A=input('nhap ma tran toan phuong A= '); [a b]=size(A); while a~=b disp('ma tran A khong vuong') A=input('nhap lai ma tran toan phuong A= '); [a b]=size(A); end while ~isequal(A,A') disp('ma tran A khong doi xung') disp(A) A=input('nhap lai ma tran toan phuong A= '); end d=eig(A);d=unique(d); disp('cac lamda cua A la:') disp(d') if isempty(d(d<=0))
13. disp(lamda') [n,~]=size(lamda); disp('co so cua cac khong gian con rieng:') for i=1:n disp(['voi lamda= ' num2str(lamda(i)) ':']) disp(null(A-lamda(i)*eye(a))') end end XÁC ĐỊNH: function xacdinh a=input('nhap ma tran A= '); [m,n]=size(a); while m~=n disp('ma tran A ko vuong, yeu cau nhap lai: ') a=input('nhap ma tran A= '); [m,n]=size(a); end for i=1:m b=a(1:i,1:i); d(i,1)=det(b); end d=double(d); disp ('cac dinh thuc duong cheo la') disp(d) b=d(d<=0); if isempty(b) disp('a xac dinh duong') else disp('a khong xac dinh duong') end end XÓA: function xoa a=input('nhap a='); b=input('nhap ma tran cac cot can xoa: '); b=unique(b);[~,m]=size(b); for i=m:-1:1 if b(i)<=length(a) a(:,b(i))=[]; end end disp(a) CƠ SỞ ẢNH: function [C DIM]= cosoanh(A, B) [n, ~]= size(A); C= rref(B, n- rank(B)); C= C(1:rank(A), :); DIM= rank(B); end CỞ SỞ NHÂN: function [C DIM]= cosonhan(A, B) [n, ~]= size(A); E= []; for i= 1: n D= rref([A, B(:, i)]); E=[E; D(:, n +1)']; end C= null(E);