Add inlet outlet bc

src/Control/Inciter/InputDeck/InputDeck.hpp

@@ -47,7 +47,6 @@ using ncomp_t = std::size_t;
 using bclist = tk::TaggedTuple< brigand::list<
   tag::dirichlet,   std::vector< std::size_t >,
   tag::symmetry,    std::vector< std::size_t >,
-  tag::inlet,       std::vector< std::size_t >,
   tag::outlet,      std::vector< std::size_t >,
   tag::farfield,    std::vector< std::size_t >,
   tag::extrapolate, std::vector< std::size_t >,
@@ -128,23 +127,29 @@ using materialList = tk::TaggedTuple< brigand::list<
 
 // Boundary conditions block
 using bcList = tk::TaggedTuple< brigand::list<
-  tag::mesh,           std::vector< std::size_t >,
-  tag::dirichlet,      std::vector< std::size_t >,
-  tag::symmetry,       std::vector< std::size_t >,
-  tag::inlet,          std::vector< std::size_t >,
-  tag::outlet,         std::vector< std::size_t >,
-  tag::farfield,       std::vector< std::size_t >,
-  tag::extrapolate,    std::vector< std::size_t >,
-  tag::noslipwall,     std::vector< std::size_t >,
-  tag::stag_point,     std::vector< tk::real >,
-  tag::radius,         tk::real,
-  tag::velocity,       std::vector< tk::real >,
-  tag::pressure,       tk::real,
-  tag::density,        tk::real,
-  tag::temperature,    tk::real,
+  tag::mesh,        std::vector< std::size_t >,
+  tag::dirichlet,   std::vector< std::size_t >,
+  tag::symmetry,    std::vector< std::size_t >,
+  tag::inlet,       std::vector<
+    tk::TaggedTuple< brigand::list<
+      tag::sideset,      std::vector< uint64_t >,
+      tag::velocity,     std::vector< tk::real >,
+      tag::materialid,   std::size_t
+    > >
+  >,
+  tag::outlet,      std::vector< std::size_t >,
+  tag::farfield,    std::vector< std::size_t >,
+  tag::extrapolate, std::vector< std::size_t >,
+  tag::noslipwall,  std::vector< std::size_t >,
+  tag::stag_point,  std::vector< tk::real >,
+  tag::radius,      tk::real,
+  tag::velocity,    std::vector< tk::real >,
+  tag::pressure,    tk::real,
+  tag::density,     tk::real,
+  tag::temperature, tk::real,
   tag::mass_fractions, std::vector< tk::real >,
-  tag::materialid,     std::size_t,
-  tag::timedep,        std::vector<
+  tag::materialid,  std::size_t,
+  tag::timedep,     std::vector<
     tk::TaggedTuple< brigand::list<
       tag::sideset,    std::vector< uint64_t >,
       tag::fn,         std::vector< tk::real >

src/Control/Inciter/InputDeck/LuaParser.cpp

@@ -1087,8 +1087,22 @@ LuaParser::storeInputDeck(
       storeVecIfSpecd< uint64_t >(sol_bc[i+1], "symmetry",
         bc_deck[i].get< tag::symmetry >(), {});
 
-      storeVecIfSpecd< uint64_t >(sol_bc[i+1], "inlet",
-        bc_deck[i].get< tag::inlet >(), {});
+      if (sol_bc[i+1]["inlet"].valid()) {
+        const sol::table& sol_inbc = sol_bc[i+1]["inlet"];
+        auto& inbc_deck = bc_deck[i].get< tag::inlet >();
+        inbc_deck.resize(sol_inbc.size());
+
+        for (std::size_t j=0; j<inbc_deck.size(); ++j) {
+          storeVecIfSpecd< uint64_t >(sol_inbc[j+1], "sideset",
+            inbc_deck[j].get< tag::sideset >(), {});
+          storeVecIfSpecd< tk::real >(sol_inbc[j+1], "velocity",
+            inbc_deck[j].get< tag::velocity >(), {0.0, 0.0, 0.0});
+          if (inbc_deck[j].get< tag::velocity >().size() != 3)
+            Throw("Inlet velocity requires 3 components.");
+          storeIfSpecd< std::size_t >(sol_inbc[j+1], "materialid",
+            inbc_deck[j].get< tag::materialid >(), 1);
+        }
+      }
 
       storeVecIfSpecd< uint64_t >(sol_bc[i+1], "outlet",
         bc_deck[i].get< tag::outlet >(), {});

src/Inciter/Transporter.cpp

@@ -475,14 +475,18 @@ Transporter::matchBCs( std::map< int, std::vector< std::size_t > >& bnd )
 // *****************************************************************************
 {
   // Query side set ids at which BCs assigned for all BC types for all PDEs
-  using bclist = ctr::bclist::Keys;
   std::unordered_set< int > usedsets;
-  brigand::for_each< bclist >( UserBC( g_inputdeck, usedsets ) );
 
-  // Query side sets of time dependent BCs (since tag::bctimedep is not a part
+  // Query side sets of time dependent and inlet BCs (since these are not a part
   // of tag::bc)
+  using bclist = ctr::bclist::Keys;
   const auto& bcs = g_inputdeck.get< tag::bc >();
+  brigand::for_each< bclist >( UserBC( g_inputdeck, usedsets ) );
   for (const auto& bci : bcs) {
+    for (const auto& b : bci.get< tag::inlet >()) {
+      for (auto i : b.get< tag::sideset >())
+        usedsets.insert(static_cast<int>(i));
+    }
     for (const auto& b : bci.get< tag::timedep >()) {
       for (auto i : b.get< tag::sideset >())
         usedsets.insert(static_cast<int>(i));

src/PDE/CompFlow/DGCompFlow.hpp

@@ -72,23 +72,21 @@ class CompFlow {
       // associate boundary condition configurations with state functions, the
       // order in which the state functions listed matters, see ctr::bc::Keys
       brigand::for_each< ctr::bclist::Keys >( ConfigBC( m_bc,
-        // BC State functions
-        { dirichlet
-        , symmetry
-        , invalidBC         // Inlet BC not implemented
-        , invalidBC         // Outlet BC not implemented
-        , farfield
-        , extrapolate
-        , invalidBC },      // No slip wall BC not implemented
-        // BC Gradient functions
-        { noOpGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad }
-        ) );
+      // BC State functions
+      { dirichlet
+      , symmetry
+      , invalidBC         // Outlet BC not implemented
+      , farfield
+      , extrapolate
+      , invalidBC },      // No slip wall BC not implemented
+      // BC Gradient functions
+      { noOpGrad
+      , noOpGrad
+      , noOpGrad
+      , noOpGrad
+      , noOpGrad
+      , noOpGrad }
+      ) );
 
       // EoS initialization
       const auto& matprop =

src/PDE/ConfigureTransport.cpp

@@ -101,8 +101,10 @@ infoTransport( std::map< ctr::PDEType, tk::ncomp_t >& cnt )
     const auto& bcinlet =
       ib.get< tag::inlet >();
     if (!bcinlet.empty())
-      nfo.emplace_back( "Inlet boundary [" + std::to_string( ncomp ) + "]",
-        parameters( bcinlet ) );
+      for (const auto& bndry : bcinlet) {
+        nfo.emplace_back( "Inlet boundary [" + std::to_string( ncomp ) + "]",
+          parameters(bndry.get< tag::sideset >()) );
+    }
 
     const auto& bcoutlet =
       ib.get< tag::outlet >();

src/PDE/MultiMat/BCFunctions.hpp

@@ -116,6 +116,85 @@ namespace inciter {
     return {{ std::move(ul), std::move(ur) }};
   }
 
+  //! \brief Boundary state function providing the left and right state of a
+  //!   face at inlet boundaries
+  //! \param[in] ncomp Number of scalar components in this PDE system
+  //! \param[in] ul Left (domain-internal) state
+  //! \return Left and right states for all scalar components in this PDE
+  //!   system
+  //! \details The inlet boundary condition specifies a velocity at a
+  //!   sideset and assumes a zero bulk pressure and density gradient
+  //! \note The function signature must follow tk::StateFn
+  static tk::StateFn::result_type
+  inlet( ncomp_t ncomp,
+            const std::vector< EOS >& mat_blk,
+            const std::vector< tk::real >& ul,
+            tk::real, tk::real, tk::real, tk::real,
+            const std::array< tk::real, 3 >&  )
+  {
+    auto nmat = g_inputdeck.get< tag::multimat, tag::nmat >();
+    const auto& solidx = g_inputdeck.get< tag::matidxmap, tag::solidx >();
+    auto& inbc = g_inputdeck.get< tag::bc >()[0].get< tag::inlet >();
+
+    // inlet velocity and material
+    auto u_in = inbc[0].get< tag::velocity >();
+    auto mat_in = inbc[0].get< tag::materialid >() - 1;
+
+    [[maybe_unused]] auto nsld = numSolids(nmat, solidx);
+
+    Assert( ul.size() == ncomp+nmat+3+nsld*6, "Incorrect size for appended "
+            "internal state vector" );
+
+    auto ur = ul;
+
+    // Internal cell velocity components
+    auto v1l = ul[ncomp+velocityIdx(nmat, 0)];
+    auto v2l = ul[ncomp+velocityIdx(nmat, 1)];
+    auto v3l = ul[ncomp+velocityIdx(nmat, 2)];
+
+    // External cell velocity, such that velocity = v_in at face
+    auto v1r = 2.0*u_in[0] - v1l;
+    auto v2r = 2.0*u_in[1] - v2l;
+    auto v3r = 2.0*u_in[2] - v3l;
+
+    // Calculate bulk pressure/density to ensure zero gradient
+    tk::real p(0.0);
+    tk::real rho(0.0);
+    for (std::size_t k=0; k<nmat; ++k) {
+      p += ul[ncomp+pressureIdx(nmat,k)];
+      rho += ul[densityIdx(nmat,k)];
+    }
+
+    tk::real alphamin = 1e-12;
+    for (std::size_t k=0; k<nmat; ++k) {
+      if (k == mat_in)
+        ur[volfracIdx(nmat,k)] = 1.0 -
+          (static_cast< tk::real >(nmat-1))*alphamin;
+      else
+        ur[volfracIdx(nmat,k)] = alphamin;
+
+      // zero bulk pressure and density gradient
+      ur[ncomp+pressureIdx(nmat, k)] = ur[volfracIdx(nmat, k)] * p;
+      ur[densityIdx(nmat,k)] = ur[volfracIdx(nmat,k)] * rho;
+      ur[energyIdx(nmat,k)] = ur[volfracIdx(nmat,k)] *
+        mat_blk[k].compute< EOS::totalenergy >(rho, v1r, v2r, v3r, p);
+    }
+
+    ur[momentumIdx(nmat, 0)] = rho * v1r;
+    ur[momentumIdx(nmat, 1)] = rho * v2r;
+    ur[momentumIdx(nmat, 2)] = rho * v3r;
+
+    // velocity
+    ur[ncomp+velocityIdx(nmat, 0)] = v1r;
+    ur[ncomp+velocityIdx(nmat, 1)] = v2r;
+    ur[ncomp+velocityIdx(nmat, 2)] = v3r;
+
+    Assert( ur.size() == ncomp+nmat+3+nsld*6, "Incorrect size for appended "
+            "boundary state vector" );
+
+    return {{ std::move(ul), std::move(ur) }};
+  }
+
   //! \brief Boundary state function providing the left and right state of a
   //!   face at farfield boundaries
   //! \param[in] ncomp Number of scalar components in this PDE system

src/PDE/MultiMat/DGMultiMat.hpp

@@ -76,23 +76,36 @@ class MultiMat {
     {
       // associate boundary condition configurations with state functions
       brigand::for_each< ctr::bclist::Keys >( ConfigBC( m_bc,
-        // BC State functions
-        { dirichlet
-        , symmetry
-        , invalidBC         // Inlet BC not implemented
-        , invalidBC         // Outlet BC not implemented
-        , farfield
-        , extrapolate
-        , noslipwall },
-        // BC Gradient functions
-        { noOpGrad
-        , symmetryGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad
-        , noOpGrad }
-        ) );
+      // BC State functions
+      { dirichlet
+      , symmetry
+      , invalidBC         // Outlet BC not implemented
+      , farfield
+      , extrapolate
+      , noslipwall },
+      // BC Gradient functions
+      { noOpGrad
+      , symmetryGrad
+      , noOpGrad
+      , noOpGrad
+      , noOpGrad
+      , noOpGrad }
+      ) );
+
+      // Inlet BC has a different structure than above BCs, so it must be 
+      // handled differently than with ConfigBC
+      const auto& bc = g_inputdeck.get< tag::bc >();
+      std::vector< std::size_t > v;
+      for (const auto& ib : bc) {
+        const auto& in = ib.get< tag::inlet >();
+        if (!in.empty()) {
+          for (const auto& bndry : in) {
+            const auto& sideset = bndry.get< tag::sideset >();
+            v.insert(v.end(), sideset.begin(), sideset.end());
+            m_bc.push_back( { v, inlet, noOpGrad } );
+          }
+        }
+      }
 
       // EoS initialization
       initializeMaterialEoS( m_mat_blk );

src/PDE/MultiMat/FVMultiMat.hpp

@@ -75,7 +75,6 @@ class MultiMat {
         // BC State functions
         { dirichlet
         , symmetry
-        , invalidBC         // Inlet BC not implemented
         , invalidBC         // Outlet BC not implemented
         , farfield
         , extrapolate
@@ -86,10 +85,25 @@ class MultiMat {
         , noOpGrad
         , noOpGrad
         , noOpGrad
-        , noOpGrad
         , noOpGrad }
         ) );
 
+      // Inlet BC has a different structure than above BCs, so it must be 
+      // handled differently than with ConfigBC
+      const auto& bc = g_inputdeck.get< tag::bc >();
+      std::vector< std::size_t > v;
+      for (const auto& ib : bc) {
+        const auto& in = ib.get< tag::inlet >();
+        if (!in.empty()) {
+          for (const auto& bndry : in) {
+            const auto& sideset = bndry.get< tag::sideset >();
+            v.clear();
+            v.insert(v.end(), sideset.begin(), sideset.end());
+            m_bc.push_back( { v, inlet, noOpGrad } );
+          }
+        }
+      }
+
       // EoS initialization
       initializeMaterialEoS( m_mat_blk );
     }

src/PDE/Transport/DGTransport.hpp

@@ -71,7 +71,6 @@ class Transport {
         // BC State functions
         { dirichlet
         , invalidBC  // Symmetry BC not implemented
-        , inlet
         , outlet
         , invalidBC  // Characteristic BC not implemented
         , extrapolate
@@ -82,9 +81,24 @@ class Transport {
         , noOpGrad
         , noOpGrad
         , noOpGrad
-        , noOpGrad
         , noOpGrad }
         ) );
+
+      // Inlet BC has a different structure than above BCs, so it must be 
+      // handled differently than with ConfigBC
+      const auto& bc = g_inputdeck.get< tag::bc >();
+      std::vector< std::size_t > v;
+      for (const auto& ib : bc) {
+        const auto& in = ib.get< tag::inlet >();
+        if (!in.empty()) {
+          for (const auto& bndry : in) {
+            const auto& sideset = bndry.get< tag::sideset >();
+            v.clear();
+            v.insert(v.end(), sideset.begin(), sideset.end());
+            m_bc.push_back( { v, inlet, noOpGrad } );
+          }
+        }
+      }
       m_problem.errchk( m_ncomp );
     }
 

tests/regression/inciter/mesh_refinement/dtref/gauss_hump.lua

@@ -18,7 +18,11 @@ inciter = {
   bc = {
     {
       extrapolate = { 1 },
-      inlet = { 2 },
+      inlet = {
+        {
+          sideset = { 2 }
+        }
+      },
       outlet = { 3 }
     }
   },

tests/regression/inciter/mesh_refinement/dtref/gauss_hump_reord.lua

@@ -19,7 +19,11 @@ inciter = {
   bc = {
     {
       extrapolate = { 1 },
-      inlet = { 2 },
+      inlet = {
+        {
+          sideset = { 2 }
+        }
+      },
       outlet = { 3 }
     }
   },