1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17
18
19
20 #include "libmesh/libmesh_config.h"
21
22 #ifdef LIBMESH_TRILINOS_HAVE_DTK
23
24 #include "libmesh/dtk_solution_transfer.h"
25 #include "libmesh/system.h"
26
27 #include "libmesh/ignore_warnings.h"
28
29 // Trilinos Includes
30 #include <Teuchos_RCP.hpp>
31 #include <Teuchos_GlobalMPISession.hpp>
32 #include <Teuchos_Ptr.hpp>
33 #include <Teuchos_DefaultMpiComm.hpp>
34 #include <Teuchos_OpaqueWrapper.hpp>
35
36 // DTK Includes
37 #include <DTK_MeshManager.hpp>
38 #include <DTK_MeshContainer.hpp>
39 #include <DTK_FieldManager.hpp>
40 #include <DTK_FieldContainer.hpp>
41 #include <DTK_FieldTools.hpp>
42 #include <DTK_CommTools.hpp>
43 #include <DTK_CommIndexer.hpp>
44
45 #include "libmesh/restore_warnings.h"
46
47 namespace libMesh
48 {
49
DTKSolutionTransfer(const libMesh::Parallel::Communicator & comm)50 DTKSolutionTransfer::DTKSolutionTransfer(const libMesh::Parallel::Communicator & comm) :
51 SolutionTransfer(comm)
52 {
53 //comm_default = Teuchos::DefaultComm<int>::getComm();
54 comm_default = Teuchos::rcp(new Teuchos::MpiComm<int>(Teuchos::rcp(new Teuchos::OpaqueWrapper<MPI_Comm>(comm.get()))));
55 }
56
~DTKSolutionTransfer()57 DTKSolutionTransfer::~DTKSolutionTransfer()
58 {
59 for (auto & pr : adapters)
60 delete pr.second;
61
62 for (auto & pr : dtk_maps)
63 delete pr.second;
64 }
65
66 void
transfer(const Variable & from_var,const Variable & to_var)67 DTKSolutionTransfer::transfer(const Variable & from_var,
68 const Variable & to_var)
69 {
70 libmesh_experimental();
71
72 EquationSystems * from_es = &from_var.system()->get_equation_systems();
73 EquationSystems * to_es = &to_var.system()->get_equation_systems();
74
75 // Possibly make an Adapter for from_es
76 if (adapters.find(from_es) == adapters.end())
77 adapters[from_es] = new DTKAdapter(comm_default, *from_es);
78
79 // Possibly make an Adapter for to_es
80 if (adapters.find(to_es) == adapters.end())
81 adapters[to_es] = new DTKAdapter(comm_default, *to_es);
82
83 DTKAdapter * from_adapter = adapters[from_es];
84 DTKAdapter * to_adapter = adapters[to_es];
85
86 std::pair<EquationSystems *, EquationSystems *> from_to(from_es, to_es);
87
88 // If we haven't created a map for this pair of EquationSystems yet, do it now
89 if (dtk_maps.find(from_to) == dtk_maps.end())
90 {
91 libmesh_assert(from_es->get_mesh().mesh_dimension() == to_es->get_mesh().mesh_dimension());
92
93 shared_domain_map_type * map = new shared_domain_map_type(comm_default, from_es->get_mesh().mesh_dimension(), true);
94 dtk_maps[from_to] = map;
95
96 // The tolerance here is for the "contains_point()" implementation in DTK. Set a larger value for a looser tolerance...
97 map->setup(from_adapter->get_mesh_manager(), to_adapter->get_target_coords(), 30*Teuchos::ScalarTraits<double>::eps());
98 }
99
100 DTKAdapter::RCP_Evaluator from_evaluator = from_adapter->get_variable_evaluator(from_var.name());
101 Teuchos::RCP<DataTransferKit::FieldManager<DTKAdapter::FieldContainerType>> to_values = to_adapter->get_values_to_fill(to_var.name());
102
103 dtk_maps[from_to]->apply(from_evaluator, to_values);
104
105 if (dtk_maps[from_to]->getMissedTargetPoints().size())
106 libMesh::out<<"Warning: Some points were missed in the transfer of "<<from_var.name()<<" to "<<to_var.name()<<"!"<<std::endl;
107
108 to_adapter->update_variable_values(to_var.name());
109 }
110
111 } // namespace libMesh
112
113 #endif // #ifdef LIBMESH_TRILINOS_HAVE_DTK
114