====================================================================== Part A, Question 1 ================== Digits:=20; f:=x^6-12*x^4+39*x^2+2*x=28; ans1:=fsolve( diff(f,x)=0, x=2..3); ans2:=subs( x=ans1, f); r1:=fsolve( f=0, x=0..1); % could just do "x=1" r2:=fsolve( f=0, x=-1.2..-1); % could just do "x=-1" ans3:=int( -f ,x=r1..r2); % the minus is not really nec ====================================================================== Part A, Question 2 ================== Make an m-file called ft.m function y=ft(x,t) y=(1+t*sin(x)+x^4)/(1+4*x^2); Now do the following (can write the commands in an m-file if you wish): t=[1:0.05:4]; m=zeros(size(t)); for i=1:size(t,2) [a m(i)] = fminsearch( @ft, 0, [], t(i) ); end plot(t,m) ====================================================================== Part A, Question 3 ================== f:=piecewise( x<2*Pi/3, exp(x)*sin(x), x>2*Pi/3, 1/(a+x^2) ); a:=solve( limit(f,x=2*Pi/3,left) = limit(f,x=2*Pi/3,right) , a); plot(f,x=0..6); ====================================================================== Part A, Question 4 ================== A is a 50 times 50 matrix, with entries random numbers in the range [0,1] x is a 50 dimensional vector, with elements spaced equally from -1 to 1 y is a 50 dimensional vector, the eigenvalues of A + A^T - ones(50,50) z is a 50 dimensional vector, a simple cubic function of x in the plot y is plotted against x using points (y is not a function) and z is plotted against x using a line (z is a function) we see the eigenvalues of the (apparently) random matrix A+A^T-1 are not at all random, but are more or less the same as the values of the cubic function z ! ====================================================================== Part A, Question 5 ================== with(LinearAlgebra); M:=<<1|r|s>,<2|-2|1>,<3|4|1>>; p:=CharacteristicPolynomial(M,lambda); subs(lambda=1,p)=0; % this is the condition M2:=<<1|r|-3*r/17>,<2|-2|1>,<3|4|1>>; % or solve for r ColumnSpace( M2 - IdentityMatrix(3,3) ); ====================================================================== Part B, Question 1 ================== a) This is easiest in Maple r1:=fsolve( 4*sin(x)=2/(1+x), x=0..1); r2:=fsolve( x^2/3=2/(1+x), x=1..2); r3:=fsolve( x^2/3=4*sin(x), x=2..3); ans:=int( 4*sin(x)-2/(1=x), x=r1..r2 ) + int( 4*sin(x)-x^2/3, x=r2..r3 ) b) This is easiest in Matlab. We do not know where f is "over" g, or g "over" f, and of course we only want the x's where 4sin(x) is over both of them. Remember that to use quad we need a function that takes a VECTOR! Assume we already have files f.m, g.m, but these might work only on scalars :-( Make an mfile fq.m: function z=fq(x) z = zeros(size(x)) for i=1:size(x,1) for j=1:size(x,2) % make it work for matrices! vectors is really enough s4=4*sin(x(i,j)); fx=f(x(i,j)); gx=g(x(i,j)); if s4>max(fx,gx) z(i,j)=s4-max(fx,gx); % otherwise zero end end end Now just need to do quad(@fq,0,pi) ====================================================================== Part B, Question 2 ================== a) The quick and dirty way. File nch.m function z=nch(X,Y,r,s) Yp=Y; Yp(r)=Y(s); Yp(s)=Y(r); if norm(X-Yp)<norm(X-Y) z=Yp else z=Y end % in fact this is a waste because the condition can be rewritten % (X(r)-Yp(r))^2 + (X(s)-Yp(s))^2 < (X(r)-Y(r))^2 + (X(s)-Y(s))^2 % and in other even easier ways b) It is _not_ enough to write d=size(X,1) % assume column vectors for r=1:(d-1) for s=(i+1):d Y=nch(X,Y,r,s); end end Y it may be necessary to repeat this procedure for a few times. so make the above a function, i.e. create a file nch2.m function Y=nch2(X,Y) d=size(X,1) % assume column vectors for r=1:(d-1) for s=(i+1):d Y=nch(X,Y,r,s); end end then do something like Yp=nch2(X,Y); while Yp~=Y Y=Yp; Yp=nch2(X,Y); end Y ====================================================================== Part B, Question 3 ================== with(LinearAlgebra) p:=proc(A,b,Ps) % Ps is a list of points (vectors) local numpoints,dmin,dmax,dtot,dav,d,i,lA; lA=norm(A,2); numpoints:=nops(Ps); dmin:=infinity; dmax:=0; dtot:=0; for i from 1 to numpoints do d:= abs( DotProduct(A,op(i,Ps))-b )/lA; if d>dmax then dmax:=d end if: if d<dmin then dmin:=d end if: dtot:=dtot+d; end: return( dmax, dmin, dtot/numpoints) end proc; ====================================================================== Part B, Question 4 ================== we need to make a function that computes (lambda(s,t)-1)^2 provided (s,t) in the unit circle. zero otherwise. it has to work when s is a vector!! M-file: findl.m function z=findl(s,t) z=zeros(size(s)); for i=1:size(s,2) if s(i)^2+t^2<1 M=[2 s(i) t; s(i) -1 4; t 4 -1]; % matrix v=eig(M)-1; % values of lambda-1 q=v(1); % pick the closest to zero if abs(v(2))<abs(q) q=v(2); end if abs(v(3))<abs(q) q=v(3); end z(i)=q^2; end end finally dblquad( @findl, -1, 1, -1, 1 ) ======================================================================
Back to course homepage
Back to my main teaching page
Back to my main page