How to calculate smoothing spline coefficients

Question

I am attempting to calculate smoothing spline coefficients based on the description in Reinsch's 1967 paper, but I'm having some trouble. The first derivative is not continuous. Here are the equations I'm attempting to follow:

And this is the MATLAB code. There are some simplifications because all h_sub_i and weight values are equal to 1.

pts = [0      0
  1      1
  2      2
  3      3
  4      2
  5      1
  6      0];

n = length(pts)-1;
p = 10; % closeness of fit

Q = -2*eye(n-1) + diag(ones(1,n-2),1) + diag(ones(1,n-2),-1);
Q(n+1, end) = 0;
Q(n, end) = 1

T = 4/3*eye(n-1) + diag(1/3*ones(1,n-2),1) + diag(1/3*ones(1,n-2),-1)
Ap = transpose(Q)*Q + p*T

c = Ap^-1 * p * transpose(Q) * pts(:, 2);
a = pts(:, 2) - (1/p)*Q*c;
c = [0; c; 0];
d = (c(2:end) - c(1:end-1)) / 3;
b = (a(2:end) - a(1:end-1)) - c(1:end-1) - d;

N = 10000;
breaks = 0:6;
pp = mkpp(breaks, [d(1:n) c(1:n) b(1:n) a(1:n)]);
x = 0:0.01:6;
y = ppval(pp, x);

figure(1)
clf
hold on
grid on
plot(pts(:,1), pts(:,2), '--bo');
plot(x, y, 'r');

The following output is produced:

Q =

    -2     1     0     0     0
     1    -2     1     0     0
     0     1    -2     1     0
     0     0     1    -2     1
     0     0     0     1    -2
     0     0     0     0     1
     0     0     0     0     0


T =

       1.3333      0.33333            0            0            0
      0.33333       1.3333      0.33333            0            0
            0      0.33333       1.3333      0.33333            0
            0            0      0.33333       1.3333      0.33333
            0            0            0      0.33333       1.3333


Ap =

       18.333     -0.66667            1            0            0
     -0.66667       19.333     -0.66667            1            0
            1     -0.66667       19.333     -0.66667            1
            0            1     -0.66667       19.333     -0.66667
            0            0            1     -0.66667       19.333

What is my error?