% The function illustrates midpoint-based second order
% accurate finite elements;
function fem_02
% Size of the problem;
N1 = 6;
% Viewing the basis functions
Nx = 200;
dx = 1.0 / (Nx - 1);
X = 0.0 : dx : 1.0;
figure(1);
plot( X, phi( N1, 1, X ), X, phi( N1, 2, X ), X, phi( N1, 3, X ), ...
X, phi( N1, N1, X ) ); hold on;
plot( X, psi( N1, 1, X ),':', X, psi( N1, 2, X ),':', X, psi( N1, 3, X ),':', ...
X, psi( N1, N1, X ),':' ); hold off;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% The basis functions based on the interval endpoints;
function [fun] = phi( N, k, x )
Nx = length( x );
h = 1.0 / N;
for i=1:Nx
if ( (x(i) > (k-1) * h) & (x(i) <= k * h) )
xi = ( k * h + (k-1) * h ) / 2.0;
fun(i) = 2.0 * ( x(i) - (k-1) * h ) * ( x(i) - xi ) / h^2;
elseif ( (x(i) > k * h) & (x(i) <= (k+1) * h) )
xi = ( (k+1) * h + k * h ) / 2.0;
fun(i) = 2.0 * ( x(i) - (k+1) * h ) * ( x(i) - xi ) / h^2;
else
fun(i) = 0.0;
end;
end;
% The basis functions based on the midpoints;
function [fun] = psi( N, k, x )
Nx = length( x );
h = 1.0 / N;
for i=1:Nx
if ( (x(i) > (k-1) * h) & (x(i) < k * h) )
fun(i) = 4.0 * ( k * h - x(i) ) * ( x(i) - (k-1)*h ) / h^2;
else
fun(i) = 0.0;
end;
end;