Advanced Command Line Postprocessing and Visualization

Banner - Advanced Command Line Postprocessing and Visualization

The following post breaks down and explains the postprocessing and visualization steps used to create the images and video for the axisymmetric stress-strain brake disk analysis model. By using the Octave or MATLAB command line and exported finite element fea data struct one can create custom graphics and visualizations that is currently not possible to do within the FEATool graphical user interface (GUI).

The first step in creating the advanced visualizations is to run the FEATool model example m-script file to generate and return the finite element analysis struct fea, together with the solution stored in the fea.sol.u field (in this case the automatic pre-defined plots can be turned off with the ‘iplot’, ‘0’ flag, since postprocessing will be done manually)

fea = ex_axistressstrain4( 'iplot', 0 );

Next, model constants and parameters used in during postprocessing are prescribed

deltad = 5.5e-3;     % Disk thickness.
rd     = 66e-3;      % Disk inner radius.
rp     = 75.5e-3;    % Pad inner radius.
Rd     = 113.5e-3;   % Disk outer radius.

Continuing, the solution vector for all time steps is saved in a temporary variable u (the solution at each time step corresponds to a column in u). The parameters cmax and cmin are also calculated which later are used to fix the colormap range (otherwise not much difference between hot and cold maxima can be seen)

u    = fea.sol.u;
cmax = max(u(:)) - 273.15;
cmin = 0;

Coordinates to plot and visualize the contours and boundaries of the brake disk must also be defined (here the assumption is that the disk is centered at the origin with the radial coordinate spanning the x-z plane)

% Points around the circumference.
npth = 72;
th   = linspace( 0, 2*pi, npth );

% Inner coordinates.
x    = rd*cos(th);
z    = rd*sin(th);

% Outer coordinates.
xx   = Rd*cos(th);
zz   = Rd*sin(th);

y    = deltad*ones(size(x));

The ringgrid function is also used to generate coordinates around the disk in the grid struct variable g and evaluation points in rz. These will be used to interpolate the temperature field from axisymmetric coordinates to three dimensions (the ring grid can be seen and visualized by itself using the plotgrid( g ) command)

npr  = 25;
g    = ringgrid( npr-1, npth-1, rd, Rd );
r    = sqrt( g.p(1,:).^2 + g.p(2,:).^2 );
rz   = [rd+sqrt(eps) rd+(Rd-rd)/2 Rd-sqrt(eps);deltad/2*[1 1 1]];

Now, the actual plotting is initialized by creating a figure window to draw in, an empty array for the stresses at the three points (inner radius, outer radius, and the mid point), and starting the postprocessing loop

st = zeros(3,1);
for i=1:size(u,2)

The use subplot MATLAB command is used to hold two horizontal images. First is the 3D temperature field, starting with drawing the black outline of the brake disk

  hold on

  plot3(x,  y, z,'k-')
  plot3(x, -y, z,'k-')
  plot3(xx, y,zz,'k-')

The i:th solution is selected as the only solution in fea.sol.u (since evalexpr per default uses the solution at the last time step), after which the temperature is evaluated along the radius at the disk center line coordinate (z=0).

  fea.sol.u = u(:,i);
  T = evalexpr( 'T-273.15', [r;zeros(size(r))], fea );

Using the assumption that the solution is axisymmetric, one can extrapolate the data all along the brake disk with the help of the annular grid struct g created earlier. The MATLAB patch command is used to plot the x-z temperature field

  h = patch( 'faces',     g.c', ...
             'vertices',  [g.p(1,:)' zeros(size(g.p,2),1) g.p(2,:)'], ...
             'facevertexcdata', T, ...
             'facecolor', 'interp', ...
             'linestyle', 'none' );

Furthermore, the temperature for the original patch in the positive and mirrored (negative) z-directions is also plotted

  postplot( fea, 'surfexpr', 'T-273.15', 'colorbar', 'off' )
  fea.grid.p(2,:) = -fea.grid.p(2,:);
  postplot( fea, 'surfexpr', 'T-273.15', 'colorbar', 'off' )
  fea.grid.p(2,:) = -fea.grid.p(2,:);

Lastly, a title is set, the color map range fixed, and a suitable viewing angle is chosen

  title( 'Temperature' )
  caxis( [cmin cmax] )
  axis off
  axis tight

For the second subplot the idea is to follow the von Mieses stress as a function of time for the three points prescribed in rz. Just as for the temperature, the stress can be evaluated with the evalexpr command where the stress expression svm has been defined earlier in the fea.const field of the model fea struct

  hold on
  grid on
  if( i>1 )   % Skip first zero solution
    st = [st evalexpr( '(svm)*1e-6', rz, fea )];

The stress is plotted as a function of time, with additional text and legends added for identification

    plot( tlist(1:i), st(1,:) )
    plot( tlist(1:i), st(2,:) )
    plot( tlist(1:i), st(3,:) )

    text( tlist(i), st(1,end), 'r = rd' )
    text( tlist(i), st(2,end), 'r = (rd+Rd)/2' )
    text( tlist(i), st(3,end), 'r = Rd' )

    axis( [0 tmax 0 max(st(:))] )
  xlabel( 'Time [s]' )
  ylabel( 'von Mieses stress [MPa]' )

Finally, before closing the loop, a drawing update is forced with the drawnow command, and optionally an jpeg image file is written out and saved for later processing

  print( '-r200', '-djpeg', sprintf('img%03i.jpg',i) )

After all the MATLAB postprocessing has finished the images can easily be converted to a video with a suitable tool, such as for example ffmpeg with the command

ffmpeg -i img%03d.jpg -c:v libx264 -vf "fps=25, format=yuv420p, crop=1600:1198:0:0" out.mp4

From the generated video one can now clearly see that the temperature is the highest on the outer two-thirds of the disk consistent with where the brake pad is acting. In contrast the maximum stress is achieved at the inner radius r = rd around t = 3.0 s, and is significantly higher than at the two points measured further out on the disk.

Category: postprocessing

Tags: cli


Twitter LinkedIn Facebook