%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% PLEASE TYPE the content of lines with >>
%%%% the lines with % are comments/instructions
%%%% Type help <command> or doc <command> for more info on, well, <command>
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  PART 1
>> xplot=linspace(0,1);
>> plot(xplot,sin(xplot))
>> xpoints=[.1,.3,.7];
>> upoints=sin(xpoints);
>> plot(xplot,sin(xplot),'k-',xpoints,upoints,'r*-');
%%%% So YOU"VE PRACTICED PLOTTING and how it depends on your nodes
%% now, make a mistake and recover form it
>> plot(xplot,sin(xplo))
%%%% 
%%% Now more on interpolation
%%% 
>> myfun=@(x)(sin(3*x));
>> upoints=myfun(xpoints);
>> plot(xplot,myfun(xplot),'k-',xpoints,upoints,'r*-');
>> uinterp=interp1(xpoints,upoints,xplot);
>> plot(xplot,myfun(xplot),'k-',xpoints,upoints,'r*-',xplot,uinterp,'bo');
%%% So you saw how to construct a piecewise linear interpolant.
%%% NEXT, how to construct a linear APPROXIMATION based on upoints
>> p=polyfit(xpoints,upoints,1);
>> plot(xplot,myfun(xplot),'k-',xpoints,upoints,'r*-',xplot,polyval(p,xplot));
>> plot(xplot,myfun(xplot),'k-',xpoints,upoints,'r*',xplot,polyval(p,xplot),'b.');
%%% Is it good ? Not really... so let us try a quadratic APPROXIMATION
>> p2=polyfit(xpoints,upoints,2);
>> plot(xplot,myfun(xplot),'k-',xpoints,upoints,'r*',xplot,polyval(p,xplot),'b.',xplot,polyval(p2,xplot),'m-');
>> legend('original function','piecewise linear interpolant','approximating linear polynomial','approximating quadratic')
%%%%%%%%%%%%%%%%%%
%%% 
%%% MORE PRACTICE: change myfun, xpoints, order of interpolation, order of approximation
%%% 
%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  PART 2
%%  TRY THE FOLLOWING Finite Element example and play with the order of the error
%%  download the code fem1d.m first and run it. Watch what it does
>> fem1d(0,1,5);
>> fem1d(0,1,10);
>> fem1d(0,1,20);
%% WATCH THE COMPARISON between FE and FD solution. 
%%%%%%%
%% OTHER THINGS to do mwith the code if you have time: 1) Look inside. 
%% 2) print the matrix of the linear system (find where to do the following
%% full(mat),pause
%% 3) make the code return the fd solution, and exact solution
%% and plot them from the outside of the code fem1d.m
%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  PART 3
%%% PRACTICE plotting error order using the following example
>> h=[0.2,0.1,0.05]
>> L2err = [0.081,0.021,0.0049];
>> plot(h,L2err,'*-');
>> loglog(h,L2err,'*-')
>> grid on
%%% WHAT order of error do we see ? WHY loglog plot ? 
%%% If you guessed right, the order is quadratic. Let us see ...
>> loglog(h,L2err,'k*-',h,h.^2,'b-'); grid on;
>> legend('Error(h)','Quadratic trend');
%%%%%%%%
%% MORE PRACTICE: come up with something of order 3 or 2.5 
%% (and maybe perturb it gently using rand() )
%% then test the order of convergence
%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  PART 4
%%%%%%%%% NOW FINISH UP and do the analysis on what the FE code produced
%% recall that it says, for example  
%% Approximate L2 norm of error is=0.0788696 h=0.2
>> fem1d(0,1,5);
%% what h did it use ? what error did it produce ?
%% Now record all the 'h' and error in your vector h, and L2err
%% and repeat the loglog thing. You're done !
%% 
%% Well, you could try other grid norms ... Linf. Ask me why/how !
%%