MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_MatlabUtils_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // MueLu: A package for multigrid based preconditioning
4 //
5 // Copyright 2012 NTESS and the MueLu contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef MUELU_MATLABUTILS_DECL_HPP
11 #define MUELU_MATLABUTILS_DECL_HPP
12 
13 #include "MueLu_ConfigDefs.hpp"
14 
15 #if !defined(HAVE_MUELU_MATLAB) || !defined(HAVE_MUELU_TPETRA)
16 #error "Muemex requires MATLAB, Epetra and Tpetra."
17 #else
18 
19 // Matlab fwd style declarations
20 struct mxArray_tag;
21 typedef struct mxArray_tag mxArray;
22 typedef size_t mwIndex;
23 
24 #include <string>
25 #include <complex>
26 #include <stdexcept>
28 #include <Teuchos_RCP.hpp>
29 #include <Teuchos_DefaultComm.hpp>
30 #include "MueLu_Factory.hpp"
31 #include "MueLu_Hierarchy_decl.hpp"
34 #include "MueLu_Utilities_decl.hpp"
35 #include "MueLu_Graph_fwd.hpp"
36 #ifdef HAVE_MUELU_EPETRA
37 #include "Epetra_MultiVector.h"
38 #include "Epetra_CrsMatrix.h"
39 #endif
41 #ifdef HAVE_MUELU_EPETRA
43 #endif
44 #include "Xpetra_MapFactory.hpp"
45 #include "Xpetra_CrsGraph.hpp"
46 #include "Xpetra_VectorFactory.hpp"
47 #include <Tpetra_Core.hpp>
48 
49 #include "Kokkos_DynRankView.hpp"
50 
51 namespace MueLu {
52 
53 enum MuemexType {
54  INT,
69 #ifdef HAVE_MUELU_EPETRA
72 #endif
76 #ifdef HAVE_MUELU_INTREPID2
77  ,
78  FIELDCONTAINER_ORDINAL
79 #endif
80 };
81 
82 typedef Tpetra::KokkosCompat::KokkosDeviceWrapperNode<Kokkos::Serial, Kokkos::HostSpace> mm_node_t;
83 typedef typename Tpetra::Map<>::local_ordinal_type mm_LocalOrd; // these are used for LocalOrdinal and GlobalOrdinal of all xpetra/tpetra templated types
85 typedef std::complex<double> complex_t;
103 
104 #ifdef HAVE_MUELU_INTREPID2
105 typedef Kokkos::DynRankView<mm_LocalOrd, typename mm_node_t::device_type> FieldContainer_ordinal;
106 #endif
107 
108 class MuemexArg {
109  public:
110  MuemexArg(MuemexType dataType) { type = dataType; }
112 };
113 
114 template <typename T>
115 MuemexType getMuemexType(const T& data);
116 
117 template <typename T>
118 class MuemexData : public MuemexArg {
119  public:
120  MuemexData(T& data); // Construct from pre-existing data, to pass to MATLAB.
121  MuemexData(T& data, MuemexType type); // Construct from pre-existing data, to pass to MATLAB.
122  MuemexData(const mxArray* mxa); // Construct from MATLAB array, to get from MATLAB.
123  mxArray* convertToMatlab(); // Create a MATLAB object and copy this data to it
124  T& getData(); // Set and get methods
125  void setData(T& data);
126 
127  private:
128  T data;
129 };
130 
131 template <typename T>
132 MuemexType getMuemexType(const T& data);
133 
134 template <typename T>
136 
137 template <typename T>
138 T loadDataFromMatlab(const mxArray* mxa);
139 
140 template <typename T>
141 mxArray* saveDataToMatlab(T& data);
142 
143 // Add data to level. Set the keep flag on the data to "user-provided" so it's not deleted.
144 template <typename T>
145 void addLevelVariable(const T& data, std::string& name, Level& lvl, const FactoryBase* fact = NoFactory::get());
146 
147 template <typename T>
148 const T& getLevelVariable(std::string& name, Level& lvl);
149 
150 // Functions used to put data through matlab factories - first arg is "this" pointer of matlab factory
151 template <typename Scalar = double, typename LocalOrdinal = mm_LocalOrd, typename GlobalOrdinal = mm_GlobalOrd, typename Node = mm_node_t>
152 std::vector<Teuchos::RCP<MuemexArg> > processNeeds(const Factory* factory, std::string& needsParam, Level& lvl);
153 
154 template <typename Scalar = double, typename LocalOrdinal = mm_LocalOrd, typename GlobalOrdinal = mm_GlobalOrd, typename Node = mm_node_t>
155 void processProvides(std::vector<Teuchos::RCP<MuemexArg> >& mexOutput, const Factory* factory, std::string& providesParam, Level& lvl);
156 
157 // create a sparse array in Matlab
158 template <typename Scalar>
159 mxArray* createMatlabSparse(int numRows, int numCols, int nnz);
160 template <typename Scalar>
161 mxArray* createMatlabMultiVector(int numRows, int numCols);
162 template <typename Scalar>
163 void fillMatlabArray(Scalar* array, const mxArray* mxa, int n);
164 int* mwIndex_to_int(int N, mwIndex* mwi_array);
165 bool isValidMatlabAggregates(const mxArray* mxa);
166 bool isValidMatlabGraph(const mxArray* mxa);
167 std::vector<std::string> tokenizeList(const std::string& param);
168 // The two callback functions that MueLu can call to run anything in MATLAB
169 void callMatlabNoArgs(std::string function);
170 std::vector<Teuchos::RCP<MuemexArg> > callMatlab(std::string function, int numOutputs, std::vector<Teuchos::RCP<MuemexArg> > args);
173 
174 // trim from start
175 static inline std::string& ltrim(std::string& s) {
176  s.erase(0, s.find_first_not_of(" "));
177  return s;
178 }
179 
180 // trim from end
181 static inline std::string& rtrim(std::string& s) {
182  s.erase(s.find_last_not_of(" "), std::string::npos);
183  return s;
184 }
185 
186 // trim from both ends
187 static inline std::string& trim(std::string& s) {
188  return ltrim(rtrim(s));
189 }
190 
191 } // namespace MueLu
192 
193 #endif // HAVE_MUELU_MATLAB error handler
194 #endif // MUELU_MATLABUTILS_DECL_HPP guard
Xpetra::MultiVector< double, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_MultiVector_double
bool isValidMatlabAggregates(const mxArray *mxa)
template mxArray * saveDataToMatlab(bool &data)
LocalOrdinal local_ordinal_type
std::vector< std::string > tokenizeList(const std::string &params)
Xpetra::MultiVector< complex_t, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_MultiVector_complex
Xpetra::CrsGraph< mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_CrsGraph
GlobalOrdinal global_ordinal_type
Container class for aggregation information.
Tpetra::CrsMatrix< complex_t, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Tpetra_CrsMatrix_complex
Tpetra::Map::local_ordinal_type mm_LocalOrd
Tpetra::Map::global_ordinal_type mm_GlobalOrd
std::vector< RCP< MuemexArg > > callMatlab(std::string function, int numOutputs, std::vector< RCP< MuemexArg > > args)
void fillMatlabArray(Scalar *array, const mxArray *mxa, int n)
T loadDataFromMatlab(const mxArray *mxa)
Tpetra::KokkosCompat::KokkosDeviceWrapperNode< Kokkos::Serial, Kokkos::HostSpace > mm_node_t
bool isValidMatlabGraph(const mxArray *mxa)
MueLu::Hierarchy< complex_t, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Hierarchy_complex
MuemexType getMuemexType(const T &data)
std::vector< Teuchos::RCP< MuemexArg > > processNeeds(const Factory *factory, std::string &needsParam, Level &lvl)
void addLevelVariable(const T &data, std::string &name, Level &lvl, const FactoryBase *fact=NoFactory::get())
static const NoFactory * get()
Teuchos::RCP< Teuchos::ParameterList > getInputParamList()
struct mxArray_tag mxArray
Base class for factories (e.g., R, P, and A_coarse).
Xpetra::Vector< mm_LocalOrd, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_ordinal_vector
MuemexArg(MuemexType dataType)
MueLu::AmalgamationInfo< mm_LocalOrd, mm_GlobalOrd, mm_node_t > MAmalInfo
MueLu::DefaultScalar Scalar
int * mwIndex_to_int(int N, mwIndex *mwi_array)
Class that holds all level-specific information.
Definition: MueLu_Level.hpp:63
MueLu::Hierarchy< double, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Hierarchy_double
Xpetra::Matrix< double, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_Matrix_double
const T & getLevelVariable(std::string &name, Level &lvl)
MueLu::LWGraph< mm_LocalOrd, mm_GlobalOrd, mm_node_t > MGraph
Xpetra::Matrix< complex_t, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_Matrix_complex
static std::string & trim(std::string &s)
void callMatlabNoArgs(std::string function)
static std::string & ltrim(std::string &s)
std::complex< double > complex_t
static std::string & rtrim(std::string &s)
MueLu::Aggregates< mm_LocalOrd, mm_GlobalOrd, mm_node_t > MAggregates
Tpetra::Map muemex_map_type
Tpetra::MultiVector< complex_t, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Tpetra_MultiVector_complex
Lightweight MueLu representation of a compressed row storage graph.
Xpetra::Map< mm_LocalOrd, mm_GlobalOrd, mm_node_t > Xpetra_map
mxArray * createMatlabSparse(int numRows, int numCols, int nnz)
Tpetra::MultiVector< double, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Tpetra_MultiVector_double
void processProvides(std::vector< Teuchos::RCP< MuemexArg > > &mexOutput, const Factory *factory, std::string &providesParam, Level &lvl)
int mwIndex
minimal container class for storing amalgamation information
mxArray * createMatlabMultiVector(int numRows, int numCols)
Tpetra::CrsMatrix< double, mm_LocalOrd, mm_GlobalOrd, mm_node_t > Tpetra_CrsMatrix_double
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
Teuchos::RCP< MuemexArg > convertMatlabVar(const mxArray *mxa)