From 05b942043b75a77839a6b7a4d78da3de28d523fc Mon Sep 17 00:00:00 2001
From: Rahil Doshi <rahil.doshi@fau.de>
Date: Mon, 3 Mar 2025 16:55:15 +0100
Subject: [PATCH] Remove pybind11 interface and associated C++ files, update
package dependencies and requirements
---
pyproject.toml | 2 +-
requirements.txt | 2 +-
setup.py | 10 +-
.../core/cpp/binary_search_interpolation.cpp | 329 ------------------
.../core/cpp/double_lookup_interpolation.cpp | 246 -------------
.../cpp/include/binary_search_interpolation.h | 20 --
.../cpp/include/double_lookup_interpolation.h | 18 -
7 files changed, 7 insertions(+), 620 deletions(-)
delete mode 100644 src/pymatlib/core/cpp/binary_search_interpolation.cpp
delete mode 100644 src/pymatlib/core/cpp/double_lookup_interpolation.cpp
delete mode 100644 src/pymatlib/core/cpp/include/binary_search_interpolation.h
delete mode 100644 src/pymatlib/core/cpp/include/double_lookup_interpolation.h
diff --git a/pyproject.toml b/pyproject.toml
index d5499a5..92725cc 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -42,7 +42,7 @@ test = [
requires = [
"setuptools>=69",
"wheel",
- "pybind11",
+ # "pybind11",
"ruamel.yaml>=0.15.70",
#"PyYAML",
"pandas",
diff --git a/requirements.txt b/requirements.txt
index 30df18c..fc39e17 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -10,6 +10,6 @@ matplotlib~=3.10.0
scipy~=1.15.1
setuptools~=75.8.0
versioneer~=0.29
-pybind11~=2.14.0.dev1
+# pybind11~=2.14.0.dev1
pandas~=2.2.3
openpyxl~=3.1.5
diff --git a/setup.py b/setup.py
index 594f94f..5e7af27 100644
--- a/setup.py
+++ b/setup.py
@@ -1,9 +1,9 @@
from setuptools import setup, find_packages, Extension
-import pybind11
+# import pybind11
# Define the extension module
-ext_modules = [
+'''ext_modules = [
Extension(
"pymatlib.core.cpp.fast_interpolation", # Module name in Python
[
@@ -16,7 +16,7 @@ ext_modules = [
extra_compile_args=['-O3', '-std=c++11'], # Enable high optimization and C++11
language='c++'
),
-]
+]'''
setup(
name='pymatlib',
@@ -44,7 +44,7 @@ setup(
'sympy>=1.7.0',
'pytest>=6.0.0',
'pystencils@git+https://i10git.cs.fau.de/pycodegen/pystencils.git@v2.0-dev',
- 'pybind11>=2.6.0',
+ # 'pybind11>=2.6.0',
],
extras_require={
'dev': [
@@ -53,6 +53,6 @@ setup(
'black', # For code formatting
],
},
- ext_modules=ext_modules,
+ # ext_modules=ext_modules,
include_package_data=True, # Include package data specified in MANIFEST.in
)
diff --git a/src/pymatlib/core/cpp/binary_search_interpolation.cpp b/src/pymatlib/core/cpp/binary_search_interpolation.cpp
deleted file mode 100644
index 7887af7..0000000
--- a/src/pymatlib/core/cpp/binary_search_interpolation.cpp
+++ /dev/null
@@ -1,329 +0,0 @@
-#include "include/binary_search_interpolation.h"
-#include <iostream>
-#include <cmath>
-
-
-double interpolate_binary_search(
- const py::array_t<double>& temperature_array,
- double h_in,
- const py::array_t<double>& energy_density_array) {
-
- static constexpr double EPSILON = 1e-6;
-
- // Get array views
- auto temp_arr = temperature_array.unchecked<1>();
- auto energy_arr = energy_density_array.unchecked<1>();
- const size_t n = temp_arr.shape(0);
-
- // Input validation
- if (temperature_array.size() != energy_density_array.size() || n < 2) {
- throw std::runtime_error("Invalid array sizes");
- }
-
- // Determine array order and cache indices
- const bool is_ascending = temp_arr(0) < temp_arr(n-1);
- //std::cout << "C++ - is_ascending: " << is_ascending << std::endl;
- //std::cout << "C++ - h_in: " << h_in << std::endl;
- //std::cout << "C++ - First few T values: ";
- //for(size_t i = 0; i < 5; i++) std::cout << temp_arr(i) << " ";
- //std::cout << std::endl;
- //std::cout << "C++ - First few E values: ";
- //for(size_t i = 0; i < 5; i++) std::cout << energy_arr(i) << " ";
- //std::cout << std::endl;
-
- const size_t start_idx = is_ascending ? 0 : n-1;
- const size_t end_idx = is_ascending ? n-1 : 0;
-
- // Boundary checks
- if (h_in <= energy_arr(start_idx)) return temp_arr(start_idx);
- if (h_in >= energy_arr(end_idx)) return temp_arr(end_idx);
-
- // Binary search with optimized memory access
- size_t left = 0;
- size_t right = n - 1;
-
- while (left <= right) {
- const size_t mid = (left + right) / 2;
- const double mid_val = energy_arr(mid);
- //std::cout << "C++ - mid: " << mid << ", mid_val: " << mid_val << std::endl;
-
- if (std::abs(mid_val - h_in) < EPSILON) {
- //std::cout << "C++ - Exact match found at " << mid << std::endl;
- return temp_arr(mid);
- }
-
- const bool go_left = (mid_val > h_in) == is_ascending;
- if (go_left) {
- right = mid - 1;
- } else {
- left = mid + 1;
- }
- //std::cout << "C++ - left: " << left << ", right: " << right << std::endl;
- }
-
- //std::cout << "C++ - Final indices - left: " << left << ", right: " << right << std::endl;
-
- // Linear interpolation
- const double x0 = energy_arr(right);
- const double x1 = energy_arr(left);
- const double y0 = temp_arr(right);
- const double y1 = temp_arr(left);
-
- //std::cout << "C++ - Interpolation points:" << std::endl;
- //std::cout << "C++ - x0, x1: " << x0 << ", " << x1 << std::endl;
- //std::cout << "C++ - y0, y1: " << y0 << ", " << y1 << std::endl;
-
- double result = y0 + (y1 - y0) * (h_in - x0) / (x1 - x0);
- //std::cout << "C++ - Result: " << result << std::endl;
- return result;
-}
-
-/*
-double interpolate_binary_search(
- const py::array_t<double>& temperature_array,
- double h_in,
- const py::array_t<double>& energy_density_array) {
-
- static constexpr double EPSILON = 1e-6;
-
- // Cache buffer requests and get pointers
- const auto temp_buf = temperature_array.request();
- const auto energy_buf = energy_density_array.request();
-
- // Use __restrict__ to inform compiler about non-aliasing
- const double* __restrict__ temp_ptr = static_cast<double*>(temp_buf.ptr);
- const double* __restrict__ energy_ptr = static_cast<double*>(energy_buf.ptr);
- const size_t n = temp_buf.size;
-
- // Input validation
- if (temp_buf.size != energy_buf.size || n < 2) {
- throw std::runtime_error("Invalid array sizes");
- }
-
- // Determine array order and cache indices
- const bool is_ascending = temp_ptr[0] < temp_ptr[n-1];
- const size_t start_idx = is_ascending ? 0 : n-1;
- const size_t end_idx = is_ascending ? n-1 : 0;
- std::cout << "C++ - is_ascending: " << is_ascending << std::endl;
- std::cout << "C++ - h_in: " << h_in << std::endl;
- std::cout << "C++ - First few T values: ";
- for(size_t i = 0; i < 5; i++) std::cout << temp_ptr[i] << " ";
- std::cout << std::endl;
- std::cout << "C++ - First few E values: ";
- for(size_t i = 0; i < 5; i++) std::cout << energy_ptr[i] << " ";
- std::cout << std::endl;
-
- if (energy_ptr[start_idx] >= energy_ptr[end_idx]) {
- throw std::runtime_error("Energy density must increase with temperature");
- }
-
- // Boundary checks
- if (h_in <= energy_ptr[start_idx]) return temp_ptr[start_idx];
- if (h_in >= energy_ptr[end_idx]) return temp_ptr[end_idx];
-
- // Binary search with optimized memory access
- size_t left = 0;
- size_t right = n - 1;
-
- while (left <= right) {
- const size_t mid = (left + right) / 2;
- const double mid_val = energy_ptr[mid];
-
- if (std::abs(mid_val - h_in) < EPSILON) {
- return temp_ptr[mid];
- }
-
- const bool go_left = (mid_val > h_in) == is_ascending;
- if (go_left) {
- right = mid - 1;
- } else {
- left = mid + 1;
- }
- }
-
- // Linear interpolation
- const size_t idx1 = is_ascending ? right : left;
- const size_t idx2 = is_ascending ? left : right;
- const double x0 = energy_ptr[idx1];
- const double dx = energy_ptr[idx2] - x0;
- const double y0 = temp_ptr[idx1];
-
- return y0 + (temp_ptr[idx2] - y0) * (h_in - x0) / dx;
-}
-*/
-/*
-double interpolate_binary_search(
- const py::array_t<double>& temperature_array,
- double h_in,
- const py::array_t<double>& energy_density_array) {
-
- static constexpr double EPSILON = 1e-6;
-
- const auto temp_buf = temperature_array.request();
- const auto energy_buf = energy_density_array.request();
-
- if (temp_buf.size != energy_buf.size) {
- throw std::runtime_error("Arrays must have same length");
- }
-
- const auto temp_ptr = static_cast<double*>(temp_buf.ptr);
- const auto energy_ptr = static_cast<double*>(energy_buf.ptr);
- const size_t n = temp_buf.size;
-
- if (n < 2) {
- throw std::runtime_error("Arrays must have at least 2 elements");
- }
-
- // Determine array order
- const bool is_ascending = temp_ptr[0] < temp_ptr[n-1];
- const size_t start_idx = is_ascending ? 0 : n-1;
- const size_t end_idx = is_ascending ? n-1 : 0;
-
- // Validate energy density increases with temperature
- if (energy_ptr[start_idx] >= energy_ptr[end_idx]) {
- throw std::runtime_error("Energy density must increase with temperature");
- }
-
- // Quick boundary checks
- if (h_in <= energy_ptr[start_idx]) return temp_ptr[start_idx];
- if (h_in >= energy_ptr[end_idx]) return temp_ptr[end_idx];
-
- // Binary search
- size_t left = 0;
- size_t right = n - 1;
-
- while (left <= right) {
- const size_t mid = (left + right) / 2;
- const double mid_val = energy_ptr[mid];
-
- if (std::abs(mid_val - h_in) < EPSILON) {
- return temp_ptr[mid];
- }
-
- if (mid_val > h_in) {
- if (is_ascending) {
- right = mid - 1;
- } else {
- left = mid + 1;
- }
- } else {
- if (is_ascending) {
- left = mid + 1;
- } else {
- right = mid - 1;
- }
- }
- }
-
- // Linear interpolation
- const size_t idx1 = is_ascending ? right : left;
- const size_t idx2 = is_ascending ? left : right;
-
- const double x0 = energy_ptr[idx1];
- const double x1 = energy_ptr[idx2];
- const double y0 = temp_ptr[idx1];
- const double y1 = temp_ptr[idx2];
-
- const double slope = (y1 - y0) / (x1 - x0);
- return y0 + slope * (h_in - x0);
-}
-*/
-/*
-namespace detail {
- static constexpr double EPSILON = 1e-6;
-
- inline size_t interpolation_guess(const size_t left, const size_t right,
- const double* energy_ptr, const double h_in) {
- const double left_val = energy_ptr[left];
- const double right_val = energy_ptr[right];
-
- // Interpolation guess
- const double pos = static_cast<double>(left) +
- ((h_in - left_val) * static_cast<double>(right - left)) /
- (right_val - left_val);
-
- return static_cast<size_t>(std::round(pos));
- }
-}
-
-double interpolate_binary_search(
- const py::array_t<double>& temperature_array,
- double h_in,
- const py::array_t<double>& energy_density_array) {
-
- const auto temp_buf = temperature_array.request();
- const auto energy_buf = energy_density_array.request();
-
- if (temp_buf.size != energy_buf.size) {
- throw std::runtime_error("Arrays must have same length");
- }
-
- const double* temp_ptr = static_cast<double*>(temp_buf.ptr);
- const double* energy_ptr = static_cast<double*>(energy_buf.ptr);
- const size_t n = temp_buf.size;
-
- if (n < 2) {
- throw std::runtime_error("Arrays must have at least 2 elements");
- }
-
- // Determine array order
- const bool is_ascending = temp_ptr[0] < temp_ptr[n-1];
- const size_t start_idx = is_ascending ? 0 : n-1;
- const size_t end_idx = is_ascending ? n-1 : 0;
-
- // Validate energy density increases with temperature
- if (energy_ptr[start_idx] >= energy_ptr[end_idx]) {
- throw std::runtime_error("Energy density must increase with temperature");
- }
-
- // Quick boundary checks
- if (h_in <= energy_ptr[start_idx]) return temp_ptr[start_idx];
- if (h_in >= energy_ptr[end_idx]) return temp_ptr[end_idx];
-
- // Search with interpolation and binary fallback
- size_t left = 0;
- size_t right = n - 1;
-
- while (left <= right) {
- // Try interpolation search first
- size_t mid = detail::interpolation_guess(left, right, energy_ptr, h_in);
-
- // Fallback to binary search if interpolation guess is out of bounds
- if (mid < left || mid > right) {
- mid = left + (right - left) / 2;
- }
-
- const double mid_val = energy_ptr[mid];
-
- if (std::abs(mid_val - h_in) < detail::EPSILON) {
- return temp_ptr[mid];
- }
-
- if (mid_val > h_in) {
- if (is_ascending) {
- right = mid - 1;
- } else {
- left = mid + 1;
- }
- } else {
- if (is_ascending) {
- left = mid + 1;
- } else {
- right = mid - 1;
- }
- }
- }
-
- // Linear interpolation
- const size_t idx1 = is_ascending ? right : left;
- const size_t idx2 = is_ascending ? left : right;
-
- const double x0 = energy_ptr[idx1];
- const double x1 = energy_ptr[idx2];
- const double y0 = temp_ptr[idx1];
- const double y1 = temp_ptr[idx2];
-
- const double slope = (y1 - y0) / (x1 - x0);
- return y0 + slope * (h_in - x0);
-}
-*/
\ No newline at end of file
diff --git a/src/pymatlib/core/cpp/double_lookup_interpolation.cpp b/src/pymatlib/core/cpp/double_lookup_interpolation.cpp
deleted file mode 100644
index cdf5b76..0000000
--- a/src/pymatlib/core/cpp/double_lookup_interpolation.cpp
+++ /dev/null
@@ -1,246 +0,0 @@
-#include "include/double_lookup_interpolation.h"
-#include <iostream>
-#include <cmath>
-
-
-double interpolate_double_lookup(
- double E_target,
- const py::array_t<double>& T_eq,
- const py::array_t<double>& E_neq,
- const py::array_t<double>& E_eq,
- double inv_delta_E_eq,
- const py::array_t<int32_t>& idx_map) {
-
- //std::cout << "\nC++ Debug:" << std::endl;
- //std::cout << "E_target: " << E_target << std::endl;
-
- // Get array views with bounds checking disabled for performance
- auto T_eq_arr = T_eq.unchecked<1>();
- auto E_neq_arr = E_neq.unchecked<1>();
- auto E_eq_arr = E_eq.unchecked<1>();
- auto idx_map_arr = idx_map.unchecked<1>();
-
- // Cache array size
- const size_t n = T_eq_arr.shape(0);
-
- //std::cout << "First few values:" << std::endl;
- //std::cout << "T_eq: ";
- //for(int i = 0; i < 5; i++) std::cout << T_eq_arr[i] << " ";
- //std::cout << std::endl;
- //std::cout << "E_neq: ";
- //for(int i = 0; i < 5; i++) std::cout << E_neq_arr[i] << " ";
- //std::cout << std::endl;
- //std::cout << "E_eq: ";
- //for(int i = 0; i < 5; i++) std::cout << E_eq_arr[i] << " ";
- //std::cout << std::endl;
- //std::cout << "idx_map: ";
- //for(int i = 0; i < 5; i++) std::cout << idx_map_arr[i] << " ";
- //std::cout << std::endl;
-
- // Cache frequently accessed values
- //const size_t last_idx = E_neq_arr.shape(0) - 1;
- //const double first_e = E_neq_arr(0);
- //const double last_e = E_neq_arr(last_idx);
-
- // Quick boundary checks with cached values
- // if (__builtin_expect(E_target <= E_neq_arr(0), 0)) {
- if (E_target <= E_neq_arr(0)) {
- //std::cout << "Lower boundary case: returning " << T_eq_arr[0] << std::endl;
- return T_eq_arr(0);
- }
-
- // if (__builtin_expect(E_target >= E_neq_arr(n-1), 0)) {
- if (E_target >= E_neq_arr(n-1)) {
- //std::cout << "Upper boundary case: returning " << T_eq_arr(last_idx) << std::endl;
- return T_eq_arr(n-1);
- }
-
- // Precompute and cache interpolation parameters
- // const double E_eq_start = E_eq_arr(0);
- // const double delta_E_eq = E_eq_arr(1) - E_eq_start;
- // const double inv_delta_E = 1.0 / delta_E_eq;
-
- /* const int idx_E_eq = std::min(
- static_cast<int>((E_target - E_eq_arr(0)) * inv_delta_E_eq),
- static_cast<int>(idx_map_arr.shape(0) - 1)
- ); */
- const int idx_E_eq = static_cast<int>((E_target - E_eq_arr(0)) * inv_delta_E_eq);
- //std::cout << "idx_E_eq: " << idx_E_eq << std::endl;
-
- int idx_E_neq = idx_map_arr(idx_E_eq);
- //std::cout << "idx_E_neq initial: " << idx_E_neq << std::endl;
- /*if (E_neq_arr(idx_E_neq + 1) < E_target) {
- ++idx_E_neq;
- //std::cout << "idx_E_neq adjusted: " << idx_E_neq << std::endl;
- }*/
- idx_E_neq += E_neq_arr(idx_E_neq + 1) < E_target;
-
- // Get interpolation index
- const double E1 = E_neq_arr(idx_E_neq);
- const double E2 = E_neq_arr(idx_E_neq + 1);
- const double T1 = T_eq_arr(idx_E_neq);
- const double T2 = T_eq_arr(idx_E_neq + 1);
-
- //std::cout << "Interpolation points:" << std::endl;
- //std::cout << "E1, E2: " << E1 << ", " << E2 << std::endl;
- //std::cout << "T1, T2: " << T1 << ", " << T2 << std::endl;
-
- // Optimized linear interpolation
- const double inv_dE = 1.0 / (E2 - E1);
- double result = T1 + (T2 - T1) * (E_target - E1) * inv_dE;
- //std::cout << "Final result: " << result << std::endl;
- return result;
-}
-
-/*
-double interpolate_double_lookup(
- double E_target,
- const py::array_t<double>& T_eq,
- const py::array_t<double>& E_neq,
- const py::array_t<double>& E_eq,
- const py::array_t<int32_t>& idx_map) {
-
- // Cache all buffer requests at start
- const auto T_eq_buf = T_eq.request();
- const auto E_neq_buf = E_neq.request();
- const auto E_eq_buf = E_eq.request();
- const auto idx_map_buf = idx_map.request();
-
- // Get all pointers
- const double* __restrict__ T_eq_ptr = static_cast<double*>(T_eq_buf.ptr);
- const double* __restrict__ E_neq_ptr = static_cast<double*>(E_neq_buf.ptr);
- const double* __restrict__ E_eq_ptr = static_cast<double*>(E_eq_buf.ptr);
- const int32_t* __restrict__ idx_map_ptr = static_cast<int32_t*>(idx_map_buf.ptr);
- std::cout << "C++ - First few T_eq_ptr values: ";
- for(size_t i = 0; i < 5; i++) std::cout << T_eq_ptr[i] << " ";
- std::cout << std::endl;
- std::cout << "C++ - First few E_neq_ptr values: ";
- for(size_t i = 0; i < 5; i++) std::cout << E_neq_ptr[i] << " ";
- std::cout << std::endl;
-
- // Quick boundary checks
- if (E_target <= E_neq_ptr[0]) {
- return T_eq_ptr[0];
- }
-
- const size_t last_idx = E_neq_buf.size - 1;
- if (E_target >= E_neq_ptr[last_idx]) {
- return T_eq_ptr[last_idx];
- }
-
- // Cache delta_E to avoid recomputation
- const double delta_E = E_eq_ptr[1] - E_eq_ptr[0];
- const double E_eq_start = E_eq_ptr[0];
-
- // Compute indices
- const int idx_E_eq = std::min(
- static_cast<int>((E_target - E_eq_start) / delta_E),
- static_cast<int>(idx_map_buf.size - 1)
- );
-
- int idx_E_neq = idx_map_ptr[idx_E_eq];
- if (E_neq_ptr[idx_E_neq + 1] < E_target) {
- ++idx_E_neq;
- }
-
- // Linear interpolation with minimal temporaries
- const double E1 = E_neq_ptr[idx_E_neq];
- const double dE = E_neq_ptr[idx_E_neq + 1] - E1;
- const double T1 = T_eq_ptr[idx_E_neq];
-
- return T1 + (T_eq_ptr[idx_E_neq + 1] - T1) * (E_target - E1) / dE;
-}
-*/
-/*
-double interpolate_double_lookup(
- double E_target,
- const py::array_t<double>& T_eq,
- const py::array_t<double>& E_neq,
- const py::array_t<double>& E_eq,
- const py::array_t<int32_t>& idx_map) {
-
- auto T_eq_buf = T_eq.request();
- auto E_neq_buf = E_neq.request();
- auto E_eq_buf = E_eq.request();
- auto idx_map_buf = idx_map.request();
-
- const double* T_eq_ptr = static_cast<double*>(T_eq_buf.ptr);
- const double* E_neq_ptr = static_cast<double*>(E_neq_buf.ptr);
- const double* E_eq_ptr = static_cast<double*>(E_eq_buf.ptr);
- const int32_t* idx_map_ptr = static_cast<int32_t*>(idx_map_buf.ptr);
-
- if (E_target <= E_neq_ptr[0]) {
- return T_eq_ptr[0];
- }
- if (E_target >= E_neq_ptr[E_neq_buf.size - 1]) {
- return T_eq_ptr[T_eq_buf.size - 1];
- }
-
- int idx_E_eq = static_cast<int>((E_target - E_eq_ptr[0]) / (E_eq_ptr[1] - E_eq_ptr[0]));
- idx_E_eq = std::min(idx_E_eq, static_cast<int>(idx_map_buf.size - 1));
-
- int idx_E_neq = idx_map_ptr[idx_E_eq];
-
- if (E_neq_ptr[idx_E_neq + 1] < E_target) {
- ++idx_E_neq;
- }
-
- double E1 = E_neq_ptr[idx_E_neq];
- double E2 = E_neq_ptr[idx_E_neq + 1];
- double T1 = T_eq_ptr[idx_E_neq];
- double T2 = T_eq_ptr[idx_E_neq + 1];
-
- return T1 + (T2 - T1) * (E_target - E1) / (E2 - E1);
-}
-*/
-/*
-double interpolate_double_lookup(
- double E_target,
- const py::array_t<double>& T_eq,
- const py::array_t<double>& E_neq,
- const py::array_t<double>& E_eq,
- const py::array_t<int32_t>& idx_map) {
-
- // Get buffer info once and cache sizes
- const auto& E_neq_buf = E_neq.request();
- const auto E_neq_size = E_neq_buf.size;
- const double* const E_neq_ptr = static_cast<double*>(E_neq_buf.ptr);
-
- // Quick boundary checks
- if (E_target <= E_neq_ptr[0]) {
- return static_cast<double*>(T_eq.request().ptr)[0];
- }
- if (E_target >= E_neq_ptr[E_neq_size - 1]) {
- return static_cast<double*>(T_eq.request().ptr)[E_neq_size - 1];
- }
-
- // Cache frequently used values
- const double* const E_eq_ptr = static_cast<double*>(E_eq.request().ptr);
- const double delta_E = E_eq_ptr[1] - E_eq_ptr[0];
-
- // Calculate index in equidistant grid
- const int idx_E_eq = std::min(
- static_cast<int>((E_target - E_eq_ptr[0]) / delta_E),
- static_cast<int>(idx_map.request().size - 1)
- );
-
- // Get mapped index and data pointers
- int idx_E_neq = static_cast<int32_t*>(idx_map.request().ptr)[idx_E_eq];
- const double* const T_eq_ptr = static_cast<double*>(T_eq.request().ptr);
-
- // Adjust index if needed
- if (E_neq_ptr[idx_E_neq + 1] < E_target) {
- ++idx_E_neq;
- }
-
- // Linear interpolation using direct array access
- const double E1 = E_neq_ptr[idx_E_neq];
- const double E2 = E_neq_ptr[idx_E_neq + 1];
- const double T1 = T_eq_ptr[idx_E_neq];
- const double T2 = T_eq_ptr[idx_E_neq + 1];
-
- return T1 + ((T2 - T1) * (E_target - E1)) / (E2 - E1);
-}
-*/
-
-
diff --git a/src/pymatlib/core/cpp/include/binary_search_interpolation.h b/src/pymatlib/core/cpp/include/binary_search_interpolation.h
deleted file mode 100644
index 2b4dc4c..0000000
--- a/src/pymatlib/core/cpp/include/binary_search_interpolation.h
+++ /dev/null
@@ -1,20 +0,0 @@
-#pragma once
-#include <pybind11/pybind11.h>
-#include <pybind11/numpy.h>
-#include "interpolate_binary_search_cpp.h"
-
-namespace py = pybind11;
-
-/**
- * Fast temperature interpolation using interpolation search with binary fallback
- *
- * @param temperature_array Array of temperature values (monotonically decreasing)
- * @param h_in Energy density value to interpolate
- * @param energy_density_array Array of energy density values (monotonically decreasing)
- * @return Interpolated temperature value
- * @throws std::runtime_error if arrays have different lengths or wrong monotonicity
- */
-double interpolate_binary_search(
- const py::array_t<double>& temperature_array,
- double h_in,
- const py::array_t<double>& energy_density_array);
diff --git a/src/pymatlib/core/cpp/include/double_lookup_interpolation.h b/src/pymatlib/core/cpp/include/double_lookup_interpolation.h
deleted file mode 100644
index 62efae5..0000000
--- a/src/pymatlib/core/cpp/include/double_lookup_interpolation.h
+++ /dev/null
@@ -1,18 +0,0 @@
-#ifndef DOUBLE_LOOKUP_INTERPOLATION_H
-#define DOUBLE_LOOKUP_INTERPOLATION_H
-
-#include <pybind11/pybind11.h>
-#include <pybind11/numpy.h>
-#include "interpolate_double_lookup_cpp.h"
-
-namespace py = pybind11;
-
-double interpolate_double_lookup(
- double E_target,
- const py::array_t<double>& T_eq,
- const py::array_t<double>& E_neq,
- const py::array_t<double>& E_eq,
- double inv_delta_E_eq,
- const py::array_t<int32_t>& idx_map);
-
-#endif // DOUBLE_LOOKUP_INTERPOLATION_H
--
GitLab