Demo: A 3-D steady state diffusion code

This code is developed to solve 3-D steady-state (stationary) heat transfer problems. The code is very similar to previous demo codes: stdc2 and ttd2.

  • We construct an instance of class IPF to manage 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. 

       IPF data("std3 - 1.0",argv[1]);
   int output_flag = data.getOutput();
   int save_flag = data.getSave();

  • After outputting some information about the code, we declare an instance of class Mesh. 

       Mesh ms(data.getMeshFile());

  • In order to read various prescriptions (boundary conditions, heat sources, fluxes), we use the class Prescription. The name of the file containing these prescriptions is given in data.getDataFile() 

       Prescription p(ms,data.getDataFile());

  • We declare vectors that will contain the solution and various data. 

       Vect<double> u(ms);
   Vect<double> bc(ms), bf(ms), sf(ms,BOUNDARY_SIDE_DOF,1);

  • Now we can read various prescriptions: we first instantiate vectors bc, bf and sf that will contain respectively Dirichlet boundary conditions and sources. Neumann boundary conditions (imposed fluxes) are stored in a vector with mentioning that this one is a SIDE_FIELD. Then we use the member function get with various arguments to fill properly these vectors. 

       p.get(BOUNDARY_CONDITION,bc);
   p.get(BODY_FORCE,bf);
   p.get(BOUNDARY_FORCE,sf,0);

  • We declare as an equation the Diffusion-Convection equation in 3-D using P1 (Tetrahedral), namely class DC3DT4: 

       DC3DT4 eq(ms,u);

  • We define as inputs of the equations various data as usual: 

       eq.setInput(BOUNDARY_CONDITION,bc);
   eq.setInput(SOURCE,bf);
   eq.setInput(FLUX,sf);

  • We next define equation terms. Here we add the diffusion term, without convection and then define the Conjugate Gradient, combined with the diagonal ILU preconditioner, as iterative solver of the resulting linear system: 

       eq.setTerms(DIFFUSION);
   eq.setSolver(CG_SOLVER,DILU_PREC);

  • We finally run the equation, print the solution and save this solution in a Gmsh file. 

       eq.run();
   if (Verbosity>4)
      cout << u;
   if (save_flag)
      saveField(u,"beam.pos",GMSH);

      

    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">
      <mesh_file>beam.m</mesh_file>
      <verbose>2</verbose>
      <save>1</save>
   </Project>
   <Prescription>
      <BoundaryCondition code="1">1.0</BoundaryCondition>
   </Prescription>
</OFELI_File>

    In summary, this file looks mainly like the one in the previous example.