Demo: A planar truss

This code models an assembly of elastic bars subject to elongation and compression. The libray contains the class Bar2DL2 that handles elastic planar with 2 degrees of freedom for each node that are x and y displacements at nodes.

  • 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("truss",argv[1]);
       int output_flag = data.getOutput();
       int save_flag = data.getSave();
       Mesh ms(data.getMeshFile());
       Prescription p(ms,data.getDataFile());

  • We declare an instance of class SkSMatrix<double> to store the stiffness matrix. Then the vectors bc and bf are declared to contain boundary conditions and loads respectively. We fill these vectors by using the class Prescription. We finally use a Vect<double> instance that contains the right-hand side and later, the solution. The bar class needs to provide the section of the bar. This is read from the project file.

       SkSMatrix<double> A(ms);
       Vect<double> bc(ms), bf(ms);
       p.get(BOUNDARY_CONDITION,bc);
       p.get(POINT_FORCE,bf);
       Vect<double> b(bf);
       double section = data.getDouble("section");

  • The global matrix and right-hand side are filled by looping over elements. For each element, an instance of class Bar2DL2 is considered. For each element, we add the stiffness contribution and then assemble the local matrix into the global one.

       MeshElements(ms) {
          Bar2DL2 eq(theElement,section);
          eq.Stiffness();
          eq.ElementAssembly(A);
       }

  • Boundary conditions are enforces by the penalty method. Then the resulting linear system is solved by factorization and backsubstitution.

       A.Prescribe(b);
       A.solve(b);

  • A last step consists here in computing the deformed structure by using the function DeformMesh. We then store the solution in appropriate file.

       IOField pf(data.getMeshFile(),data.getPlotFile(),ms,IOField::OUT);
       pf.put(b);
       DeformMesh(ms,b,40);

  

An example

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


<OFELI_File>
   <info>
      <title />
      <date />
      <author />
   </info>
   <Project name="tower">
      <mesh_file value="tower.dat"/>
      <verbose value="1"/>
      <output value="1"/>
      <save value="1"/>
      <plot value="1"/>
      <init value="1"/>
      <parameter label="section" value="0.005"/>
      <plot_file value="tower1.d"/>
   </Project>
                              
   <Mesh dim="2" nb_dof="2">
      <Nodes>
         0.00000000e+00      0.00000000e+00    11         6.00000000e+00      0.00000000e+00    11
         3.00000000e+00      1.50000000e+00     0         0.00000000e+00      3.00000000e+00     0
         6.00000000e+00      3.00000000e+00     0         3.00000000e+00      4.50000000e+00     0
         0.00000000e+00      6.00000000e+00     0         6.00000000e+00      6.00000000e+00     0
         3.00000000e+00      7.50000000e+00     0         0.00000000e+00      9.00000000e+00     0	  
         6.00000000e+00      9.00000000e+00     0         3.00000000e+00      1.05000000e+01     0
         0.00000000e+00      1.20000000e+01     0         6.00000000e+00      1.20000000e+01     0
         3.00000000e+00      1.50000000e+01     0   
      </Nodes>
      <Elements shape="line"  nodes="2">
         1       2       1          1       3       1          2       3       1
         3       4       1          3       5       1          1       4       1
         2       5       1          4       5       1          4       6       1
         5       6       1          6       7       1          6       8       1
         4       7       1          5       8       1          7       8       1
         7       9       1          8       9       1          9      10       1
         9      11       1          7      10       1          8      11       1
        10      11       1         10      12       1         11      12       1
        12      13       1         12      14       1         10      13       1
        11      14       1         13      14       1         13      15       1   
        14      15       1   
      </Elements>
      <Material>
         1   Iron
      </Material>
   </Mesh>

   <Prescription>
      <PointForce y="15.0" dof="1">5.0</PointForce>
   </Prescription>
</OFELI_File>

The mesh here above describes a tower. Note that the node with coordinate y=15 is submitted to a concentrated load equal to 15.0.