function u = u(x)
    u = -(pi^2/16)*(
	cos(pi*x).^4 - 5*cos(pi*x).^2 + 4
    ) + (3/16)*(x.^2 - x);
end
