From 08b00416e23a5a83fa38a01cb34688db87651ae3 Mon Sep 17 00:00:00 2001 From: Jrenaud-Desk Date: Tue, 1 Aug 2023 19:10:47 -0400 Subject: [PATCH 1/4] added new array funcs and cython hooks --- CHANGES.md | 14 + CyRK/_test.py | 2 +- CyRK/array/__init__.pxd | 0 CyRK/array/__init__.py | 2 +- CyRK/array/__init__.pyx | 0 CyRK/array/interp.pxd | 12 +- CyRK/array/interp.pyx | 242 +++++++++++++++++- CyRK/cy/cysolver.pxd | 10 + Tests/B_Other_Tests/test_array.py | 105 ++++++++ .../test_helpers.py | 0 .../test_a_cython.py | 0 .../test_b_cy_extra_output.py | 0 .../test_c_cy_readonly.py | 0 .../test_a_numba.py | 0 .../test_b_nb_extra_output.py | 0 Tests/D_Other_Tests/test_array.py | 35 --- pyproject.toml | 2 +- 17 files changed, 380 insertions(+), 44 deletions(-) create mode 100644 CyRK/array/__init__.pxd create mode 100644 CyRK/array/__init__.pyx create mode 100644 Tests/B_Other_Tests/test_array.py rename Tests/{D_Other_Tests => B_Other_Tests}/test_helpers.py (100%) rename Tests/{B_Cython_Tests => C_Cython_Tests}/test_a_cython.py (100%) rename Tests/{B_Cython_Tests => C_Cython_Tests}/test_b_cy_extra_output.py (100%) rename Tests/{B_Cython_Tests => C_Cython_Tests}/test_c_cy_readonly.py (100%) rename Tests/{C_Numba_Tests => D_Numba_Tests}/test_a_numba.py (100%) rename Tests/{C_Numba_Tests => D_Numba_Tests}/test_b_nb_extra_output.py (100%) delete mode 100644 Tests/D_Other_Tests/test_array.py diff --git a/CHANGES.md b/CHANGES.md index aa02014..10460dc 100644 --- a/CHANGES.md +++ b/CHANGES.md @@ -2,6 +2,20 @@ ## 2023 +#### v0.6.1 + +New Features +- Added top level parameters (like `MAX_STEP`) used in `CySolver` to `cysolver.pxd` so they can be cimported. +- Added new argument to `array.interp` and `array.interp_complex`: `provided_j` the user can provide a `j` index, +allowing the functions to skip the binary search. +- Added new function `interpj` to array module that outputs the interpolation result as well as the `j` index that was found. + +Bug Fixes +- Fixed issue with array tests not actually checking values. + +Other Changes +- Reordered tests since numba tests take the longest. + ### v0.6.0 New Features diff --git a/CyRK/_test.py b/CyRK/_test.py index a344121..f2fb711 100644 --- a/CyRK/_test.py +++ b/CyRK/_test.py @@ -54,7 +54,7 @@ def test_cysolver(): from CyRK.cy.cysolvertest import CySolverTester # TODO: Currently CySolver only works with floats not complex - CySolverTesterInst = CySolverTester(time_span, np.asarray(initial_conds, dtype=np.float64)) + CySolverTesterInst = CySolverTester(time_span, np.asarray(np.real(initial_conds), dtype=np.float64)) CySolverTesterInst.solve() assert CySolverTesterInst.success diff --git a/CyRK/array/__init__.pxd b/CyRK/array/__init__.pxd new file mode 100644 index 0000000..e69de29 diff --git a/CyRK/array/__init__.py b/CyRK/array/__init__.py index 90da19b..45da24d 100644 --- a/CyRK/array/__init__.py +++ b/CyRK/array/__init__.py @@ -1 +1 @@ -from CyRK.array.interp import interp, interp_complex, interp_array, interp_complex_array \ No newline at end of file +from CyRK.array.interp import interpj, interp_complexj, interp, interp_complex, interp_array, interp_complex_array \ No newline at end of file diff --git a/CyRK/array/__init__.pyx b/CyRK/array/__init__.pyx new file mode 100644 index 0000000..e69de29 diff --git a/CyRK/array/interp.pxd b/CyRK/array/interp.pxd index a167c7b..78eb5f4 100644 --- a/CyRK/array/interp.pxd +++ b/CyRK/array/interp.pxd @@ -2,9 +2,17 @@ cdef unsigned int binary_search_with_guess(double key, double[:] array, unsigned int length, unsigned int guess) nogil -cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_values) nogil +cpdef (double, int) interpj(double desired_x, double[:] x_domain, double[:] dependent_values, + int provided_j = *) nogil -cpdef double complex interp_complex(double desired_x, double[:] x_domain, double complex[:] dependent_values) nogil +cpdef (double complex, int) interp_complexj(double desired_x, double[:] x_domain, + double complex[:] dependent_values, int provided_j = *) nogil + +cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_values, + int provided_j = *) nogil + +cpdef double complex interp_complex(double desired_x, double[:] x_domain, double complex[:] dependent_values, + int provided_j = *) nogil cpdef void interp_array(double[:] desired_x_array, double[:] x_domain, double[:] dependent_values, double[:] desired_dependent_array) nogil diff --git a/CyRK/array/interp.pyx b/CyRK/array/interp.pyx index 68bcd25..a9d9e38 100644 --- a/CyRK/array/interp.pyx +++ b/CyRK/array/interp.pyx @@ -93,7 +93,219 @@ cdef unsigned int binary_search_with_guess(double key, double[:] array, unsigned return imin - 1 -cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_values) nogil: +cpdef (double, int) interpj(double desired_x, double[:] x_domain, double[:] dependent_values, + int provided_j = -1) nogil: + """ Interpolation function for floats. This function will return the index that it found during interpolation. + + Provided a domain, `x_domain` and a dependent array `dependent_values` search domain for value closest to + `desired_x` and return the value of `dependent_values` at that location if it is defined. Otherwise, use local + slopes of `x_domain` and `dependent_values` to interpolate a value of `dependent_values` at `desired_x`. + + Based on `numpy`'s `interp` function. + + Parameters + ---------- + desired_x : float + Location where `dependent_variables` is desired. + x_domain : np.ndarray[float] + Domain to search for the correct location. + dependent_values : np.ndarray[float] + Dependent values that are to be returned after search and interpolation. + provided_j : int + Give a j index from a previous interpolation to improve performance. + + Returns + ------- + result : float + Desired value of `dependent_values`. + j_out : int + The index that was found during binary_search_with_guess which can be used by other interpolation functions + to improve performance. + + """ + + cdef unsigned int lenx + lenx = len(x_domain) + # TODO: Needs to be at least 3 item long array. Add exception here? + + cdef double left_value + left_value = dependent_values[0] + cdef double right_value + right_value = dependent_values[lenx - 1] + + # Binary Search with Guess + cdef unsigned int j + j = 0 + cdef double slope + + cdef double result + cdef double fp_at_j + cdef double xp_at_j + cdef double fp_at_jp1 + cdef double xp_at_jp1 + + # Perform binary search with guess + if provided_j == -1: + # No j provided; search for it instead. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + elif provided_j < -1: + # Error + # TODO: How to handle exception handling in a cdef function... For now just repeat the search. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + else: + # provided_j is strictly positive and smaller size in this block so the conversion is safe. + j = provided_j + + cdef int j_out + j_out = j + + if j == 0: + result = left_value + elif j == lenx: + result = right_value + else: + fp_at_j = dependent_values[j] + xp_at_j = x_domain[j] + if j == lenx - 1: + result = fp_at_j + elif xp_at_j == desired_x: + result = fp_at_j + else: + fp_at_jp1 = dependent_values[j + 1] + xp_at_jp1 = x_domain[j + 1] + slope = (fp_at_jp1 - fp_at_j) / (xp_at_jp1 - xp_at_j) + + # If we get nan in one direction, try the other + result = slope * (desired_x - xp_at_j) + fp_at_j + if isnan(result): + result = slope * (desired_x - xp_at_jp1) + fp_at_jp1 + if isnan(result) and (fp_at_jp1 == fp_at_j): + result = fp_at_j + + return result, j_out + + +cpdef (double complex, int) interp_complexj(double desired_x, double[:] x_domain, + double complex[:] dependent_values, int provided_j = -1) nogil: + """ Interpolation function for complex numbers. + + Provided a domain, `desired_x` and a dependent array `dependent_values` search domain for value closest to + `desired_x` and return the value of `dependent_values` at that location if it is defined. Otherwise, use local + slopes of `desired_x` and `dependent_values` to interpolate a value of `dependent_values` at `desired_x`. + + Based on `numpy`'s `interp` function. + + Parameters + ---------- + desired_x : float + Location where `dependent_variables` is desired. + x_domain : np.ndarray[float] + Domain to search for the correct location. + dependent_values : np.ndarray[complex] + Dependent values that are to be returned after search and interpolation. + provided_j : int + Give a j index from a previous interpolation to improve performance. + + Returns + ------- + result : complex + Desired value of `dependent_values`. + j_out : int + The index that was found during binary_search_with_guess which can be used by other interpolation functions + to improve performance. + + """ + + cdef unsigned int lenx + lenx = len(x_domain) + # Note: Needs to be at least 3 item long array. Add exception here? + + cdef double complex left_value + left_value = dependent_values[0] + cdef double complex right_value + right_value = dependent_values[lenx - 1] + + # Binary Search with Guess + cdef unsigned int j + j = 0 + cdef double slope_real + cdef double slope_imag + cdef double x_slope_inverse + + cdef double result_real + cdef double result_imag + cdef double complex fp_at_j + cdef double fp_at_j_real + cdef double fp_at_j_imag + cdef double xp_at_j + cdef double complex fp_at_jp1 + cdef double fp_at_jp1_real + cdef double fp_at_jp1_imag + cdef double xp_at_jp1 + + # Perform binary search with guess + if provided_j == -1: + # No j provided; search for it instead. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + elif provided_j < -1: + # Error + # TODO: How to handle exception handling in a cdef function... For now just repeat the search. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + else: + # provided_j is strictly positive and smaller size in this block so the conversion is safe. + j = provided_j + + cdef int j_out + j_out = j + + if j == 0: + result_real = left_value.real + result_imag = left_value.imag + elif j == lenx: + result_real = right_value.real + result_imag = right_value.imag + else: + fp_at_j = dependent_values[j] + fp_at_j_real = fp_at_j.real + fp_at_j_imag = fp_at_j.imag + xp_at_j = x_domain[j] + if j == lenx - 1: + result_real = fp_at_j_real + result_imag = fp_at_j_imag + elif xp_at_j == desired_x: + result_real = fp_at_j_real + result_imag = fp_at_j_imag + else: + fp_at_jp1 = dependent_values[j + 1] + fp_at_jp1_real = fp_at_jp1.real + fp_at_jp1_imag = fp_at_jp1.imag + xp_at_jp1 = x_domain[j + 1] + x_slope_inverse = 1.0 / (xp_at_jp1 - xp_at_j) + slope_real = (fp_at_jp1_real - fp_at_j_real) * x_slope_inverse + slope_imag = (fp_at_jp1_imag - fp_at_j_imag) * x_slope_inverse + + # If we get nan in one direction try the other + # Real Part + result_real = slope_real * (desired_x - xp_at_j) + fp_at_j_real + if isnan(result_real): + result_real = slope_real * (desired_x - xp_at_jp1) + fp_at_jp1_real + if isnan(result_real) and (fp_at_jp1_real == fp_at_j_real): + result_real = fp_at_j_real + + # Imaginary Part + result_imag = slope_imag * (desired_x - xp_at_j) + fp_at_j_imag + if isnan(result_imag): + result_imag = slope_imag * (desired_x - xp_at_jp1) + fp_at_jp1_imag + if isnan(result_imag) and (fp_at_jp1_imag == fp_at_j_imag): + result_imag = fp_at_j_imag + + cdef double complex result + result = result_real + 1.0j * result_imag + + return result, j_out + + +cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_values, int provided_j = -1) nogil: """ Interpolation function for floats. Provided a domain, `x_domain` and a dependent array `dependent_values` search domain for value closest to @@ -110,6 +322,8 @@ cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_va Domain to search for the correct location. dependent_values : np.ndarray[float] Dependent values that are to be returned after search and interpolation. + provided_j : int + Give a j index from a previous interpolation to improve performance. Returns ------- @@ -139,7 +353,16 @@ cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_va cdef double xp_at_jp1 # Perform binary search with guess - j = binary_search_with_guess(desired_x, x_domain, lenx, j) + if provided_j == -1: + # No j provided; search for it instead. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + elif provided_j < -1: + # Error + # TODO: How to handle exception handling in a cdef function... For now just repeat the search. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + else: + # provided_j is strictly positive and smaller size in this block so the conversion is safe. + j = provided_j if j == 0: result = left_value @@ -168,7 +391,7 @@ cpdef double interp(double desired_x, double[:] x_domain, double[:] dependent_va cpdef double complex interp_complex(double desired_x, double[:] x_domain, - double complex[:] dependent_values) nogil: + double complex[:] dependent_values, int provided_j = -1) nogil: """ Interpolation function for complex numbers. Provided a domain, `desired_x` and a dependent array `dependent_values` search domain for value closest to @@ -185,6 +408,8 @@ cpdef double complex interp_complex(double desired_x, double[:] x_domain, Domain to search for the correct location. dependent_values : np.ndarray[complex] Dependent values that are to be returned after search and interpolation. + provided_j : int + Give a j index from a previous interpolation to improve performance. Returns ------- @@ -221,7 +446,16 @@ cpdef double complex interp_complex(double desired_x, double[:] x_domain, cdef double xp_at_jp1 # Perform binary search with guess - j = binary_search_with_guess(desired_x, x_domain, lenx, j) + if provided_j == -1: + # No j provided; search for it instead. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + elif provided_j < -1: + # Error + # TODO: How to handle exception handling in a cdef function... For now just repeat the search. + j = binary_search_with_guess(desired_x, x_domain, lenx, j) + else: + # provided_j is strictly positive and smaller size in this block so the conversion is safe. + j = provided_j if j == 0: result_real = left_value.real diff --git a/CyRK/cy/cysolver.pxd b/CyRK/cy/cysolver.pxd index 5925eef..4ca957e 100644 --- a/CyRK/cy/cysolver.pxd +++ b/CyRK/cy/cysolver.pxd @@ -1,5 +1,15 @@ cimport numpy as np from libcpp cimport bool as bool_cpp_t + +cdef double SAFETY +cdef double MIN_FACTOR +cdef double MAX_FACTOR +cdef double MAX_STEP +cdef double INF +cdef double EPS +cdef double EPS_10 +cdef double EPS_100 + cdef class CySolver: # Class attributes diff --git a/Tests/B_Other_Tests/test_array.py b/Tests/B_Other_Tests/test_array.py new file mode 100644 index 0000000..41f5bf5 --- /dev/null +++ b/Tests/B_Other_Tests/test_array.py @@ -0,0 +1,105 @@ +import numpy as np + +def test_array_module(): + """ Test that the functions load correctly. """ + from CyRK.array import interpj, interp_complexj, interp, interp_complex, interp_array, interp_complex_array + + +def test_interp(): + """ Test custom interpolation function (floats). """ + from CyRK.array import interp + + t = np.linspace(0., 1000., 1001, dtype=np.float64) + x = 100 * np.cos(t)**2 + t_test = 250.2235 + + numpy_interp = np.interp(t_test, t, x) + + cyrk_interp = interp(t_test, t, x) + + assert np.allclose([numpy_interp], [cyrk_interp]) + + +def test_interp_complex(): + """ Test custom interpolation function (complex). """ + from CyRK.array import interp_complex + + t = np.linspace(0., 1000., 1001, dtype=np.float64) + x = 100 * np.cos(t)**2 - 1.j * np.sin(t) + t_test = 250.2235 + + numpy_interp = np.interp(t_test, t, x) + + cyrk_interp = interp_complex(t_test, t, x) + + assert np.allclose([numpy_interp], [cyrk_interp]) + +def test_interp_with_provided_j(): + """ Test custom interpolation function with provided j. """ + from CyRK.array import interpj, interp_complexj, interp, interp_complex + + # Test float function + t = np.linspace(0., 1000., 1001, dtype=np.float64) + x = 100 * np.cos(t)**2 + t_test = 250.2235 + + numpy_interp = np.interp(t_test, t, x) + + # Get j + cyrk_interp_1, provided_j = interpj(t_test, t, x) + # Use j + cyrk_interp_2 = interp(t_test, t, x, provided_j=provided_j) + + # Two cyrk versions should be exactly the same. + assert cyrk_interp_1 == cyrk_interp_2 + + assert np.allclose([numpy_interp], [cyrk_interp_1]) + assert np.allclose([numpy_interp], [cyrk_interp_2]) + + # Test complex function + t_cmp = np.linspace(0., 1000., 1001, dtype=np.float64) + x_cmp = 100 * np.cos(t_cmp)**2 - 1.j * np.sin(t_cmp) + t_test_cmp = 250.2235 + + numpy_interp_cmp = np.interp(t_test_cmp, t_cmp, x_cmp) + + # Get j + cyrk_interp_1_cmp, provided_j_cmp = interp_complexj(t_test_cmp, t_cmp, x_cmp) + # Use j + cyrk_interp_2_cmp = interp_complex(t_test_cmp, t_cmp, x_cmp, provided_j=provided_j_cmp) + + # Two cyrk versions should be exactly the same. + assert cyrk_interp_1_cmp == cyrk_interp_2_cmp + + assert np.allclose([numpy_interp_cmp], [cyrk_interp_1_cmp]) + assert np.allclose([numpy_interp_cmp], [cyrk_interp_2_cmp]) + +def test_interp_array(): + """ Test custom interpolation function (float arrays). """ + from CyRK.array import interp_array + + t = np.linspace(0., 1000., 500, dtype=np.float64) + t_small = np.linspace(0., 1000., 50, dtype=np.float64) + x = 100 * np.cos(t)**2 + + numpy_interp = np.interp(t_small, t, x) + + cyrk_interp = np.empty(t_small.size, dtype=np.float64) + interp_array(t_small, t, x, cyrk_interp) + + assert np.allclose(numpy_interp, cyrk_interp) + +def test_interp_complex_array(): + """ Test custom interpolation function (complex arrays). """ + from CyRK.array import interp_complex_array + + t = np.linspace(0., 1000., 500, dtype=np.float64) + t_small = np.linspace(0., 1000., 50, dtype=np.float64) + x = 100 * np.cos(t)**2 - 1.j * np.sin(t) + + numpy_interp = np.interp(t_small, t, x) + + cyrk_interp = np.empty(t_small.size, dtype=np.complex128) + interp_complex_array(t_small, t, x, cyrk_interp) + + assert np.allclose(numpy_interp, cyrk_interp) diff --git a/Tests/D_Other_Tests/test_helpers.py b/Tests/B_Other_Tests/test_helpers.py similarity index 100% rename from Tests/D_Other_Tests/test_helpers.py rename to Tests/B_Other_Tests/test_helpers.py diff --git a/Tests/B_Cython_Tests/test_a_cython.py b/Tests/C_Cython_Tests/test_a_cython.py similarity index 100% rename from Tests/B_Cython_Tests/test_a_cython.py rename to Tests/C_Cython_Tests/test_a_cython.py diff --git a/Tests/B_Cython_Tests/test_b_cy_extra_output.py b/Tests/C_Cython_Tests/test_b_cy_extra_output.py similarity index 100% rename from Tests/B_Cython_Tests/test_b_cy_extra_output.py rename to Tests/C_Cython_Tests/test_b_cy_extra_output.py diff --git a/Tests/B_Cython_Tests/test_c_cy_readonly.py b/Tests/C_Cython_Tests/test_c_cy_readonly.py similarity index 100% rename from Tests/B_Cython_Tests/test_c_cy_readonly.py rename to Tests/C_Cython_Tests/test_c_cy_readonly.py diff --git a/Tests/C_Numba_Tests/test_a_numba.py b/Tests/D_Numba_Tests/test_a_numba.py similarity index 100% rename from Tests/C_Numba_Tests/test_a_numba.py rename to Tests/D_Numba_Tests/test_a_numba.py diff --git a/Tests/C_Numba_Tests/test_b_nb_extra_output.py b/Tests/D_Numba_Tests/test_b_nb_extra_output.py similarity index 100% rename from Tests/C_Numba_Tests/test_b_nb_extra_output.py rename to Tests/D_Numba_Tests/test_b_nb_extra_output.py diff --git a/Tests/D_Other_Tests/test_array.py b/Tests/D_Other_Tests/test_array.py deleted file mode 100644 index 3115204..0000000 --- a/Tests/D_Other_Tests/test_array.py +++ /dev/null @@ -1,35 +0,0 @@ -import numpy as np - -def test_array_module(): - """ Test that the functions load correctly. """ - from CyRK.array import interp, interp_complex, interp_array, interp_complex_array - -def test_interp_array(): - """ Test custom interpolation function (float arrays). """ - from CyRK.array import interp_array - - t = np.linspace(0., 1000., 500, dtype=np.float64) - t_small = np.linspace(0., 1000., 50, dtype=np.float64) - x = 100 * np.cos(t)**2 - - numpy_interp = np.interp(t_small, t, x) - - cyrk_interp = np.empty(t_small.size, dtype=np.float64) - interp_array(t_small, t, x, cyrk_interp) - - np.allclose(numpy_interp, cyrk_interp) - -def test_interp_complex_array(): - """ Test custom interpolation function (complex arrays). """ - from CyRK.array import interp_complex_array - - t = np.linspace(0., 1000., 500, dtype=np.float64) - t_small = np.linspace(0., 1000., 50, dtype=np.float64) - x = 100 * np.cos(t)**2 - 1.j * np.sin(t) - - numpy_interp = np.interp(t_small, t, x) - - cyrk_interp = np.empty(t_small.size, dtype=np.complex128) - interp_complex_array(t_small, t, x, cyrk_interp) - - np.allclose(numpy_interp, cyrk_interp) diff --git a/pyproject.toml b/pyproject.toml index bc655ed..377542f 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [project] name='CyRK' -version = '0.6.0' +version = '0.6.1dev1' description='Runge-Kutta ODE Integrator Implemented in Cython and Numba.' authors= [ {name = 'Joe P. Renaud', email = 'joe.p.renaud@gmail.com'} From 88fdecb07e2dfd5ff1fa8f288d0ae5a6356579c1 Mon Sep 17 00:00:00 2001 From: Jrenaud-Desk Date: Tue, 1 Aug 2023 19:18:35 -0400 Subject: [PATCH 2/4] changed array compile args --- CHANGES.md | 1 + CyRK/array/interp.pyx | 2 +- pyproject.toml | 2 +- 3 files changed, 3 insertions(+), 2 deletions(-) diff --git a/CHANGES.md b/CHANGES.md index 10460dc..6217c87 100644 --- a/CHANGES.md +++ b/CHANGES.md @@ -15,6 +15,7 @@ Bug Fixes Other Changes - Reordered tests since numba tests take the longest. +- Added `initializedcheck=False` to the array module compile arguments. ### v0.6.0 diff --git a/CyRK/array/interp.pyx b/CyRK/array/interp.pyx index a9d9e38..c5d639b 100644 --- a/CyRK/array/interp.pyx +++ b/CyRK/array/interp.pyx @@ -1,5 +1,5 @@ # distutils: language = c++ -# cython: boundscheck=False, wraparound=False, nonecheck=False, cdivision=True +# cython: boundscheck=False, wraparound=False, nonecheck=False, cdivision=True, initializedcheck=False import cython from cython.parallel import parallel, prange cimport cython diff --git a/pyproject.toml b/pyproject.toml index 377542f..bc0d18f 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [project] name='CyRK' -version = '0.6.1dev1' +version = '0.6.1dev2' description='Runge-Kutta ODE Integrator Implemented in Cython and Numba.' authors= [ {name = 'Joe P. Renaud', email = 'joe.p.renaud@gmail.com'} From abe5a8657527b7293105130e357c405b739f3747 Mon Sep 17 00:00:00 2001 From: Joe Renaud Date: Tue, 1 Aug 2023 21:13:18 -0400 Subject: [PATCH 3/4] Update README.md --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index a007877..7e2e141 100644 --- a/README.md +++ b/README.md @@ -11,7 +11,7 @@ --- -CyRK Version 0.6.0 Alpha +CyRK Version 0.6.1 Alpha **Runge-Kutta ODE Integrator Implemented in Cython and Numba** From 457c054ba5ae97f521fd8c07144454a963347d57 Mon Sep 17 00:00:00 2001 From: Joe Renaud Date: Tue, 1 Aug 2023 21:13:49 -0400 Subject: [PATCH 4/4] Update version --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index bc0d18f..9a169fd 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [project] name='CyRK' -version = '0.6.1dev2' +version = '0.6.1' description='Runge-Kutta ODE Integrator Implemented in Cython and Numba.' authors= [ {name = 'Joe P. Renaud', email = 'joe.p.renaud@gmail.com'}