ASPECT
python_helper.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2025 - 2026 by the authors of the ASPECT code.
3 
4  This file is part of ASPECT.
5 
6  ASPECT is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2, or (at your option)
9  any later version.
10 
11  ASPECT is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with ASPECT; see the file LICENSE. If not see
18  <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #ifndef _aspect_python_helper_h
23 #define _aspect_python_helper_h
24 
25 #include <aspect/global.h>
26 #include <deal.II/base/array_view.h>
27 
28 
29 #ifdef ASPECT_WITH_PYTHON
30 
31 // Python does not like it if this macro is already defined. This happens at
32 // least in some versions of Trilinos and can trigger only with certain unity
33 // build options:
34 #ifdef HAVE_SYS_TIME_H
35 # undef HAVE_SYS_TIME_H
36 #endif
37 #define PY_SSIZE_T_CLEAN
38 #include <Python.h>
39 
40 // Declare that we are compatible with numpy 1.7 API and don't want any
41 // compile warnings for deprecated API calls:
42 #define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION
43 
44 // Define a unique symbol name for the numpy Python API for used when building
45 // the ASPECT executable. This will be used by numpy to declare a static member
46 // variable that will be initialized using the import_array() function (see below).
47 #ifndef PY_ARRAY_UNIQUE_SYMBOL
48 # define PY_ARRAY_UNIQUE_SYMBOL ASPECT_ARRAY_API
49 #endif
50 
51 // NumPY API requires a call to import_array() to be made before usage
52 // but only once per executable. We will define
53 // ASPECT_NUMPY_DEFINE_API in exactly one translation unit per
54 // executable/shared library (in main.cc and any plugin that uses
55 // Python). All other translation units including this header will
56 // therefore define NO_IMPORT_ARRAY that asks numpy to skip defining
57 // import_array().
58 #ifndef ASPECT_NUMPY_DEFINE_API
59 # define NO_IMPORT_ARRAY
60 #endif
61 
62 #include <numpy/arrayobject.h>
63 
64 // Clean up any of the macros defined above. This is important if we
65 // are doing a unity build:
66 #ifdef ASPECT_NUMPY_DEFINE_API
67 # undef ASPECT_NUMPY_DEFINE_API
68 #endif
69 #ifdef NO_IMPORT_ARRAY
70 # undef NO_IMPORT_ARRAY
71 #endif
72 
73 namespace aspect
74 {
75  namespace PythonHelper
76  {
85  inline
86  std::unique_ptr<PyObject, void(*)(PyObject *)> vector_to_numpy_object(const std::vector<double> &vec)
87  {
88  const npy_intp size = static_cast<npy_intp>(vec.size());
89 
90  auto deleter = [](PyObject *obj)
91  {
92  Py_DECREF(obj);
93  };
94  if (size == 0)
95  return std::unique_ptr<PyObject, void(*)(PyObject *)>(PyArray_SimpleNew(1, &size, NPY_DOUBLE), deleter);
96  double *data = const_cast<double *>(vec.data());
97  return std::unique_ptr<PyObject, void(*)(PyObject *)>(PyArray_SimpleNewFromData(1, &size, NPY_DOUBLE, data), deleter);
98  }
99 
100 
101 
106  inline ::ArrayView<double> numpy_to_array_view(PyObject *obj)
107  {
108  AssertThrow(PyArray_Check(obj), ::ExcMessage("Expected a numpy array"));
109  PyArrayObject *arr = reinterpret_cast<PyArrayObject *>(obj);
110  return ::ArrayView<double>(static_cast<double *>(PyArray_DATA(arr)),
111  static_cast<size_t>(PyArray_SIZE(arr)));
112  }
113 
114  }
115 }
116 #endif
117 
118 
119 #endif