Simple Heat Transfer Tutorial¶
This tutorial provides an introduction to running simulations with Serac and demonstrates the setup of a simple steady-state heat transfer problem.
The full source code for this tutorial is available in examples/simple_conduction/without_input_file.cpp
, which demonstrates C++ configuration
of a heat transfer physics module.
The heat transfer modeled in this section is based on the formulation discussed in Heat Transfer.
Setting Up Includes and Initializing¶
The most important parts of Serac are its physics modules, each of which corresponds to a particular discretization
of a partial differential equation (e.g., continuous Galerkin finite element method for heat transfer).
In this example, we are building a heat transfer simulation, so we include Serac's HeatTransfer
module and
thermal material models:
#include "serac/physics/heat_transfer.hpp"
#include "serac/physics/materials/thermal_material.hpp"
The following header provides access to the StateManager class which manages the individual finite element states and the mesh:
#include "serac/physics/state/state_manager.hpp"
Serac also provides a set of setup/teardown functions that encapsulate much of the boilerplate setup required for each simulation, e.g., MPI initialization/finalization and logger setup/teardown, so we include their headers:
#include "serac/infrastructure/initialize.hpp"
#include "serac/infrastructure/terminator.hpp"
Finally, we include the header for Serac's mesh utilities, which includes support for reading meshes from a file and for generating meshes of common solids, like cuboids, rectangles, disks, and cylinders:
#include "serac/mesh/mesh_utils.hpp"
We're now ready to start our main()
function by initializing Serac, which performs the setup described above:
int main(int argc, char* argv[])
{
/*auto [num_procs, rank] = */ serac::initialize(argc, argv);
axom::sidre::DataStore datastore;
serac::StateManager::initialize(datastore, "without_input_file_example");
To simplify saving to an output file and restarting a simulation, Serac stores critical state information, like
the mesh and fields, in a single StateManager
object, which is initialized here.
Warning
Since Serac's initialization helper initializes MPI, you should not call MPI_Init
directly.
Constructing the Mesh¶
In this introductory example, we will use a simple square mesh with 10 quadrilateral elements in each space dimension
for 100 elements total. Once created, the primary mesh must always be registered with the StateManager
:
auto mesh = serac::mesh::refineAndDistribute(serac::buildRectangleMesh(10, 10));
std::string mesh_tag{"mesh"};
serac::StateManager::setMesh(std::move(mesh), mesh_tag);
After constructing the serial mesh, we call refineAndDistribute
to distribute it into a parallel mesh.
Constructing the Physics Module¶
// Create a Heat Transfer class instance with Order 1 and Dimensions of 2
constexpr int order = 1;
constexpr int dim = 2;
serac::HeatTransfer<order, dim> heat_transfer(
serac::heat_transfer::default_nonlinear_options, serac::heat_transfer::default_linear_options,
serac::heat_transfer::default_static_options, "thermal_solver", mesh_tag);
When using the C++ API, the HeatTransfer
constructor requires the polynomial order of the elements and the dimension
of the mesh at compile time, i.e. they are template parameters. We also need to pass the options for solving the nonlinear
system of equations and ordinary differential equations arising from the discretization. In this example, we use the default
static heat transfer options.
Configuring Material Conductivity¶
We define a material model that includes information needed for the constitutive response needed by the physics solver.
In this example, we define a linear isotropic conductor with uniform density, heat capacity, and conductivity (kappa).
That material model is then passed to the HeatTransfer
object. Note that this material model could be user-defined.
constexpr double kappa = 0.5;
serac::heat_transfer::LinearIsotropicConductor mat(1.0, 1.0, kappa);
heat_transfer.setMaterial(mat);
Setting Thermal (Dirichlet) Boundary Conditions¶
The following snippets add two Dirichlet boundary conditions:
One that constrains the temperature to 1.0 at boundary attribute 1, which for this mesh corresponds to the side of the square mesh along the x-axis, i.e., the "bottom" of the mesh.
One that constrains the temperature to \(x^2 + y - 1\) at boundary attributes 2 and 3, which for this mesh correspond to the right side and top of the mesh, respectively.
const std::set<int> boundary_constant_attributes = {1};
constexpr double boundary_constant = 1.0;
auto ebc_func = [boundary_constant](const auto&, auto) { return boundary_constant; };
heat_transfer.setTemperatureBCs(boundary_constant_attributes, ebc_func);
const std::set<int> boundary_function_attributes = {2, 3};
auto boundary_function_coef = [](const auto& vec, auto) { return vec[0] * vec[0] + vec[1] - 1; };
heat_transfer.setTemperatureBCs(boundary_function_attributes, boundary_function_coef);
Running the Simulation¶
Now that we've configured the HeatTransfer
instance using a few of its configuration
options, we're ready to run the simulation. We call completeSetup
to "finalize" the simulation
configuration, and then save off the initial state of the simulation. This also allocates and builds
all of the internal finite element data structures.
We can then perform the steady-state solve and save the end result:
heat_transfer.completeSetup();
heat_transfer.outputStateToDisk();
heat_transfer.advanceTimestep(1.0);
heat_transfer.outputStateToDisk();
Note
The dt
variable does not actually get used in a quasistatic simulation.
This should produce the following initial state:
Note the areas of the mesh where the boundary conditions were applied.
The final (steady-state) solution should look like the following:
The end result is not particularly impressive, but should be fairly intuitive.
Cleaning Up¶
To make sure Serac terminates properly, don't forget to call its exit function at the very end of your program:
serac::exitGracefully();
}