126int main(
int argc,
char* argv[])
129 "Creates boundary conditions for mesh nodes along polylines."
130 "The documentation is available at "
131 "https://docs.opengeosys.org/docs/tools/model-preparation/"
132 "create-boundary-conditions-along-a-polyline.\n\n"
133 "OpenGeoSys-6 software, version " +
136 "Copyright (c) 2012-2025, OpenGeoSys Community "
137 "(http://www.opengeosys.org)",
139 TCLAP::SwitchArg gml_arg(
"",
"gml",
"Write found nodes to gml file.");
142 TCLAP::ValueArg<std::string> output_base_fname(
143 "o",
"output-base-file-name",
144 "Output. The base name of the file the output (geometry (gli) and "
146 "condition (bc)) will be written to",
147 true,
"",
"BASE_FILENAME_OUTPUT");
148 cmd.add(output_base_fname);
150 std::vector<std::string> allowed_types_vector{
"LIQUID_FLOW",
152 TCLAP::ValuesConstraint<std::string> allowed_types(allowed_types_vector);
153 TCLAP::ValueArg<std::string> bc_type(
155 "the process type the boundary condition will be written for currently "
156 "LIQUID_FLOW (primary variable PRESSURE1) and GROUNDWATER_FLOW "
157 "(primary variable HEAD, default) are supported, ",
158 true,
"", &allowed_types);
161 TCLAP::ValueArg<double> search_length_arg(
162 "s",
"search-length",
"The size of the search length ",
false,
163 std::numeric_limits<double>::epsilon(),
"SEARCH_LENGTH");
164 cmd.add(search_length_arg);
166 TCLAP::ValueArg<std::string> geometry_fname(
167 "i",
"input-geometry",
168 "Input (.gml | .gli). The name of the input file containing the "
170 true,
"",
"INPUT_FILE");
171 cmd.add(geometry_fname);
173 TCLAP::ValueArg<std::string> mesh_arg(
175 "Input (.vtu). The name of the input file containing the mesh",
true,
179 TCLAP::ValueArg<std::string> gmsh_path_arg(
180 "g",
"gmsh-path",
"Input (.msh). The path to the gmsh binary",
false,
182 cmd.add(gmsh_path_arg);
184 cmd.parse(argc, argv);
189 INFO(
"Reading mesh '{:s}' ... ", mesh_arg.getValue());
190 std::unique_ptr<MeshLib::Mesh> subsurface_mesh(
193 INFO(
"Extracting top surface of mesh '{:s}' ... ", mesh_arg.getValue());
194 Eigen::Vector3d
const dir({0, 0, -1});
195 double const angle(90);
196 std::unique_ptr<MeshLib::Mesh> surface_mesh(
200 subsurface_mesh.reset(
nullptr);
205 gmsh_path_arg.getValue());
211 std::vector<GeoLib::Polyline*>
const* plys(
215 ERR(
"Could not get vector of polylines out of geometry '{:s}'.",
220 auto search_length_strategy =
221 std::make_unique<MeshGeoToolsLib::SearchLength>();
222 if (search_length_arg.isSet())
224 search_length_strategy.reset(
230 *surface_mesh, std::move(search_length_strategy),
232 for (std::size_t k(0); k < plys->size(); k++)
239 std::string polyline_name(
"Polyline-" + std::to_string(k));
246 if (geo_names.empty())
248 ERR(
"Did not find mesh nodes along polylines.");
252 std::string merge_name(
"AllMeshNodesAlongPolylines");
255 merge_name = geo_names[0];
259 auto const& merged_pnts(pnt_vec->
getVector());
261 std::vector<GeoLib::Point> pnts_with_id;
262 const std::size_t n_merged_pnts(merged_pnts.size());
263 for (std::size_t k(0); k < n_merged_pnts; ++k)
265 pnts_with_id.emplace_back(*(merged_pnts[k]), k);
268 std::sort(pnts_with_id.begin(), pnts_with_id.end(),
270 { return p0 < p1; });
272 double const eps(std::numeric_limits<double>::epsilon());
273 std::vector<GeoLib::Point*> surface_pnts;
277 surface_pnts.push_back(
280 std::string element_name;
282 name_id_map[element_name] = 0;
284 for (std::size_t k(1); k < n_merged_pnts; ++k)
288 if (std::abs(p0[0] - p1[0]) > eps || std::abs(p0[1] - p1[1]) > eps)
290 surface_pnts.push_back(
292 std::string element_name;
294 name_id_map[element_name] = surface_pnts.size() - 1;
299 "-MeshNodesAlongPolylines");
300 geometry_sets.
addPointVec(std::move(surface_pnts), surface_name,
301 std::move(name_id_map), 1e-6);
304 std::string
const base_fname(
307 bc_type.getValue(), gml_arg.getValue());