Five days of hard work

galpropcc: d9829b55373b3520e10653ab27a90cd49835ce44

.gitignore

+# =====================================
+# C++
+# =====================================
+# Prerequisites
+*.d
+
+# Compiled Object files
+*.slo
+*.lo
+*.o
+*.obj
+
+# Precompiled Headers
+*.gch
+*.pch
+
+# Compiled Dynamic libraries
+*.so
+*.dylib
+*.dll
+
+# Fortran module files
+*.mod
+*.smod
+
+# Compiled Static libraries
+*.lai
+*.la
+*.a
+*.lib
+
+# Executables
+*.exe
+*.out
+*.app
+
+# =====================================
+# cmake
+# =====================================
+CMakeLists.txt.user
+CMakeCache.txt
+CMakeFiles
+CMakeScripts
+Testing
+Makefile
+cmake_install.cmake
+install_manifest.txt
+compile_commands.json
+CTestTestfile.cmake
+_deps
+
+# =====================================
+# Python
+# =====================================
+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+# C extensions
+*.so
+
+# Distribution / packaging
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+share/python-wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+MANIFEST
+
+# PyInstaller
+#  Usually these files are written by a python script from a template
+#  before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.nox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*.cover
+*.py,cover
+.hypothesis/
+.pytest_cache/
+cover/
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+db.sqlite3
+db.sqlite3-journal
+
+# Flask stuff:
+instance/
+.webassets-cache
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+.pybuilder/
+target/
+
+# Jupyter Notebook
+.ipynb_checkpoints
+
+# IPython
+profile_default/
+ipython_config.py
+
+# pyenv
+#   For a library or package, you might want to ignore these files since the code is
+#   intended to run in multiple environments; otherwise, check them in:
+# .python-version
+
+# pipenv
+#   According to pypa/pipenv#598, it is recommended to include Pipfile.lock in version control.
+#   However, in case of collaboration, if having platform-specific dependencies or dependencies
+#   having no cross-platform support, pipenv may install dependencies that don't work, or not
+#   install all needed dependencies.
+#Pipfile.lock
+
+# poetry
+#   Similar to Pipfile.lock, it is generally recommended to include poetry.lock in version control.
+#   This is especially recommended for binary packages to ensure reproducibility, and is more
+#   commonly ignored for libraries.
+#   https://python-poetry.org/docs/basic-usage/#commit-your-poetrylock-file-to-version-control
+#poetry.lock
+
+# PEP 582; used by e.g. github.com/David-OConnor/pyflow
+__pypackages__/
+
+# Celery stuff
+celerybeat-schedule
+celerybeat.pid
+
+# SageMath parsed files
+*.sage.py
+
+# Environments
+.env
+.venv
+env/
+venv/
+ENV/
+env.bak/
+venv.bak/
+
+# Spyder project settings
+.spyderproject
+.spyproject
+
+# Rope project settings
+.ropeproject
+
+# mkdocs documentation
+/site
+
+# mypy
+.mypy_cache/
+.dmypy.json
+dmypy.json
+
+# Pyre type checker
+.pyre/
+
+# pytype static type analyzer
+.pytype/
+
+# Cython debug symbols
+cython_debug/
+
+# scikit-build
+_skbuild/

CMakeLists.txt

+cmake_minimum_required(VERSION 3.15...3.20)
+project(galpropy LANGUAGES CXX)
+
+set(CMAKE_CXX_FLAGS_DEBUG "-g -Wall -Wextra")
+set(CMAKE_CXX_FLAGS_RELEASE "-O2 -DNDEBUG")
+if(NOT CMAKE_BUILD_TYPE)
+    set(CMAKE_BUILD_TYPE Release)
+endif()
+message(">> CMAKE_BUILD_TYPE: ${CMAKE_BUILD_TYPE}")
+
+set(CMAKE_CXX_STANDARD 17)
+
+# set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}/cmake")
+
+if(SKBUILD)
+    set(THE_PYTHON ${PYTHON_EXECUTABLE})
+else()
+    find_package(Python COMPONENTS Interpreter Development REQUIRED)
+    set(THE_PYTHON ${Python_EXECUTABLE})
+endif()
+
+# find_package(GALPROPCC REQUIRED)
+set(GALPROPCC_INC_DIR "/usr/include/cfitsio")
+set(GALPROPCC_LIB "galpropcc")
+
+execute_process(
+    COMMAND "${THE_PYTHON}" -c
+                "import pybind11; print(pybind11.get_cmake_dir())"
+    OUTPUT_VARIABLE _tmp_dir
+    OUTPUT_STRIP_TRAILING_WHITESPACE COMMAND_ECHO STDOUT)
+
+list(APPEND CMAKE_PREFIX_PATH "${_tmp_dir}")
+
+find_package(OpenMP)
+if (OPENMP_FOUND)
+set (CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
+set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}")
+set (CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS} ${OpenMP_Fortran_FLAGS}")
+set (CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} ${OpenMP_EXE_LINKER_FLAGS}")
+endif()
+
+find_package(pybind11 CONFIG REQUIRED)
+
+FILE(GLOB TARGET_SRCS src/*.[hc]pp)
+
+set(CORE_TARGET "_bind_galpropcc")
+
+pybind11_add_module(${CORE_TARGET} MODULE ${TARGET_SRCS})
+
+target_include_directories(${CORE_TARGET} PUBLIC ${GALPROPCC_INC_DIR})
+target_link_libraries(${CORE_TARGET} PUBLIC ${GALPROPCC_LIB})
+
+install(TARGETS ${CORE_TARGET} LIBRARY DESTINATION .)

LICENSE

+MIT License
+
+Copyright (c) 2021 Chen Xia
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
\ No newline at end of file

MANIFEST.in

+include README.md
+include LICENSE
+include setup.py
+include pyproject.toml
+include CMakeLists.txt
+recursive-include src *.cpp *.hpp

README.md

+# Galpropy
+
+Python wrapper for [GALPROPCC](https://code.itp.ac.cn/ycding/galpropcc), a
+modified [GALPROP](https://galprop.stanford.edu) with redesigned C++ interface.
+The seamless experience in Python comes from the powerful
+[pybind11](https://github.com/pybind/pybind11).
+
+# Installation
+
+Install GALPROPCC and its dependencies first and then install with pip
+
+```bash
+$ pip install .
+```

pyproject.toml

+[build-system]
+requires = [
+    "setuptools>=59.0",
+    "wheel",
+    "pybind11>=2.8.0",
+    "scikit-build>=0.12",
+    "cmake>=3.21",
+    "ninja"
+]
+build-backend = "setuptools.build_meta"

setup.py

+from setuptools import find_packages
+from skbuild import setup
+
+
+def readme():
+    with open("README.md", "r") as f:
+        return f.read()
+
+
+meta_data = {
+    "name": "galpropy",
+    "version": "0.0.1",
+    "author": "Chen Xia",
+    "author_email": "xiachen@itp.ac.cn",
+    "maintainer": "Chen Xia",
+    "maintainer_email": "xiachen1996@outlook.com",
+    "description": "Python wrapper for GALPROPCC",
+    "long_description": readme(),
+    "long_description_content_type": "text/markdown",
+    "keywords": "Cosmic Ray, GALPROPCC",
+    "url": "https://code.itp.ac.cn/physcxia/galpropy",
+    "project_urls": {"Source Code":
+                     "https://code.itp.ac.cn/physcxia/galpropy"},
+    "platforms": "Linux",
+    "license": "MIT",
+    "classifiers": [
+        "Programming Language :: Python :: 3 :: Only",
+        "Programming Language :: Python :: 3.6",
+        "Programming Language :: Python :: 3.7",
+        "Programming Language :: Python :: 3.8",
+        "Programming Language :: Python :: 3.9",
+        "Programming Language :: Python :: 3.10",
+        "Programming Language :: C++",
+        "License :: OSI Approved :: MIT License",
+        "Operating System :: POSIX :: Linux",
+        "Topic :: Scientific/Engineering :: Physics",
+        "Development Status :: 3 - Alpha",
+    ]
+}
+
+
+if __name__ == '__main__':
+    setup(**meta_data,
+          packages=find_packages(where="src"),
+          package_dir={"": "src"},
+          python_requires=">=3.6",
+          zip_safe=False,
+          cmake_install_dir="src/galpropy")

src/Distribution.cpp

+#include "galpropy.hpp"
+
+void bind_Distribution(py::module_ &m)
+{
+    py::class_<Dimensions>(m, "Dimensions").def(py::init<>())
+        .def("is1", &Dimensions::is1)
+        .def("is2", &Dimensions::is2)
+        .def("is3", &Dimensions::is3)
+        .def("ndim", &Dimensions::ndim)
+        .def("npos", &Dimensions::npos);
+
+    py::bind_vector<std::vector<Distribution>>(m, "VecDistribution");
+    py::class_<Distribution>(m, "Distribution")
+        .def(py::init<>())
+        .def("init", py::overload_cast<size_t, size_t, size_t, size_t>
+             (&Distribution::init))
+        .def("init", py::overload_cast<size_t, size_t, size_t>
+             (&Distribution::init))
+        .def("init", py::overload_cast<const Dimensions&, size_t>
+             (&Distribution::init), "dims"_a, "np"_a=0)
+        .def("max", &Distribution::max)
+        .def("delete_array", &Distribution::delete_array);
+}

src/ErrorLogger.cpp

+#include <galprop/ErrorLogger.h>
+#include "galpropy.hpp"
+
+void bind_ErrorLogger(py::module_ &m)
+{
+    py::enum_<utl::ErrorLogger::SeverityLevel>(m, "SeverityLevel")
+        .value("eDebug", utl::ErrorLogger::SeverityLevel::eDebug)
+        .value("eInfo", utl::ErrorLogger::SeverityLevel::eInfo)
+        .value("eWarning", utl::ErrorLogger::SeverityLevel::eWarning)
+        .value("eError", utl::ErrorLogger::SeverityLevel::eError)
+        .value("eFatal", utl::ErrorLogger::SeverityLevel::eFatal);
+    py::enum_<utl::ErrorLogger::VerbosityLevel>(m, "VerbosityLevel")
+        .value("eDefaultVerbosity", utl::ErrorLogger::VerbosityLevel::eDefaultVerbosity)
+        .value("eTerse", utl::ErrorLogger::VerbosityLevel::eTerse)
+        .value("eVerbose", utl::ErrorLogger::VerbosityLevel::eVerbose)
+        .value("eSilent", utl::ErrorLogger::VerbosityLevel::eSilent);
+}

src/GCR.cpp

+#include <galprop/GCR.h>
+#include <galprop/DarkMatter.h>
+#include "galpropy.hpp"
+
+
+class PyPopulator : public GCR::Populator {
+public:
+    using GCR::Populator::Populator;
+
+    using ptr_t = std::shared_ptr<Particle>;
+    using map_t = std::map<std::string, ptr_t>;
+    map_t populate(const Galprop& g) const override {
+        PYBIND11_OVERRIDE_PURE(map_t, GCR::Populator, populate, g);
+    }
+};
+
+
+void bind_GCR(py::module_ &m)
+{
+    using vecptrpop = std::vector<std::shared_ptr<GCR::Populator>>;
+    py::bind_vector<vecptrpop>(m, "VecPtrPopulator");
+    py::implicitly_convertible<py::list, vecptrpop>();
+    py::implicitly_convertible<py::tuple, vecptrpop>();
+
+    py::class_<GCR> gcr(m, "GCR");
+    gcr.def(py::init<>())
+        .def("populate", &GCR::populate)
+        .def("__getitem__",
+             [](const GCR& self, const std::string& name) {
+                return self[name]; },
+             "get by name")
+        .def("__getitem__",
+             [](const GCR &self, std::size_t i) {
+                 if (i >= self.size()) throw py::index_error("out of range");
+                 return self[i]; },
+            "get by index")
+        .def("__len__", [](const GCR& self) { return self.size(); })
+        .def("__iter__", [](const GCR& self)
+             { return py::make_iterator(self.begin(), self.end()); },
+             py::keep_alive<0, 1>())
+        .def("push_back", [](GCR& self, std::shared_ptr<Particle> p) {
+             self.push_back(p); })
+        .def_readwrite("populators", &GCR::populators,
+                       py::return_value_policy::reference_internal);
+
+    py::class_<GCR::Populator, PyPopulator, std::shared_ptr<GCR::Populator>>(
+        gcr, "Populator")
+        .def(py::init_alias<>())
+        .def("populate", &GCR::Populator::populate);
+
+    py::class_<PopPrimaryNucV54, GCR::Populator,
+               std::shared_ptr<PopPrimaryNucV54>>(m, "PopPrimaryNucV54")
+        .def(py::init<>());
+    py::class_<PopPrimaryEleV54, GCR::Populator,
+               std::shared_ptr<PopPrimaryEleV54>>(m, "PopPrimaryEleV54")
+        .def(py::init<>());
+    py::class_<PopSecondaryV54, GCR::Populator,
+               std::shared_ptr<PopSecondaryV54>>(m, "PopSecondaryV54")
+        .def(py::init<>());
+
+    // dark matter TODO
+    py::class_<PopPrimaryDMSP, GCR::Populator,
+               std::shared_ptr<PopPrimaryDMSP>>(m, "PopPrimaryDMSP")
+        .def(py::init<>()).def_readwrite("dm", &PopPrimaryDMSP::dm);
+    py::class_<PopSecondaryAntinucleus, GCR::Populator,
+               std::shared_ptr<PopSecondaryAntinucleus>>(
+        m, "PopSecondaryAntinucleus").def(py::init<>());
+}

src/Galaxy.cpp

+#include <galprop/Galaxy.h>
+#include <galprop/Galdef.h>
+#include "galpropy.hpp"
+
+void bind_Galaxy(py::module_ &m)
+{
+    py::class_<Galaxy, std::shared_ptr<Galaxy>>(m, "Galaxy")
+        .def(py::init<>());
+}

src/Galdef.cpp

+#include <galprop/Galdef.h>
+#include "galpropy.hpp"
+
+void bind_Galdef(py::module_ &m)
+{
+    py::class_<Galdef>(m, "Galdef").def(py::init<>())
+        .def_readwrite("z_min", &Galdef::z_min)
+        .def_readwrite("z_max", &Galdef::z_max)
+        .def_readwrite("dz", &Galdef::dz)
+        .def_readwrite("r_min", &Galdef::r_min)
+        .def_readwrite("r_max", &Galdef::r_max)
+        .def_readwrite("primary_electrons", &Galdef::primary_electrons)
+        .def_readwrite("electron_g_0", &Galdef::electron_g_0)
+        .def_readwrite("electron_g_1", &Galdef::electron_g_1)
+        .def_readwrite("electron_g_2", &Galdef::electron_g_2)
+        .def_readwrite("electron_rigid_cutoff", &Galdef::electron_rigid_cutoff)
+        .def_readwrite("electron_rigid_br0", &Galdef::electron_rigid_br0)
+        .def_readwrite("electron_rigid_br", &Galdef::electron_rigid_br)
+        .def_readwrite("electron_rigid_br2", &Galdef::electron_rigid_br2)
+        .def_readwrite("max_Z", &Galdef::max_Z);
+}

src/Galprop.cpp

+#include "galpropy.hpp"
+
+void bind_Galprop(py::module_ &m)
+{
+    py::class_<Galprop, std::unique_ptr<Galprop, py::nodelete>>(m, "Galprop")
+        .def(py::init<>())
+        .def("create_gcr", &Galprop::create_gcr)
+        .def("create_galaxy", &Galprop::create_galaxy)
+        // .def("cr_luminosity", &Galprop::cr_luminosity)
+        .def("propagate_particles", &Galprop::propagate_particles)
+        .def("D_pp", &Galprop::D_pp)
+        .def("D_xx", &Galprop::D_xx)
+        .def("nuclei_normalize", &Galprop::nuclei_normalize)
+        .def("electron_normalize", &Galprop::electrons_normalize)
+        .def("run_DM", &Galprop::run_DM)
+        .def("propel", &Galprop::propel)
+        .def("set_dpp_limit", &Galprop::set_dpp_limit)
+        .def("load_config", &Galprop::load_config)
+        .def("run", py::overload_cast<>(&Galprop::run), "run without arg")
+        .def("run", py::overload_cast<const std::string &>(&Galprop::run),
+             "run with name")
+        .def("GetEkin", &Galprop::GetEkin)
+        .def("GetFlux", &Galprop::GetFlux, "Get the flux by ParticleNumber",
+             "ZA"_a, "ekin_toa"_a=std::valarray<double>{},
+             "mod"_a=NoModulation())
+        .def("GetFluxByName", &Galprop::GetFluxByName, "Get the flux by name",
+             "name"_a, "ekin_toa"_a=std::valarray<double>{},
+             "mod"_a=NoModulation())
+        .def_readwrite("n_species", &Galprop::n_species)
+        .def_readwrite("galdef", &Galprop::galdef)
+        .def_readwrite("cfg", &Galprop::cfg)
+        .def_readwrite("particle_lib", &Galprop::particle_lib)
+        .def_readwrite("gcr", &Galprop::gcr)
+        .def_readwrite("galaxy", &Galprop::galaxy)
+        .def_readwrite("energy_losses", &Galprop::energy_losses)
+        .def_readwrite("propel_operators", &Galprop::propel_operators)
+        .def_readwrite("energy_loss_fallback", &Galprop::energy_loss_fallback);
+        // .def_readwrite("propel_operator_fallback",
+        //                &Galprop::propel_operator_fallback);
+}

src/Particle.cpp

+#include <galprop/Particle.h>
+#include "galpropy.hpp"
+
+
+class PyParticle : public Particle {
+public:
+    using Particle::Particle;
+
+    bool depends_on(const Particle& p, Galprop& g) const override {
+        PYBIND11_OVERRIDE_PURE(bool, Particle, depends_on, p, g);
+    }
+    void fill_transport_arrays(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, Particle, fill_transport_arrays, g);
+    }
+    void fill_source_functions(Galprop& g) override {
+        PYBIND11_OVERRIDE_PURE(void, Particle, fill_source_functions, g);
+    }
+    void fill_diffusion_coeffs(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, Particle, fill_diffusion_coeffs, g);
+    }
+    void fill_energy_loss_rate(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, Particle, fill_energy_loss_rate, g);
+    }
+    void fill_destruction_rate(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, Particle, fill_destruction_rate, g);
+    }
+    void propagate(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, Particle, propagate, g);
+    }
+    double abundance(void) const override {
+        PYBIND11_OVERRIDE(double, Particle, abundance);
+    }
+};
+
+void bind_Particle(py::module_ &m)
+{
+    py::bind_vector<std::vector<std::shared_ptr<Particle>>>(
+        m, "VecPtrParticle");
+    py::implicitly_convertible<
+        py::list, std::vector<std::shared_ptr<Particle>>>();
+    py::implicitly_convertible<
+        py::tuple, std::vector<std::shared_ptr<Particle>>>();
+
+    py::bind_map<std::map<std::string, std::shared_ptr<Particle>>>(
+        m, "MapPtrParticle");
+
+    py::class_<Particle, PyParticle, std::shared_ptr<Particle>>(m, "Particle")
+        .def(py::init<std::shared_ptr<const ParticleInfo>,
+                      std::shared_ptr<const Galaxy>,
+                      const std::string &>())
+        .def("depends_on", &Particle::depends_on)
+        .def("fill_transport_arrays", &Particle::fill_transport_arrays)
+        .def("fill_source_functions", &Particle::fill_source_functions)
+        .def("fill_diffusion_coeffs", &Particle::fill_diffusion_coeffs)
+        .def("fill_energy_loss_rate", &Particle::fill_energy_loss_rate)
+        .def("fill_destruction_rate", &Particle::fill_destruction_rate)
+        .def("local_flux", py::overload_cast<double, double>(
+             &Particle::local_flux, py::const_),
+             "2-D local flux", "r"_a=Rsun, "z"_a=0)
+        .def("local_flux", py::overload_cast<double, double, double>(
+             &Particle::local_flux, py::const_),
+             "3-D local flux", "x"_a, "y"_a, "z"_a)
+        .def_readwrite("name", &Particle::name)
+        .def_readwrite("K_electron", &Particle::K_electron);
+}

src/ParticleLibrary.cpp

+#include <galprop/ParticleLibrary.h>
+#include "galpropy.hpp"
+
+void bind_ParticleLibrary(py::module_ &m)
+{
+    py::class_<ParticleInfo, std::shared_ptr<ParticleInfo>>
+        (m, "ParticleInfo")
+        .def(py::init<short, short, double, const std::string&>())
+        .def_readwrite("name", &ParticleInfo::name)
+        .def_readwrite("Ekin", &ParticleInfo::Ekin)
+        .def_readwrite("Etot", &ParticleInfo::Etot)
+        .def_readwrite("P", &ParticleInfo::P)
+        .def_readwrite("Rig", &ParticleInfo::Rig)
+        .def_readwrite("gamma", &ParticleInfo::gamma)
+        .def_readwrite("beta", &ParticleInfo::beta)
+        .def_readwrite("mnuc", &ParticleInfo::mnuc)
+        .def_readwrite("mass", &ParticleInfo::mass)
+        // .def_readwrite("t_half", &ParticleInfo::t_half) // TODO
+        .def_readwrite("beta", &ParticleInfo::beta)
+        .def_readwrite("abundance", &ParticleInfo::abundance)
+        .def_readwrite("Z", &ParticleInfo::Z)
+        .def_readwrite("A", &ParticleInfo::A)
+        .def("is_electron", &ParticleInfo::is_electron)
+        .def("is_positron", &ParticleInfo::is_positron)
+        .def("is_lepton", &ParticleInfo::is_lepton)
+        .def("is", &ParticleInfo::is)
+        .def("set_ekn", &ParticleInfo::set_ekn);
+
+    py::class_<ParticleLibrary, std::shared_ptr<ParticleLibrary>>
+        (m, "ParticleLibrary").def(py::init<>())
+        .def("add_record", &ParticleLibrary::add_record)
+        .def("get", &ParticleLibrary::get)
+        .def("set_ekn", &ParticleLibrary::set_ekn)
+        .def("init", &ParticleLibrary::init)
+        .def_static("element_name", &ParticleLibrary::element_name)
+        .def_static("isotope_name", &ParticleLibrary::isotope_name);
+}

src/PrimaryParticle.cpp

+#include <galprop/PrimaryParticle.h>
+#include "galpropy.hpp"
+
+class PyPrimaryParticle : PrimaryParticle {
+public:
+    using PrimaryParticle::PrimaryParticle;
+    bool depends_on(const Particle& p, Galprop& g) const override {
+        PYBIND11_OVERRIDE(bool, PrimaryParticle, depends_on, p, g);
+    }
+    void fill_source_functions(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, PrimaryParticle, fill_source_functions, g);
+    }
+};
+
+class PyPrimaryNucleus : PrimaryNucleus {
+public:
+    using PrimaryNucleus::PrimaryNucleus;
+    bool depends_on(const Particle& p, Galprop& g) const override {
+        PYBIND11_OVERRIDE(bool, PrimaryNucleus, depends_on, p, g);
+    }
+    void fill_source_functions(Galprop& g) override {
+        PYBIND11_OVERRIDE(void, PrimaryNucleus, fill_source_functions, g);
+    }
+    double abundance(void) const override {
+        PYBIND11_OVERRIDE(double, PrimaryNucleus, abundance,);
+    }
+};
+
+void bind_PrimaryParticle(py::module_ &m)
+{
+    py::class_<PrimaryParticle, PyPrimaryParticle, Particle,
+               std::shared_ptr<PrimaryParticle>>(m, "PrimaryParticle")
+        .def(py::init<std::shared_ptr<const ParticleInfo>,
+                      std::shared_ptr<const Galaxy>,
+                      const std::string &>())
+        .def_readwrite("spec_ptr", &PrimaryParticle::spec_ptr)
+        .def_readwrite("dist_ptr", &PrimaryParticle::dist_ptr);
+
+    py::class_<PrimaryNucleus, PyPrimaryNucleus, PrimaryParticle,
+               std::shared_ptr<PrimaryNucleus>>(m, "PrimaryNucleus")
+        .def(py::init<std::shared_ptr<const ParticleInfo>,
+                      std::shared_ptr<const Galaxy>,
+                      const std::string &>());
+}

src/PropelBase.cpp

+#include <galprop/PropelBase.h>
+#include "galpropy.hpp"
+
+class PyPropelBase : PropelBase {
+public:
+    using PropelBase::PropelBase;
+    void propel(Particle& p) override {
+        PYBIND11_OVERRIDE_PURE(void, PropelBase, propel, p);
+    }
+};
+
+class PyPropelOperatorFallback : PropelOperatorFallback {
+public:
+    using PropelOperatorFallback::PropelOperatorFallback;
+    void propel(Particle& p) override {
+        PYBIND11_OVERRIDE(void, PropelOperatorFallback, propel, p);
+    }
+};
+
+void bind_PropelBase(py::module_ &m)
+{
+    py::class_<PropelBase, PyPropelBase>(m, "PropelBase")
+        .def(py::init<const Galdef&>())
+        .def("propel", &PropelBase::propel);
+
+    py::class_<PropelOperatorFallback, PyPropelOperatorFallback>(
+        m, "PropelOperatorFallback")
+        .def(py::init<const Galdef&>())
+        .def("propel", &PropelOperatorFallback::propel);
+
+}

src/SolarModulation.cpp

+#include <pybind11/stl.h>
+#include <galprop/SolarModulation.h>
+
+#include "galpropy.hpp"
+
+class PySolarModulation : public SolarModulation {
+public:
+    using SolarModulation::SolarModulation;
+
+    std::valarray<double> operator()(
+        const std::valarray<double>& ekin_toa,
+        const std::valarray<double>& ekin_is,
+        const std::valarray<double>& flux_is,
+        const ParticleNumber& ZA) const override {
+        PYBIND11_OVERRIDE_PURE_NAME(std::valarray<double>, SolarModulation,
+                                    "__call__", operator(),
+                                    ekin_toa, ekin_is, flux_is, ZA);
+    }
+};
+
+void bind_SolarModulation(py::module_ &m)
+{
+    py::class_<SolarModulation, PySolarModulation>(m, "SolarModulation")
+        .def(py::init<>())
+        .def("__call__", &SolarModulation::operator());
+
+    py::class_<NoModulation, SolarModulation>(m, "NoModulation")
+        .def(py::init<>());
+
+    py::class_<ForceField, SolarModulation>(m, "ForceField")
+        .def(py::init<double>())
+        .def_readwrite("phi", &ForceField::phi);
+}

src/SourceDistribution.cpp

+#include <galprop/SourceDistribution.h>
+#include "galpropy.hpp"
+
+class PySourceDistributionVirtual : public SourceDistributionVirtual {
+public:
+    using SourceDistributionVirtual::SourceDistributionVirtual;
+    const Galprop* galprop_ptr(void) const override {
+        PYBIND11_OVERRIDE(const Galprop*, SourceDistributionVirtual,
+                          galprop_ptr,);
+    }
+    double density_2d(double r, double z) override {
+        PYBIND11_OVERRIDE(double, SourceDistributionVirtual, density_2d, r, z);
+    }
+    double density_3d(double x, double y, double z) override {
+        PYBIND11_OVERRIDE_PURE(double, SourceDistributionVirtual,
+                               density_3d, x, y, z);
+    }
+    std::valarray<double> distribution_2d(
+        const std::valarray<double>& r,
+        const std::valarray<double>& z) override {
+        PYBIND11_OVERRIDE(std::valarray<double>, SourceDistributionVirtual,
+                          distribution_2d, r, z);
+    }
+    std::valarray<double> distribution_3d(
+        const std::valarray<double>& x,
+        const std::valarray<double>& y,
+        const std::valarray<double>& z) override {
+        PYBIND11_OVERRIDE(std::valarray<double>, SourceDistributionVirtual,
+                          distribution_3d, x, y, z);
+    }
+};
+
+
+class PySrcDistParameterized : public SrcDistParameterized {
+public:
+    using SrcDistParameterized::SrcDistParameterized;
+    double density_3d(double x, double y, double z) override {
+        PYBIND11_OVERRIDE(double, SrcDistParameterized,
+                          density_3d, x, y, z);
+    }
+    bool is_electron(void) const override {
+        PYBIND11_OVERRIDE_PURE(bool, SrcDistParameterized, is_electron,);
+    }
+};
+
+class PyEleSrcDistParameterized : public EleSrcDistParameterized {
+public:
+    using EleSrcDistParameterized::EleSrcDistParameterized;
+    bool is_electron(void) const override {
+        PYBIND11_OVERRIDE(bool, EleSrcDistParameterized, is_electron,);
+    }
+};
+
+class PyNucSrcDistParameterized : public NucSrcDistParameterized {
+public:
+    using NucSrcDistParameterized::NucSrcDistParameterized;
+    bool is_electron(void) const override {
+        PYBIND11_OVERRIDE(bool, NucSrcDistParameterized, is_electron,);
+    }
+};
+
+void bind_SourceDistribution(py::module_ &m)
+{
+    py::class_<SourceDistributionVirtual, PySourceDistributionVirtual,
+               std::shared_ptr<SourceDistributionVirtual>>(
+        m, "SourceDistributionVirtual")
+        .def(py::init_alias<>())
+        .def("galprop_ptr", &SourceDistributionVirtual::galprop_ptr)
+        .def("density_2d", &SourceDistributionVirtual::density_2d)
+        .def("density_3d", &SourceDistributionVirtual::density_3d)
+        .def("distribution_2d", &SourceDistributionVirtual::distribution_2d)
+        .def("distribution_3d", &SourceDistributionVirtual::distribution_3d)
+        .def_readwrite("pG", &SourceDistributionVirtual::pG);
+
+    // TODO use trampoline template
+    py::class_<SrcDistParameterized, PySrcDistParameterized,
+               SourceDistributionVirtual,
+               std::shared_ptr<SrcDistParameterized>>(
+        m, "SrcDistParameterized")
+        .def(py::init_alias<>())
+        .def("is_electron", &SrcDistParameterized::is_electron)
+        .def("density_3d", &SrcDistParameterized::density_3d);
+
+    py::class_<EleSrcDistParameterized, PyEleSrcDistParameterized,
+               SrcDistParameterized, std::shared_ptr<EleSrcDistParameterized>>(
+        m, "EleSrcDistParameterized")
+        .def(py::init_alias<>())
+        .def("is_electron", &EleSrcDistParameterized::is_electron);
+
+    py::class_<NucSrcDistParameterized, PyNucSrcDistParameterized,
+               SrcDistParameterized, std::shared_ptr<NucSrcDistParameterized>>(
+        m, "NucSrcDistParameterized")
+        .def(py::init_alias<>())
+        .def("is_electron", &NucSrcDistParameterized::is_electron);
+}

src/SourceSpectrum.cpp

+#include <galprop/SourceSpectrum.h>
+#include "galpropy.hpp"
+
+class PySourceSpectrumVirtual : public SourceSpectrumVirtual {
+public:
+    using SourceSpectrumVirtual::SourceSpectrumVirtual;
+    std::valarray<double> spectrum(const Particle& p) const override {
+        PYBIND11_OVERRIDE_PURE(std::valarray<double>, SourceSpectrumVirtual,
+                               spectrum, p);
+    }
+};
+
+void bind_SourceSpectrum(py::module_ &m)
+{
+    using srcseg = SrcSpecSBPL::SmoothlyBrokenSegment;
+    using vecseg = std::vector<srcseg>;
+    py::bind_vector<vecseg>(m, "VecSmoothlyBrokenSegment");
+    py::implicitly_convertible<py::list, vecseg>();
+    py::implicitly_convertible<py::tuple, vecseg>();
+
+    py::class_<SourceSpectrumVirtual, PySourceSpectrumVirtual,
+               std::shared_ptr<SourceSpectrumVirtual>>(
+        m, "SourceSpectrumVirtual")
+        .def(py::init_alias<>());
+
+    py::class_<SrcSpecSBPL, SourceSpectrumVirtual,
+               std::shared_ptr<SrcSpecSBPL>> src(m, "SrcSpecSBPL");
+    src.def(py::init<>())
+        .def("set_norm_ekn_GeV", &SrcSpecSBPL::set_norm_ekn_GeV)
+        .def("set_norm_rig_GV", &SrcSpecSBPL::set_norm_rig_GV)
+        .def("spectrum", &SrcSpecSBPL::spectrum)
+        .def_readwrite("abundance", &SrcSpecSBPL::abundance)
+        .def_readwrite("norm_ekn", &SrcSpecSBPL::norm_ekn)
+        .def_readwrite("norm_rig", &SrcSpecSBPL::norm_rig)
+        .def_readwrite("cutoff_rig", &SrcSpecSBPL::cutoff_rig)
+        .def_readwrite("segments", &SrcSpecSBPL::segments);
+
+    py::class_<srcseg>(src, "SmoothlyBrokenSegment")
+        .def(py::init<>())
+        .def(py::init<double, double, double>(),
+             "gamma"_a, "rig"_a=HUGE_VAL, "s"_a=0)
+        .def(py::init([](const py::tuple& args) {
+                std::size_t arg_len = py::len(args);
+                if (arg_len < 1 || arg_len > 3) {
+                    throw py::value_error("tuple length incorrect");
+                }
+                double gamma = args[0].cast<double>(), rig = HUGE_VAL, s = 0;
+                if (arg_len > 1) {
+                    rig = args[1].cast<double>();
+                    if (arg_len > 2) {
+                        s = args[2].cast<double>();
+                    }
+                }
+                return srcseg {gamma, rig, s};
+            }))
+        .def_readwrite("gamma", &srcseg::gamma)
+        .def_readwrite("rig", &srcseg::rig)
+        .def_readwrite("s", &srcseg::s);
+
+    py::implicitly_convertible<py::tuple, srcseg>();
+}