LCOV - code coverage report
Current view: top level - src/transfers - MultiAppGeneralFieldTransfer.C (source / functions) Hit Total Coverage
Test: idaholab/moose framework: 6f668f Lines: 810 918 88.2 %
Date: 2025-09-22 20:01:15 Functions: 40 41 97.6 %
Legend: Lines: hit not hit

          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 "MultiAppGeneralFieldTransfer.h"
      11             : 
      12             : // MOOSE includes
      13             : #include "DisplacedProblem.h"
      14             : #include "FEProblem.h"
      15             : #include "MooseMesh.h"
      16             : #include "MooseTypes.h"
      17             : #include "MooseVariableFE.h"
      18             : #include "MeshDivision.h"
      19             : #include "Positions.h"
      20             : #include "Factory.h"
      21             : #include "MooseAppCoordTransform.h"
      22             : 
      23             : // libmesh includes
      24             : #include "libmesh/point_locator_base.h"
      25             : #include "libmesh/enum_point_locator_type.h"
      26             : 
      27             : // TIMPI includes
      28             : #include "timpi/communicator.h"
      29             : #include "timpi/parallel_sync.h"
      30             : 
      31             : using namespace libMesh;
      32             : 
      33             : namespace GeneralFieldTransfer
      34             : {
      35             : Number BetterOutOfMeshValue = std::numeric_limits<Real>::infinity();
      36             : }
      37             : 
      38             : InputParameters
      39       73462 : MultiAppGeneralFieldTransfer::validParams()
      40             : {
      41       73462 :   InputParameters params = MultiAppConservativeTransfer::validParams();
      42             :   // Expansion default to make a node in the target mesh overlapping with a node in the origin
      43             :   // mesh always register as being inside the origin application bounding box. The contains_point
      44             :   // bounding box checks uses exact comparisons
      45      367310 :   params.addRangeCheckedParam<Real>("bbox_factor",
      46      146924 :                                     1.00000001,
      47             :                                     "bbox_factor>0",
      48             :                                     "Factor to inflate or deflate the source app bounding boxes");
      49      440772 :   params.addRangeCheckedParam<std::vector<Real>>(
      50             :       "fixed_bounding_box_size",
      51             :       "fixed_bounding_box_size >= 0",
      52             :       "Override source app bounding box size(s) for searches. App bounding boxes will still be  "
      53             :       "centered on the same coordinates. Only non-zero components passed will override.");
      54      220386 :   params.addParam<Real>(
      55             :       "extrapolation_constant",
      56      146924 :       0,
      57             :       "Constant to use when no source app can provide a valid value for a target location.");
      58             : 
      59             :   // Block restrictions
      60      293848 :   params.addParam<std::vector<SubdomainName>>(
      61             :       "from_blocks",
      62             :       "Subdomain restriction to transfer from (defaults to all the origin app domain)");
      63      293848 :   params.addParam<std::vector<SubdomainName>>(
      64             :       "to_blocks", "Subdomain restriction to transfer to, (defaults to all the target app domain)");
      65             : 
      66             :   // Boundary restrictions
      67      293848 :   params.addParam<std::vector<BoundaryName>>(
      68             :       "from_boundaries",
      69             :       "The boundary we are transferring from (if not specified, whole domain is used).");
      70      293848 :   params.addParam<std::vector<BoundaryName>>(
      71             :       "to_boundaries",
      72             :       "The boundary we are transferring to (if not specified, whole domain is used).");
      73      293848 :   MooseEnum nodes_or_sides("nodes sides", "sides");
      74      293848 :   params.addParam<MooseEnum>("elemental_boundary_restriction",
      75             :                              nodes_or_sides,
      76             :                              "Whether elemental variable boundary restriction is considered by "
      77             :                              "element side or element nodes");
      78             : 
      79             :   // Mesh division restriction
      80      293848 :   params.addParam<MeshDivisionName>("from_mesh_division",
      81             :                                     "Mesh division object on the origin application");
      82      293848 :   params.addParam<MeshDivisionName>("to_mesh_division",
      83             :                                     "Mesh division object on the target application");
      84             :   MooseEnum mesh_division_uses("spatial_restriction matching_division matching_subapp_index none",
      85      293848 :                                "none");
      86      293848 :   params.addParam<MooseEnum>("from_mesh_division_usage",
      87             :                              mesh_division_uses,
      88             :                              "How to use the source mesh division in the transfer. See object "
      89             :                              "documentation for description of each option");
      90      293848 :   params.addParam<MooseEnum>("to_mesh_division_usage",
      91             :                              mesh_division_uses,
      92             :                              "How to use the target mesh division in the transfer. See object "
      93             :                              "documentation for description of each option");
      94             : 
      95             :   // Array and vector variables
      96      220386 :   params.addParam<std::vector<unsigned int>>("source_variable_components",
      97      146924 :                                              std::vector<unsigned int>(),
      98             :                                              "The source array or vector variable component(s).");
      99      220386 :   params.addParam<std::vector<unsigned int>>("target_variable_components",
     100      146924 :                                              std::vector<unsigned int>(),
     101             :                                              "The target array or vector variable component(s).");
     102             : 
     103             :   // Search options
     104      220386 :   params.addParam<bool>(
     105             :       "greedy_search",
     106      146924 :       false,
     107             :       "Whether or not to send a point to all the domains. If true, all the processors will be "
     108             :       "checked for a given point."
     109             :       "The code will be slow if this flag is on but it will give a better solution.");
     110      220386 :   params.addParam<bool>(
     111             :       "error_on_miss",
     112      146924 :       true,
     113             :       "Whether or not to error in the case that a target point is not found in the source domain.");
     114      220386 :   params.addParam<bool>("use_bounding_boxes",
     115      146924 :                         true,
     116             :                         "When set to false, bounding boxes will not be used to restrict the source "
     117             :                         "of the transfer. Either source applications must be set using the "
     118             :                         "from_mesh_division parameter, or a greedy search must be used.");
     119      220386 :   params.addParam<bool>(
     120             :       "use_nearest_app",
     121      146924 :       false,
     122             :       "When True, transfers from a child application will work by finding the nearest (using "
     123             :       "the `position` + mesh centroid) sub-app and query that app for the value to transfer.");
     124      293848 :   params.addParam<PositionsName>(
     125             :       "use_nearest_position",
     126             :       "Name of the the Positions object (in main app) such that transfers to/from a child "
     127             :       "application will work by finding the nearest position to a target and query only the "
     128             :       "app / points closer to this position than to any other position for the value to transfer.");
     129      220386 :   params.addParam<bool>(
     130             :       "from_app_must_contain_point",
     131      146924 :       false,
     132             :       "Wether on not the origin mesh must contain the point to evaluate data at. If false, this "
     133             :       "allows for interpolation between origin app meshes. Origin app bounding boxes are still "
     134             :       "considered so you may want to increase them with 'fixed_bounding_box_size'");
     135      220386 :   params.addParam<bool>("search_value_conflicts",
     136      146924 :                         false,
     137             :                         "Whether to look for potential conflicts between two valid and different "
     138             :                         "source values for any target point");
     139      220386 :   params.addParam<unsigned int>(
     140             :       "value_conflicts_output",
     141      146924 :       10,
     142             :       "Maximum number of conflicts to output if value-conflicts, from equidistant sources to a "
     143             :       "given transfer target location, search is turned on");
     144             : 
     145      293848 :   params.addParamNamesToGroup(
     146             :       "to_blocks from_blocks to_boundaries from_boundaries elemental_boundary_restriction "
     147             :       "from_mesh_division from_mesh_division_usage to_mesh_division to_mesh_division_usage",
     148             :       "Transfer spatial restriction");
     149      293848 :   params.addParamNamesToGroup(
     150             :       "greedy_search use_bounding_boxes use_nearest_app use_nearest_position "
     151             :       "search_value_conflicts",
     152             :       "Search algorithm");
     153      293848 :   params.addParamNamesToGroup("error_on_miss from_app_must_contain_point extrapolation_constant",
     154             :                               "Extrapolation behavior");
     155      220386 :   params.addParamNamesToGroup("bbox_factor fixed_bounding_box_size", "Source app bounding box");
     156      146924 :   return params;
     157       73462 : }
     158             : 
     159        7243 : MultiAppGeneralFieldTransfer::MultiAppGeneralFieldTransfer(const InputParameters & parameters)
     160             :   : MultiAppConservativeTransfer(parameters),
     161        7243 :     _from_var_components(getParam<std::vector<unsigned int>>("source_variable_components")),
     162       14486 :     _to_var_components(getParam<std::vector<unsigned int>>("target_variable_components")),
     163       14486 :     _use_bounding_boxes(getParam<bool>("use_bounding_boxes")),
     164       14486 :     _use_nearest_app(getParam<bool>("use_nearest_app")),
     165        7243 :     _nearest_positions_obj(
     166       14486 :         isParamValid("use_nearest_position")
     167        7735 :             ? &_fe_problem.getPositionsObject(getParam<PositionsName>("use_nearest_position"))
     168             :             : nullptr),
     169       14486 :     _source_app_must_contain_point(getParam<bool>("from_app_must_contain_point")),
     170       14486 :     _from_mesh_division_behavior(getParam<MooseEnum>("from_mesh_division_usage")),
     171       14486 :     _to_mesh_division_behavior(getParam<MooseEnum>("to_mesh_division_usage")),
     172        7243 :     _elemental_boundary_restriction_on_sides(
     173       14486 :         getParam<MooseEnum>("elemental_boundary_restriction") == "sides"),
     174       14486 :     _greedy_search(getParam<bool>("greedy_search")),
     175       14486 :     _search_value_conflicts(getParam<bool>("search_value_conflicts")),
     176        7243 :     _already_output_search_value_conflicts(false),
     177       14486 :     _search_value_conflicts_max_log(getParam<unsigned int>("value_conflicts_output")),
     178       14486 :     _error_on_miss(getParam<bool>("error_on_miss")),
     179       14486 :     _default_extrapolation_value(getParam<Real>("extrapolation_constant")),
     180       14486 :     _bbox_factor(getParam<Real>("bbox_factor")),
     181       36215 :     _fixed_bbox_size(isParamValid("fixed_bounding_box_size")
     182        7243 :                          ? getParam<std::vector<Real>>("fixed_bounding_box_size")
     183       36215 :                          : std::vector<Real>(3, 0))
     184             : {
     185        7243 :   _var_size = _to_var_names.size();
     186        9913 :   if (_to_var_names.size() != _from_var_names.size() && !parameters.isPrivate("source_variable"))
     187           8 :     paramError("variable", "The number of variables to transfer to and from should be equal");
     188             : 
     189             :   // Check the parameters of the components of the array / vector variable
     190        7239 :   if (_from_var_names.size() != _from_var_components.size() && _from_var_components.size() > 0)
     191           0 :     paramError("source_variable_components",
     192             :                "This parameter must be equal to the number of source variables");
     193        7239 :   if (_to_var_names.size() != _to_var_components.size() && _to_var_components.size() > 0)
     194           0 :     paramError("target_variable_components",
     195             :                "This parameter must be equal to the number of target variables");
     196             : 
     197             :   // Make simple 'use_nearest_app' parameter rely on Positions
     198        7239 :   if (_use_nearest_app)
     199             :   {
     200         158 :     if (_nearest_positions_obj)
     201           8 :       paramError("use_nearest_app", "Cannot use nearest-app and nearest-position together");
     202         154 :     if (!hasFromMultiApp())
     203           8 :       paramError("use_nearest_app",
     204             :                  "Should have a 'from_multiapp' when using the nearest-app informed search");
     205         300 :     auto pos_params = _app.getFactory().getValidParams("MultiAppPositions");
     206         600 :     pos_params.set<std::vector<MultiAppName>>("multiapps") = {getFromMultiApp()->name()};
     207         450 :     _fe_problem.addReporter("MultiAppPositions", "_created_for_" + name(), pos_params);
     208         150 :     _nearest_positions_obj = &_fe_problem.getPositionsObject("_created_for_" + name());
     209         150 :   }
     210             : 
     211             :   // Dont let users get wrecked by bounding boxes if it looks like they are trying to extrapolate
     212       14372 :   if (!_source_app_must_contain_point && _use_bounding_boxes &&
     213       27878 :       (_nearest_positions_obj || isParamSetByUser("from_app_must_contain_point")))
     214        1967 :     if (!isParamSetByUser("bbox_factor") && !isParamSetByUser("fixed_bounding_box_size"))
     215           4 :       mooseWarning(
     216             :           "Extrapolation (nearest-source options, outside-app source) parameters have "
     217             :           "been passed, but no subapp bounding box expansion parameters have been passed.");
     218             : 
     219        7309 :   if (!_use_bounding_boxes &&
     220        7637 :       (isParamValid("fixed_bounding_box_size") || isParamSetByUser("bbox_factor")))
     221           0 :     paramError("use_bounding_boxes",
     222             :                "Cannot pass additional bounding box parameters (sizes, expansion, etc) if we are "
     223             :                "not using bounding boxes");
     224        7527 : }
     225             : 
     226             : void
     227        7007 : MultiAppGeneralFieldTransfer::initialSetup()
     228             : {
     229        7007 :   MultiAppConservativeTransfer::initialSetup();
     230             : 
     231             :   // Use IDs for block and boundary restriction
     232             :   // Loop over all source problems
     233       16544 :   for (const auto i_from : index_range(_from_problems))
     234             :   {
     235        9573 :     const auto & from_moose_mesh = _from_problems[i_from]->mesh(_displaced_source_mesh);
     236       28719 :     if (isParamValid("from_blocks"))
     237             :     {
     238        2932 :       const auto & block_names = getParam<std::vector<SubdomainName>>("from_blocks");
     239             : 
     240        3186 :       for (const auto & b : block_names)
     241        1724 :         if (!MooseMeshUtils::hasSubdomainName(from_moose_mesh.getMesh(), b))
     242           8 :           paramError("from_blocks", "The block '", b, "' was not found in the mesh");
     243             : 
     244        1462 :       if (!block_names.empty())
     245             :       {
     246        1462 :         const auto ids = from_moose_mesh.getSubdomainIDs(block_names);
     247        1462 :         _from_blocks.insert(ids.begin(), ids.end());
     248        1462 :       }
     249             :     }
     250             : 
     251       28707 :     if (isParamValid("from_boundaries"))
     252             :     {
     253         390 :       const auto & boundary_names = getParam<std::vector<BoundaryName>>("from_boundaries");
     254         558 :       for (const auto & bn : boundary_names)
     255         367 :         if (!MooseMeshUtils::hasBoundaryName(from_moose_mesh.getMesh(), bn))
     256           8 :           paramError("from_boundaries", "The boundary '", bn, "' was not found in the mesh");
     257             : 
     258         191 :       if (!boundary_names.empty())
     259             :       {
     260         191 :         const auto boundary_ids = from_moose_mesh.getBoundaryIDs(boundary_names);
     261         191 :         _from_boundaries.insert(boundary_ids.begin(), boundary_ids.end());
     262         191 :       }
     263             :     }
     264             : 
     265       28695 :     if (isParamValid("from_mesh_division"))
     266             :     {
     267         716 :       const auto & mesh_div_name = getParam<MeshDivisionName>("from_mesh_division");
     268         358 :       _from_mesh_divisions.push_back(&_from_problems[i_from]->getMeshDivision(mesh_div_name));
     269             :       // Check that the behavior set makes sense
     270         358 :       if (_from_mesh_division_behavior == MeshDivisionTransferUse::RESTRICTION)
     271             :       {
     272         158 :         if (_from_mesh_divisions[i_from]->coversEntireMesh())
     273           0 :           mooseInfo("'from_mesh_division_usage' is set to use a spatial restriction but the "
     274           0 :                     "'from_mesh_division' for source app of global index " +
     275           0 :                     std::to_string(getGlobalSourceAppIndex(i_from)) +
     276             :                     " covers the entire mesh. Do not expect any restriction from a mesh "
     277             :                     "division that covers the entire mesh");
     278             :       }
     279         310 :       else if (_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX &&
     280         530 :                !isParamValid("to_mesh_division"))
     281           0 :         paramError("to_mesh_division_usage",
     282             :                    "Source mesh division cannot match target mesh division if no target mesh "
     283             :                    "division is specified");
     284         200 :       else if (_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX)
     285             :       {
     286          90 :         if (!hasToMultiApp())
     287           8 :           paramError("from_mesh_division_usage",
     288             :                      "Cannot match source mesh division index to target subapp index if there is "
     289             :                      "only one target: the parent app (not a subapp)");
     290         172 :         else if (getToMultiApp()->numGlobalApps() !=
     291          86 :                  _from_mesh_divisions[i_from]->getNumDivisions())
     292           4 :           mooseWarning("Attempting to match target subapp index with the number of source mesh "
     293           4 :                        "divisions, which is " +
     294           8 :                        std::to_string(_from_mesh_divisions[i_from]->getNumDivisions()) +
     295           8 :                        " while there are " + std::to_string(getToMultiApp()->numGlobalApps()) +
     296             :                        " target subapps");
     297          82 :         if (_to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX)
     298             :           // We do not support it because it would require sending the point + target app index +
     299             :           // target app division index, and we only send the Point + one number
     300           8 :           paramError("from_mesh_division_usage",
     301             :                      "We do not support using target subapp index for source division behavior and "
     302             :                      "matching the division index for the target mesh division behavior.");
     303             :       }
     304         110 :       else if (_from_mesh_division_behavior == "none")
     305           0 :         paramError("from_mesh_division_usage", "User must specify a 'from_mesh_division_usage'");
     306             :     }
     307        9207 :     else if (_from_mesh_division_behavior != "none")
     308           8 :       paramError("from_mesh_division",
     309             :                  "'from_mesh_division' must be specified if the usage method is specified");
     310             :   }
     311             : 
     312             :   // Loop over all target problems
     313       16486 :   for (const auto i_to : index_range(_to_problems))
     314             :   {
     315        9539 :     const auto & to_moose_mesh = _to_problems[i_to]->mesh(_displaced_target_mesh);
     316       28617 :     if (isParamValid("to_blocks"))
     317             :     {
     318        2520 :       const auto & block_names = getParam<std::vector<SubdomainName>>("to_blocks");
     319        2774 :       for (const auto & b : block_names)
     320        1518 :         if (!MooseMeshUtils::hasSubdomainName(to_moose_mesh.getMesh(), b))
     321           8 :           paramError("to_blocks", "The block '", b, "' was not found in the mesh");
     322             : 
     323        1256 :       if (!block_names.empty())
     324             :       {
     325        1256 :         const auto ids = to_moose_mesh.getSubdomainIDs(block_names);
     326        1256 :         _to_blocks.insert(ids.begin(), ids.end());
     327        1256 :       }
     328             :     }
     329             : 
     330       28605 :     if (isParamValid("to_boundaries"))
     331             :     {
     332        2844 :       const auto & boundary_names = getParam<std::vector<BoundaryName>>("to_boundaries");
     333        4766 :       for (const auto & bn : boundary_names)
     334        3348 :         if (!MooseMeshUtils::hasBoundaryName(to_moose_mesh.getMesh(), bn))
     335           8 :           paramError("to_boundaries", "The boundary '", bn, "' was not found in the mesh");
     336             : 
     337        1418 :       if (!boundary_names.empty())
     338             :       {
     339        1418 :         const auto boundary_ids = to_moose_mesh.getBoundaryIDs(boundary_names);
     340        1418 :         _to_boundaries.insert(boundary_ids.begin(), boundary_ids.end());
     341        1418 :       }
     342             :     }
     343             : 
     344       28593 :     if (isParamValid("to_mesh_division"))
     345             :     {
     346         708 :       const auto & mesh_div_name = getParam<MeshDivisionName>("to_mesh_division");
     347         354 :       _to_mesh_divisions.push_back(&_to_problems[i_to]->getMeshDivision(mesh_div_name));
     348             :       // Check that the behavior set makes sense
     349         354 :       if (_to_mesh_division_behavior == MeshDivisionTransferUse::RESTRICTION)
     350             :       {
     351         158 :         if (_to_mesh_divisions[i_to]->coversEntireMesh())
     352           0 :           mooseInfo("'to_mesh_division_usage' is set to use a spatial restriction but the "
     353           0 :                     "'to_mesh_division' for target application of global index " +
     354           0 :                     std::to_string(getGlobalSourceAppIndex(i_to)) +
     355             :                     " covers the entire mesh. Do not expect any restriction from a mesh "
     356             :                     "division that covers the entire mesh");
     357             :       }
     358         196 :       else if (_to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX)
     359             :       {
     360         330 :         if (!isParamValid("from_mesh_division"))
     361           0 :           paramError("to_mesh_division_usage",
     362             :                      "Target mesh division cannot match source mesh division if no source mesh "
     363             :                      "division is specified");
     364         220 :         else if ((*_from_mesh_divisions.begin())->getNumDivisions() !=
     365         110 :                  _to_mesh_divisions[i_to]->getNumDivisions())
     366           4 :           mooseWarning("Source and target mesh divisions do not have the same number of bins. If "
     367             :                        "this is what you expect, please reach out to a MOOSE or app developer to "
     368             :                        "ensure appropriate use");
     369             :       }
     370          86 :       else if (_to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX)
     371             :       {
     372          86 :         if (!hasFromMultiApp())
     373           8 :           paramError(
     374             :               "to_mesh_division_usage",
     375             :               "Cannot match target mesh division index to source subapp index if there is only one "
     376             :               "source: the parent app (not a subapp)");
     377          82 :         else if (getFromMultiApp()->numGlobalApps() != _to_mesh_divisions[i_to]->getNumDivisions())
     378           4 :           mooseWarning("Attempting to match source subapp index with the number of target mesh "
     379           4 :                        "divisions, which is " +
     380           8 :                        std::to_string(_to_mesh_divisions[i_to]->getNumDivisions()) +
     381           8 :                        " while there are " + std::to_string(getFromMultiApp()->numGlobalApps()) +
     382             :                        " source subapps");
     383          78 :         if (_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX)
     384           0 :           paramError(
     385             :               "from_mesh_division_usage",
     386             :               "We do not support using source subapp index for the target division behavior and "
     387             :               "matching the division index for the source mesh division behavior.");
     388             :       }
     389           0 :       else if (_to_mesh_division_behavior == "none")
     390           0 :         paramError("to_mesh_division_usage", "User must specify a 'to_mesh_division_usage'");
     391             :     }
     392        9177 :     else if (_to_mesh_division_behavior != "none")
     393           4 :       paramError("to_mesh_division",
     394           4 :                  "'to_mesh_division' must be specified if usage method '" +
     395           8 :                      Moose::stringify(_to_mesh_division_behavior) + "' is specified");
     396             :   }
     397             : 
     398             :   // Check if components are set correctly if using an array variable
     399       16448 :   for (const auto i_from : index_range(_from_problems))
     400             :   {
     401       17564 :     for (const auto var_index : make_range(_from_var_names.size()))
     402             :     {
     403             :       MooseVariableFieldBase & from_var =
     404       16126 :           _from_problems[i_from]->getVariable(0,
     405        8063 :                                               _from_var_names[var_index],
     406             :                                               Moose::VarKindType::VAR_ANY,
     407             :                                               Moose::VarFieldType::VAR_FIELD_ANY);
     408        8063 :       if (from_var.count() > 1 && _from_var_components.empty())
     409           8 :         paramError("source_variable_components", "Component must be passed for an array variable");
     410        8059 :       if (_from_var_components.size() && from_var.count() < _from_var_components[var_index])
     411           8 :         paramError("source_variable_components",
     412             :                    "Component passed is larger than size of variable");
     413             :     }
     414             :   }
     415       16438 :   for (const auto i_to : index_range(_to_problems))
     416             :   {
     417       19343 :     for (const auto var_index : make_range(_to_var_names.size()))
     418             :     {
     419             :       MooseVariableFieldBase & to_var =
     420       19688 :           _to_problems[i_to]->getVariable(0,
     421        9844 :                                           _to_var_names[var_index],
     422             :                                           Moose::VarKindType::VAR_ANY,
     423             :                                           Moose::VarFieldType::VAR_FIELD_ANY);
     424        9844 :       if (to_var.count() > 1 && _to_var_components.empty())
     425           8 :         paramError("target_variable_components", "Component must be passed for an array variable");
     426        9840 :       if (_to_var_components.size() && to_var.count() < _to_var_components[var_index])
     427           8 :         paramError("target_variable_components",
     428             :                    "Component passed is larger than size of variable");
     429             :     }
     430             :   }
     431             : 
     432             :   // Cache some quantities to avoid having to get them on every transferred point
     433        6931 :   if (_to_problems.size())
     434             :   {
     435        6916 :     _to_variables.resize(_to_var_names.size());
     436       14020 :     for (const auto i_var : index_range(_to_var_names))
     437       14208 :       _to_variables[i_var] = &_to_problems[0]->getVariable(
     438        7104 :           0, _to_var_names[i_var], Moose::VarKindType::VAR_ANY, Moose::VarFieldType::VAR_FIELD_ANY);
     439             :   }
     440        6931 : }
     441             : 
     442             : void
     443       63715 : MultiAppGeneralFieldTransfer::getAppInfo()
     444             : {
     445       63715 :   MultiAppFieldTransfer::getAppInfo();
     446             : 
     447             :   // Create the point locators to locate evaluation points in the origin mesh(es)
     448       63715 :   _from_point_locators.resize(_from_problems.size());
     449      132323 :   for (const auto i_from : index_range(_from_problems))
     450             :   {
     451       68608 :     const auto & from_moose_mesh = _from_problems[i_from]->mesh(_displaced_source_mesh);
     452       68608 :     _from_point_locators[i_from] =
     453      137216 :         PointLocatorBase::build(TREE_LOCAL_ELEMENTS, from_moose_mesh.getMesh());
     454       68608 :     _from_point_locators[i_from]->enable_out_of_mesh_mode();
     455             :   }
     456       63715 : }
     457             : 
     458             : void
     459       56708 : MultiAppGeneralFieldTransfer::execute()
     460             : {
     461       56708 :   TIME_SECTION(
     462             :       "MultiAppGeneralFieldTransfer::execute()_" + name(), 5, "Transfer execution " + name());
     463       56708 :   getAppInfo();
     464             : 
     465             :   // Set up bounding boxes, etc
     466       56708 :   prepareToTransfer();
     467             : 
     468             :   // loop over the vector of variables and make the transfer one by one
     469      113515 :   for (const auto i : make_range(_var_size))
     470       56854 :     transferVariable(i);
     471             : 
     472       56661 :   postExecute();
     473       56661 : }
     474             : 
     475             : void
     476       56708 : MultiAppGeneralFieldTransfer::prepareToTransfer()
     477             : {
     478             :   // Get the bounding boxes for the "from" domains.
     479             :   // Clean up _from_bboxes from the previous transfer execution
     480       56708 :   _from_bboxes.clear();
     481             : 
     482             :   // NOTE: This ignores the app's bounding box inflation and padding
     483       56708 :   _from_bboxes = getRestrictedFromBoundingBoxes();
     484             : 
     485             :   // Expand bounding boxes. Some desired points might be excluded
     486             :   // without an expansion
     487       56708 :   extendBoundingBoxes(_bbox_factor, _from_bboxes);
     488             : 
     489             :   // Figure out how many "from" domains each processor owns.
     490       56708 :   _froms_per_proc.clear();
     491       56708 :   _froms_per_proc = getFromsPerProc();
     492             : 
     493             :   // Get the index for the first source app every processor owns
     494       56708 :   _global_app_start_per_proc = getGlobalStartAppPerProc();
     495             : 
     496             :   // No need to keep searching for conflicts if the mesh has not changed
     497       56708 :   if (_already_output_search_value_conflicts && !_displaced_source_mesh && !_displaced_target_mesh)
     498           0 :     _search_value_conflicts = false;
     499       56708 : }
     500             : 
     501             : void
     502       56661 : MultiAppGeneralFieldTransfer::postExecute()
     503             : {
     504       56661 :   MultiAppConservativeTransfer::postExecute();
     505       56661 :   if (_search_value_conflicts)
     506         180 :     _already_output_search_value_conflicts = true;
     507       56661 : }
     508             : 
     509             : void
     510       56854 : MultiAppGeneralFieldTransfer::transferVariable(unsigned int i)
     511             : {
     512             :   mooseAssert(i < _var_size, "The variable of index " << i << " does not exist");
     513             : 
     514             :   // Find outgoing target points
     515             :   // We need to know what points we need to send which processors
     516             :   // One processor will receive many points from many processors
     517             :   // One point may go to different processors
     518       56854 :   ProcessorToPointVec outgoing_points;
     519       56854 :   extractOutgoingPoints(i, outgoing_points);
     520             : 
     521       56846 :   if (_from_var_names.size())
     522       55523 :     prepareEvaluationOfInterpValues(i);
     523             :   else
     524        1323 :     prepareEvaluationOfInterpValues(-1);
     525             : 
     526             :   // Fill values and app ids for incoming points
     527             :   // We are responsible to compute values for these incoming points
     528             :   auto gather_functor =
     529       85064 :       [this](processor_id_type /*pid*/,
     530             :              const std::vector<std::pair<Point, unsigned int>> & incoming_locations,
     531             :              std::vector<std::pair<Real, Real>> & outgoing_vals)
     532             :   {
     533       85064 :     outgoing_vals.resize(
     534             :         incoming_locations.size(),
     535             :         {GeneralFieldTransfer::BetterOutOfMeshValue, GeneralFieldTransfer::BetterOutOfMeshValue});
     536             :     // Evaluate interpolation values for these incoming points
     537       85064 :     evaluateInterpValues(incoming_locations, outgoing_vals);
     538      141910 :   };
     539             : 
     540      113692 :   DofobjectToInterpValVec dofobject_to_valsvec(_to_problems.size());
     541       56846 :   InterpCaches interp_caches(_to_problems.size(), getMaxToProblemsBBoxDimensions());
     542       56846 :   InterpCaches distance_caches(_to_problems.size(), getMaxToProblemsBBoxDimensions());
     543             : 
     544             :   // Copy data out to incoming_vals_ids
     545       85064 :   auto action_functor = [this, &i, &dofobject_to_valsvec, &interp_caches, &distance_caches](
     546             :                             processor_id_type pid,
     547             :                             const std::vector<std::pair<Point, unsigned int>> & my_outgoing_points,
     548             :                             const std::vector<std::pair<Real, Real>> & incoming_vals)
     549             :   {
     550       85064 :     auto & pointInfoVec = _processor_to_pointInfoVec[pid];
     551             : 
     552             :     // Cache interpolation values for each dof object / points
     553       85064 :     cacheIncomingInterpVals(pid,
     554             :                             i,
     555             :                             pointInfoVec,
     556             :                             my_outgoing_points,
     557             :                             incoming_vals,
     558             :                             dofobject_to_valsvec,
     559             :                             interp_caches,
     560             :                             distance_caches);
     561      141910 :   };
     562             : 
     563             :   // We assume incoming_vals_ids is ordered in the same way as outgoing_points
     564             :   // Hopefully, pull_parallel_vector_data will not mess up this
     565       56846 :   const std::pair<Real, Real> * ex = nullptr;
     566       56846 :   libMesh::Parallel::pull_parallel_vector_data(
     567       56846 :       comm(), outgoing_points, gather_functor, action_functor, ex);
     568             : 
     569             :   // Check for conflicts and overlaps from the maps that were built during the transfer
     570       56846 :   if (_search_value_conflicts)
     571         211 :     outputValueConflicts(i, dofobject_to_valsvec, distance_caches);
     572             : 
     573             :   // Set cached values into solution vector
     574       56815 :   setSolutionVectorValues(i, dofobject_to_valsvec, interp_caches);
     575       56807 : }
     576             : 
     577             : void
     578     5948561 : MultiAppGeneralFieldTransfer::locatePointReceivers(const Point point,
     579             :                                                    std::set<processor_id_type> & processors)
     580             : {
     581             :   // Check which processors have apps that may include or be near this point
     582             :   // A point may be close enough to several problems, hosted on several processes
     583     5948561 :   bool found = false;
     584             : 
     585             :   // Additional process-restriction techniques we could use (TODOs):
     586             :   // - create a heuristic for using nearest-positions
     587             :   // - from_mesh_divisions could be polled for which divisions they possess on each
     588             :   //   process, depending on the behavior chosen. This could limit potential senders.
     589             :   //   This should be done ahead of this function call, for all points at once
     590             : 
     591             :   // Determine the apps which will be receiving points (then sending values) using various
     592             :   // heuristics
     593     5948561 :   if (_use_nearest_app)
     594             :   {
     595             :     // Find the nearest position for the point
     596        9047 :     const bool initial = _fe_problem.getCurrentExecuteOnFlag() == EXEC_INITIAL;
     597             :     // The apps form the nearest positions here, this is the index of the nearest app
     598        9047 :     const auto nearest_index = _nearest_positions_obj->getNearestPositionIndex(point, initial);
     599             : 
     600             :     // Find the apps that are nearest to the same position
     601             :     // Global search over all applications
     602       20815 :     for (processor_id_type i_proc = 0; i_proc < n_processors(); ++i_proc)
     603             :     {
     604             :       // We need i_from to correspond to the global app index
     605       11768 :       unsigned int from0 = _global_app_start_per_proc[i_proc];
     606       33840 :       for (unsigned int i_from = from0; i_from < from0 + _froms_per_proc[i_proc]; ++i_from)
     607             :       {
     608       22072 :         if (_greedy_search || _search_value_conflicts || i_from == nearest_index)
     609             :         {
     610       12127 :           processors.insert(i_proc);
     611       12127 :           found = true;
     612             :         }
     613             :         mooseAssert(i_from < getFromMultiApp()->numGlobalApps(), "We should not reach this");
     614             :       }
     615             :     }
     616             :     mooseAssert((getFromMultiApp()->numGlobalApps() < n_processors() || processors.size() == 1) ||
     617             :                     _greedy_search || _search_value_conflicts,
     618             :                 "Should only be one source processor when using more processors than source apps");
     619             :   }
     620     5939514 :   else if (_use_bounding_boxes)
     621             :   {
     622             :     // We examine all (global) bounding boxes and find the minimum of the maximum distances within a
     623             :     // bounding box from the point. This creates a sphere around the point of interest. Any app
     624             :     // with a bounding box that intersects this sphere (with a bboxMinDistance <
     625             :     // nearest_max_distance) will be considered a potential source
     626             :     // NOTE: This is a heuristic. We could try others
     627             :     // NOTE: from_bboxes are in the reference space, as is the point.
     628     5936945 :     Real nearest_max_distance = std::numeric_limits<Real>::max();
     629    14529802 :     for (const auto & bbox : _from_bboxes)
     630             :     {
     631     8592857 :       Real distance = bboxMaxDistance(point, bbox);
     632     8592857 :       if (distance < nearest_max_distance)
     633     7187956 :         nearest_max_distance = distance;
     634             :     }
     635             : 
     636     5936945 :     unsigned int from0 = 0;
     637    13935828 :     for (processor_id_type i_proc = 0; i_proc < n_processors();
     638     7998883 :          from0 += _froms_per_proc[i_proc], ++i_proc)
     639             :       // i_from here is a hybrid index based on the cumulative sum of the apps per processor
     640    16591740 :       for (unsigned int i_from = from0; i_from < from0 + _froms_per_proc[i_proc]; ++i_from)
     641             :       {
     642     8592857 :         Real distance = bboxMinDistance(point, _from_bboxes[i_from]);
     643             :         // We will not break here because we want to send a point to all possible source domains
     644     8764994 :         if (_greedy_search || distance <= nearest_max_distance ||
     645      172137 :             _from_bboxes[i_from].contains_point(point))
     646             :         {
     647     8421438 :           processors.insert(i_proc);
     648     8421438 :           found = true;
     649             :         }
     650             :       }
     651             :   }
     652             :   // Greedy search will contact every single processor. It's not scalable, but if there's valid data
     653             :   // on any subapp on any process, it will find it
     654        2569 :   else if (_greedy_search)
     655             :   {
     656        2340 :     found = true;
     657        5460 :     for (const auto i_proc : make_range(n_processors()))
     658        3120 :       processors.insert(i_proc);
     659             :   }
     660             :   // Since we indicated that we only wanted values from a subapp with the same global index as the
     661             :   // target mesh division, we might as well only communicate with the process that owns this app
     662         454 :   else if (!_to_mesh_divisions.empty() &&
     663         225 :            _to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX)
     664             :   {
     665             :     // The target point could have a different index in each target mesh division. So on paper, we
     666             :     // would need to check all of them.
     667         225 :     auto saved_target_div = MooseMeshDivision::INVALID_DIVISION_INDEX;
     668         450 :     for (const auto i_to : index_range(_to_meshes))
     669             :     {
     670         450 :       const auto target_div = _to_mesh_divisions[i_to]->divisionIndex(
     671         225 :           _to_transforms[getGlobalTargetAppIndex(i_to)]->mapBack(point));
     672             :       // If it's the same division index, do not redo the search
     673         225 :       if (target_div == saved_target_div)
     674           0 :         continue;
     675             :       else
     676         225 :         saved_target_div = target_div;
     677             : 
     678             :       // Look for the processors owning a source-app with an index equal to the target mesh division
     679         525 :       for (const auto i_proc : make_range(n_processors()))
     680        1200 :         for (const auto i_from : make_range(_froms_per_proc[i_proc]))
     681         900 :           if (target_div == _global_app_start_per_proc[i_proc] + i_from)
     682             :           {
     683         225 :             processors.insert(i_proc);
     684         225 :             found = true;
     685             :           }
     686             :     }
     687             :   }
     688             :   else
     689           4 :     mooseError("No algorithm were selected to find which processes may send value data "
     690             :                "for a each target point. Please either specify using bounding boxes, "
     691             :                "greedy search, or to_mesh_division-based parameters");
     692             : 
     693             :   // Error out if we could not find this point when ask us to do so
     694     5948557 :   if (!found && _error_on_miss)
     695           0 :     mooseError(
     696             :         "Cannot find a source application to provide a value at point: ",
     697             :         point,
     698             :         " \n ",
     699             :         "It must be that mismatched meshes, between the source and target application, are being "
     700             :         "used.\nIf you are using the bounding boxes or nearest-app heuristics, or mesh-divisions, "
     701             :         "please consider using the greedy_search to confirm. Then consider choosing a different "
     702             :         "transfer type.\nThis check can be turned off by setting 'error_on_miss' to false. The "
     703             :         "'extrapolation_constant' parameter will be used to set the local value at missed points.");
     704     5948557 : }
     705             : 
     706             : void
     707     5948561 : MultiAppGeneralFieldTransfer::cacheOutgoingPointInfo(const Point point,
     708             :                                                      const dof_id_type dof_object_id,
     709             :                                                      const unsigned int problem_id,
     710             :                                                      ProcessorToPointVec & outgoing_points)
     711             : {
     712     5948561 :   std::set<processor_id_type> processors;
     713             :   // Find which processors will receive point data so they can send back value data
     714             :   // The list can be larger than needed, depending on the heuristic / algorithm used to make
     715             :   // the call on whether a processor (and the apps it runs) should be involved
     716     5948561 :   processors.clear();
     717     5948561 :   locatePointReceivers(point, processors);
     718             : 
     719             :   // We need to send this location data to these processors so they can send back values
     720    13924735 :   for (const auto pid : processors)
     721             :   {
     722             :     // Select which from_mesh_division the source data must come from for this point
     723     7976178 :     unsigned int required_source_division = 0;
     724     7976178 :     if (_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX)
     725        8784 :       required_source_division = getGlobalTargetAppIndex(problem_id);
     726     7967394 :     else if (_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX ||
     727    15931668 :              _to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX ||
     728     7964274 :              _to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX)
     729        7890 :       required_source_division = _to_mesh_divisions[problem_id]->divisionIndex(
     730        7890 :           _to_transforms[getGlobalTargetAppIndex(problem_id)]->mapBack(point));
     731             : 
     732             :     // Skip if we already know we don't want the point
     733     7976178 :     if (required_source_division == MooseMeshDivision::INVALID_DIVISION_INDEX)
     734           0 :       continue;
     735             : 
     736             :     // Store outgoing information for every source process
     737     7976178 :     outgoing_points[pid].push_back(std::pair<Point, unsigned int>(point, required_source_division));
     738             : 
     739             :     // Store point information locally for processing received data
     740             :     // We can use these information when inserting values into the solution vector
     741             :     PointInfo pointinfo;
     742     7976178 :     pointinfo.problem_id = problem_id;
     743     7976178 :     pointinfo.dof_object_id = dof_object_id;
     744     7976178 :     _processor_to_pointInfoVec[pid].push_back(pointinfo);
     745             :   }
     746     5948557 : }
     747             : 
     748             : void
     749       56854 : MultiAppGeneralFieldTransfer::extractOutgoingPoints(const unsigned int var_index,
     750             :                                                     ProcessorToPointVec & outgoing_points)
     751             : {
     752             :   // Get the variable name, with the accommodation for array/vector names
     753       56854 :   const auto & var_name = getToVarName(var_index);
     754             : 
     755             :   // Clean up the map from processor to pointInfo vector
     756             :   // This map should be consistent with outgoing_points
     757       56854 :   _processor_to_pointInfoVec.clear();
     758             : 
     759             :   // Loop over all problems
     760      115967 :   for (const auto i_to : index_range(_to_problems))
     761             :   {
     762       59121 :     const auto global_i_to = getGlobalTargetAppIndex(i_to);
     763             : 
     764             :     // libMesh EquationSystems
     765       59121 :     auto & es = getEquationSystem(*_to_problems[i_to], _displaced_target_mesh);
     766             :     // libMesh system that has this variable
     767       59121 :     System * to_sys = find_sys(es, var_name);
     768       59121 :     auto sys_num = to_sys->number();
     769       59121 :     auto var_num = _to_variables[var_index]->number();
     770       59121 :     auto & fe_type = _to_variables[var_index]->feType();
     771       59121 :     bool is_nodal = _to_variables[var_index]->isNodal();
     772             : 
     773             :     // Moose mesh
     774       59121 :     const auto & to_moose_mesh = _to_problems[i_to]->mesh(_displaced_target_mesh);
     775       59121 :     const auto & to_mesh = to_moose_mesh.getMesh();
     776             : 
     777             :     // We support more general variables via libMesh GenericProjector
     778       59121 :     if (fe_type.order > CONSTANT && !is_nodal)
     779             :     {
     780         370 :       GeneralFieldTransfer::RecordRequests<Number> f;
     781         370 :       GeneralFieldTransfer::RecordRequests<Gradient> g;
     782         370 :       GeneralFieldTransfer::NullAction<Number> nullsetter;
     783         370 :       const std::vector<unsigned int> varvec(1, var_num);
     784             : 
     785             :       libMesh::GenericProjector<GeneralFieldTransfer::RecordRequests<Number>,
     786             :                                 GeneralFieldTransfer::RecordRequests<Gradient>,
     787             :                                 Number,
     788             :                                 GeneralFieldTransfer::NullAction<Number>>
     789         370 :           request_gather(*to_sys, f, &g, nullsetter, varvec);
     790             : 
     791             :       // Defining only boundary values will not be enough to describe the variable, disallow it
     792         370 :       if (_to_boundaries.size() && (_to_variables[var_index]->getContinuity() == DISCONTINUOUS))
     793           4 :         mooseError("Higher order discontinuous elemental variables are not supported for "
     794             :                    "target-boundary "
     795             :                    "restricted transfers");
     796             : 
     797             :       // Not implemented as the target mesh division could similarly be cutting elements in an
     798             :       // arbitrary way with not enough requested points to describe the target variable
     799         366 :       if (!_to_mesh_divisions.empty() && !_to_mesh_divisions[i_to]->coversEntireMesh())
     800           0 :         mooseError("Higher order variable support not implemented for target mesh division "
     801             :                    "unless the mesh is fully covered / indexed in the mesh division. This must be "
     802             :                    "set programmatically in the MeshDivision object used.");
     803             : 
     804             :       // We dont look at boundary restriction, not supported for higher order target variables
     805             :       // Same for mesh divisions
     806         366 :       const auto & to_begin = _to_blocks.empty()
     807         615 :                                   ? to_mesh.active_local_elements_begin()
     808         366 :                                   : to_mesh.active_local_subdomain_set_elements_begin(_to_blocks);
     809             : 
     810         366 :       const auto & to_end = _to_blocks.empty()
     811         615 :                                 ? to_mesh.active_local_elements_end()
     812         366 :                                 : to_mesh.active_local_subdomain_set_elements_end(_to_blocks);
     813             : 
     814         366 :       ConstElemRange to_elem_range(to_begin, to_end);
     815             : 
     816         366 :       request_gather.project(to_elem_range);
     817             : 
     818         366 :       dof_id_type point_id = 0;
     819      298716 :       for (const Point & p : f.points_requested())
     820             :         // using the point number as a "dof_object_id" will serve to identify the point if we ever
     821             :         // rework interp/distance_cache into the dof_id_to_value maps
     822      298350 :         this->cacheOutgoingPointInfo(
     823      298350 :             (*_to_transforms[global_i_to])(p), point_id++, i_to, outgoing_points);
     824             : 
     825             :       // This is going to require more complicated transfer work
     826         366 :       if (!g.points_requested().empty())
     827           0 :         mooseError("We don't currently support variables with gradient degrees of freedom");
     828         366 :     }
     829       58751 :     else if (is_nodal)
     830             :     {
     831     5790305 :       for (const auto & node : to_mesh.local_node_ptr_range())
     832             :       {
     833             :         // Skip this node if the variable has no dofs at it.
     834     5735080 :         if (node->n_dofs(sys_num, var_num) < 1)
     835      295207 :           continue;
     836             : 
     837             :         // Skip if it is a block restricted transfer and current node does not have
     838             :         // specified blocks
     839     5439873 :         if (!_to_blocks.empty() && !inBlocks(_to_blocks, to_moose_mesh, node))
     840       14818 :           continue;
     841             : 
     842     5425055 :         if (!_to_boundaries.empty() && !onBoundaries(_to_boundaries, to_moose_mesh, node))
     843      190960 :           continue;
     844             : 
     845             :         // Skip if the node does not meet the target mesh division behavior
     846             :         // We cannot know from which app the data will come from so we cannot know
     847             :         // the source mesh division index and the source app global index
     848     5234095 :         if (!_to_mesh_divisions.empty() && _to_mesh_divisions[i_to]->divisionIndex(*node) ==
     849             :                                                MooseMeshDivision::INVALID_DIVISION_INDEX)
     850        4950 :           continue;
     851             : 
     852             :         // Cache point information
     853             :         // We will use this information later for setting values back to solution vectors
     854     5229145 :         cacheOutgoingPointInfo(
     855     5229145 :             (*_to_transforms[global_i_to])(*node), node->id(), i_to, outgoing_points);
     856       55225 :       }
     857             :     }
     858             :     else // Elemental, constant monomial
     859             :     {
     860        3522 :       for (const auto & elem :
     861      481365 :            as_range(to_mesh.local_elements_begin(), to_mesh.local_elements_end()))
     862             :       {
     863             :         // Skip this element if the variable has no dofs at it.
     864      474321 :         if (elem->n_dofs(sys_num, var_num) < 1)
     865        2458 :           continue;
     866             : 
     867             :         // Skip if the element is not inside the block restriction
     868      471863 :         if (!_to_blocks.empty() && !inBlocks(_to_blocks, elem))
     869        5116 :           continue;
     870             : 
     871             :         // Skip if the element does not have a side on the boundary
     872      466747 :         if (!_to_boundaries.empty() && !onBoundaries(_to_boundaries, to_moose_mesh, elem))
     873       42315 :           continue;
     874             : 
     875             :         // Skip if the element is not indexed within the mesh division
     876      424432 :         if (!_to_mesh_divisions.empty() && _to_mesh_divisions[i_to]->divisionIndex(*elem) ==
     877             :                                                MooseMeshDivision::INVALID_DIVISION_INDEX)
     878        3366 :           continue;
     879             : 
     880             :         // Cache point information
     881             :         // We will use this information later for setting values back to solution vectors
     882      421066 :         cacheOutgoingPointInfo((*_to_transforms[global_i_to])(elem->vertex_average()),
     883      421066 :                                elem->id(),
     884             :                                i_to,
     885             :                                outgoing_points);
     886        3522 :       } // for
     887             :     } // else
     888             :   } // for
     889       56846 : }
     890             : 
     891             : void
     892        4326 : MultiAppGeneralFieldTransfer::extractLocalFromBoundingBoxes(std::vector<BoundingBox> & local_bboxes)
     893             : {
     894        4326 :   local_bboxes.resize(_froms_per_proc[processor_id()]);
     895             :   // Find the index to the first of this processor's local bounding boxes.
     896        4326 :   unsigned int local_start = 0;
     897        5376 :   for (processor_id_type i_proc = 0; i_proc < n_processors() && i_proc != processor_id(); ++i_proc)
     898        1050 :     local_start += _froms_per_proc[i_proc];
     899             : 
     900             :   // Extract the local bounding boxes.
     901       10073 :   for (const auto i_from : make_range(_froms_per_proc[processor_id()]))
     902        5747 :     local_bboxes[i_from] = _from_bboxes[local_start + i_from];
     903        4326 : }
     904             : 
     905             : void
     906       85064 : MultiAppGeneralFieldTransfer::cacheIncomingInterpVals(
     907             :     processor_id_type pid,
     908             :     const unsigned int var_index,
     909             :     std::vector<PointInfo> & pointInfoVec,
     910             :     const std::vector<std::pair<Point, unsigned int>> & point_requests,
     911             :     const std::vector<std::pair<Real, Real>> & incoming_vals,
     912             :     DofobjectToInterpValVec & dofobject_to_valsvec,
     913             :     InterpCaches & interp_caches,
     914             :     InterpCaches & distance_caches)
     915             : {
     916             :   mooseAssert(pointInfoVec.size() == incoming_vals.size(),
     917             :               "Number of dof objects does not equal to the number of incoming values");
     918             : 
     919       85064 :   dof_id_type val_offset = 0;
     920     8061242 :   for (const auto & pointinfo : pointInfoVec)
     921             :   {
     922             :     // Retrieve target information from cached point infos
     923     7976178 :     const auto problem_id = pointinfo.problem_id;
     924     7976178 :     const auto dof_object_id = pointinfo.dof_object_id;
     925             : 
     926     7976178 :     auto & fe_type = _to_variables[var_index]->feType();
     927     7976178 :     bool is_nodal = _to_variables[var_index]->isNodal();
     928             : 
     929             :     // In the higher order elemental variable case, we receive point values, not nodal or
     930             :     // elemental. We use an InterpCache to store the values. The distance_cache is necessary to
     931             :     // choose between multiple origin problems sending values. This code could be unified with the
     932             :     // lower order order case by using the dofobject_to_valsvec
     933     7976178 :     if (fe_type.order > CONSTANT && !is_nodal)
     934             :     {
     935             :       // Cache solution on target mesh in its local frame of reference
     936      393825 :       InterpCache & value_cache = interp_caches[problem_id];
     937      393825 :       InterpCache & distance_cache = distance_caches[problem_id];
     938      393825 :       Point p = _to_transforms[getGlobalTargetAppIndex(problem_id)]->mapBack(
     939      393825 :           point_requests[val_offset].first);
     940      393825 :       const Number val = incoming_vals[val_offset].first;
     941             : 
     942             :       // Initialize distance to be able to compare
     943      393825 :       if (!distance_cache.hasKey(p))
     944      298350 :         distance_cache[p] = std::numeric_limits<Real>::max();
     945             : 
     946             :       // We should only have one closest value for each variable at any given point.
     947             :       // While there are shared Qps, on vertices for higher order variables usually,
     948             :       // the generic projector only queries each point once
     949      393825 :       if (_search_value_conflicts && !GeneralFieldTransfer::isBetterOutOfMeshValue(val) &&
     950      393825 :           value_cache.hasKey(p) != 0 && !MooseUtils::absoluteFuzzyEqual(value_cache[p], val) &&
     951           0 :           MooseUtils::absoluteFuzzyEqual(distance_cache[p], incoming_vals[val_offset].second))
     952           0 :         registerConflict(problem_id, dof_object_id, p, incoming_vals[val_offset].second, false);
     953             : 
     954             :       // if we use the nearest app, even if the value is bad we want to save the distance because
     955             :       // it's the distance to the app, if that's the closest app then so be it with the bad value
     956      469760 :       if ((!GeneralFieldTransfer::isBetterOutOfMeshValue(val) || _use_nearest_app) &&
     957       75935 :           MooseUtils::absoluteFuzzyGreaterThan(distance_cache[p], incoming_vals[val_offset].second))
     958             :       {
     959             :         // NOTE: We store the distance as well as the value. We really only need the
     960             :         // value to construct the variable, but the distance is used to make decisions in nearest
     961             :         // node schemes on which value to use
     962       72358 :         value_cache[p] = val;
     963       72358 :         distance_cache[p] = incoming_vals[val_offset].second;
     964             :       }
     965             :     }
     966             :     else
     967             :     {
     968             :       // Using the dof object pointer, so we can handle
     969             :       // both element and node using the same code
     970             : #ifndef NDEBUG
     971             :       auto var_num = _to_variables[var_index]->number();
     972             :       auto & to_sys = _to_variables[var_index]->sys();
     973             : 
     974             :       const MeshBase & to_mesh = _to_problems[problem_id]->mesh(_displaced_target_mesh).getMesh();
     975             :       const DofObject * dof_object_ptr = nullptr;
     976             :       const auto sys_num = to_sys.number();
     977             :       // It is a node
     978             :       if (is_nodal)
     979             :         dof_object_ptr = to_mesh.node_ptr(dof_object_id);
     980             :       // It is an element
     981             :       else
     982             :         dof_object_ptr = to_mesh.elem_ptr(dof_object_id);
     983             : 
     984             :       // We should only be supporting nodal and constant elemental
     985             :       // variables in this code path; if we see multiple DoFs on one
     986             :       // object we should have been using GenericProjector
     987             :       mooseAssert(dof_object_ptr->n_dofs(sys_num, var_num) == 1,
     988             :                   "Unexpectedly found " << dof_object_ptr->n_dofs(sys_num, var_num)
     989             :                                         << "dofs instead of 1");
     990             : #endif
     991             : 
     992     7582353 :       auto & dofobject_to_val = dofobject_to_valsvec[problem_id];
     993             : 
     994             :       // Check if we visited this dof object earlier
     995     7582353 :       auto values_ptr = dofobject_to_val.find(dof_object_id);
     996             :       // We did not visit this
     997     7582353 :       if (values_ptr == dofobject_to_val.end())
     998             :       {
     999             :         // Values for this dof object
    1000     5650207 :         auto & val = dofobject_to_val[dof_object_id];
    1001             :         // Interpolation value
    1002     5650207 :         val.interp = incoming_vals[val_offset].first;
    1003             :         // Where this value came from
    1004     5650207 :         val.pid = pid;
    1005             :         // Distance
    1006     5650207 :         val.distance = incoming_vals[val_offset].second;
    1007             :       }
    1008             :       else
    1009             :       {
    1010     1932146 :         auto & val = values_ptr->second;
    1011             : 
    1012             :         // Look for value conflicts
    1013     1932146 :         if (detectConflict(val.interp,
    1014     1932146 :                            incoming_vals[val_offset].first,
    1015             :                            val.distance,
    1016     1932146 :                            incoming_vals[val_offset].second))
    1017             :         {
    1018             :           // Keep track of distance and value
    1019             :           const Point p =
    1020           0 :               getPointInTargetAppFrame(point_requests[val_offset].first,
    1021             :                                        problem_id,
    1022             :                                        "Registration of received equi-distant value conflict");
    1023           0 :           registerConflict(problem_id, dof_object_id, p, incoming_vals[val_offset].second, false);
    1024             :         }
    1025             : 
    1026             :         // We adopt values that are, in order of priority
    1027             :         // - valid (or from nearest app)
    1028             :         // - closest distance
    1029             :         // - the smallest rank with the same distance
    1030             :         // It is debatable whether we want invalid values from the nearest app. It could just be
    1031             :         // that the app position was closer but the extent of another child app was large enough
    1032     3864292 :         if ((!GeneralFieldTransfer::isBetterOutOfMeshValue(incoming_vals[val_offset].first) ||
    1033     3646220 :              _use_nearest_app) &&
    1034     1714074 :             (MooseUtils::absoluteFuzzyGreaterThan(val.distance, incoming_vals[val_offset].second) ||
    1035     1363995 :              ((val.pid > pid) &&
    1036      590231 :               MooseUtils::absoluteFuzzyEqual(val.distance, incoming_vals[val_offset].second))))
    1037             :         {
    1038      350492 :           val.interp = incoming_vals[val_offset].first;
    1039      350492 :           val.pid = pid;
    1040      350492 :           val.distance = incoming_vals[val_offset].second;
    1041             :         }
    1042             :       }
    1043             :     }
    1044             : 
    1045             :     // Move it to next position
    1046     7976178 :     val_offset++;
    1047             :   }
    1048       85064 : }
    1049             : 
    1050             : void
    1051         438 : MultiAppGeneralFieldTransfer::registerConflict(
    1052             :     unsigned int problem, dof_id_type dof_id, Point p, Real dist, bool local)
    1053             : {
    1054             :   // NOTE We could be registering the same conflict several times, we could count them instead
    1055         438 :   if (local)
    1056         438 :     _local_conflicts.push_back(std::make_tuple(problem, dof_id, p, dist));
    1057             :   else
    1058           0 :     _received_conflicts.push_back(std::make_tuple(problem, dof_id, p, dist));
    1059         438 : }
    1060             : 
    1061             : void
    1062         211 : MultiAppGeneralFieldTransfer::examineReceivedValueConflicts(
    1063             :     const unsigned int var_index,
    1064             :     const DofobjectToInterpValVec & dofobject_to_valsvec,
    1065             :     const InterpCaches & distance_caches)
    1066             : {
    1067         211 :   const auto var_name = getToVarName(var_index);
    1068             :   // We must check a posteriori because we could have two
    1069             :   // equidistant points with different values from two different problems, but a third point from
    1070             :   // another problem is actually closer, so there is no conflict because only that last one
    1071             :   // matters We check here whether the potential conflicts actually were the nearest points Loop
    1072             :   // over potential conflicts
    1073         211 :   for (auto conflict_it = _received_conflicts.begin(); conflict_it != _received_conflicts.end();)
    1074             :   {
    1075           0 :     const auto potential_conflict = *conflict_it;
    1076           0 :     bool overlap_found = false;
    1077             : 
    1078             :     // Extract info for the potential conflict
    1079           0 :     const unsigned int problem_id = std::get<0>(potential_conflict);
    1080           0 :     const dof_id_type dof_object_id = std::get<1>(potential_conflict);
    1081           0 :     const Point p = std::get<2>(potential_conflict);
    1082           0 :     const Real distance = std::get<3>(potential_conflict);
    1083             : 
    1084             :     // Extract target variable info
    1085           0 :     auto & es = getEquationSystem(*_to_problems[problem_id], _displaced_target_mesh);
    1086           0 :     System * to_sys = find_sys(es, var_name);
    1087           0 :     auto var_num = to_sys->variable_number(var_name);
    1088           0 :     auto & fe_type = to_sys->variable_type(var_num);
    1089           0 :     bool is_nodal = _to_variables[var_index]->isNodal();
    1090             : 
    1091             :     // Higher order elemental
    1092           0 :     if (fe_type.order > CONSTANT && !is_nodal)
    1093             :     {
    1094           0 :       auto cached_distance = distance_caches[problem_id].find(p);
    1095           0 :       if (cached_distance == distance_caches[problem_id].end())
    1096           0 :         mooseError("Conflict point was not found in the map of all origin-target distances");
    1097             :       // Distance is still the distance when we detected a potential overlap
    1098           0 :       if (MooseUtils::absoluteFuzzyEqual(cached_distance->second, distance))
    1099           0 :         overlap_found = true;
    1100             :     }
    1101             :     // Nodal and const monomial variable
    1102           0 :     else if (MooseUtils::absoluteFuzzyEqual(
    1103           0 :                  dofobject_to_valsvec[problem_id].find(dof_object_id)->second.distance, distance))
    1104           0 :       overlap_found = true;
    1105             : 
    1106             :     // Map will only keep the actual overlaps
    1107           0 :     if (!overlap_found)
    1108           0 :       _received_conflicts.erase(conflict_it);
    1109             :     else
    1110           0 :       ++conflict_it;
    1111             :   }
    1112         211 : }
    1113             : 
    1114             : void
    1115         211 : MultiAppGeneralFieldTransfer::examineLocalValueConflicts(
    1116             :     const unsigned int var_index,
    1117             :     const DofobjectToInterpValVec & dofobject_to_valsvec,
    1118             :     const InterpCaches & distance_caches)
    1119             : {
    1120         211 :   const auto var_name = getToVarName(var_index);
    1121             :   // We must check a posteriori because we could have:
    1122             :   // - two equidistant points with different values from two different problems
    1123             :   // - two (or more) equidistant points with different values from the same problem
    1124             :   // but a third point/value couple from another problem is actually closer, so there is no
    1125             :   // conflict because only that last one matters. We check here whether the potential conflicts
    1126             :   // actually were the nearest points. We use several global reductions. If there are not too many
    1127             :   // potential conflicts (and there should not be in a well-posed problem) it should be manageably
    1128             :   // expensive
    1129             : 
    1130             :   // Move relevant conflict info (location, distance) to a smaller data structure
    1131         211 :   std::vector<std::tuple<Point, Real>> potential_conflicts;
    1132         211 :   potential_conflicts.reserve(_local_conflicts.size());
    1133             : 
    1134             :   // Loop over potential conflicts to broadcast all the conflicts
    1135         649 :   for (auto conflict_it = _local_conflicts.begin(); conflict_it != _local_conflicts.end();
    1136         438 :        ++conflict_it)
    1137             :   {
    1138             :     // Extract info for the potential conflict
    1139         438 :     const auto potential_conflict = *conflict_it;
    1140         438 :     const unsigned int i_from = std::get<0>(potential_conflict);
    1141         438 :     Point p = std::get<2>(potential_conflict);
    1142         438 :     const Real distance = std::get<3>(potential_conflict);
    1143             :     // If not using nearest-positions: potential conflict was saved in the source frame
    1144             :     // If using nearest-positions: potential conflict was saved in the reference frame
    1145         438 :     if (!_nearest_positions_obj)
    1146             :     {
    1147         390 :       const auto from_global_num = getGlobalSourceAppIndex(i_from);
    1148         390 :       p = (*_from_transforms[from_global_num])(p);
    1149             :     }
    1150             : 
    1151             :     // Send data in the global frame of reference
    1152         438 :     potential_conflicts.push_back(std::make_tuple(p, distance));
    1153             :   }
    1154         211 :   _communicator.allgather(potential_conflicts, false);
    1155             :   // conflicts could have been reported multiple times within a tolerance
    1156         211 :   std::sort(potential_conflicts.begin(), potential_conflicts.end());
    1157         211 :   potential_conflicts.erase(unique(potential_conflicts.begin(),
    1158             :                                    potential_conflicts.end(),
    1159         407 :                                    [](auto l, auto r)
    1160             :                                    {
    1161         451 :                                      return std::get<0>(l).absolute_fuzzy_equals(std::get<0>(r)) &&
    1162         451 :                                             std::abs(std::get<1>(l) - std::get<1>(r)) < TOLERANCE;
    1163             :                                    }),
    1164         211 :                             potential_conflicts.end());
    1165             : 
    1166         211 :   std::vector<std::tuple<Point, Real>> real_conflicts;
    1167         211 :   real_conflicts.reserve(potential_conflicts.size());
    1168             : 
    1169             :   // For each potential conflict, we need to identify what problem asked for that value
    1170         605 :   for (auto conflict_it = potential_conflicts.begin(); conflict_it != potential_conflicts.end();
    1171         394 :        ++conflict_it)
    1172             :   {
    1173             :     // Extract info for the potential conflict
    1174         394 :     auto potential_conflict = *conflict_it;
    1175         394 :     const Point p = std::get<0>(potential_conflict);
    1176         394 :     const Real distance = std::get<1>(potential_conflict);
    1177             : 
    1178             :     // Check all the problems to try to find this requested point in the data structures filled
    1179             :     // with the received information
    1180         394 :     bool target_found = false;
    1181         394 :     bool conflict_real = false;
    1182         788 :     for (const auto i_to : index_range(_to_problems))
    1183             :     {
    1184             :       // Extract variable info
    1185         394 :       auto & es = getEquationSystem(*_to_problems[i_to], _displaced_target_mesh);
    1186         394 :       System * to_sys = find_sys(es, var_name);
    1187         394 :       auto var_num = to_sys->variable_number(var_name);
    1188         394 :       auto & fe_type = to_sys->variable_type(var_num);
    1189         394 :       bool is_nodal = _to_variables[var_index]->isNodal();
    1190             : 
    1191             :       // Move to the local frame of reference for the target problem
    1192             :       Point local_p =
    1193         788 :           getPointInTargetAppFrame(p, i_to, "Resolution of local value conflicts detected");
    1194             : 
    1195             :       // Higher order elemental
    1196         394 :       if (fe_type.order > CONSTANT && !is_nodal)
    1197             :       {
    1198             :         // distance_caches finds use a binned floating point search
    1199           0 :         auto cached_distance = distance_caches[i_to].find(local_p);
    1200           0 :         if (cached_distance != distance_caches[i_to].end())
    1201             :         {
    1202           0 :           target_found = true;
    1203             :           // Distance between source & target is still the distance we found in the sending
    1204             :           // process when we detected a potential overlap while gathering values to send
    1205           0 :           if (MooseUtils::absoluteFuzzyEqual(cached_distance->second, distance))
    1206           0 :             conflict_real = true;
    1207             :         }
    1208             :       }
    1209             :       // Nodal-value-dof-only and const monomial variable
    1210             :       else
    1211             :       {
    1212             :         // Find the dof id for the variable to be set
    1213         394 :         dof_id_type dof_object_id = std::numeric_limits<dof_id_type>::max();
    1214         394 :         auto pl = _to_problems[i_to]->mesh().getPointLocator();
    1215         394 :         pl->enable_out_of_mesh_mode();
    1216         394 :         if (is_nodal)
    1217             :         {
    1218         366 :           auto node = pl->locate_node(local_p);
    1219         366 :           if (node)
    1220             :             // this is not the dof_id for the variable, but the dof_object_id
    1221         366 :             dof_object_id = node->id();
    1222             :         }
    1223             :         else
    1224             :         {
    1225          28 :           auto elem = (*pl)(local_p);
    1226          28 :           if (elem)
    1227          28 :             dof_object_id = elem->id();
    1228             :         }
    1229         394 :         pl->disable_out_of_mesh_mode();
    1230             : 
    1231             :         // point isn't even in mesh
    1232         394 :         if (dof_object_id == std::numeric_limits<dof_id_type>::max())
    1233           0 :           continue;
    1234             : 
    1235             :         // this dof was not requested by this problem on this process
    1236         394 :         if (dofobject_to_valsvec[i_to].find(dof_object_id) == dofobject_to_valsvec[i_to].end())
    1237           0 :           continue;
    1238             : 
    1239         394 :         target_found = true;
    1240             :         // Check the saved distance in the vector of saved results. If the same, then the local
    1241             :         // conflict we detected with that distance is still an issue after receiving all values
    1242         788 :         if (MooseUtils::absoluteFuzzyEqual(
    1243         394 :                 dofobject_to_valsvec[i_to].find(dof_object_id)->second.distance, distance))
    1244         394 :           conflict_real = true;
    1245         394 :       }
    1246             :     }
    1247             :     // Only keep the actual conflicts / overlaps
    1248         394 :     if (target_found && conflict_real)
    1249         394 :       real_conflicts.push_back(potential_conflict);
    1250             :   }
    1251             : 
    1252             :   // Communicate real conflicts to all so they can be checked by every process
    1253         211 :   _communicator.allgather(real_conflicts, false);
    1254             : 
    1255             :   // Delete potential conflicts that were resolved
    1256             :   // Each local list of conflicts will now be updated. It's important to keep conflict lists local
    1257             :   // so we can give more context like the sending processor id (the domain of which can be
    1258             :   // inspected by the user)
    1259         649 :   for (auto conflict_it = _local_conflicts.begin(); conflict_it != _local_conflicts.end();)
    1260             :   {
    1261             :     // Extract info for the potential conflict
    1262         438 :     const auto potential_conflict = *conflict_it;
    1263         438 :     const unsigned int i_from = std::get<0>(potential_conflict);
    1264         438 :     Point p = std::get<2>(potential_conflict);
    1265         438 :     const Real distance = std::get<3>(potential_conflict);
    1266         438 :     if (!_nearest_positions_obj)
    1267             :     {
    1268         390 :       const auto from_global_num = getGlobalSourceAppIndex(i_from);
    1269         390 :       p = (*_from_transforms[from_global_num])(p);
    1270             :     }
    1271             : 
    1272             :     // If not in the vector of real conflicts, was not real so delete it
    1273         438 :     if (std::find_if(real_conflicts.begin(),
    1274             :                      real_conflicts.end(),
    1275        6599 :                      [p, distance](const auto & item)
    1276             :                      {
    1277        7037 :                        return std::get<0>(item).absolute_fuzzy_equals(p) &&
    1278        7037 :                               std::abs(std::get<1>(item) - distance) < TOLERANCE;
    1279         876 :                      }) == real_conflicts.end())
    1280           0 :       _local_conflicts.erase(conflict_it);
    1281             :     else
    1282         438 :       ++conflict_it;
    1283             :   }
    1284         211 : }
    1285             : 
    1286             : void
    1287         211 : MultiAppGeneralFieldTransfer::outputValueConflicts(
    1288             :     const unsigned int var_index,
    1289             :     const DofobjectToInterpValVec & dofobject_to_valsvec,
    1290             :     const InterpCaches & distance_caches)
    1291             : {
    1292             :   // Remove potential conflicts that did not materialize, the value did not end up being used
    1293         211 :   examineReceivedValueConflicts(var_index, dofobject_to_valsvec, distance_caches);
    1294         211 :   examineLocalValueConflicts(var_index, dofobject_to_valsvec, distance_caches);
    1295             : 
    1296             :   // Output the conflicts from the selection of local values (evaluateInterpValues-type routines)
    1297             :   // to send in response to value requests at target points
    1298         211 :   const std::string rank_str = std::to_string(_communicator.rank());
    1299         211 :   if (_local_conflicts.size())
    1300             :   {
    1301          31 :     unsigned int num_outputs = 0;
    1302          62 :     std::string local_conflicts_string = "";
    1303             :     std::string potential_reasons =
    1304             :         "Are some points in target mesh equidistant from the sources "
    1305          31 :         "(nodes/centroids/apps/positions, depending on transfer) in origin mesh(es)?\n";
    1306          31 :     if (hasFromMultiApp() && _from_problems.size() > 1)
    1307          23 :       potential_reasons += "Are multiple subapps overlapping?\n";
    1308         469 :     for (const auto & conflict : _local_conflicts)
    1309             :     {
    1310         438 :       const unsigned int problem_id = std::get<0>(conflict);
    1311         438 :       Point p = std::get<2>(conflict);
    1312         438 :       num_outputs++;
    1313             : 
    1314         438 :       std::string origin_domain_message;
    1315         438 :       if (hasFromMultiApp() && !_nearest_positions_obj)
    1316             :       {
    1317             :         // NOTES:
    1318             :         // - The origin app for a conflict may not be unique.
    1319             :         // - The conflicts vectors only store the conflictual points, not the original one
    1320             :         //   The original value found with a given distance could be retrieved from the main
    1321             :         //   caches
    1322         390 :         const auto app_id = _from_local2global_map[problem_id];
    1323         390 :         origin_domain_message = "In source child app " + std::to_string(app_id) + " mesh,";
    1324             :       }
    1325             :       // We can't locate the source app when considering nearest positions, so we saved the data
    1326             :       // in the reference space. So we return the conflict location in the target app (parent or
    1327             :       // sibling) instead
    1328          48 :       else if (hasFromMultiApp() && _nearest_positions_obj)
    1329             :       {
    1330          48 :         if (_to_problems.size() == 1 || _skip_coordinate_collapsing)
    1331             :         {
    1332          48 :           p = (*_to_transforms[0])(p);
    1333          48 :           origin_domain_message = "In target app mesh,";
    1334             :         }
    1335             :         else
    1336           0 :           origin_domain_message = "In reference (post-coordinate collapse) mesh,";
    1337             :       }
    1338             :       else
    1339           0 :         origin_domain_message = "In source parent app mesh,";
    1340             : 
    1341         438 :       if (num_outputs < _search_value_conflicts_max_log)
    1342         348 :         local_conflicts_string += origin_domain_message + " point: (" + std::to_string(p(0)) +
    1343         696 :                                   ", " + std::to_string(p(1)) + ", " + std::to_string(p(2)) +
    1344         696 :                                   "), equi-distance: " + std::to_string(std::get<3>(conflict)) +
    1345         174 :                                   "\n";
    1346         264 :       else if (num_outputs == _search_value_conflicts_max_log)
    1347             :         local_conflicts_string +=
    1348             :             "Maximum output of the search for value conflicts has been reached. Further conflicts "
    1349           8 :             "will not be output.\nIncrease 'search_value_conflicts_max_log' to output more.";
    1350         438 :     }
    1351             :     // Explicitly name source to give more context
    1352          31 :     std::string source_str = "unknown";
    1353          31 :     if (_from_var_names.size())
    1354          27 :       source_str = "variable '" + getFromVarName(var_index) + "'";
    1355          12 :     else if (isParamValid("source_user_object"))
    1356          12 :       source_str = "user object '" + getParam<UserObjectName>("source_user_object") + "'";
    1357             : 
    1358          62 :     mooseWarning("On rank " + rank_str +
    1359             :                  ", multiple valid values from equidistant points were "
    1360          31 :                  "found in the origin mesh for source " +
    1361          62 :                  source_str + " for " + std::to_string(_local_conflicts.size()) +
    1362          31 :                  " target points.\n" + potential_reasons + "Conflicts detected at :\n" +
    1363             :                  local_conflicts_string);
    1364           0 :   }
    1365             : 
    1366             :   // Output the conflicts discovered when receiving values from multiple origin problems
    1367         180 :   if (_received_conflicts.size())
    1368             :   {
    1369           0 :     unsigned int num_outputs = 0;
    1370           0 :     std::string received_conflicts_string = "";
    1371             :     std::string potential_reasons =
    1372             :         "Are some points in target mesh equidistant from the sources "
    1373           0 :         "(nodes/centroids/apps/positions, depending on transfer) in origin mesh(es)?\n";
    1374           0 :     if (hasToMultiApp() && _to_problems.size() > 1)
    1375           0 :       potential_reasons += "Are multiple subapps overlapping?\n";
    1376           0 :     for (const auto & conflict : _received_conflicts)
    1377             :     {
    1378             :       // Extract info for the potential overlap
    1379           0 :       const unsigned int problem_id = std::get<0>(conflict);
    1380           0 :       const Point p = std::get<2>(conflict);
    1381           0 :       num_outputs++;
    1382             : 
    1383           0 :       std::string target_domain_message;
    1384           0 :       if (hasToMultiApp())
    1385             :       {
    1386           0 :         const auto app_id = _to_local2global_map[problem_id];
    1387           0 :         target_domain_message = "In target child app " + std::to_string(app_id) + " mesh,";
    1388             :       }
    1389             :       else
    1390           0 :         target_domain_message = "In target parent app mesh,";
    1391             : 
    1392           0 :       if (num_outputs < _search_value_conflicts_max_log)
    1393           0 :         received_conflicts_string += target_domain_message + " point: (" + std::to_string(p(0)) +
    1394           0 :                                      ", " + std::to_string(p(1)) + ", " + std::to_string(p(2)) +
    1395           0 :                                      "), equi-distance: " + std::to_string(std::get<3>(conflict)) +
    1396           0 :                                      "\n";
    1397           0 :       else if (num_outputs == _search_value_conflicts_max_log)
    1398             :         received_conflicts_string +=
    1399             :             "Maximum output of the search for value conflicts has been reached. Further conflicts "
    1400           0 :             "will not be output.\nIncrease 'search_value_conflicts_max_log' to output more.";
    1401           0 :     }
    1402           0 :     mooseWarning("On rank " + rank_str +
    1403             :                  ", multiple valid values from equidistant points were "
    1404           0 :                  "received for target variable '" +
    1405           0 :                  getToVarName(var_index) + "' for " + std::to_string(_received_conflicts.size()) +
    1406           0 :                  " target points.\n" + potential_reasons + "Conflicts detected at :\n" +
    1407             :                  received_conflicts_string);
    1408           0 :   }
    1409             : 
    1410         180 :   if (_local_conflicts.empty() && _received_conflicts.empty())
    1411             :   {
    1412         540 :     if (isParamSetByUser("search_value_conflict"))
    1413           0 :       mooseInfo("Automated diagnosis did not detect floating point indetermination in transfer");
    1414         180 :     else if (_to_problems.size() > 10 || _from_problems.size() > 10 || _communicator.size() > 10)
    1415           0 :       mooseInfo(
    1416             :           "Automated diagnosis did not detect any floating point indetermination in "
    1417             :           "the transfer. You may consider turning it off using `search_value_conflicts=false` "
    1418             :           "to improve performance/scalability.");
    1419             :   }
    1420             : 
    1421             :   // Reset the conflicts vectors, to be used for checking conflicts when transferring the next
    1422             :   // variable
    1423         180 :   _local_conflicts.clear();
    1424         180 :   _received_conflicts.clear();
    1425         180 : }
    1426             : 
    1427             : void
    1428       56815 : MultiAppGeneralFieldTransfer::setSolutionVectorValues(
    1429             :     const unsigned int var_index,
    1430             :     const DofobjectToInterpValVec & dofobject_to_valsvec,
    1431             :     const InterpCaches & interp_caches)
    1432             : {
    1433             :   // Get the variable name, with the accommodation for array/vector names
    1434       56815 :   const auto & var_name = getToVarName(var_index);
    1435             : 
    1436      115889 :   for (const auto problem_id : index_range(_to_problems))
    1437             :   {
    1438       59082 :     auto & dofobject_to_val = dofobject_to_valsvec[problem_id];
    1439             : 
    1440             :     // libMesh EquationSystems
    1441             :     // NOTE: we would expect to set variables from the displaced equation system here
    1442       59082 :     auto & es = getEquationSystem(*_to_problems[problem_id], false);
    1443             : 
    1444             :     // libMesh system
    1445       59082 :     System * to_sys = find_sys(es, var_name);
    1446             : 
    1447             :     // libMesh mesh
    1448       59082 :     const MeshBase & to_mesh = _to_problems[problem_id]->mesh(_displaced_target_mesh).getMesh();
    1449       59082 :     auto var_num = to_sys->variable_number(var_name);
    1450       59082 :     auto sys_num = to_sys->number();
    1451             : 
    1452       59082 :     auto & fe_type = _to_variables[var_index]->feType();
    1453       59082 :     bool is_nodal = _to_variables[var_index]->isNodal();
    1454             : 
    1455       59082 :     if (fe_type.order > CONSTANT && !is_nodal)
    1456             :     {
    1457             :       // We may need to use existing data values in places where the
    1458             :       // from app domain doesn't overlap
    1459         366 :       MeshFunction to_func(es, *to_sys->current_local_solution, to_sys->get_dof_map(), var_num);
    1460         366 :       to_func.init();
    1461             : 
    1462             :       GeneralFieldTransfer::CachedData<Number> f(
    1463         366 :           interp_caches[problem_id], to_func, _default_extrapolation_value);
    1464         366 :       libMesh::VectorSetAction<Number> setter(*to_sys->solution);
    1465         366 :       const std::vector<unsigned int> varvec(1, var_num);
    1466             : 
    1467             :       libMesh::GenericProjector<GeneralFieldTransfer::CachedData<Number>,
    1468             :                                 GeneralFieldTransfer::CachedData<Gradient>,
    1469             :                                 Number,
    1470             :                                 libMesh::VectorSetAction<Number>>
    1471         366 :           set_solution(*to_sys, f, nullptr, setter, varvec);
    1472             : 
    1473             :       // We dont look at boundary restriction, not supported for higher order target variables
    1474         366 :       const auto & to_begin = _to_blocks.empty()
    1475         615 :                                   ? to_mesh.active_local_elements_begin()
    1476         366 :                                   : to_mesh.active_local_subdomain_set_elements_begin(_to_blocks);
    1477             : 
    1478         366 :       const auto & to_end = _to_blocks.empty()
    1479         615 :                                 ? to_mesh.active_local_elements_end()
    1480         366 :                                 : to_mesh.active_local_subdomain_set_elements_end(_to_blocks);
    1481             : 
    1482         366 :       ConstElemRange active_local_elem_range(to_begin, to_end);
    1483             : 
    1484         366 :       set_solution.project(active_local_elem_range);
    1485         366 :     }
    1486             :     else
    1487             :     {
    1488     5704132 :       for (const auto & val_pair : dofobject_to_val)
    1489             :       {
    1490     5645424 :         const auto dof_object_id = val_pair.first;
    1491             : 
    1492     5645424 :         const DofObject * dof_object = nullptr;
    1493     5645424 :         if (is_nodal)
    1494     5225394 :           dof_object = to_mesh.node_ptr(dof_object_id);
    1495             :         else
    1496      420030 :           dof_object = to_mesh.elem_ptr(dof_object_id);
    1497             : 
    1498     5645424 :         const auto dof = dof_object->dof_number(sys_num, var_num, 0);
    1499     5645424 :         const auto val = val_pair.second.interp;
    1500             : 
    1501             :         // This will happen if meshes are mismatched
    1502     5645424 :         if (_error_on_miss && GeneralFieldTransfer::isBetterOutOfMeshValue(val))
    1503             :         {
    1504             :           const auto target_location =
    1505           8 :               hasToMultiApp()
    1506           8 :                   ? " on target app " + std::to_string(getGlobalTargetAppIndex(problem_id))
    1507          20 :                   : " on parent app";
    1508           8 :           const auto info_msg = "\nThis check can be turned off by setting 'error_on_miss' to "
    1509             :                                 "false. The 'extrapolation_constant' parameter will be used to set "
    1510             :                                 "the local value at missed points.";
    1511           8 :           if (is_nodal)
    1512           4 :             mooseError("No source value for node ",
    1513             :                        dof_object_id,
    1514             :                        target_location,
    1515             :                        " could be located. Node details:\n",
    1516           4 :                        _to_meshes[problem_id]->nodePtr(dof_object_id)->get_info(),
    1517             :                        "\n",
    1518             :                        info_msg);
    1519             :           else
    1520           4 :             mooseError("No source value for element ",
    1521             :                        dof_object_id,
    1522             :                        target_location,
    1523             :                        " could be located. Element details:\n",
    1524           4 :                        _to_meshes[problem_id]->elemPtr(dof_object_id)->get_info(),
    1525             :                        "\n",
    1526             :                        info_msg);
    1527           0 :         }
    1528             : 
    1529             :         // We should not put garbage into our solution vector
    1530             :         // but it can be that we want to set it to a different value than what was already there
    1531             :         // for example: the source app has been displaced and was sending an indicator of its
    1532             :         // position
    1533     5645416 :         if (GeneralFieldTransfer::isBetterOutOfMeshValue(val))
    1534             :         {
    1535      608586 :           if (!GeneralFieldTransfer::isBetterOutOfMeshValue(_default_extrapolation_value))
    1536      608082 :             to_sys->solution->set(dof, _default_extrapolation_value);
    1537      608586 :           continue;
    1538             :         }
    1539             : 
    1540     5036830 :         to_sys->solution->set(dof, val);
    1541             :       }
    1542             :     }
    1543             : 
    1544       59074 :     to_sys->solution->close();
    1545             :     // Sync local solutions
    1546       59074 :     to_sys->update();
    1547             :   }
    1548       56807 : }
    1549             : 
    1550             : bool
    1551     2027204 : MultiAppGeneralFieldTransfer::acceptPointInOriginMesh(unsigned int i_from,
    1552             :                                                       const std::vector<BoundingBox> & local_bboxes,
    1553             :                                                       const Point & pt,
    1554             :                                                       const unsigned int only_from_mesh_div,
    1555             :                                                       Real & distance) const
    1556             : {
    1557     2027204 :   if (_use_bounding_boxes && !local_bboxes[i_from].contains_point(pt))
    1558     1643277 :     return false;
    1559             :   else
    1560             :   {
    1561      383927 :     auto * pl = _from_point_locators[i_from].get();
    1562      383927 :     const auto from_global_num = getGlobalSourceAppIndex(i_from);
    1563      383927 :     const auto transformed_pt = _from_transforms[from_global_num]->mapBack(pt);
    1564             : 
    1565             :     // Check point against source block restriction
    1566      383927 :     if (!_from_blocks.empty() && !inBlocks(_from_blocks, pl, transformed_pt))
    1567       13937 :       return false;
    1568             : 
    1569             :     // Check point against source boundary restriction. Block restriction will speed up the search
    1570      383725 :     if (!_from_boundaries.empty() &&
    1571           0 :         !onBoundaries(_from_boundaries, _from_blocks, *_from_meshes[i_from], pl, transformed_pt))
    1572           0 :       return false;
    1573             : 
    1574             :     // Check point against the source mesh division
    1575      391205 :     if ((!_from_mesh_divisions.empty() || !_to_mesh_divisions.empty()) &&
    1576        7480 :         !acceptPointMeshDivision(transformed_pt, i_from, only_from_mesh_div))
    1577        4180 :       return false;
    1578             : 
    1579             :     // Get nearest position (often a subapp position) for the target point
    1580             :     // We want values from the child app that is closest to the same position as the target
    1581      379545 :     Point nearest_position_source;
    1582      379545 :     if (_nearest_positions_obj)
    1583             :     {
    1584       16444 :       const bool initial = _fe_problem.getCurrentExecuteOnFlag() == EXEC_INITIAL;
    1585             :       // The search for the nearest position is done in the reference frame
    1586       16444 :       const Point nearest_position = _nearest_positions_obj->getNearestPosition(pt, initial);
    1587       16444 :       nearest_position_source = _nearest_positions_obj->getNearestPosition(
    1588       16444 :           (*_from_transforms[from_global_num])(Point(0, 0, 0)), initial);
    1589             : 
    1590       16444 :       if (!_skip_coordinate_collapsing &&
    1591           0 :           _from_transforms[from_global_num]->hasNonTranslationTransformation())
    1592           0 :         mooseError("Rotation and scaling currently unsupported with nearest positions transfer.");
    1593             : 
    1594             :       // Compute distance to nearest position and nearest position source
    1595       16444 :       const Real distance_to_position_nearest_source = (pt - nearest_position_source).norm();
    1596       16444 :       const Real distance_to_nearest_position = (pt - nearest_position).norm();
    1597             : 
    1598             :       // Source (usually app position) is not closest to the same positions as the target, dont
    1599             :       // send values. We check the distance instead of the positions because if they are the same
    1600             :       // that means there's two equidistant positions and we would want to capture that as a "value
    1601             :       // conflict"
    1602       16444 :       if (!MooseUtils::absoluteFuzzyEqual(distance_to_position_nearest_source,
    1603             :                                           distance_to_nearest_position))
    1604        9395 :         return false;
    1605             : 
    1606             :       // Set the distance as the distance from the nearest position to the target point
    1607        7049 :       distance = distance_to_position_nearest_source;
    1608             :     }
    1609             : 
    1610             :     // Check that the app actually contains the origin point
    1611             :     // We dont need to check if we already found it in a block or a boundary
    1612      370310 :     if (_from_blocks.empty() && _from_boundaries.empty() && _source_app_must_contain_point &&
    1613         160 :         !inMesh(pl, transformed_pt))
    1614         160 :       return false;
    1615             :   }
    1616      369990 :   return true;
    1617             : }
    1618             : 
    1619             : bool
    1620         160 : MultiAppGeneralFieldTransfer::inMesh(const PointLocatorBase * const pl, const Point & point) const
    1621             : {
    1622             :   // Note: we do not take advantage of a potential block restriction of the mesh here. This is
    1623             :   // because we can avoid this routine by calling inBlocks() instead
    1624         160 :   const Elem * elem = (*pl)(point);
    1625         160 :   return (elem != nullptr);
    1626             : }
    1627             : 
    1628             : bool
    1629       53252 : MultiAppGeneralFieldTransfer::inBlocks(const std::set<SubdomainID> & blocks,
    1630             :                                        const Elem * elem) const
    1631             : {
    1632       53252 :   return blocks.find(elem->subdomain_id()) != blocks.end();
    1633             : }
    1634             : 
    1635             : bool
    1636       34200 : MultiAppGeneralFieldTransfer::inBlocks(const std::set<SubdomainID> & blocks,
    1637             :                                        const MooseMesh & /* mesh */,
    1638             :                                        const Elem * elem) const
    1639             : {
    1640       34200 :   return inBlocks(blocks, elem);
    1641             : }
    1642             : 
    1643             : bool
    1644       42573 : MultiAppGeneralFieldTransfer::inBlocks(const std::set<SubdomainID> & blocks,
    1645             :                                        const MooseMesh & mesh,
    1646             :                                        const Node * node) const
    1647             : {
    1648       42573 :   const auto & node_blocks = mesh.getNodeBlockIds(*node);
    1649       42573 :   std::set<SubdomainID> u;
    1650       42573 :   std::set_intersection(blocks.begin(),
    1651             :                         blocks.end(),
    1652             :                         node_blocks.begin(),
    1653             :                         node_blocks.end(),
    1654             :                         std::inserter(u, u.begin()));
    1655       85146 :   return !u.empty();
    1656       42573 : }
    1657             : 
    1658             : bool
    1659       10471 : MultiAppGeneralFieldTransfer::inBlocks(const std::set<SubdomainID> & blocks,
    1660             :                                        const PointLocatorBase * const pl,
    1661             :                                        const Point & point) const
    1662             : {
    1663       10471 :   const Elem * elem = (*pl)(point, &blocks);
    1664       10471 :   return (elem != nullptr);
    1665             : }
    1666             : 
    1667             : bool
    1668      505427 : MultiAppGeneralFieldTransfer::onBoundaries(const std::set<BoundaryID> & boundaries,
    1669             :                                            const MooseMesh & mesh,
    1670             :                                            const Node * node) const
    1671             : {
    1672      505427 :   const BoundaryInfo & bnd_info = mesh.getMesh().get_boundary_info();
    1673      505427 :   std::vector<BoundaryID> vec_to_fill;
    1674      505427 :   bnd_info.boundary_ids(node, vec_to_fill);
    1675      505427 :   std::set<BoundaryID> vec_to_fill_set(vec_to_fill.begin(), vec_to_fill.end());
    1676      505427 :   std::set<BoundaryID> u;
    1677      505427 :   std::set_intersection(boundaries.begin(),
    1678             :                         boundaries.end(),
    1679             :                         vec_to_fill_set.begin(),
    1680             :                         vec_to_fill_set.end(),
    1681             :                         std::inserter(u, u.begin()));
    1682     1010854 :   return !u.empty();
    1683      505427 : }
    1684             : 
    1685             : bool
    1686       77184 : MultiAppGeneralFieldTransfer::onBoundaries(const std::set<BoundaryID> & boundaries,
    1687             :                                            const MooseMesh & mesh,
    1688             :                                            const Elem * elem) const
    1689             : {
    1690             :   // Get all boundaries each side of the element is part of
    1691       77184 :   const BoundaryInfo & bnd_info = mesh.getMesh().get_boundary_info();
    1692       77184 :   std::vector<BoundaryID> vec_to_fill;
    1693       77184 :   std::vector<BoundaryID> vec_to_fill_temp;
    1694       77184 :   if (_elemental_boundary_restriction_on_sides)
    1695      477288 :     for (const auto side : make_range(elem->n_sides()))
    1696             :     {
    1697      409104 :       bnd_info.boundary_ids(elem, side, vec_to_fill_temp);
    1698      409104 :       vec_to_fill.insert(vec_to_fill.end(), vec_to_fill_temp.begin(), vec_to_fill_temp.end());
    1699             :     }
    1700             :   else
    1701       81000 :     for (const auto node_index : make_range(elem->n_nodes()))
    1702             :     {
    1703       72000 :       bnd_info.boundary_ids(elem->node_ptr(node_index), vec_to_fill_temp);
    1704       72000 :       vec_to_fill.insert(vec_to_fill.end(), vec_to_fill_temp.begin(), vec_to_fill_temp.end());
    1705             :     }
    1706       77184 :   std::set<BoundaryID> vec_to_fill_set(vec_to_fill.begin(), vec_to_fill.end());
    1707             : 
    1708             :   // Look for a match between the boundaries from the restriction and those near the element
    1709       77184 :   std::set<BoundaryID> u;
    1710       77184 :   std::set_intersection(boundaries.begin(),
    1711             :                         boundaries.end(),
    1712             :                         vec_to_fill_set.begin(),
    1713             :                         vec_to_fill_set.end(),
    1714             :                         std::inserter(u, u.begin()));
    1715      154368 :   return !u.empty();
    1716       77184 : }
    1717             : 
    1718             : bool
    1719           0 : MultiAppGeneralFieldTransfer::onBoundaries(const std::set<BoundaryID> & boundaries,
    1720             :                                            const std::set<SubdomainID> & block_restriction,
    1721             :                                            const MooseMesh & mesh,
    1722             :                                            const PointLocatorBase * const pl,
    1723             :                                            const Point & point) const
    1724             : {
    1725             :   // Find the element containing the point and use the block restriction if known for speed
    1726             :   const Elem * elem;
    1727           0 :   if (block_restriction.empty())
    1728           0 :     elem = (*pl)(point);
    1729             :   else
    1730           0 :     elem = (*pl)(point, &block_restriction);
    1731             : 
    1732           0 :   if (!elem)
    1733           0 :     return false;
    1734           0 :   return onBoundaries(boundaries, mesh, elem);
    1735             : }
    1736             : 
    1737             : bool
    1738        7480 : MultiAppGeneralFieldTransfer::acceptPointMeshDivision(
    1739             :     const Point & pt, const unsigned int i_local, const unsigned int only_from_this_mesh_div) const
    1740             : {
    1741             :   // This routine can also be called to examine if the to_mesh_division index matches the current
    1742             :   // source subapp index
    1743        7480 :   unsigned int source_mesh_div = MooseMeshDivision::INVALID_DIVISION_INDEX - 1;
    1744        7480 :   if (!_from_mesh_divisions.empty())
    1745        5318 :     source_mesh_div = _from_mesh_divisions[i_local]->divisionIndex(pt);
    1746             : 
    1747             :   // If the point is not indexed in the source division
    1748        7480 :   if (!_from_mesh_divisions.empty() && source_mesh_div == MooseMeshDivision::INVALID_DIVISION_INDEX)
    1749        3190 :     return false;
    1750             :   // If the point is not the at the same index in the target and the origin meshes, reject
    1751        4290 :   else if ((_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX ||
    1752        4290 :             _to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_DIVISION_INDEX) &&
    1753             :            source_mesh_div != only_from_this_mesh_div)
    1754           0 :     return false;
    1755             :   // If the point is at a certain division index that is not the same as the index of the subapp
    1756             :   // we wanted the information to be from for that point, reject
    1757        4290 :   else if (_from_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX &&
    1758             :            source_mesh_div != only_from_this_mesh_div)
    1759         504 :     return false;
    1760        4497 :   else if (_to_mesh_division_behavior == MeshDivisionTransferUse::MATCH_SUBAPP_INDEX &&
    1761         711 :            only_from_this_mesh_div != getGlobalSourceAppIndex(i_local))
    1762         486 :     return false;
    1763             :   else
    1764        3300 :     return true;
    1765             : }
    1766             : 
    1767             : bool
    1768      153475 : MultiAppGeneralFieldTransfer::closestToPosition(unsigned int pos_index, const Point & pt) const
    1769             : {
    1770             :   mooseAssert(_nearest_positions_obj, "Should not be here without a positions object");
    1771      153475 :   if (!_skip_coordinate_collapsing)
    1772           0 :     paramError("skip_coordinate_collapsing", "Coordinate collapsing not implemented");
    1773      153475 :   bool initial = _fe_problem.getCurrentExecuteOnFlag() == EXEC_INITIAL;
    1774      153475 :   if (!_search_value_conflicts)
    1775             :     // Faster to just compare the index
    1776       19227 :     return pos_index == _nearest_positions_obj->getNearestPositionIndex(pt, initial);
    1777             :   else
    1778             :   {
    1779             :     // Get the distance to the position and see if we are missing a value just because the position
    1780             :     // is not officially the closest, but it is actually at the same distance
    1781      134248 :     const auto nearest_position = _nearest_positions_obj->getNearestPosition(pt, initial);
    1782      134248 :     const auto nearest_position_at_index = _nearest_positions_obj->getPosition(pos_index, initial);
    1783      134248 :     Real distance_to_position_at_index = (pt - nearest_position_at_index).norm();
    1784      134248 :     const Real distance_to_nearest_position = (pt - nearest_position).norm();
    1785             : 
    1786      134248 :     if (!MooseUtils::absoluteFuzzyEqual(distance_to_position_at_index,
    1787             :                                         distance_to_nearest_position))
    1788       89553 :       return false;
    1789             :     // Actually the same position (point)
    1790       44695 :     else if (nearest_position == nearest_position_at_index)
    1791       44695 :       return true;
    1792             :     else
    1793             :     {
    1794           0 :       mooseWarning("Two equidistant positions ",
    1795             :                    nearest_position,
    1796             :                    " and ",
    1797             :                    nearest_position_at_index,
    1798             :                    " detected near point ",
    1799             :                    pt);
    1800           0 :       return true;
    1801             :     }
    1802             :   }
    1803             : }
    1804             : 
    1805             : Real
    1806     8592857 : MultiAppGeneralFieldTransfer::bboxMaxDistance(const Point & p, const BoundingBox & bbox) const
    1807             : {
    1808     8592857 :   std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
    1809             : 
    1810     8592857 :   std::array<Point, 8> all_points;
    1811    25778571 :   for (unsigned int x = 0; x < 2; x++)
    1812    51557142 :     for (unsigned int y = 0; y < 2; y++)
    1813   103114284 :       for (unsigned int z = 0; z < 2; z++)
    1814    68742856 :         all_points[x + 2 * y + 4 * z] =
    1815   137485712 :             Point(source_points[x](0), source_points[y](1), source_points[z](2));
    1816             : 
    1817     8592857 :   Real max_distance = 0.;
    1818             : 
    1819    77335713 :   for (unsigned int i = 0; i < 8; i++)
    1820             :   {
    1821    68742856 :     Real distance = (p - all_points[i]).norm();
    1822    68742856 :     if (distance > max_distance)
    1823    18293854 :       max_distance = distance;
    1824             :   }
    1825             : 
    1826     8592857 :   return max_distance;
    1827             : }
    1828             : 
    1829             : Real
    1830     8592857 : MultiAppGeneralFieldTransfer::bboxMinDistance(const Point & p, const BoundingBox & bbox) const
    1831             : {
    1832     8592857 :   std::array<Point, 2> source_points = {{bbox.first, bbox.second}};
    1833             : 
    1834     8592857 :   std::array<Point, 8> all_points;
    1835    25778571 :   for (unsigned int x = 0; x < 2; x++)
    1836    51557142 :     for (unsigned int y = 0; y < 2; y++)
    1837   103114284 :       for (unsigned int z = 0; z < 2; z++)
    1838    68742856 :         all_points[x + 2 * y + 4 * z] =
    1839   137485712 :             Point(source_points[x](0), source_points[y](1), source_points[z](2));
    1840             : 
    1841     8592857 :   Real min_distance = std::numeric_limits<Real>::max();
    1842             : 
    1843    77335713 :   for (unsigned int i = 0; i < 8; i++)
    1844             :   {
    1845    68742856 :     Real distance = (p - all_points[i]).norm();
    1846    68742856 :     if (distance < min_distance)
    1847    18958532 :       min_distance = distance;
    1848             :   }
    1849             : 
    1850     8592857 :   return min_distance;
    1851             : }
    1852             : 
    1853             : std::vector<BoundingBox>
    1854       56708 : MultiAppGeneralFieldTransfer::getRestrictedFromBoundingBoxes() const
    1855             : {
    1856       56708 :   std::vector<std::pair<Point, Point>> bb_points(_from_meshes.size());
    1857       56708 :   const Real min_r = std::numeric_limits<Real>::lowest();
    1858       56708 :   const Real max_r = std::numeric_limits<Real>::max();
    1859             : 
    1860      115711 :   for (const auto j : make_range(_from_meshes.size()))
    1861             :   {
    1862       59003 :     Point min(max_r, max_r, max_r);
    1863       59003 :     Point max(min_r, min_r, min_r);
    1864       59003 :     bool at_least_one = false;
    1865       59003 :     const auto & from_mesh = _from_problems[j]->mesh(_displaced_source_mesh);
    1866             : 
    1867      118006 :     for (const auto & elem : as_range(from_mesh.getMesh().local_elements_begin(),
    1868     5288215 :                                       from_mesh.getMesh().local_elements_end()))
    1869             :     {
    1870     5111206 :       if (!_from_blocks.empty() && !inBlocks(_from_blocks, from_mesh, elem))
    1871       19296 :         continue;
    1872             : 
    1873    28942466 :       for (const auto & node : elem->node_ref_range())
    1874             :       {
    1875    23850556 :         if (!_from_boundaries.empty() && !onBoundaries(_from_boundaries, from_mesh, &node))
    1876      181768 :           continue;
    1877             : 
    1878    23668788 :         at_least_one = true;
    1879    94675152 :         for (const auto i : make_range(LIBMESH_DIM))
    1880             :         {
    1881    71006364 :           min(i) = std::min(min(i), node(i));
    1882    71006364 :           max(i) = std::max(max(i), node(i));
    1883             :         }
    1884             :       }
    1885       59003 :     }
    1886             : 
    1887             :     // For 2D RZ problems, we need to amend the bounding box to cover the whole XYZ projection
    1888             :     // - The XYZ-Y axis is assumed aligned with the RZ-Z axis
    1889             :     // - RZ systems also cover negative coordinates hence the use of the maximum R
    1890             :     // NOTE: We will only support the case where there is only one coordinate system
    1891       59003 :     if ((from_mesh.getUniqueCoordSystem() == Moose::COORD_RZ) && (LIBMESH_DIM == 3))
    1892             :     {
    1893         132 :       min(0) = -max(0);
    1894         132 :       min(2) = -max(0);
    1895         132 :       max(2) = max(0);
    1896             :     }
    1897             : 
    1898       59003 :     BoundingBox bbox(min, max);
    1899       59003 :     if (!at_least_one)
    1900          31 :       bbox.min() = max; // If we didn't hit any nodes, this will be _the_ minimum bbox
    1901             :     else
    1902             :     {
    1903             :       // Translate the bounding box to the from domain's position. We may have rotations so we
    1904             :       // must be careful in constructing the new min and max (first and second)
    1905       58972 :       const auto from_global_num = getGlobalSourceAppIndex(j);
    1906       58972 :       transformBoundingBox(bbox, *_from_transforms[from_global_num]);
    1907             :     }
    1908             : 
    1909             :     // Cast the bounding box into a pair of points (so it can be put through
    1910             :     // MPI communication).
    1911       59003 :     bb_points[j] = static_cast<std::pair<Point, Point>>(bbox);
    1912             :   }
    1913             : 
    1914             :   // Serialize the bounding box points.
    1915       56708 :   _communicator.allgather(bb_points);
    1916             : 
    1917             :   // Recast the points back into bounding boxes and return.
    1918       56708 :   std::vector<BoundingBox> bboxes(bb_points.size());
    1919      144545 :   for (const auto i : make_range(bb_points.size()))
    1920       87837 :     bboxes[i] = static_cast<BoundingBox>(bb_points[i]);
    1921             : 
    1922             :   // TODO move up
    1923             :   // Check for a user-set fixed bounding box size and modify the sizes as appropriate
    1924       56708 :   if (_fixed_bbox_size != std::vector<Real>(3, 0))
    1925        1380 :     for (const auto i : make_range(LIBMESH_DIM))
    1926        1035 :       if (!MooseUtils::absoluteFuzzyEqual(_fixed_bbox_size[i], 0))
    1927        2142 :         for (const auto j : make_range(bboxes.size()))
    1928             :         {
    1929        1331 :           const auto current_width = (bboxes[j].second - bboxes[j].first)(i);
    1930        1331 :           bboxes[j].first(i) -= (_fixed_bbox_size[i] - current_width) / 2;
    1931        1331 :           bboxes[j].second(i) += (_fixed_bbox_size[i] - current_width) / 2;
    1932             :         }
    1933             : 
    1934      113416 :   return bboxes;
    1935       56708 : }
    1936             : 
    1937             : std::vector<unsigned int>
    1938       56708 : MultiAppGeneralFieldTransfer::getGlobalStartAppPerProc() const
    1939             : {
    1940       56708 :   std::vector<unsigned int> global_app_start_per_proc(1, -1);
    1941       56708 :   if (_from_local2global_map.size())
    1942       29705 :     global_app_start_per_proc[0] = _from_local2global_map[0];
    1943       56708 :   _communicator.allgather(global_app_start_per_proc, true);
    1944       56708 :   return global_app_start_per_proc;
    1945           0 : }
    1946             : 
    1947             : VariableName
    1948       60492 : MultiAppGeneralFieldTransfer::getFromVarName(unsigned int var_index)
    1949             : {
    1950             :   mooseAssert(var_index < _from_var_names.size(), "No source variable at this index");
    1951       60492 :   VariableName var_name = _from_var_names[var_index];
    1952       60492 :   if (_from_var_components.size())
    1953         312 :     var_name += "_" + std::to_string(_from_var_components[var_index]);
    1954       60492 :   return var_name;
    1955           0 : }
    1956             : 
    1957             : VariableName
    1958      114091 : MultiAppGeneralFieldTransfer::getToVarName(unsigned int var_index)
    1959             : {
    1960             :   mooseAssert(var_index < _to_var_names.size(), "No target variable at this index");
    1961      114091 :   VariableName var_name = _to_var_names[var_index];
    1962      114091 :   if (_to_var_components.size())
    1963         384 :     var_name += "_" + std::to_string(_to_var_components[var_index]);
    1964      114091 :   return var_name;
    1965           0 : }
    1966             : 
    1967             : Point
    1968      113692 : MultiAppGeneralFieldTransfer::getMaxToProblemsBBoxDimensions() const
    1969             : {
    1970             :   Point max_dimension = {std::numeric_limits<Real>::min(),
    1971             :                          std::numeric_limits<Real>::min(),
    1972      113692 :                          std::numeric_limits<Real>::min()};
    1973             : 
    1974      231918 :   for (const auto & to_mesh : _to_meshes)
    1975             :   {
    1976      118226 :     const auto bbox = to_mesh->getInflatedProcessorBoundingBox();
    1977      472904 :     for (const auto dim : make_range(LIBMESH_DIM))
    1978      354678 :       max_dimension(dim) = std::max(
    1979      709356 :           max_dimension(dim), std::max(std::abs(bbox.first(dim)), std::abs(bbox.second(dim))));
    1980             :   }
    1981             : 
    1982      113692 :   return max_dimension;
    1983             : }
    1984             : 
    1985             : bool
    1986     2302136 : MultiAppGeneralFieldTransfer::detectConflict(Real current_value,
    1987             :                                              Real new_value,
    1988             :                                              Real current_distance,
    1989             :                                              Real new_distance) const
    1990             : {
    1991             :   // No conflict if we're not looking for them
    1992     2302136 :   if (_search_value_conflicts)
    1993             :     // Only consider conflicts if the values are valid and different
    1994        3680 :     if (current_value != GeneralFieldTransfer::BetterOutOfMeshValue &&
    1995        5654 :         new_value != GeneralFieldTransfer::BetterOutOfMeshValue &&
    1996        1974 :         !MooseUtils::absoluteFuzzyEqual(current_value, new_value))
    1997             :       // Conflict only occurs if the origin points are equidistant
    1998        1767 :       if (MooseUtils::absoluteFuzzyEqual(current_distance, new_distance))
    1999          86 :         return true;
    2000     2302050 :   return false;
    2001             : }

Generated by: LCOV version 1.14