Exam solutions

Question 1

for r from 1 to 10 do 
   Digits:=2*r+2;
   fsolve( sin(x)=exp(x), x=-r*Pi );
end do;

Question 2

make m-file mb2.m:
function z=mb2(x)
   z=x(1)^2+x(2)^2+sin(x(1)+x(2))/2+sin(x(1)+2*x(2))/3;

in command window:
format long
fminsearch( @mb2, [0,0], [optimset('tolx',1e-8)])

Question 3

for i from 0 to 100 do
   c:=-5+i/10;
   p[i]:=[c,fsolve(2*x^3+3*x+c=0)];
end do:
plot( [seq(p[i],i=0..100)], style=point );

Question 4

c=[-5:0.1:5]; 
for i=0:100 
   v=roots(2,0,3,c(i));
   for j=1:size(v,1)
      if im(j)==0 
         d(i)=v(j);
      end
   end
end
plot(c,d); 

Question 5

with(LinearAlgebra);
assume(t,real);
v:=<sqrt(t)*cos(t),sqrt(t)*sin(t),t>;
w:=Vector(3);
for i from 1 to 3 do w[i]:=diff(v[i],t) end do;
int( w.w, t=0..1);

Question 6

make m-file mb6.m:
function y=mb6(x,p)
   y=sin(x)./(1+p*x.*x);

in command window (or another m-file):
p=[0:0.1:2];
q=zeros(size(p));
for i=1:size(p,2)
   q(i)=quad(@mb6,0,p(i),[],[],p(i));
end
plot(p,q)

Question 7

make m-file mb7.m:
function y=mb7(s)
M=[1.9 1.1 1.4 1.4; 1.1 1.5 s 0.7; 1.4 s 1.2 0.9; 1.4 0.7 0.9 0.8];
v=eigs(M);
vmin=v(1); 
vmax=v(1); 
for i=1:size(v)
   if v(i)<vmin
     vmin=v(i)
   end
   if v(i)>vmax
     vmax=v(i)
   end
end
y=(vmax-vmin)/2/s;

At this stage it is tempting to write 
quad( @mb7, 10, 20) 
but this does not work as mb7 is not vectorized
So make a better m-file mb7a.m:
function y=mb7a(s)

y=zeros(size(s));
for k=1:size(s,1)
for l=1:size(s,2)

M=[1.9 1.1 1.4 1.4; 1.1 1.5 s(k,l) 0.7; 1.4 s(k,l) 1.2 0.9; 1.4 0.7 0.9 0.8];
v=eigs(M);
vmin=v(1); 
vmax=v(1); 
for i=1:size(v)
   if v(i)<vmin
     vmin=v(i);
   end
   if v(i)>vmax
     vmax=v(i);
   end
end
y(k,l)=(vmax-vmin)/2/s(k,l);

end
end

and now you can run quad( @mb7a, 10, 20) 

Question 8

So long as the m-file f.m works on a pair of variables x,y where
x is a vector and y scalar (giving a vector answer) you can write
dblquad(@f,0,1,0,1)

Implement the 2d trap rule by the following m-file:

function S=trap2(f,N)

h=1/N;
S1=feval(f,[0,0])+feval(f,[0,1])+feval(f,[1,0])+feval(f,[1,1]);
S2=0;
for i=1:(N-1)
   S2=S2+feval(f,[i*h,0])+feval(f,[0,i*h])+feval(f,[i*h,1])+feval(f,[1,i*h]);
end
S3=0;
for i=1:(N-1)
for j=1:(N-1)
   S3=S3+feval(f,[i*h,j*h]);
end
end
S=(S1/4+S2/2+S3)*h^2;

Here f is a function that takes as input a 2-vector [x,y]. 

Question 9

a)

q9a:=proc(p1,p2,l)        
local s,t;
s:=op(1,l)*op(1,p1)+op(2,l)*op(2,p1)+op(3,l);
t:=op(1,l)*op(1,p2)+op(2,l)*op(2,p2)+op(3,l);
if (s*t<0) then return(1) else return(0) end if;
end proc;

Here p1 is the point [x1,y1] (a list with 2 elements)
     p2 is the point [x2,y2] (a list with 2 elements)
     l  is the line  [a,b,c] (a list with 3 elements)
the function returns 1 if the points are on opposite sides of the line
   and 0 if they are on the same side

b) 

q9b:=proc( ps , ls)     % ps are the points, ls are the lines
local i,j,k,s,t;
s:=1;                   % this will change to zero if there are 
                        % two points without a line between them
for i from 1 to nops(ps) do
for j from (i+1) to nops(ps) do
     t:=0; 
     for k from 1 to nops(ls) do
         if (q9a(op(i,ps),op(j,ps),op(k,ls))=1) 
         then 
            t:=t+1;      % t becomes 1 if there is a line between the points
            break;
         end if;
     end do:
s:=s*t;
end do:
end do:                    
return(s);   % 1 if there is a line between every pair of points
end proc;

c)

q9c:=proc(p1,p2,m)        
local s,t;
s:=op(1,m)*op(1,p1)+op(2,m)*op(2,p1)+op(3,m)*op(3,p1)+op(4,m);
t:=op(1,m)*op(1,p2)+op(2,m)*op(2,p2)+op(3,m)*op(3,p2)+op(4,m);
if (s*t<0) then return(1) else return(0) end if;
end proc;

Here p1 is the point [x1,y1,z1] (a list with 3 elements)
     p2 is the point [x2,y2,z2] (a list with 3 elements)
     m  is the plane  [a,b,c,d] (a list with 4 elements)
the function returns 1 if the points are on opposite sides of the plane
   and 0 if they are on the same side

d) 

q9d:=proc( ps , ms)     % ps are the points, ms are the planes
local i,j,k,s,t;
s:=1;                   % this will change to zero if there are 
                        % two points without a plane between them
for i from 1 to nops(ps) do
for j from (i+1) to nops(ps) do
     t:=0; 
     for k from 1 to nops(ms) do
         if (q9a(op(i,ps),op(j,ps),op(k,ms))=1) 
         then 
            t:=t+1;      % t becomes 1 if there is a plane between the points
            break;
         end if;
     end do:
s:=s*t;
end do:
end do:                    
return(s);   % 1 if there ia a plane between every pair of points
end proc;

Back to course homepage
Back to my main teaching page
Back to my main page