Back to tutorial index

# Lab Solutions

## Introduction

clear all;
close all;

Back to the top

## Problem 1

% get coefficients for the polynomial
clf;
c = [5, -2, 4, -1];
xs = linspace(-5,5,100);
ys = polyval(c,xs);
p1 = plot(xs,ys,'k');
hold on;

c_deriv = polyder(c);
ys_deriv = polyval(c_deriv,xs);

p2 = plot(xs,ys_deriv,'r');
legend([p1 p2],{'p(x)','p''(x)'});
snapnow;

Back to the top

## Problem 2

clear all;
clf;

c = [5, -2, 4, -1];
xs = linspace(-5,5,100);
ys = polyval(c,xs-1);
p1 = plot(xs,ys,'k');
hold on;

c_deriv = polyder(c);
ys_deriv = polyval(c_deriv,xs-1);

p2 = plot(xs,ys_deriv,'r');
legend([p1 p2],{'p(x)','p''(x)'});
snapnow;

Back to the top

## Problem 3

clf;
clear all;
v = [-1, 2, 5];
format long;
disp(norm(v,2))
disp(sqrt(dot(v,v)));
disp(sqrt(v(1)^2 + v(2)^2 + v(3)^2));
   5.477225575051660

5.477225575051661

5.477225575051661



Back to the top

## Problem 4

clf;
clear all;

xdata = [-2,-1,0,1,2];
ydata = xdata.*(xdata.^2 - 1);
pp = spline(xdata,ydata);
xs = linspace(-2,2,100);
ys = ppval(pp,xs);

p1 = plot(xs,ys,'b','linewidth',2);
hold on;
plot(xs,xs.*(xs.^2 - 1),'r','linewidth',1);

d2_coefs = zeros(pp.pieces,3);

for i = 1:pp.pieces
d2_coefs(i,:) = polyder(pp.coefs(i,:));
end

dpp = mkpp(pp.breaks,d2_coefs);

ys_deriv = ppval(dpp,xs);

plot(xs,ys_deriv,'k','linewidth',2);
plot(xs,3*xs.^2 - 1,'r','linewidth',1);
snapnow;

Back to the top

## Problem 5

clear all;
clf;

N = 8;
xdata = linspace(-5,5,N);
ydata = -2 + 4*rand(1,N);

pp = spline(xdata,ydata);
xs = linspace(-5,5,100);
ys = ppval(pp,xs);
plot(xs,ys);
hold on;
plot(xdata,ydata,'.','markersize',30);

Back to the top

## Problem 6

clear all;
clf;

N = 10;
xdata = linspace(-5,5,N);
ydata = -2 + 4*rand(1,N);
pp = spline(xdata,ydata);
xs = linspace(-5,5,200);
ys = ppval(pp,xs);
plot(xs,ys,'r','linewidth',2);
hold on;

d_coefs = zeros(pp.pieces,3);
d2_coefs = zeros(pp.pieces,2);

for j = 1:pp.pieces,
c = pp.coefs(j,:);
d_coefs(j,:) = polyder(c);
d2_coefs(j,:) = polyder(d_coefs(j,:));
end

% First and second derivatives of the spline.
dpp = mkpp(pp.breaks,d_coefs);
d2pp = mkpp(pp.breaks,d2_coefs);

% Compute the tangents and normals at the data points.  We will make
% r'(x) and r''(x) 3-component vectors so that the norm and cross product
% functions in matlab can be used. The third component will just be zero
% since the curve is restricted to the x-y plane.
fprime = ppval(dpp,xdata);
rprime = [ones(1,N); fprime; zeros(1,N)];
f2prime = ppval(d2pp,xdata);
r2prime = [zeros(1,N); f2prime; zeros(1,N)]; % Cross product requires

tangents = zeros(3,N);
normals = zeros(3,N);
% Loop over the values in xdata to compute the tangents and normals.
for j = 1:N
g = norm(rprime(:,j),2);
tangents(:,j) = rprime(:,j)/g;
gprime = dot(rprime(:,j),r2prime(:,j))/g;
kappa = norm(cross(r2prime(:,j),rprime(:,j)),2)/g^3;
normals(:,j) = (g*r2prime(:,j)-gprime*rprime(:,j))/(g^3*kappa);
end

% Plot the original data
plot(xdata,ydata,'k.','markersize',30);

% Plot the tangents, use only the first two components since this is a
% 2D curve
quiver(xdata,ydata,tangents(1,:),tangents(2,:),0.25,'k-','LineWidth',1);

% Plot the normals, use only the first two components since this is a
% 2D curve
quiver(xdata,ydata,normals(1,:),normals(2,:),0.25,'b-','LineWidth',1);

title('Tangents and Normals (spline)','fontsize',18);
xlim([-5.5 5.5]);
ylim([-3 3]);
daspect([1 1 1]);
snapnow;

Notice how the normal always points towards the radius of curvature.

Back to the top

## Problem 7

for i = 1:N-1,
fprintf('T(i) dot N(i) = %12.4e\n',dot(tangents(:,i),normals(:,i)));
end
T(i) dot N(i) =  -1.1102e-16
T(i) dot N(i) =   1.6653e-16
T(i) dot N(i) =   1.1102e-16
T(i) dot N(i) =   2.2204e-16
T(i) dot N(i) =  -1.6653e-16
T(i) dot N(i) =   5.5511e-17
T(i) dot N(i) =   1.1102e-16
T(i) dot N(i) =  -5.5511e-17
T(i) dot N(i) =   8.3267e-16


Back to the top

## Problem 8

binormals = zeros(3,N);
for i = 1:N-1,
b = cross(tangents(:,i),normals(:,i));
binormals(:,i) = b(:)';
fprintf('%24.16f %24.16f %24.16f\n',b(1),b(2),b(3));

end
     -0.0000000000000000       0.0000000000000000       1.0000000000000000
0.0000000000000000      -0.0000000000000000       1.0000000000000000
0.0000000000000000       0.0000000000000000      -0.9999999999999999
0.0000000000000000      -0.0000000000000000      -1.0000000000000000
0.0000000000000000      -0.0000000000000000       0.9999999999999999
0.0000000000000000       0.0000000000000000      -0.9999999999999999
-0.0000000000000000       0.0000000000000000       1.0000000000000000
0.0000000000000000       0.0000000000000000      -1.0000000000000000
-0.0000000000000000       0.0000000000000000       1.0000000000000000


Back to the top

## Problem 9

for i = 1:N-1,
fprintf('B(i) dot N(i) = %12.4e\n',dot(binormals(:,i),normals(:,i)));
end
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00
B(i) dot N(i) =   0.0000e+00

fprintf('\n');
for i = 1:N-1,
fprintf('B(i) dot T(i) = %12.4e\n',dot(binormals(:,i),tangents(:,i)));
end
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00
B(i) dot T(i) =   0.0000e+00


Back to the top

## Problem 10

for i = 1:N-1,
d1 = norm(binormals(:,i),2);
d2 = norm(tangents(:,i),2);
d3 = norm(normals(:,i),2);
fprintf('%24.16f %24.16f %24.16f\n',d1,d2,d3);
end
      1.0000000000000000       0.9999999999999999       1.0000000000000002
1.0000000000000000       1.0000000000000000       1.0000000000000000
0.9999999999999999       1.0000000000000000       0.9999999999999997
1.0000000000000000       0.9999999999999999       1.0000000000000002
0.9999999999999999       1.0000000000000000       0.9999999999999999
0.9999999999999999       1.0000000000000000       0.9999999999999998
1.0000000000000000       1.0000000000000000       0.9999999999999998
1.0000000000000000       1.0000000000000000       1.0000000000000000
1.0000000000000000       1.0000000000000000       0.9999999999999998


Back to lab exercises.

Back to the top