FEATool Multiphysics  v1.16.1 Finite Element Analysis Toolbox
ex_spanner.m File Reference

## Description

EX_SPANNER STL geometry import and stress calculation of a spanner.

[ FEA, OUT ] = EX_SPANNER( VARARGIN ) Example to import a STL CAD geometry and calculate displacements on a fixed spanner. The load force may be distributed in the tangential load direction with the force fraction parameter FRAC.

Accepts the following property/value pairs.

Input       Value/{Default}        Description
-----------------------------------------------------------------------------------
E           scalar {190e3}         Modulus of elasticity [N/mm^2]
nu          scalar {0.29}          Poissons ratio
force       scalar {1000}          Load force [N]
frac        scalar {0}             Fraction of stress against pulling direction
sfun        string {sflag1}        Shape function for displacements
iplot       scalar 0/{1}           Plot solution (=1)
.
Output      Value/(Size)           Description
-----------------------------------------------------------------------------------
fea         struct                 Problem definition struct
out         struct                 Output struct


# Code listing

 cOptDef = { 'E',        190e3;
'nu',       0.29;
'force',    1000;
'frac',     0;
'sfun',     'sflag1';
'iplot',    1;
'fid',      1 };
[got,opt] = parseopt(cOptDef,varargin{:});

% Define scaling factor m to mm.
USE_METERS = true;
s = double( ~USE_METERS + USE_METERS*1e-3 );

% Import (and scale) grid.
fea.sdim = {'x','y','z'};
fea.geom = impexp_stl( 'spanner1.stl', 'import', [], 'solid', 2, 'extend', 0 );
fea.grid = gridgen( fea, 'hmax', 5, 'gridgen', 'robust', 'fid', opt.fid );
fea.grid.p = fea.grid.p*s;

% Add linear elasticity physics mode and define material parameters.
fea = addphys( fea, @linearelasticity );
fea.phys.el.sfun            = { opt.sfun opt.sfun opt.sfun };
fea.phys.el.eqn.coef{1,end} = { opt.nu };
fea.phys.el.eqn.coef{2,end} = { opt.E/s^2 };

% Set all boundaries to no load per default.
n_bdr  = max(fea.grid.b(3,:));
bc_sel = cell(3,n_bdr);
[bc_sel{:}] = deal(0);

% Fix all displacements on mandible boundaries.
i_fix   = [13 14];
i_force = [3 1];
faxis   = [1 3];
[bc_sel{:,i_fix}] = deal(1);
fea.phys.el.bdr.coef{5} = bc_sel;

% Apply force for x > 140 mm.
force = opt.force/(6*s*80*s);
fea.phys.el.bdr.coef{7}{faxis(1),i_force(1)} = ['-',num2str((1-opt.frac)*force),'*(y>140*',num2str(s),')'];
fea.phys.el.bdr.coef{7}{faxis(2),i_force(2)} = [num2str(opt.frac*force),'*(y>140*',num2str(s),')'];

% Parse and solve problem.
fea = parsephys(fea);
fea = parseprob(fea);
fea.sol.u = solvestat( fea, 'fid', opt.fid, 'icub', 1+str2num(strrep(opt.sfun,'sflag','')) );

% Postprocessing.
if( opt.iplot>0 )
subplot(1,2,1)
postplot( fea, 'surfexpr', ['sqrt(u^2+v^2+w^2)/',num2str(s)] )
view([30 20])
title('Total displacement (mm)')

subplot(1,2,2)
DSCALE = 5;
dp = zeros(size(fea.grid.p));
for i=1:3
dp(i,:) = DSCALE*evalexpr( fea.dvar{i}, fea.grid.p, fea );
end
fea_disp.grid   = fea.grid;
fea_disp.grid.p = fea_disp.grid.p + dp;
plotgrid( fea_disp )
title(['Displacement plot'])
view([30 20])
end

% Error checking.
u = fea.sol.u(unique(fea.eqn.dofm{1}(:)))/s;
v = fea.sol.u(unique(fea.eqn.ndof(1)+fea.eqn.dofm{2}(:)))/s;
w = fea.sol.u(unique(sum(fea.eqn.ndof(1:2))+fea.eqn.dofm{3}(:)))/s;
out.disp = sqrt( u.^2 + v.^2 + w.^2 );
out.pass = nan;
if( ~(got.frac || got.E || got.nu || got.force) )
out.pass = abs( max(out.disp) - 2.5 )/2.5 < 0.1;
end

if( nargout==0 )
clear fea out
else
fea.grid.p = fea.grid.p/s;
end