Demo: A linear elasticity problem with planar deformations

This code solves linearized elastostatics problems with planar deformations.

  • We start by defining an instance of class IPF to deal with project data. We extract useful information from this instance: An output flag to deal with output level and a save flag to say if the solution is to be saved in file. We then construct the Mesh and the Prescription instances.

       IPF data("lelas2d - 1.1 ",argv[1]);
       int output_flag = data.getOutput();
       int save_flag = data.getSave();
       Mesh ms(data.getMeshFile());
       Prescription p(ms,data.getDataFile());

  • We now declare vectors to contain the solution, boundary conditions, body and boundary forces and then initialize these vectors from prescribed data.

       Vect<double> u(ms), bc(ms), body_f(ms), bound_f(ms,2,BOUNDARY_SIDE_DOF);
       p.get(BOUNDARY_CONDITION,bc,0);
       p.get(BODY_FORCE,body_f);
       p.get(TRACTION,bound_f,0);

  • We can now instantiate the equation to solve: Here we consider the linearized elasticity equations in 2-D using bilinear quadrilateral elements. Namely, class Elas2DQ4.

       Elas2DQ4 eq(ms,u);
       eq.setInput(BOUNDARY_CONDITION,bc);
       eq.setInput(BODY_FORCE,body_f);
       eq.setInput(TRACTION,bound_f);
       eq.run();

  • We can save the resulting displacement vector for plotting or any other purpose. For this, the class IOField is used.

       if (save_flag) {
          IOField pl_file(data.getPlotFile(),IOField::OUT);
          pl_file.put(u);
       }

  • Let us illustrate how to do some postprocessing. We want here to compute principal and Von Mises stresses. For this, we have to use the member function Stress. The element principal stresses are stored in the instance st while the vector vm will contain the Von-Mises stresses. These vectors are instantiated here and computed using the member function Stress:

       Vect<double> st, vm;
       eq.Stress(st,vm);

  

An example

Let us run this program with the data presented in the following project file:

<?xml version="1.0" encoding="ISO-8859-1" ?>
<OFELI_File>
   <info>
      <title></title>
      <date></date>
      <author></author>
   </info>
   <Project name="beam">
      <verbose>1</verbose>
      <output>1</output>
      <save>1</save>
      <plot>1</plot>
      <init>1</init>
      <mesh_file>beam.m</mesh_file>
      <plot_file>beam.d</plot_file>
      <prescription_file>beam.pr</prescription_file>
   </Project>
</OFELI_File>

The only thing that may differ from other seen project files is that prescription data is given in a separate file: beam.pr. Let us give this one:

<?xml version="1.0" encoding="ISO-8859-1" ?>
<OFELI_File>
   <info>
      <title></title>
      <date></date>
      <author></author>
   </info>
   <Prescription>
      <BoundaryCondition code="1">0.0</BoundaryCondition>
      <Traction code="1" dof="2">-1000.0</Traction>
   </Prescription>
</OFELI_File>

Clearly, we have prescribed by this file a (Dirichlet) boundary condition that enforces a value of 0.0 for degrees of freedom with code 1. In addition, we have prescribed a traction of value -1000 for the second degrees of freedom (y-traction) having code 1.