Skip to content

add PythonHelper and tests#6847

Open
tjhei wants to merge 2 commits intogeodynamics:mainfrom
tjhei:python-helper
Open

add PythonHelper and tests#6847
tjhei wants to merge 2 commits intogeodynamics:mainfrom
tjhei:python-helper

Conversation

@tjhei
Copy link
Member

@tjhei tjhei commented Feb 2, 2026

Python and the numpy C API make it very difficult to use in several translation units due to the import_array() facility that internally defines some static functions/members. It took me some time to get this right, but it now works correctly for me. I decided to move the shared functionality into a new header python_helper.h.
@bangerth will hopefully like the idea that we can put helper functions there and use them from several places.

FYI @danieldouglas92

@tjhei
Copy link
Member Author

tjhei commented Feb 2, 2026

Copy link
Contributor

@danieldouglas92 danieldouglas92 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This looks great to me!

@danieldouglas92
Copy link
Contributor

I'm not sure why the periodic_box_mpi4 test would be failing with these changes

Copy link
Member

@gassmoeller gassmoeller left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I like the direction to collect python parts in python_helper.h. I have a few questions, maybe I am just not understanding the purpose correctly.


#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION

//
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This needs a comment. What is the purpose of PY_ARRAY_UNIQUE_SYMBOL, and what is ASPECT_ARRAY_API?

// ASPECT_NUMPY_DEFINE_API in exactly one translation unit per
// executable/shared library (in main.cc and any plugin that uses
// Python). All other translation units including this header will
// therefor define NO_IMPORT_ARRAY that asks numpy to skip defining
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
// therefor define NO_IMPORT_ARRAY that asks numpy to skip defining
// therefore define NO_IMPORT_ARRAY that asks numpy to skip defining

#endif

#define ASPECT_NUMPY_DEFINE_API
#include <aspect/python_helper.h>
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

does ASPECT_NUMPY_DEFINE_API need to be undefined after the include (for unity builds where other translation units may be combined with main.cc)?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good point, thank you! I am now #undefing both

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you still need the #define here, though?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes! in exactly one translation unit I need to define this macro (see description in the header file)

return PyArray_SimpleNewFromData(1, &size, NPY_DOUBLE, data);
}
#define ASPECT_NUMPY_DEFINE_API
#include <aspect/python_helper.h>
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same question here, shouldnt we undefine ASPECT_NUMPY_DEFINE_API after the include?

@gassmoeller
Copy link
Member

Hm, yes and I dont know why the test fails either. Maybe it is unrelated to this PR, we could have picked up that test failure in the past weeks while the deal.II development tester didnt work for us.

@gassmoeller
Copy link
Member

I checked and periodic_box_mpi4 also fails on ASPECTs main branch with deal.II master. So there must be something going on unrelated to this PR.

@tjhei
Copy link
Member Author

tjhei commented Feb 6, 2026

I checked and periodic_box_mpi4 also fails on ASPECTs main branch with deal.II master. So there must be something going on unrelated to this PR.

Is this a fallout from the "ghosted vector" changes?

@tjhei
Copy link
Member Author

tjhei commented Feb 6, 2026

I think I addressed all your points.

Comment on lines 66 to 71
#ifdef ASPECT_NUMPY_DEFINE_API
#undef ASPECT_NUMPY_DEFINE_API
#endif
#ifdef NO_IMPORT_ARRAY
#undef NO_IMPORT_ARRAY
#endif
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
#ifdef ASPECT_NUMPY_DEFINE_API
#undef ASPECT_NUMPY_DEFINE_API
#endif
#ifdef NO_IMPORT_ARRAY
#undef NO_IMPORT_ARRAY
#endif
#ifdef ASPECT_NUMPY_DEFINE_API
# undef ASPECT_NUMPY_DEFINE_API
#endif
#ifdef NO_IMPORT_ARRAY
# undef NO_IMPORT_ARRAY
#endif

Comment on lines 75 to 86
/**
* Convert a std::vector<double> to a numpy array.
* @param[in] vec The vector to convert.
* @return A PyObject* to the numpy array.
*/
inline PyObject *vector_to_numpy_object(const std::vector<double> &vec)
{
npy_intp size = static_cast<npy_intp>(vec.size());
if (size == 0)
return PyArray_SimpleNew(1, &size, NPY_DOUBLE);
double *data = const_cast<double *>(vec.data());
return PyArray_SimpleNewFromData(1, &size, NPY_DOUBLE, data);
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could you comment on whether the object has to be explicitly destroyed? I imagine so, and perhaps in that case it would be useful to return not just a plain pointer, but a std::unique_ptr<PyObject, void(*)(PyObject*)> where when you create the object, you add a lambda that calls the destructor on the PyObject.

#endif

#define ASPECT_NUMPY_DEFINE_API
#include <aspect/python_helper.h>
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you still need the #define here, though?

@tjhei
Copy link
Member Author

tjhei commented Feb 10, 2026

I am not a big fan of the unique_ptr approach, but it works.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants