Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/Inciter/NodeBC.cpp
4 : : \copyright 2012-2015 J. Bakosi,
5 : : 2016-2018 Los Alamos National Security, LLC.,
6 : : 2019-2021 Triad National Security, LLC.
7 : : All rights reserved. See the LICENSE file for details.
8 : : \brief Boundary conditions for nodal discretizations
9 : : \details Boundary conditions for nodal discretizations, such as continuous
10 : : Galerkin finite elements.
11 : : */
12 : : // *****************************************************************************
13 : :
14 : : #include <vector>
15 : : #include <map>
16 : : #include <unordered_map>
17 : : #include <algorithm>
18 : :
19 : : #include "NodeBC.hpp"
20 : : #include "CGPDE.hpp"
21 : : #include "Fields.hpp"
22 : : #include "Vector.hpp"
23 : :
24 : : namespace inciter {
25 : :
26 : : extern std::vector< CGPDE > g_cgpde;
27 : :
28 : : std::unordered_map< std::size_t, std::vector< std::pair< bool, tk::real > > >
29 : 45114 : match( std::size_t meshid,
30 : : [[maybe_unused]] tk::ncomp_t ncomp,
31 : : tk::real t,
32 : : tk::real dt,
33 : : const std::vector< tk::real >& tp,
34 : : const std::vector< tk::real >& dtp,
35 : : const tk::UnsMesh::Coords& coord,
36 : : const std::unordered_map< std::size_t, std::size_t >& lid,
37 : : const std::map< int, std::vector< std::size_t > >& bnode,
38 : : bool increment )
39 : : // *****************************************************************************
40 : : // Match user-specified boundary conditions at nodes for side sets
41 : : //! \param[in] ncomp Number of scalar components in PDE system
42 : : //! \param[in] t Physical time at which to query boundary conditions
43 : : //! \param[in] dt Time step size (for querying BC increments in time)
44 : : //! \param[in] tp Physical time for each mesh node
45 : : //! \param[in] dtp Time step size for each mesh node
46 : : //! \param[in] coord Mesh node coordinates
47 : : //! \param[in] lid Local node IDs associated to local node IDs
48 : : //! \param[in] bnode Map storing global mesh node IDs mapped to side set ids
49 : : //! \param[in] increment If true, evaluate the solution increment between
50 : : //! t and t+dt for Dirichlet BCs. If false, evlauate the solution instead.
51 : : //! \return Vector of pairs of bool and boundary condition value associated to
52 : : //! local mesh node IDs at which the user has set Dirichlet boundary
53 : : //! conditions for all systems of PDEs integrated. The bool indicates whether
54 : : //! the BC is set at the node for that component: if true, the real value is
55 : : //! the increment (from t to dt) in (or the value of) the BC specified for a
56 : : //! component.
57 : : //! \details Boundary conditions (BC), mathematically speaking, are applied on
58 : : //! finite surfaces. These finite surfaces are given by element sets (i.e., a
59 : : //! list of elements). This function queries Dirichlet boundary condition
60 : : //! values from all PDEs in the multiple systems of PDEs integrated at the
61 : : //! node lists associated to side set IDs, given by bnode. Each
62 : : //! PDE system returns a BC data structure. Note that the BC mesh nodes that
63 : : //! this function results in (stored in dirbc) only contains those nodes that
64 : : //! are supplied via bnode. i.e., in parallel only a part of the mesh is
65 : : //! worked on.
66 : : // *****************************************************************************
67 : : {
68 : : using inciter::g_cgpde;
69 : :
70 : : // Vector of pairs of bool and boundary condition value associated to mesh
71 : : // node IDs at which the user has set Dirichlet BCs for all PDEs integrated.
72 : : std::unordered_map< std::size_t,
73 : 45114 : std::vector< std::pair< bool, tk::real > > > dirbc;
74 : :
75 : : // Details for the algorithm below: PDE::dirbc() returns a new map that
76 : : // associates a vector of pairs associated to local node IDs. (The pair is a
77 : : // pair of bool and real value, the former is the fact that the BC is to be
78 : : // set while the latter is the value if it is to be set). The length of this
79 : : // NodeBC vector, returning from each system of PDEs equals to the number of
80 : : // scalar components the given PDE integrates. Here we contatenate this map
81 : : // for all PDEs being integrated. If there are multiple BCs set at a mesh node
82 : : // (dirbc::key), either because (1) in the same PDE system the user prescribed
83 : : // BCs on side sets that share nodes or (2) because more than a single PDE
84 : : // system assigns BCs to a given node (on different variables), the NodeBC
85 : : // vector must be correctly stored. "Correctly" here means that the size of
86 : : // the NodeBC vectors must all be the same and equal to the sum of all scalar
87 : : // components integrated by all PDE systems integrated. Example: single-phase
88 : : // compressible flow (density, momentum, energy = 5) + transported scalars of
89 : : // 10 variables -> NodeBC vector length = 15. Note that in case (1) above a
90 : : // new node encountered must "overwrite" the already existing space for the
91 : : // NodeBC vector. "Overwrite" here means that it should keep the existing BCs
92 : : // and add the new ones yielding the union the two prescription for BCs but in
93 : : // the same space that already exist in the NodeBC vector. In case (2),
94 : : // however, the NodeBC pairs must go to the location in the vector assigned to
95 : : // the given PDE system, i.e., using the above example BCs for the 10 (or
96 : : // less) scalars should go in the positions starting at 5, leaving the first 5
97 : : // false, indicating no BCs for the flow variables.
98 : : //
99 : : // Note that the logic described above is only partially implemented at this
100 : : // point. What works is the correct insertion of multiple BCs for nodes shared
101 : : // among multiple side sets, e.g., corners, originating from the same PDE
102 : : // system. What is not yet implemented is the case when there are no BCs set
103 : : // for flow variables but there are BCs for transport, the else branch below
104 : : // will incorrectly NOT skip the space for the flow variables. In other words,
105 : : // this only works for a single PDE system and a sytem of systems. This
106 : : // machinery is only tested with a single system of PDEs at this point.
107 : : //
108 : : // When a particular node belongs to two or more side sets with different BCs,
109 : : // there is an ambiguity as to which of the multiple BCs should be applied to
110 : : // the node. This issue is described in case (1) above. In the current
111 : : // implementation, every side set applies the BC to the common node in
112 : : // question, successively overwriting the BC applied by the previous side set.
113 : : // Effectively, the BC corresponding to the last side set ID is applied to the
114 : : // common node. Since bnode is an ordered map, the side set with a larger
115 : : // id wins if a node belongs to multiple side sets.
116 : :
117 : : // Lambda to convert global to local node ids of a list of nodes
118 : 9055643 : auto local = [ &lid ]( const std::vector< std::size_t >& gnodes ){
119 [ + - ]: 73684 : std::vector< std::size_t > lnodes( gnodes.size() );
120 [ + + ]: 9055643 : for (std::size_t i=0; i<gnodes.size(); ++i)
121 [ + - ]: 8981959 : lnodes[i] = tk::cref_find( lid, gnodes[i] );
122 : 73684 : return lnodes;
123 : 45114 : };
124 : :
125 : : // Query Dirichlet BCs for all PDEs integrated and assign to nodes
126 [ + + ]: 118798 : for (const auto& s : bnode) { // for all side sets passed in
127 : 73684 : std::size_t c = 0;
128 [ + - ]: 147368 : auto l = local(s.second); // generate local node ids on side set
129 : : // query Dirichlet BCs at nodes of this side set
130 : : auto eqbc =
131 [ + - ][ + - ]: 147368 : g_cgpde[meshid].dirbc( t, dt, tp, dtp, {s.first,l}, coord, increment );
132 [ + + ]: 856881 : for (const auto& n : eqbc) {
133 : 783197 : auto id = n.first; // BC node ID
134 : 783197 : const auto& bcs = n.second; // BCs
135 [ + - ]: 783197 : auto& nodebc = dirbc[ id ]; // BCs to be set for node
136 [ + + ]: 783197 : if (nodebc.size() > c) { // node already has BCs from this PDE
137 [ - + ][ - - ]: 153836 : Assert( nodebc.size() == c+bcs.size(), "Size mismatch" );
[ - - ][ - - ]
138 [ + + ]: 923016 : for (std::size_t i=0; i<bcs.size(); i++) {
139 [ + - ]: 769180 : if (bcs[i].first) nodebc[c+i] = bcs[i];
140 : : }
141 : : } else { // node does not yet have BCs from this PDE
142 : : // This branch needs to be completed for system of systems of PDEs.
143 : : // See note above.
144 [ + - ]: 629361 : nodebc.insert( end(nodebc), begin(bcs), end(bcs) );
145 : : }
146 : : }
147 [ + + ]: 73684 : if (!eqbc.empty()) c += eqbc.cbegin()->second.size();
148 : : }
149 : :
150 : : // Verify the size of each NodeBC vectors. They must have the same lengths and
151 : : // equal to the total number of scalar components for all systems of PDEs
152 : : // integrated.
153 [ + - ][ - + ]: 674475 : Assert( std::all_of( begin(dirbc), end(dirbc),
[ - - ][ - - ]
[ - - ]
154 : : [ ncomp ]( const auto& n ){ return n.second.size() == ncomp; } ),
155 : : "Size of NodeBC vector incorrect" );
156 : :
157 : 90228 : return dirbc;
158 : : }
159 : :
160 : : } // inciter::
|