
%\documentstyle{siam}
\documentclass[final,psfrag,twoside,epsfrag]{article}
\input{zmacros}
\textwidth 6.4in
\textheight 9in
\oddsidemargin -0.3in
\evensidemargin -0.3in
\topmargin -.5in

%\input{zmacros.tex}
\def\be{\begin{enumerate}}
\def\ee{\end{enumerate}}
%\newcommand{\vthin}{\vskip 5pt}

\newcommand{\bftau}{\mbox{\boldmath $\tau$}}

\pagestyle{empty}
\begin{document}
\hfill\vbox{\hbox{MA 780}\hbox{Homework \#2}
\hbox{\bf Due  Feb. 17}}
\vskip .2in
\baselineskip 14pt


\be
\vthin
\item  Given $\lp x_i,f(x_i),f'(x_i) \rp$, $i=0,1,\cdots, n$, $n>0$.
Let $p_{2n+1}(x)$ be the Hermite interpolation polynomial satisfying
\eqmno
p_{2n+1}(x_i)=f(x_i), \quad p'_{2n+1}(x_i)= f'(x_i), \qquad i=0,1,\cdots, n.
\enmno
Show that if $f(x)\in C^{2n+2}$, then we have the following error estimate
\eqmno
R(x) = f(x) - p_{2n+1}(x)=\fr{f^{(2n+2)}(\xi(x))}{(2n+2)!} \omega^2_{n+1}(x).
\enmno
Give a {\em least upper bound} for the case $n=1$.

\ignore{
{\bf Hint:} Consider the auxiliary function
\eqmno
\phi(t) = f(t)-p_{2n+1}(t) - \fr{\omega^2_{n+1}(t)}{\omega^2_{n+1}(x)} \lp
f(x)-p_{2n+1}(x) \rp.
\enmno }

\vthin
\item % Modified from problem 9 on page 377.
\be
\item Show that there is no polynomial $p_3(x)$ that interpolates
\eqmno
f(-1)=1, \quad f'(-1)=1, \quad f'(1)=2, \quad f(2)=1.
\enmno
{\bf Hint:} Use the un-determined coefficient method and show that the
system of equations for the coefficients
does not have a solution. Explain as much as you can.
\item Now we add two conditions
\eqmno
f(1)=0, \quad f''(1)=2.
\enmno
Use the divided difference table to find the polynomial interpolation and
verify your result.
\ee

\vthin

\item Let $\{ x_i\} _0 ^n$, $x_0=a$, $x_n=b$ be a equally spaced nodal points for the interval $[a, b]$, that is,  $x_{i+1}-x_i=h= (b-a)/n$ is a constant. 
 Assume also that $u(x)\in C^5[a,b]$ (can this be relaxed?). 
\be
\item   In the interval $[x_{i-1}, x_{i+1}]$, we can approximate $u(x)$ by a quadratic interpolation using nodal points $x_{i-1}, x_i, x_{i+1}$. From your interpolation function, find an approximation to  $u''(x_i)$ and an error estimate. {\bf Hint:} Use the second order divided difference.

\item Using the formula we can get a linear system of equations to approximate the boundary value problem
\eqmno
  u''(x) = f(x), \qquad a< x < b, \qquad u(a) = 0, \quad u(b) = 0,
\enmno 
as 
\eqmno
  \frac{U_{i-1} - 2 U_i + U_{i-1} }{h^2} = f(x_i),  \qquad U_0 = 0, \quad U_n = 0, \quad i=1,2,\cdots n-1,
\enmno
where $U_i \approx u(x_i)$. Write down the matrix-vector form ( $ A U = F$) of the linear system of equations above. 

\item Write down the component form of the SOR($\omega$) method for the system of equations above. Do the Jacobi, Gauss-Seidel,  SOR($\omega)$ converge?
\ee

{\bf Hint:}  The SOR($\omega)$ method is ${\bf x}^{k+1}= (1- \omega)  {\bf x}^{k} + \omega {\bf x}^{k+1}_{Gauss-Seidel}$, which converges for $0< \omega < 2$ for this problem. 


\vthin

\item Write down the linear system of equations for the moments
$S_3''(x_j)=M_j$ for the periodic cubic spline and simplify it for the
equally spaced nodal points by setting $h=x_{j}-x_{j-1}$. Do some research about how to solve the linear system of equations. 


\vthin
\item 
\be 
\ignore{
\item Derive a quadratic spline $S_2(x) \in C^0 $ that satisfies 
\eqm \label{p2}
S_2(x_i)=f(x_i), \quad i=0,1,\cdots,n, \quad \mbox{and} \quad
S_2'(x_0)=f'(x_0).
\enm
and $S_2'(x_n)=f'(x_n)$ that uses local basis functions, that is, $S_2(x) = \dsp \sum_{i=0}^{n+1} \alpha_i  \phi_i(x) + \sum_{j=1}^{n} \alpha_j  \bar{\phi}_j(x) $, 
where $\phi_i(x)$'s and $\bar{\phi}_j(x)$ are piecewise quadratic function in $C^0$  satisfying
\eqmno
   \phi_i(x_i) = \delta_i^j,  \quad  \phi_i\left ( \frac{x_j+ x_{j+1}}{2} \right) =0; \qquad   \bar{\phi}_j(x_i)=0, \quad \bar{\phi}_j \left ( \frac{x_j+ x_{j+1}}{2} \right) =1.
\enmno }

\item  Derive a quadratic spline $S_2(x) \in C^1$ that satisfies % \eqn{p2}. 
\eqm \label{p2}
S_2(x_i)=f(x_i), \quad i=0,1,\cdots,n, \quad \mbox{and} \quad
S_2'(x_0)=f'(x_0).
\enm
{\bf Hint:} Write the quadratic spline $S_2(x)$ as
\eqmno
S_2(x)=a_i + b_i(x-x_i) + c_i(x-x_i)(x-x_{i+1}), \quad x\in \left [x_i, \,
x_{i+1} \right ],
\enmno
and derive a system of equations for $a_i$, $b_i$ and $c_i$.

\item
\fbox{\bf Extra credit:} Let $R(x)=f(x)-S_2(x)$ and $h=\max \{ |x_{i+1}-x_i|\}$.
Show the following error estimates if $f(x)\in C^3$.
\eqmno
&& \|R''\|_{\infty} \le h \|f'''\|_{\infty}, \quad \mbox{excluding the nodes,}\\
&& \|R'\|_{\infty} \le h^2\|f'''\|_{\infty}/8, \\
&& \|R\|_{\infty} \le (x_n-x_0) h^3\|f'''\|_{\infty}/8.
\enmno
\ee


\item {\bf Programming Part: } % Can be done as group work up to 4 people. 

 Use  {\em a cubic spline} (either the cardinal or B-splines, see text  and pseudo-code there) assuming $S_3'(x_0)=f'(x_0)$ and $S_3'(x_n)=f'(x_n)$ to approximate functions $f(x) =  \cos(10 \pi x)$
and  $f(x) =\fr{1}{ 1 + 25 x^2}$ in HW\#1. 

Plot the function and  your approximation  on the same plot. Plot also  of the errors and compare  with the results from HW 1. 

% Note that: More points are awarded if  you choose the spline approach, which can also be done as group work up to 4 people.   the {\em  piecewise linear}, or

{\em The project part can be done as group  up to 4 people.} You need to state your group in you HW. Analysis can be done individually or in group.  If you can not get your code working, you can use Hermite cubic interpolation, or even piecewise linear interpolations. A couple of points may be deducted though. 

\vthin

\item \textbf{An extra credit project for the semester.} Group work is OK. 

This is example of how to simulate  mean-curvature flows using a cubic spline. Assume that we have a two-dimensional smooth and closed curve ${\bf X}(s,t)=(X(s,t), Y(s,t))$ whose motion satisfies $\frac {d {\bf X}}{dt}  = {\bf n} f(\kappa)$, where $\kappa$ is the curvature of the curve, ${\bf n}$ is the normal direction, for example $f(\kappa)=1$, or $f(\kappa)=-1$, or $f(\kappa)=C \kappa$.

\be 
 \item Initially, we can use a set of discrete points $(X_k, Y_k)$ to represent the curve at $t=0$, $k=1,2, \cdots N$, $X_1= X_{N+1}$, $Y_1= Y_{N+1}$. We can treat the points $(X_k, Y_k)$ as the function of the arc-length $s$, such that $s_0=0$, $s_{k+1}= s_k + \sqrt{ (X_{k+1}-X_k )^2 +  (Y_{k+1}-Y_k )^2} $. Thus, we have $X_k^0 \approx X(s_k, t^0)$, $Y_k \approx Y(s_k,t^0)$, where the superscript is the time level, 
 $t^{m+1} = t^m + \Delta t$ with $t^0 = 0$. 
 
 \item Use a cubic spline with the period boundary condition to approximate the curve to obtain $(X(s, t^m), Y(s, t^m))$ so that we can get the tangential ($\bftau$)  and normal directions ${ \bf n}$, and the curvature information, for example, $\bftau = (\frac {d X(s_k, t^k)}{ d s},  \frac {d Y(s_k, t^m)}{ d s})/|\bftau |$, where  $|\bftau | = \sqrt{  (\frac {d X(s_k, t^m)}{ d s})^2 +  (\frac {d Y(s_k, t^m)}{ d s} )^2} $.
 
 \item Evolve the boundary using the Euler's method ${\bf X}^{m+1} = {\bf X}^{m} + \Delta t \,  {\bf n}^m  f(\kappa^m)$. 
  
\ee

\vthin

 {\bf Hint:} In continuous case, this is an initial value problem
\eqm
  \frac {d {\bf X(s,t)}}{dt}  = {\bf n} f(\kappa), \qquad t>0, \qquad    {\bf X(s,0)} = {\bf X}_0(s). 
\enm
We can use forward Euler method to solve the problem
\eqm
 \frac { {\bf X}^{k+1} - {\bf X}^{k} }{\Delta t}  = {\bf n^k} f(\kappa^k), \quad k=0,1,\cdots, 
 \enm
where ${\bf X}^{k} $ is the initial position of the curve. A pseudo code is given below:
\begin{verbatim}
 % Set-up
     clear all; close all;
     n=80; h=2/n; x=-1:h:1; y=x;     % Domain: [-1, 1] x [-1, 1]
     n1 = n; tfinal=0.5; ds = 2*pi/n1; theta=0: ds:2*pi-ds;
     for m=1:n1;   % An initial ellipse
        X(m)=0.6*cos(theta(m)); Y(m) = 0.4*sin(theta(m));
      end
      plot(X,Y)
     axis([-1 1 -1 1]); axis('square'); hold % Visualize the initial curve
     
% Time loop
     
     dt = h*h/2; t=0;
     while t < tfinal       
        % Spline interpolation to get 
        for m=1: n1-1.   % Need to do it for m=n1 separately.
           dx = X(m+1)-X(m); dy = Y(m+1)-X(m);  ds = sqrt( dx*dx + dy*dy + 1e-15);
           dx = dx/ds dy = dy/ds;  % The tangent direction
           nx = -dy;  ny=dx;  % The normal direction
           X1(m) = X(m) - dt*nx; Y1(m) = Y(m) - dt*ny;   % f=-1;
        end
        m = n1; 
           %%%%% Fill in
        
        t = t + dt; 
        plot(X1,Y1)
        X=X1; Y= Y1;   % Overwrite the old boundary data
     end
     
     
\end{verbatim}




\ignore{
Given a quadrilateral mesh (rectangles), can we find a piecewise bi-polynomial interpolation
 \eqm
  Q_h^k(x,y) =  \sum_{i=0,j=0}^{i+j\le k} a_{ij}x^i y^j
 \enm
 such that:
 \be
 \item $Q_h^k(x,y)\in C^2$;
 %\item $Q_h^k(x_i,y_j)= u(x_i,y_j)$ at all $(x_i,y_j)$, where $u(x_i,y_j)$'s are given.
 \item $D^{\alpha}(x_i,y_j)$, $0\le \|\alpha\|_1\le 2$ are given, that is, the function values, and up to all the second partial derivatives are given at $(x_i,y_j)$.
 \ee
 Find the minimum $k$; derive,  programming, and validate the algorithm of your derivation;  and carry out necessary analysis. \\
 
 Finally, if we only impose some of $D^{\alpha}(x_i,y_j)$, $0\le \|\alpha\|_1\le 2$ values, say, its function and Laplacian values at the four corners, can we derive such basis functions?}
 
\ee 

\end{document}
