FEATool Multiphysics  v1.16.5 Finite Element Analysis Toolbox
ex_periodic2.m File Reference

## Description

EX_PERIODIC2 2D Periodic Poisson equation example.

[ FEA, OUT ] = EX_PERIODIC2( VARARGIN ) 2D Periodic Poisson equation example. Accepts the following property/value pairs.

Input       Value/{Default}        Description
-----------------------------------------------------------------------------------
isol        scalar 0/{1}           Stationary/Time dependent solution
iplot       scalar 0/{1}           Plot solution and error (=1)
.
Output      Value/(Size)           Description
-----------------------------------------------------------------------------------
fea         struct                 Problem definition struct
out         struct                 Output struct


# Code listing

 cOptDef = { 'isol',     0;
'iplot',    1;
'fid',      1 };
[got,opt] = parseopt(cOptDef,varargin{:});
fid       = opt.fid;

% Unit square grid.
fea.sdim = {'x' 'y'};
fea.grid = rectgrid( 40 );

% Define Poisson equation with custom physics mode.
fea = addphys( fea, @customeqn, { 'u' } );
fea.phys.ce.eqn.seqn = '- ux_x - uy_y = f + u*eps';

% Boundary conditions.
DirBC = true;
NeuBC = ~DirBC;
NeuBC_coef = @periodic_bc;   % Assign periodic_bc solve hook function as boundary coefficient.
fea.phys.ce.bdr.coef = { 'bcnd_ce', '', '', {}, ...
{ DirBC, NeuBC,      DirBC, NeuBC }, [], ...
{ 0    , NeuBC_coef, 0,     0     } };

% Source term.
fea.expr = { 'f' 'x*sin(5*pi*y) + exp(-((x-0.5)^2 + (y-0.5)^2)/0.02)' };

% Parse and solve problem.
fea = parsephys( fea );
fea = parseprob( fea );
if( opt.isol==0 )
fea.sol.u = solvestat( fea, 'fid', opt.fid );
else
fea.sol.u = solvetime( fea, 'tstop', -eps, 'fid', opt.fid );
end

% Postprocessing.
y  = linspace( 0, 1, 100 );
ul = evalexpr( 'u', [zeros(1,100);y], fea );
ur = evalexpr( 'u', [ones(1,100);y],  fea );
if( opt.iplot>0 )
subplot(1,2,1)
postplot(fea, 'surfexpr', 'u', 'surfhexpr', 'u' )

% Plot solution on left and right sides.
hold on
fea.grid.p(1,:) = fea.grid.p(1,:) - 1;
postplot(fea, 'surfexpr', 'u', 'surfhexpr', 'u' )
fea.grid.p(1,:) = fea.grid.p(1,:) + 2;
postplot(fea, 'surfexpr', 'u', 'surfhexpr', 'u' )

subplot(1,2,2)
plot( y, ul, 'k-' )
hold on
plot( y, ur, 'r.' )
grid on
xlabel( 'y' )
ylabel( 'u' )
legend( 'u(x=0)', 'u(x=1)', 'Location', 'South' )
end

% Error checking.
err = norm( ur - ul );

out.err  = err;
out.pass = err<1e-6;
if ( nargout==0 )
clear fea out
end

%