Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "SideSetsBetweenSubdomainsGenerator.h"
11 : #include "InputParameters.h"
12 : #include "MooseTypes.h"
13 : #include "MooseMeshUtils.h"
14 : #include "CastUniquePointer.h"
15 :
16 : #include "libmesh/remote_elem.h"
17 :
18 : registerMooseObject("MooseApp", SideSetsBetweenSubdomainsGenerator);
19 :
20 : InputParameters
21 21099 : SideSetsBetweenSubdomainsGenerator::validParams()
22 : {
23 21099 : InputParameters params = SideSetsGeneratorBase::validParams();
24 :
25 21099 : params.renameParam("included_subdomains",
26 : "primary_block",
27 : "The primary set of blocks for which to draw a sideset between");
28 21099 : params.makeParamRequired<std::vector<SubdomainName>>("primary_block");
29 21099 : params.renameParam("included_neighbors",
30 : "paired_block",
31 : "The paired set of blocks for which to draw a sideset between");
32 21099 : params.makeParamRequired<std::vector<SubdomainName>>("paired_block");
33 21099 : params.addClassDescription("MeshGenerator that creates a sideset composed of the nodes located "
34 : "between two or more subdomains.");
35 :
36 : // TODO: Implement each of these in the generate() routine using utilities in SidesetGeneratorBase
37 21099 : params.suppressParameter<bool>("fixed_normal");
38 21099 : params.suppressParameter<bool>("include_only_external_sides");
39 :
40 21099 : return params;
41 0 : }
42 :
43 3413 : SideSetsBetweenSubdomainsGenerator::SideSetsBetweenSubdomainsGenerator(
44 3413 : const InputParameters & parameters)
45 3413 : : SideSetsGeneratorBase(parameters)
46 : {
47 3413 : }
48 :
49 : std::unique_ptr<MeshBase>
50 3270 : SideSetsBetweenSubdomainsGenerator::generate()
51 : {
52 3270 : std::unique_ptr<MeshBase> mesh = std::move(_input);
53 :
54 : // construct the FE object so we can compute normals of faces
55 3270 : setup(*mesh);
56 :
57 : std::vector<boundary_id_type> boundary_ids =
58 3262 : MooseMeshUtils::getBoundaryIDs(*mesh, _boundary_names, true);
59 :
60 : // Get a reference to our BoundaryInfo object for later use
61 3262 : BoundaryInfo & boundary_info = mesh->get_boundary_info();
62 :
63 : // Prepare to query about sides adjacent to remote elements if we're
64 : // on a distributed mesh
65 3262 : const processor_id_type my_n_proc = mesh->n_processors();
66 3262 : const processor_id_type my_proc_id = mesh->processor_id();
67 : typedef std::vector<std::pair<dof_id_type, unsigned int>> vec_type;
68 3262 : std::vector<vec_type> queries(my_n_proc);
69 :
70 : // Request to compute normal vectors
71 3262 : const std::vector<Point> & face_normals = _fe_face->get_normals();
72 :
73 577446 : for (const auto & elem : mesh->active_element_ptr_range())
74 : {
75 : // We only need to loop over elements in the primary subdomain
76 287092 : if (_check_subdomains && !elementSubdomainIdInList(elem, _included_subdomain_ids))
77 153944 : continue;
78 :
79 794486 : for (const auto & side : make_range(elem->n_sides()))
80 : {
81 661338 : const Elem * neighbor = elem->neighbor_ptr(side);
82 :
83 : // On a replicated mesh, we add all subdomain sides ourselves.
84 : // On a distributed mesh, we may have missed sides which
85 : // neighbor remote elements. We should query any such cases.
86 661338 : if (neighbor == remote_elem)
87 : {
88 1199 : queries[elem->processor_id()].push_back(std::make_pair(elem->id(), side));
89 : }
90 660139 : else if (neighbor != NULL)
91 : {
92 564813 : _fe_face->reinit(elem, side);
93 : // We'll just use the normal of the first qp
94 564813 : const Point & face_normal = face_normals[0];
95 : // Add the boundaries, if appropriate
96 564813 : if (elemSideSatisfiesRequirements(elem, side, *mesh, _normal, face_normal))
97 : {
98 : // Add the boundaries
99 34974 : if (_replace)
100 0 : boundary_info.remove_side(elem, side);
101 69948 : for (const auto & boundary_id : boundary_ids)
102 34974 : boundary_info.add_side(elem, side, boundary_id);
103 : }
104 : }
105 : }
106 3262 : }
107 :
108 3262 : if (!mesh->is_serial())
109 : {
110 299 : const auto queries_tag = mesh->comm().get_unique_tag(),
111 299 : replies_tag = mesh->comm().get_unique_tag();
112 :
113 299 : std::vector<Parallel::Request> side_requests(my_n_proc - 1), reply_requests(my_n_proc - 1);
114 :
115 : // Make all requests
116 840 : for (const auto & p : make_range(my_n_proc))
117 : {
118 541 : if (p == my_proc_id)
119 299 : continue;
120 :
121 242 : Parallel::Request & request = side_requests[p - (p > my_proc_id)];
122 :
123 242 : mesh->comm().send(p, queries[p], request, queries_tag);
124 : }
125 :
126 : // Reply to all requests
127 299 : std::vector<vec_type> responses(my_n_proc - 1);
128 :
129 541 : for (const auto & p : make_range(uint(1), my_n_proc))
130 : {
131 242 : vec_type query;
132 :
133 242 : Parallel::Status status(mesh->comm().probe(Parallel::any_source, queries_tag));
134 242 : const processor_id_type source_pid = cast_int<processor_id_type>(status.source());
135 :
136 242 : mesh->comm().receive(source_pid, query, queries_tag);
137 :
138 242 : Parallel::Request & request = reply_requests[p - 1];
139 :
140 1441 : for (const auto & q : query)
141 : {
142 1199 : const Elem * elem = mesh->elem_ptr(q.first);
143 1199 : const unsigned int side = q.second;
144 1199 : const Elem * neighbor = elem->neighbor_ptr(side);
145 :
146 1199 : if (neighbor != NULL)
147 : {
148 1199 : _fe_face->reinit(elem, side);
149 : // We'll just use the normal of the first qp
150 1199 : const Point & face_normal = _fe_face->get_normals()[0];
151 : // Add the boundaries, if appropriate
152 1199 : if (elemSideSatisfiesRequirements(elem, side, *mesh, _normal, face_normal))
153 94 : responses[p - 1].push_back(std::make_pair(elem->id(), side));
154 : }
155 : }
156 :
157 242 : mesh->comm().send(source_pid, responses[p - 1], request, replies_tag);
158 242 : }
159 :
160 : // Process all incoming replies
161 541 : for (processor_id_type p = 1; p != my_n_proc; ++p)
162 : {
163 242 : Parallel::Status status(this->comm().probe(Parallel::any_source, replies_tag));
164 242 : const processor_id_type source_pid = cast_int<processor_id_type>(status.source());
165 :
166 242 : vec_type response;
167 :
168 242 : this->comm().receive(source_pid, response, replies_tag);
169 :
170 336 : for (const auto & r : response)
171 : {
172 94 : const Elem * elem = mesh->elem_ptr(r.first);
173 94 : const unsigned int side = r.second;
174 :
175 94 : if (_replace)
176 0 : boundary_info.remove_side(elem, side);
177 188 : for (const auto & boundary_id : boundary_ids)
178 94 : boundary_info.add_side(elem, side, boundary_id);
179 : }
180 242 : }
181 :
182 299 : Parallel::wait(side_requests);
183 299 : Parallel::wait(reply_requests);
184 299 : }
185 :
186 6524 : for (const auto & i : make_range(boundary_ids.size()))
187 3262 : boundary_info.sideset_name(boundary_ids[i]) = _boundary_names[i];
188 :
189 3262 : mesh->set_isnt_prepared();
190 6524 : return dynamic_pointer_cast<MeshBase>(mesh);
191 3262 : }
|