diff --git a/.github/workflows/cd.yml b/.github/workflows/cd.yml
index 6fe2590..f766947 100644
--- a/.github/workflows/cd.yml
+++ b/.github/workflows/cd.yml
@@ -43,7 +43,7 @@ jobs:
- name: Install uv
uses: astral-sh/setup-uv@v5
- name: Set up Python
- run: uv python install 3.13
+ run: uv python install 3.14
- name: Bump version
run: |
uv run bumpversion.py ${{ github.event.inputs.bumptype }}
diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml
index 47a7b41..876ebff 100644
--- a/.github/workflows/ci.yml
+++ b/.github/workflows/ci.yml
@@ -10,7 +10,7 @@ jobs:
strategy:
fail-fast: false
matrix:
- python-version: ["3.10", "3.11", "3.12", "3.13"]
+ python-version: ["3.10", "3.11", "3.12", "3.13", "3.14"]
os: [ubuntu-latest]
runs-on: ${{ matrix.os }}
steps:
diff --git a/.github/workflows/pages.yml b/.github/workflows/pages.yml
index d3a2909..e20d566 100644
--- a/.github/workflows/pages.yml
+++ b/.github/workflows/pages.yml
@@ -16,7 +16,7 @@ jobs:
- name: Install uv
uses: astral-sh/setup-uv@v5
with:
- python-version: "3.13"
+ python-version: "3.14"
enable-cache: true
- name: Set up Python
run: uv python install
diff --git a/.gitignore b/.gitignore
index 27f83ae..3c8fb2f 100644
--- a/.gitignore
+++ b/.gitignore
@@ -105,4 +105,6 @@ uv.lock
**data/
-pkg_graph.*
\ No newline at end of file
+pkg_graph.*
+
+/notebooks/tutorials-repository/
diff --git a/CHANGELOG.rst b/CHANGELOG.rst
index a82eb4c..926d54b 100644
--- a/CHANGELOG.rst
+++ b/CHANGELOG.rst
@@ -1,3 +1,14 @@
+Unreleased
+----------
+- **perf(numba)**: cache the generic Numba root-finder specializations, add
+ opt-in `warmup_numba()`, and enable disk caching on the direct Numba kernels.
+- **build(numba)**: widen Python < 3.14 support to Numba 0.61 through 0.65
+ without gaps, while keeping Python 3.14 on Numba 0.65.
+- **fix(compatibility)**: add Pandas 3.x compatibility fixes in the
+ cycle-count aggregation helpers and notebook examples.
+- **docs(numba)**: refresh the Numba speedup report, repair broken docs links,
+ and update the tutorial notebooks so they execute on Python 3.13 and 3.14.
+
2.0.4
------
- **bump(patch)**: Bump version to 2.0.4 [skip ci]
diff --git a/notebooks/0.0.0 - PD.Generic.development-guide.ipynb b/notebooks/0.0.0 - PD.Generic.development-guide.ipynb
index 2176078..5bd1a12 100644
--- a/notebooks/0.0.0 - PD.Generic.development-guide.ipynb
+++ b/notebooks/0.0.0 - PD.Generic.development-guide.ipynb
@@ -1,222 +1,373 @@
{
- "cells": [
- {
- "cell_type": "markdown",
- "metadata": {
- "toc": true
- },
- "source": [
- "Table of Contents
\n",
- ""
- ]
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Python 3.13.13 (py-fatigue; uv)",
+ "name": "Python 3.13.13 (py-fatigue; uv)",
+ "language": "python"
+ },
+ "language_info": {
+ "name": "python",
+ "version": "3.13.1",
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "file_extension": ".py",
+ "mimetype": "text/x-python"
+ },
+ "toc": {
+ "base_numbering": 1,
+ "nav_menu": {},
+ "number_sections": true,
+ "sideBar": true,
+ "skip_h1_title": true,
+ "title_cell": "Table of Contents",
+ "title_sidebar": "Contents",
+ "toc_cell": true,
+ "toc_position": {},
+ "toc_section_display": true,
+ "toc_window_display": false
+ },
+ "latex_envs": {
+ "LaTeX_envs_menu_present": true,
+ "autoclose": false,
+ "autocomplete": true,
+ "bibliofile": "biblio.bib",
+ "cite_by": "apalike",
+ "current_citInitial": 1,
+ "eqLabelWithNumbers": true,
+ "eqNumInitial": 1,
+ "hotkeys": {
+ "equation": "Ctrl-E",
+ "itemize": "Ctrl-I"
+ },
+ "labels_anchors": false,
+ "latex_user_defs": false,
+ "report_style_numbering": false,
+ "user_envs_cfg": false
+ },
+ "hide_input": false
},
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "# Development with the notebook\n",
- "This notebook server can be used as a scratchpad. Work in progress code, algorithms, visualizations, tests, ...\n",
- "\n",
- "The ultimate goal here is to refactor these snippets in well documented and tested python classes and methods so they can be imported easily. \n",
- "\n",
- "The suggested workflow is as follows:\n",
- "1. Use the notebook to prototype and test drive your development\n",
- "2. Migrate your developed code into the project. Create a python file, organize it in classes,...\n",
- "3. Document your migrated code (docstring) so it's functionality is documented in the sphinx docs.\n",
- "4. Import the migrated code back into this notebook and develop tests around it.\n",
- "5. Include these tests in the project.\n",
- "6. Cleanup your notebook as the functionality is now in the python project and can be easily imported in the future.\n"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "## Import the source code of the project\n",
- "Updates to the project's source need to be reflected here so it's possible to add new functionality and fix bugs which are directly reflected here. Jupyter supports this behaviour with the `autoreload` magic. \n",
- "\n",
- "The 2 cells below will setup this autoreload functionality (so the original import is not cached) and add the project path to the system path. Now the project source can be imported (without installing it as a python package)."
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {
- "ExecuteTime": {
- "end_time": "2022-06-16T14:47:08.852940Z",
- "start_time": "2022-06-16T14:47:08.805966Z"
- }
- },
- "outputs": [],
- "source": [
- "%load_ext autoreload\n",
- "%autoreload 2\n",
- "# %matplotlib inline"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {
- "ExecuteTime": {
- "end_time": "2022-06-16T14:47:10.411668Z",
- "start_time": "2022-06-16T14:47:09.429114Z"
- }
- },
- "outputs": [],
- "source": [
- "import datetime as dt\n",
- "import os\n",
- "import platform\n",
- "import sys\n",
- "import numpy as np\n",
- "import pandas as pd\n",
- "import matplotlib.pyplot as plt\n",
- "import plotly.graph_objs as go\n",
- "\n",
- "PROJECT_PATH = os.path.dirname(os.getcwd())\n",
- "print(f'PROJECT_PATH = {PROJECT_PATH}')\n",
- "\n",
- "if not PROJECT_PATH in sys.path:\n",
- " sys.path.append(PROJECT_PATH)\n",
- " \n",
- "print('Platform:', platform.platform())\n",
- "print('Python version:', sys.version)\n",
- "print('numpy version:', np.__version__)\n",
- "np.set_printoptions(threshold=sys.maxsize)"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {
- "ExecuteTime": {
- "end_time": "2022-06-16T14:47:10.475159Z",
- "start_time": "2022-06-16T14:47:10.411668Z"
- }
- },
- "outputs": [],
- "source": [
- "#matplotlib tweaking\n",
- "plt.rcParams['figure.figsize'] = (10.5/2.514/0.7, 4.5/2.514/0.7)\n",
- "plt.rcParams[\"font.family\"] = \"Sans\"\n",
- "plt.rcParams[\"font.size\"] = 10\n",
- "# plt.rcParams[\"axes.grid\"] = True\n",
- "# plt.rcParams['grid.color'] = \"#CCCCCC\"\n",
- "# plt.rcParams['grid.linestyle'] = \"-\"\n",
- "# plt.rcParams['grid.color'] = \"#DDDDDD\"\n",
- "# plt.rcParams['grid.linestyle'] = \"-\"\n",
- "plt.rcParams['axes.spines.right'] = False\n",
- "plt.rcParams['axes.spines.top'] = False\n",
- "plt.rcParams['lines.markersize'] = 3\n",
- "plt.rcParams['xtick.bottom'] = False\n",
- "plt.rcParams['xtick.labelbottom'] = True\n",
- "plt.rcParams['ytick.left'] = False\n",
- "plt.rcParams['ytick.labelleft'] = True"
- ]
- },
- {
- "cell_type": "markdown",
- "metadata": {},
- "source": [
- "## Import and use the functionality"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {
- "ExecuteTime": {
- "end_time": "2022-06-16T14:47:12.468692Z",
- "start_time": "2022-06-16T14:47:11.748517Z"
+ "nbformat": 4,
+ "nbformat_minor": 5,
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "b638a54f-cae2-4232-8d56-863cec9784c8",
+ "metadata": {
+ "toc": true
+ },
+ "source": [
+ "Table of Contents
\n",
+ "\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "bcec5489-4436-4c2a-b49b-6c22e639887d",
+ "metadata": {},
+ "source": [
+ "# Development with the notebook\n",
+ "This notebook server can be used as a scratchpad. Work in progress code, algorithms, visualizations, tests, ...\n",
+ "\n",
+ "The ultimate goal here is to refactor these snippets in well documented and tested python classes and methods so they can be imported easily. \n",
+ "\n",
+ "The suggested workflow is as follows:\n",
+ "1. Use the notebook to prototype and test drive your development\n",
+ "2. Migrate your developed code into the project. Create a python file, organize it in classes,...\n",
+ "3. Document your migrated code (docstring) so it's functionality is documented in the sphinx docs.\n",
+ "4. Import the migrated code back into this notebook and develop tests around it.\n",
+ "5. Include these tests in the project.\n",
+ "6. Cleanup your notebook as the functionality is now in the python project and can be easily imported in the future.\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "38b73e6c-e8f2-473b-b4f6-79a35d3aa1ef",
+ "metadata": {},
+ "source": [
+ "## Import the source code of the project\n",
+ "Updates to the project's source need to be reflected here so it's possible to add new functionality and fix bugs which are directly reflected here. Jupyter supports this behaviour with the `autoreload` magic. \n",
+ "\n",
+ "The 2 cells below will setup this autoreload functionality (so the original import is not cached) and add the project path to the system path. Now the project source can be imported (without installing it as a python package).\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "id": "7503552a-3170-45d2-8424-bee080872ebe",
+ "metadata": {
+ "ExecuteTime": {
+ "end_time": "2022-06-16T14:47:08.852940Z",
+ "start_time": "2022-06-16T14:47:08.805966Z"
+ }
+ },
+ "execution_count": 1,
+ "source": [
+ "%load_ext autoreload\n",
+ "%autoreload 2\n",
+ "# %matplotlib inline\n"
+ ],
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "id": "01a8057c-ccb2-45ee-98ee-24e66b1aaf9c",
+ "metadata": {
+ "ExecuteTime": {
+ "end_time": "2022-06-16T14:47:10.411668Z",
+ "start_time": "2022-06-16T14:47:09.429114Z"
+ }
+ },
+ "execution_count": 2,
+ "source": [
+ "import datetime as dt\n",
+ "import os\n",
+ "import platform\n",
+ "import sys\n",
+ "import numpy as np\n",
+ "import pandas as pd\n",
+ "import matplotlib.pyplot as plt\n",
+ "import plotly.graph_objs as go\n",
+ "\n",
+ "PROJECT_PATH = os.path.dirname(os.getcwd())\n",
+ "print(f'PROJECT_PATH = {PROJECT_PATH}')\n",
+ "\n",
+ "if not PROJECT_PATH in sys.path:\n",
+ " sys.path.append(PROJECT_PATH)\n",
+ " \n",
+ "print('Platform:', platform.platform())\n",
+ "print('Python version:', sys.version)\n",
+ "print('numpy version:', np.__version__)\n",
+ "np.set_printoptions(threshold=sys.maxsize)\n"
+ ],
+ "outputs": [
+ {
+ "output_type": "stream",
+ "name": "stdout",
+ "text": [
+ "PROJECT_PATH = /home/pietro/Code\n",
+ "Platform: Linux-6.17.0-35-generic-x86_64-with-glibc2.39\n",
+ "Python version: 3.13.13 (main, May 10 2026, 19:26:54) [Clang 22.1.3 ]\n",
+ "numpy version: 2.2.6\n"
+ ]
+ }
+ ]
+ },
+ {
+ "cell_type": "code",
+ "id": "7b33381d-3393-4d19-9370-e0b6a5cf6709",
+ "metadata": {
+ "ExecuteTime": {
+ "end_time": "2022-06-16T14:47:10.475159Z",
+ "start_time": "2022-06-16T14:47:10.411668Z"
+ }
+ },
+ "execution_count": 3,
+ "source": [
+ "#matplotlib tweaking\n",
+ "plt.rcParams['figure.figsize'] = (10.5/2.514/0.7, 4.5/2.514/0.7)\n",
+ "plt.rcParams[\"font.family\"] = \"Sans\"\n",
+ "plt.rcParams[\"font.size\"] = 10\n",
+ "# plt.rcParams[\"axes.grid\"] = True\n",
+ "# plt.rcParams['grid.color'] = \"#CCCCCC\"\n",
+ "# plt.rcParams['grid.linestyle'] = \"-\"\n",
+ "# plt.rcParams['grid.color'] = \"#DDDDDD\"\n",
+ "# plt.rcParams['grid.linestyle'] = \"-\"\n",
+ "plt.rcParams['axes.spines.right'] = False\n",
+ "plt.rcParams['axes.spines.top'] = False\n",
+ "plt.rcParams['lines.markersize'] = 3\n",
+ "plt.rcParams['xtick.bottom'] = False\n",
+ "plt.rcParams['xtick.labelbottom'] = True\n",
+ "plt.rcParams['ytick.left'] = False\n",
+ "plt.rcParams['ytick.labelleft'] = True\n"
+ ],
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "id": "887c6147-f370-4d86-bac9-c064d36693c0",
+ "metadata": {},
+ "source": [
+ "## Import and use the functionality\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "id": "7a121539-3276-4659-8476-50fa1f4cdd5e",
+ "metadata": {
+ "ExecuteTime": {
+ "end_time": "2022-06-16T14:47:12.468692Z",
+ "start_time": "2022-06-16T14:47:11.748517Z"
+ }
+ },
+ "execution_count": 4,
+ "source": [
+ "# Use the package\n",
+ "import py_fatigue as pf\n",
+ "import py_fatigue.cycle_count.rainflow as rf\n",
+ "import py_fatigue.cycle_count.histogram as ht\n",
+ "from py_fatigue.version import parse_version, __version__\n",
+ "v = parse_version(__version__)\n"
+ ],
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "id": "e77007e0-be2e-4d6a-940a-973696d0e664",
+ "metadata": {},
+ "execution_count": 5,
+ "source": [
+ "print(v)\n"
+ ],
+ "outputs": [
+ {
+ "output_type": "stream",
+ "name": "stdout",
+ "text": [
+ "Version(major=2, minor=1, patch=0, release=None, num=None)\n"
+ ]
+ }
+ ]
+ },
+ {
+ "cell_type": "code",
+ "id": "fa6c3c46-48fb-46c9-8bc5-d79995a97bbf",
+ "metadata": {},
+ "execution_count": 6,
+ "source": [
+ "\n",
+ "s = np.array([4,7,2,10,5,9,3,4,2,12,5,11,1,4,3,10,6,12,4,8,1,9,4,6], dtype=float)\n",
+ "# s = [0, 5, 0, 5, 0, 5, 0, 0 , 5, 0, 5, 0, 5, 0]\n",
+ "cc_4p = rf.rainflow(s, method='fourpoint', extended_output=True)\n",
+ "cc_astm = rf.rainflow(s, method='astm', extended_output=True)\n",
+ "display(cc_4p)\n",
+ "display(cc_astm)\n"
+ ],
+ "outputs": [
+ {
+ "output_type": "display_data",
+ "data": {
+ "text/plain": [
+ "(array([[ 2. , 7. , 1. , 4. , 2. ],\n",
+ " [ 0.5, 3.5, 1. , 6. , 2. ],\n",
+ " [ 4. , 6. , 1. , 3. , 10. ],\n",
+ " [ 3. , 8. , 1. , 10. , 2. ],\n",
+ " [ 0.5, 3.5, 1. , 13. , 2. ],\n",
+ " [ 2. , 8. , 1. , 15. , 2. ],\n",
+ " [ 2. , 6. , 1. , 18. , 2. ],\n",
+ " [ 5.5, 6.5, 1. , 12. , 10. ],\n",
+ " [ 1.5, 5.5, 0.5, 0. , 2. ],\n",
+ " [ 2.5, 4.5, 0.5, 1. , 2. ],\n",
+ " [ 5. , 7. , 0.5, 2. , 14. ],\n",
+ " [ 5.5, 6.5, 0.5, 9. , 22. ],\n",
+ " [ 4. , 5. , 0.5, 20. , 2. ],\n",
+ " [ 2.5, 6.5, 0.5, 21. , 2. ],\n",
+ " [ 1. , 5. , 0.5, 22. , 2. ]]),\n",
+ " array([ 4., 7., 2., 12., 1., 9., 4., 6.]),\n",
+ " array([ 0, 1, 2, 9, 20, 21, 22, 23]))\n"
+ ]
+ },
+ "metadata": {}
+ },
+ {
+ "output_type": "display_data",
+ "data": {
+ "text/plain": [
+ "(array([[ 1.5, 5.5, 0.5, 0. , 2. ],\n",
+ " [ 2.5, 4.5, 0.5, 1. , 2. ],\n",
+ " [ 2. , 7. , 1. , 4. , 2. ],\n",
+ " [ 0.5, 3.5, 1. , 6. , 2. ],\n",
+ " [ 4. , 6. , 0.5, 2. , 2. ],\n",
+ " [ 4. , 6. , 0.5, 3. , 10. ],\n",
+ " [ 3. , 8. , 1. , 10. , 2. ],\n",
+ " [ 5. , 7. , 0.5, 8. , 2. ],\n",
+ " [ 0.5, 3.5, 1. , 13. , 2. ],\n",
+ " [ 2. , 8. , 1. , 15. , 2. ],\n",
+ " [ 5.5, 6.5, 0.5, 9. , 6. ],\n",
+ " [ 2. , 6. , 1. , 18. , 2. ],\n",
+ " [ 5.5, 6.5, 0.5, 12. , 10. ],\n",
+ " [ 5.5, 6.5, 0.5, 17. , 6. ],\n",
+ " [ 4. , 5. , 0.5, 20. , 2. ],\n",
+ " [ 2.5, 6.5, 0.5, 21. , 2. ],\n",
+ " [ 1. , 5. , 0.5, 22. , 2. ]]),\n",
+ " array([ 4., 7., 2., 10., 2., 12., 1., 12., 1., 9., 4., 6.]),\n",
+ " array([ 0, 1, 2, 3, 8, 9, 12, 17, 20, 21, 22, 23]))\n"
+ ]
+ },
+ "metadata": {}
+ }
+ ]
+ },
+ {
+ "cell_type": "code",
+ "id": "d1fab177-b5c4-4b47-9959-c51449545052",
+ "metadata": {},
+ "execution_count": 14,
+ "source": [
+ "from py_fatigue.material import sn_curve as sn\n",
+ "#SN curves definition\n",
+ "b1_c = sn.SNCurve(\n",
+ " slope=3, # SNCurve can handle multiple slope and\n",
+ " intercept=12.436, # intercept types, as long as their sizes \n",
+ " norm='DNVGL-RP-C203', # are compatible. The slope and intercept \n",
+ " environment='Free corrosion', # attribute will be stored as numpy arrays.\n",
+ " curve='B1'\n",
+ ")\n",
+ "c_c = sn.SNCurve(\n",
+ " [3],\n",
+ " 12.115,\n",
+ " norm='DNVGL-RP-C203',\n",
+ " environment='Free corrosion',\n",
+ " curve='C'\n",
+ ")\n",
+ "e_c = sn.SNCurve(\n",
+ " 3,\n",
+ " intercept=[11.533,],\n",
+ " norm='DNVGL-RP-C203',\n",
+ " environment='Free corrosion',\n",
+ " curve='E'\n",
+ ")\n",
+ "w3_c = sn.SNCurve(\n",
+ " [3],\n",
+ " intercept=[10.493,],\n",
+ " norm='DNVGL-RP-C203',\n",
+ " environment='Free corrosion',\n",
+ " curve=\"W3-C\"\n",
+ ")\n",
+ "\n",
+ "#%%\n"
+ ],
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "id": "ab3977ff-4419-4126-b776-37d43bf86ea6",
+ "metadata": {},
+ "execution_count": 12,
+ "source": [
+ "c_c\n"
+ ],
+ "outputs": [
+ {
+ "output_type": "display_data",
+ "data": {
+ "text/plain": [
+ "\n"
+ ]
+ },
+ "metadata": {}
+ }
+ ]
+ },
+ {
+ "cell_type": "code",
+ "id": "e9d1caef-3e70-4a10-9045-8a2254bb0b8c",
+ "metadata": {},
+ "execution_count": null,
+ "source": [],
+ "outputs": []
}
- },
- "outputs": [],
- "source": [
- "# Use the package\n",
- "import py_fatigue as pf\n",
- "import py_fatigue.cycle_count.rainflow as rf\n",
- "import py_fatigue.cycle_count.histogram as ht\n",
- "from py_fatigue.version import parse_version, __version__\n",
- "v = parse_version(__version__)"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {},
- "outputs": [],
- "source": [
- "print(v)"
- ]
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "metadata": {},
- "outputs": [],
- "source": [
- "\n",
- "s = np.array([4,7,2,10,5,9,3,4,2,12,5,11,1,4,3,10,6,12,4,8,1,9,4,6], dtype=float)\n",
- "# s = [0, 5, 0, 5, 0, 5, 0, 0 , 5, 0, 5, 0, 5, 0]\n",
- "cc_4p = rf.rainflow(s, method='fourpoint', extended_output=True)\n",
- "cc_astm = rf.rainflow(s, method='astm', extended_output=True)\n",
- "display(cc_4p)\n",
- "display(cc_astm)"
- ]
- }
- ],
- "metadata": {
- "hide_input": false,
- "kernelspec": {
- "display_name": "py-fatigue (3.13.1)",
- "language": "python",
- "name": "python3"
- },
- "language_info": {
- "codemirror_mode": {
- "name": "ipython",
- "version": 3
- },
- "file_extension": ".py",
- "mimetype": "text/x-python",
- "name": "python",
- "nbconvert_exporter": "python",
- "pygments_lexer": "ipython3",
- "version": "3.13.1"
- },
- "latex_envs": {
- "LaTeX_envs_menu_present": true,
- "autoclose": false,
- "autocomplete": true,
- "bibliofile": "biblio.bib",
- "cite_by": "apalike",
- "current_citInitial": 1,
- "eqLabelWithNumbers": true,
- "eqNumInitial": 1,
- "hotkeys": {
- "equation": "Ctrl-E",
- "itemize": "Ctrl-I"
- },
- "labels_anchors": false,
- "latex_user_defs": false,
- "report_style_numbering": false,
- "user_envs_cfg": false
- },
- "toc": {
- "base_numbering": 1,
- "nav_menu": {},
- "number_sections": true,
- "sideBar": true,
- "skip_h1_title": true,
- "title_cell": "Table of Contents",
- "title_sidebar": "Contents",
- "toc_cell": true,
- "toc_position": {},
- "toc_section_display": true,
- "toc_window_display": false
- }
- },
- "nbformat": 4,
- "nbformat_minor": 2
-}
+ ]
+}
\ No newline at end of file
diff --git a/py_fatigue/__init__.py b/py_fatigue/__init__.py
index 131a614..9700e13 100644
--- a/py_fatigue/__init__.py
+++ b/py_fatigue/__init__.py
@@ -1,17 +1,26 @@
# -*- coding: utf-8 -*-
"""Py-fatigue bundles the main functionality for performing cyclic stress (fatigue) analysis and cycle-counting.""" # noqa: E501 # pylint: disable=C0301
+
+from typing import TYPE_CHECKING
+
from .cycle_count.cycle_count import CycleCount
from .material.sn_curve import SNCurve
from .material.crack_growth_curve import ParisCurve, WalkerCurve
from .version import __version__
from . import cycle_count, geometry, material, damage, styling, testing
+if TYPE_CHECKING:
+ from .damage import crack_growth
+ from .utils import warmup_numba
+
__all__ = [
"CycleCount",
"SNCurve",
"ParisCurve",
"WalkerCurve",
+ "warmup_numba",
"cycle_count",
+ "crack_growth",
"damage",
"geometry",
"material",
@@ -19,3 +28,17 @@
"testing",
"__version__",
]
+
+
+def __getattr__(name: str):
+ """Lazily expose optional package-level helpers."""
+
+ if name == "warmup_numba":
+ from .utils import warmup_numba # pylint: disable=C0415
+
+ return warmup_numba
+ if name == "crack_growth":
+ from .damage import crack_growth # pylint: disable=C0415
+
+ return crack_growth
+ raise AttributeError(f"module {__name__!r} has no attribute {name!r}")
diff --git a/py_fatigue/cycle_count/__init__.py b/py_fatigue/cycle_count/__init__.py
index f621e79..a9ca616 100644
--- a/py_fatigue/cycle_count/__init__.py
+++ b/py_fatigue/cycle_count/__init__.py
@@ -1,7 +1,10 @@
# -*- coding: utf-8 -*-
"""Cycle counting algorithms for fatigue analysis."""
+
from .cycle_count import CycleCount
+from . import utils
__all__ = [
"CycleCount",
+ "utils",
]
diff --git a/py_fatigue/cycle_count/cycle_count.py b/py_fatigue/cycle_count/cycle_count.py
index d947d7e..0902286 100644
--- a/py_fatigue/cycle_count/cycle_count.py
+++ b/py_fatigue/cycle_count/cycle_count.py
@@ -126,9 +126,9 @@ def _build_input_data_from_json( # noqa: C901
) -> dict:
"""Function returning the same data it receives"""
data = _assess_json_keys(data)
- the_hist = np.empty(0)
- range_bin_centers = np.empty(0)
- mean_bin_centers = np.empty(0)
+ the_hist: Any = np.empty(0)
+ range_bin_centers: Any = np.empty(0)
+ mean_bin_centers: Any = np.empty(0)
if "hist" not in data:
data["hist"] = []
if (
@@ -2004,6 +2004,8 @@ def pbar_sum(cc_list: Sequence[CycleCount]) -> CycleCount:
"""
if len(cc_list) == 1:
+ if hasattr(cc_list, "iloc"):
+ return cc_list.iloc[0] # type: ignore[return-value]
return cc_list[0]
chunk_size = int(np.sqrt(len(cc_list)))
diff --git a/py_fatigue/cycle_count/rainflow.py b/py_fatigue/cycle_count/rainflow.py
index 3f0bdec..cb4b73b 100644
--- a/py_fatigue/cycle_count/rainflow.py
+++ b/py_fatigue/cycle_count/rainflow.py
@@ -205,7 +205,7 @@ def findcross_indices(xn: Any) -> np.ndarray:
return ind[:idx]
-@njit(int64(int64[:], int8[:]))
+@njit(int64(int64[:], int8[:]), cache=True)
def _findcross(ind, y):
"""Return indices to zero up and downcrossings of a vector
@@ -422,7 +422,7 @@ def findrfc_astm(tp: np.ndarray, t: Optional[np.ndarray] = None) -> np.ndarray:
def _extract_cycles_from_sequence(
- seq: List[Tuple[float, Optional[float], int]]
+ seq: List[Tuple[float, Optional[float], int]],
) -> Tuple[List[Dict[str, Any]], List[Tuple[float, Optional[float], int]]]:
"""
Internal 4-point extraction following Amzallag et al. (1994).
@@ -474,7 +474,7 @@ def _extract_cycles_from_sequence(
def _duplicate_and_select_crossing_cycles(
- residue_seq: List[Tuple[float, Optional[float], int]]
+ residue_seq: List[Tuple[float, Optional[float], int]],
) -> List[Dict[str, Any]]:
"""
Decompose residue by duplication and select cycles crossing boundary.
@@ -507,7 +507,7 @@ def _duplicate_and_select_crossing_cycles(
def _process_residue(
- residue_seq: List[Tuple[float, Optional[float], int]]
+ residue_seq: List[Tuple[float, Optional[float], int]],
) -> Tuple[List[Dict[str, Any]], List[Tuple[float, Optional[float], int]]]:
"""
Fully decompose residue by iterative duplication until stable.
diff --git a/py_fatigue/cycle_count/utils.py b/py_fatigue/cycle_count/utils.py
index 743ffd3..f7a2c1e 100644
--- a/py_fatigue/cycle_count/utils.py
+++ b/py_fatigue/cycle_count/utils.py
@@ -1,5 +1,6 @@
# -*- coding: utf-8 -*-
"""Utilities for the cycle count module."""
+
# pragma: no cover
from __future__ import annotations
@@ -105,7 +106,8 @@ def aggregate_cc( # pragma: no cover
The time window to cluster the dataframe by. It must be a valid pandas
date offset frequency string or 'all'.
For all the frequency string aliases offered by pandas, see:
- `pandas-timeseries.html#dateoffset-objects `_.
+ `pandas timeseries offsets
+ `_.
save_residuals : bool, optional
If True, the residuals sequences of each aggregated CycleCount are
saved
@@ -154,9 +156,15 @@ def aggregate_cc( # pragma: no cover
# Retrieving the low-frequency fatigue dynamics on the aggregated dataframe
print("\33[36m4. Retrieving LFFD on aggregated \33[1mdf\33[22m.\33[0m")
- df_agg_rr = df_agg.applymap(
- partial(solve_lffd, rainflow_method=rainflow_method)
- )
+ map_cycles = getattr(df_agg, "map", None)
+ if map_cycles is None:
+ df_agg_rr = df_agg.applymap(
+ partial(solve_lffd, rainflow_method=rainflow_method)
+ )
+ else:
+ df_agg_rr = map_cycles(
+ partial(solve_lffd, rainflow_method=rainflow_method)
+ )
cc_cols: list[str] = [
col for col in df_agg_rr.columns if col.startswith("CC_")
@@ -287,13 +295,23 @@ def calc_aggregated_damage( # pragma: no cover
cc_cols: list[str] = [col for col in df.columns if col.startswith("CC_")]
damages = pd.DataFrame()
for _, sn_curve in sn.items():
- df_1 = df[cc_cols].applymap(
- lambda x, sk=sn_curve: np.sum(
- get_pm(cycle_count=x, sn_curve=sk)
- if isinstance(x, CycleCount)
- else 0
+ map_cycles = getattr(df[cc_cols], "map", None)
+ if map_cycles is None:
+ df_1 = df[cc_cols].applymap(
+ lambda x, sk=sn_curve: np.sum(
+ get_pm(cycle_count=x, sn_curve=sk)
+ if isinstance(x, CycleCount)
+ else 0
+ )
+ )
+ else:
+ df_1 = map_cycles(
+ lambda x, sk=sn_curve: np.sum(
+ get_pm(cycle_count=x, sn_curve=sk)
+ if isinstance(x, CycleCount)
+ else 0
+ )
)
- )
df_1["sn_curve"] = f"m={sn_curve.slope}"
damages = pd.concat([damages, df_1])
del df_1
diff --git a/py_fatigue/damage/crack_growth.py b/py_fatigue/damage/crack_growth.py
index b26a173..234a3df 100644
--- a/py_fatigue/damage/crack_growth.py
+++ b/py_fatigue/damage/crack_growth.py
@@ -18,7 +18,6 @@
from ..geometry.generic import AbstractCrackGeometry
from ..geometry.cylinder import f_hol_cyl_01
-
try:
# delete the accessor to avoid warning
del pd.DataFrame.cg # type: ignore
@@ -538,9 +537,9 @@ class CrackGrowth:
def __init__(self, pandas_obj):
# self._validate(pandas_obj)
self._obj = pandas_obj
- self.cg_curve = None
- self.crack_geometry = None
- self.final_cycles = None
+ self.cg_curve = pandas_obj.attrs.get("cg_curve")
+ self.crack_geometry = pandas_obj.attrs.get("crack_geometry")
+ self.final_cycles = pandas_obj.attrs.get("final_cycles")
@staticmethod
def _validate(obj):
@@ -629,6 +628,11 @@ def calc_growth(
geometry,
)
+ self._obj.attrs["cg_curve"] = cg_curve
+ self._obj.attrs["crack_geometry"] = crack_geometry
+ self._obj.attrs["final_cycles"] = cg_.final_cycles
+ self.cg_curve = cg_curve
+ self.crack_geometry = crack_geometry
self._obj.final_cycles = cg_.final_cycles
self.final_cycles = cg_.final_cycles
diff --git a/py_fatigue/damage/stress_life.py b/py_fatigue/damage/stress_life.py
index 86c128c..aa8918b 100644
--- a/py_fatigue/damage/stress_life.py
+++ b/py_fatigue/damage/stress_life.py
@@ -27,7 +27,7 @@
SNCurve,
)
from ..styling import py_fatigue_formatwarning
-from ..utils import make_axes, numba_bisect, _plot_damage_accumulation
+from ..utils import make_axes, _plot_damage_accumulation
try:
# delete the accessor to avoid warning
@@ -55,12 +55,6 @@ def _validate(obj):
"""Validate the input DataFrame. Raise an error if the input
DataFrame does not contain the right columns.
"""
- if {
- "cycles_to_failure",
- }.issubset(obj.columns):
- e_msg = "'cycles_to_failure' already calculated"
- raise AttributeError(e_msg)
-
if not {"count_cycle", "mean_stress", "stress_range"}.issubset(
obj.columns
):
@@ -86,6 +80,9 @@ def damage(self, sn_curve: SNCurve):
f"sn_curve ({sn_curve.unit}) do not match."
)
raise ValueError(e_msg)
+ if {"cycles_to_failure"}.issubset(self._obj.columns):
+ e_msg = "'cycles_to_failure' already calculated"
+ raise AttributeError(e_msg)
self._obj.sn_curve = sn_curve
self.sn_curve = sn_curve
self._obj["cycles_to_failure"] = sn_curve.get_cycles(
@@ -1184,8 +1181,8 @@ def find_sn_curve_intersection(
endurance: float,
weight: float,
res_stress: float,
- n_min: float,
- n_max: float,
+ n_min: float = 1e0,
+ n_max: float = 1e10,
tol=1e-6,
):
"""
@@ -1236,8 +1233,7 @@ def find_sn_curve_intersection(
ValueError
If the bisection method fails to find a solution.
"""
- return numba_bisect(
- __jit_sn_curve_residuals,
+ return __jit_sn_curve_residuals(
n_min,
n_max,
tol,
@@ -1260,6 +1256,12 @@ def find_sn_curve_intersection(
# )
+def find_sn_curve_intersection_2(*args, **kwargs):
+ """Backwards-compatible alias for the current SN intersection solver."""
+
+ return find_sn_curve_intersection(*args, **kwargs)
+
+
@nb.njit(
fastmath=True,
cache=True,
@@ -1497,3 +1499,19 @@ def calc_theil_sn_damage(
return tuple((label, *tup) for tup, label in zip(hist_dict.values(),
hist_dict.keys()))
# fmt: on
+
+
+def calc_theil_cycles_to_failure(
+ stress_range,
+ count_cycle,
+ sn_curve: SNCurve,
+ to_failure: bool = False,
+):
+ """Backwards-compatible alias for the current Theil SN helper."""
+
+ return calc_theil_sn_damage(
+ stress_range,
+ count_cycle,
+ sn_curve,
+ to_failure=to_failure,
+ )
diff --git a/py_fatigue/material/crack_growth_curve.py b/py_fatigue/material/crack_growth_curve.py
index d3be8cd..baa20f9 100755
--- a/py_fatigue/material/crack_growth_curve.py
+++ b/py_fatigue/material/crack_growth_curve.py
@@ -11,8 +11,8 @@
# Standard imports
import abc
-
import io
+import os
# import itertools
# import warnings
@@ -39,6 +39,8 @@
ensure_array,
)
+NUMBA_VERSION_INFO = tuple(int(part) for part in nb.__version__.split(".")[:2])
+
class AbstractCrackGrowthCurve(metaclass=abc.ABCMeta):
"""Abstract Paris curve. It also implements mean stress effect according
@@ -388,6 +390,8 @@ def __init__(
intercept: Union[int, float, list, np.ndarray],
threshold: Union[int, float] = 0,
critical: Union[int, float] = np.inf,
+ load_ratio: Union[int, float] = 0,
+ walker_exponent: Union[int, float] = 0,
environment: Optional[str] = None,
curve: Optional[str] = None,
norm: Optional[str] = None,
@@ -422,16 +426,17 @@ def __init__(
color : str, optional
RGBS or HEX string for color, by default None
"""
+ _ = (load_ratio, walker_exponent)
super().__init__(
- slope,
- intercept,
- threshold,
- critical,
- environment,
- curve,
- norm,
- unit_string,
- color,
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ environment=environment,
+ curve=curve,
+ norm=norm,
+ unit_string=unit_string,
+ color=color,
)
self.__threshold = threshold
self.__critical = critical
@@ -755,7 +760,7 @@ def get_knee_sif(
np.ndarray
knee SIF
"""
- knee_sif = np.empty(self.walker_intercept.size - 1, dtype=np.float64)
+ knee_sif = np.empty(self.walker_intercept.size - 1)
if not self.linear:
if self.walker_intercept.size > 1:
for i in range(self.walker_intercept.size - 1):
@@ -988,15 +993,15 @@ def __init__( # pylint: disable=R0913, R0917
RGBS or HEX string for color, by default None
"""
super().__init__(
- slope,
- intercept,
- threshold,
- critical,
- environment,
- curve,
- norm,
- unit_string,
- color,
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ environment=environment,
+ curve=curve,
+ norm=norm,
+ unit_string=unit_string,
+ color=color,
)
self.__slope, self.__intercept = _check_param_couple_types(
@@ -1057,41 +1062,44 @@ def walker_intercept(self) -> np.ndarray:
return self.__intercept * self.walker_correction
+@nb.njit(cache=True)
+def get_array_min(values): # pragma: no cover
+ """Return the minimum value of a one-dimensional array."""
+
+ min_value = values[0]
+ for i in range(1, values.size):
+ if np.isnan(values[i]) or values[i] < min_value:
+ min_value = values[i]
+ return min_value
+
+
@nb.njit(
# 'float64[::1](float64[::1], float64[::1], float64[::1])',
fastmath=False,
+ cache=True,
# parallel=True,
)
-def _calc_growth_rate(
+def _calc_growth_rate_numba(
sif, slope, intercept, threshold, critical
): # pragma: no cover # noqa: E501 # pylint: disable=C0301
# pylint: disable=not-an-iterable
assert intercept.size > 0 and intercept.size == slope.size
- assert np.nanmin(sif) >= 0
+ assert get_array_min(sif) >= 0
assert 0 <= threshold < critical <= np.inf
- knees_sif = np.empty(intercept.size - 1, dtype=np.float64)
- if intercept.size > 1:
- for i in nb.prange(intercept.size - 1):
- knees_sif[i] = (intercept[i] / intercept[i + 1]) ** (
- 1 / (slope[i + 1] - slope[i])
- )
- knees_sif = np.hstack(
- (
- np.array([0.9999999999 * threshold]),
- knees_sif,
- np.array([critical / 0.9999999999]),
- )
- )
e_msg = (
- "Knee(s) not in between threshold and critical SIF."
+ "Knee(s) are not ordered."
+ "\nCheck the definitions of slope, intercept, threshold, critical."
)
- assert np.all(np.diff(knees_sif) > 0), e_msg
- # print("knees_sif:", knees_sif)
- idx = np.digitize(sif, knees_sif, right=False) - 1
- # print("idx:", idx)
- the_growth_rate = np.empty(sif.size, dtype=np.float64)
+ if intercept.size > 1:
+ prev_knee = (intercept[0] / intercept[1]) ** (1 / (slope[1] - slope[0]))
+ for i in range(1, intercept.size - 1):
+ knee = (intercept[i] / intercept[i + 1]) ** (
+ 1 / (slope[i + 1] - slope[i])
+ )
+ assert knee > prev_knee, e_msg
+ prev_knee = knee
+ the_growth_rate = sif.copy()
for i in nb.prange(sif.size):
if sif[i] < threshold: # below threshold
the_growth_rate[i] = 0
@@ -1099,7 +1107,14 @@ def _calc_growth_rate(
if sif[i] > critical or sif[i] == np.inf:
the_growth_rate[i] = np.inf
continue
- the_growth_rate[i] = intercept[idx[i]] * sif[i] ** slope[idx[i]]
+ idx = 0
+ for j in range(intercept.size - 1):
+ knee = (intercept[j] / intercept[j + 1]) ** (
+ 1 / (slope[j + 1] - slope[j])
+ )
+ if sif[i] >= knee:
+ idx = j + 1
+ the_growth_rate[i] = intercept[idx] * sif[i] ** slope[idx]
return the_growth_rate
@@ -1107,58 +1122,105 @@ def _calc_growth_rate(
@nb.njit(
# 'float64[::1](float64[::1], float64[::1], float64[::1])',
fastmath=False,
+ cache=True,
# parallel=True,
)
-def _calc_sif(
+def _calc_sif_numba(
growth_rate, slope, intercept, threshold, critical
): # pragma: no cover # noqa: E501 # pylint: disable=C0301
# pylint: disable=not-an-iterable
assert intercept.size > 0 and intercept.size == slope.size
- assert np.nanmin(growth_rate) >= 0
+ assert get_array_min(growth_rate) >= 0
assert 0 <= threshold < critical <= np.inf
- knees_growth_rate = np.empty(intercept.size - 1, dtype=np.float64)
- if intercept.size > 1:
- for i in nb.prange(intercept.size - 1):
- m_i = slope[i + 1] / (slope[i + 1] - slope[i])
- m_i_p_1 = slope[i] / (slope[i + 1] - slope[i])
- knees_growth_rate[i] = (
- intercept[i] ** m_i / intercept[i + 1] ** m_i_p_1
- )
-
- knees_growth_rate = np.hstack(
- (
- np.array(
- [intercept[0] * (0.9999999999999999 * threshold) ** slope[0]]
- ),
- knees_growth_rate,
- np.array(
- [intercept[-1] * (critical / 0.9999999999999999) ** slope[-1]]
- ),
- )
- )
+ lower_bound = intercept[0] * (0.9999999999999999 * threshold) ** slope[0]
+ upper_bound = intercept[-1] * (critical / 0.9999999999999999) ** slope[-1]
e_msg = (
"Knee(s) not in between threshold and critical SIF."
+ "\nCheck the definitions of slope, intercept, threshold, critical."
)
- assert np.all(np.diff(knees_growth_rate) > 0), e_msg
- # print("knees_growth_rate:", knees_growth_rate)
- idx = np.digitize(growth_rate, knees_growth_rate, right=False) - 1
- # print("idx:", idx)
- the_sif = np.empty(growth_rate.size, dtype=np.float64)
- nr_knees = intercept.size - 1
+ if intercept.size > 1:
+ prev_growth_rate = lower_bound
+ for i in range(intercept.size - 1):
+ m_i = slope[i + 1] / (slope[i + 1] - slope[i])
+ m_i_p_1 = slope[i] / (slope[i + 1] - slope[i])
+ knee_growth_rate = intercept[i] ** m_i / intercept[i + 1] ** m_i_p_1
+ assert knee_growth_rate > prev_growth_rate, e_msg
+ prev_growth_rate = knee_growth_rate
+ assert upper_bound > prev_growth_rate, e_msg
+ else:
+ assert upper_bound > lower_bound, e_msg
+
+ the_sif = growth_rate.copy()
for i in nb.prange(growth_rate.size):
- if idx[i] <= 0: # below threshold
+ if growth_rate[i] < lower_bound: # below threshold
the_sif[i] = threshold
continue
- if idx[i] > nr_knees: # above threshold
+ if growth_rate[i] >= upper_bound: # above threshold
the_sif[i] = critical
continue
- the_sif[i] = (growth_rate[i] / intercept[idx[i]]) ** (1 / slope[idx[i]])
+ idx = 0
+ for j in range(intercept.size - 1):
+ m_i = slope[j + 1] / (slope[j + 1] - slope[j])
+ m_i_p_1 = slope[j] / (slope[j + 1] - slope[j])
+ knee_growth_rate = intercept[j] ** m_i / intercept[j + 1] ** m_i_p_1
+ if growth_rate[i] >= knee_growth_rate:
+ idx = j + 1
+ the_sif[i] = (growth_rate[i] / intercept[idx]) ** (1 / slope[idx])
return the_sif
+def use_python_crack_growth_kernels() -> bool:
+ """Return whether crack-growth kernels should bypass Numba dispatch."""
+
+ return (
+ NUMBA_VERSION_INFO <= (0, 61)
+ or getattr(nb.config, "DISABLE_JIT", False)
+ or os.environ.get("NUMBA_DISABLE_JIT") == "1"
+ )
+
+
+def _calc_growth_rate(
+ sif, slope, intercept, threshold, critical
+): # pragma: no cover # noqa: E501 # pylint: disable=C0301
+ if use_python_crack_growth_kernels():
+ return _calc_growth_rate_numba.py_func(
+ sif,
+ slope,
+ intercept,
+ threshold,
+ critical,
+ )
+ return _calc_growth_rate_numba(
+ sif,
+ slope,
+ intercept,
+ threshold,
+ critical,
+ )
+
+
+def _calc_sif(
+ growth_rate, slope, intercept, threshold, critical
+): # pragma: no cover # noqa: E501 # pylint: disable=C0301
+ if use_python_crack_growth_kernels():
+ return _calc_sif_numba.py_func(
+ growth_rate,
+ slope,
+ intercept,
+ threshold,
+ critical,
+ )
+ return _calc_sif_numba(
+ growth_rate,
+ slope,
+ intercept,
+ threshold,
+ critical,
+ )
+
+
def _paris_curve_data_points(pc: ParisCurve) -> tuple:
"""Draw Paris' curve data points, including possible knees and
threshold/critical SIF.
@@ -1192,7 +1254,12 @@ def _paris_curve_data_points(pc: ParisCurve) -> tuple:
) # 0.999999999999 to avoid rounding errors
if not pc.linear:
- sif_plot = np.sort(np.append(sif_plot, np.array(pc.get_knee_sif())))
+ knee_sif = np.array(pc.get_knee_sif())
+ knee_sif = knee_sif[
+ (knee_sif > pc.threshold)
+ & (knee_sif < 0.999999999999 * pc.critical)
+ ]
+ sif_plot = np.sort(np.append(sif_plot, knee_sif))
growth_rate_plot = pc.get_growth_rate(sif_plot)
if pc.threshold > 0:
growth_rate_plot = np.hstack(
diff --git a/py_fatigue/material/sn_curve.py b/py_fatigue/material/sn_curve.py
index 2ce2b30..650c7b2 100644
--- a/py_fatigue/material/sn_curve.py
+++ b/py_fatigue/material/sn_curve.py
@@ -18,6 +18,7 @@
import abc
import io
import itertools
+import os
import warnings
# Non-standard "from" imports
@@ -689,7 +690,13 @@ def get_cycles(
# axis=1
# ), self.endurance)
- return _calc_cycles_2(
+ calc_cycles = _calc_cycles
+ if (
+ getattr(nb.config, "DISABLE_JIT", False)
+ or os.environ.get("NUMBA_DISABLE_JIT") == "1"
+ ):
+ calc_cycles = getattr(_calc_cycles, "py_func", _calc_cycles)
+ return calc_cycles(
stress_range, self.slope, self.intercept, self.endurance
)
@@ -711,9 +718,13 @@ def get_stress(
# cycles.shape[0],-1),
# axis=1
# ), endurance_stress)
- return _calc_stress_2(
- cycles, self.slope, self.intercept, self.endurance
- )
+ calc_stress = _calc_stress
+ if (
+ getattr(nb.config, "DISABLE_JIT", False)
+ or os.environ.get("NUMBA_DISABLE_JIT") == "1"
+ ):
+ calc_stress = getattr(_calc_stress, "py_func", _calc_stress)
+ return calc_stress(cycles, self.slope, self.intercept, self.endurance)
def n( # pylint: disable=invalid-name
self, sigma: int | float | list | np.ndarray
@@ -941,35 +952,24 @@ def plot(
return fig, ax
-@nb.njit(
- # 'float64[::1](float64[::1], float64[::1], float64[::1])',
- fastmath=False,
- # parallel=True,
-)
-def _calc_cycles(stress, slope, intercept, endurance): # pragma: no cover
- # pylint: disable=not-an-iterable
- assert intercept.size > 0 and intercept.size == slope.size
- assert np.min(stress) >= 0
- log10 = np.log(10)
- the_cycles = np.empty(stress.size, dtype=np.float64)
- for i in nb.prange(stress.size):
- log_stress = np.log(stress[i])
- max_i = intercept[0] * log10 - slope[0] * log_stress
- for j in range(1, len(intercept)):
- value = intercept[j] * log10 - slope[j] * log_stress
- max_i = max(max_i, value)
- the_cycles[i] = np.exp(max_i)
- if endurance < np.inf:
- the_cycles[the_cycles > endurance] = np.inf
- return the_cycles
+@nb.njit(cache=True)
+def get_sn_array_min(values): # pragma: no cover
+ """Return the minimum value of a one-dimensional array."""
+
+ min_value = values[0]
+ for i in range(1, values.size):
+ if np.isnan(values[i]) or values[i] < min_value:
+ min_value = values[i]
+ return min_value
@nb.njit(
# 'float64[::1](float64[::1], float64[::1], float64[::1])',
fastmath=False,
+ cache=True,
# parallel=True,
)
-def _calc_cycles_2(stress, slope, intercept, endurance): # pragma: no cover
+def _calc_cycles(stress, slope, intercept, endurance): # pragma: no cover
"""
Calculate the number of cycles to failure for given stress levels.
@@ -983,11 +983,11 @@ def _calc_cycles_2(stress, slope, intercept, endurance): # pragma: no cover
numpy.ndarray: Array of calculated cycles to failure.
"""
assert intercept.size > 0 and intercept.size == slope.size
- assert np.min(stress) >= 0
+ assert get_sn_array_min(stress) >= 0
log_stress = np.log10(stress)
log_endurance = np.log10(endurance)
- log_knee_stress = np.empty(intercept.size - 1, dtype=np.float64)
+ log_knee_stress = np.empty(intercept.size - 1)
if intercept.size > 1:
for i in nb.prange(intercept.size - 1): # pylint: disable=E1133
log_knee_stress[i] = (intercept[i + 1] - intercept[i]) / (
@@ -1002,7 +1002,7 @@ def _calc_cycles_2(stress, slope, intercept, endurance): # pragma: no cover
)
)
idx = np.digitize(log_stress, log_knee_stress, right=False) - 1
- the_cycles = np.empty(stress.size, dtype=np.float64)
+ the_cycles = np.empty(stress.size)
nr_knees = intercept.size - 1
for i in nb.prange(stress.size): # pylint: disable=E1133
if idx[i] <= 0:
@@ -1021,36 +1021,10 @@ def _calc_cycles_2(stress, slope, intercept, endurance): # pragma: no cover
@nb.njit(
# 'float64[::1](float64[::1], float64[::1], float64[::1])',
fastmath=False,
+ cache=True,
# parallel=True,
)
def _calc_stress(cycles, slope, intercept, endurance): # pragma: no cover
- # pylint: disable=not-an-iterable
- assert intercept.size > 0 and intercept.size == slope.size
- assert np.min(cycles) > 0
- log10 = np.log(10)
- the_stress = np.empty(cycles.size, dtype=np.float64)
- for i in nb.prange(cycles.size):
- log_cycles = np.log(cycles[i])
- max_i = (intercept[0] * log10 - log_cycles) / slope[0]
- for j in range(1, len(intercept)):
- value = (intercept[j] * log10 - log_cycles) / slope[j]
- max_i = max(max_i, value)
- the_stress[i] = np.exp(max_i)
- if endurance < np.inf:
- endurance_stress = np.exp(
- (intercept[-1] * log10 - np.log(endurance)) / slope[-1]
- )
- the_stress[the_stress < endurance_stress] = endurance_stress
- return the_stress
-
-
-@nb.njit(
- # 'float64[::1](float64[::1], float64[::1], float64[::1])',
- fastmath=False,
- # parallel=True,
- cache=True,
-)
-def _calc_stress_2(cycles, slope, intercept, endurance): # pragma: no cover
"""
Calculate the number of cycles to failure for given stress levels.
@@ -1064,12 +1038,12 @@ def _calc_stress_2(cycles, slope, intercept, endurance): # pragma: no cover
numpy.ndarray: Array of calculated cycles to failure.
"""
assert intercept.size > 0 and intercept.size == slope.size
- assert np.nanmin(cycles) >= 0
+ assert get_sn_array_min(cycles) > 0
log_cycles = np.log10(cycles)
log_endurance = np.log10(endurance)
- log_knee_stress = np.empty(intercept.size - 1, dtype=np.float64)
- log_knee_cycles = np.empty(intercept.size - 1, dtype=np.float64)
+ log_knee_stress = np.empty(intercept.size - 1)
+ log_knee_cycles = np.empty(intercept.size - 1)
if intercept.size > 1:
for i in nb.prange(intercept.size - 1): # pylint: disable=E1133
log_knee_stress[i] = (intercept[i + 1] - intercept[i]) / (
@@ -1089,7 +1063,7 @@ def _calc_stress_2(cycles, slope, intercept, endurance): # pragma: no cover
)
endurance_stress = 10 ** log_knee_stress[-1]
idx = np.digitize(log_cycles, log_knee_cycles, right=False) - 1
- the_stress = np.empty(cycles.size, dtype=np.float64)
+ the_stress = np.empty(cycles.size)
nr_knees = intercept.size - 1
for i in nb.prange(cycles.size): # pylint: disable=E1133
if idx[i] <= 0:
@@ -1150,9 +1124,13 @@ def __sn_curve_residuals( # pragma: no cover
# res_stress = 0.
# print(f"np.array([cycles]) = {cycles}, type = {type(cycles)}")
- fail = _calc_stress_2(np.array([cycles]), slope, intercept, endurance)[0]
+ fail = _calc_stress(np.array([cycles]), slope, intercept, endurance)[0]
return fail - weight * cycles - res_stress
# Create a jitted bisection function specialized for root_func
__jit_sn_curve_residuals = compile_specialized_bisect(__sn_curve_residuals)
+
+# Backwards-compatible aliases used by older notebooks.
+sn_curve_residuals = __sn_curve_residuals
+jit_sn_curve_residuals = __jit_sn_curve_residuals
diff --git a/py_fatigue/mean_stress/corrections.py b/py_fatigue/mean_stress/corrections.py
index 88e364c..d796d37 100644
--- a/py_fatigue/mean_stress/corrections.py
+++ b/py_fatigue/mean_stress/corrections.py
@@ -5,6 +5,7 @@
"""
from __future__ import annotations
+from typing import cast
import logging
import warnings
@@ -26,7 +27,7 @@ def dnvgl_mean_stress_correction(
plot: bool = False,
) -> np.ndarray:
"""Calculates the mean stress correction according to par.
- 2.5 of`DNVGL-RP-C203 `_ which includes
+ 2.5 of DNVGL-RP-C203, which includes
an attenuation factor :math:`p` for the stress ranges if the
following cases:
@@ -320,8 +321,7 @@ def goodman_haigh_mean_stress_correction( # pylint: disable=R0912 # noqa: C901,
See Also
--------
:func:`py_fatigue.utils.numba_newton`,
- :func:`py_fatigue.utils.compile_specialized_newton`,
- :func:`py_fatigue.mean_stress.corrections.__goodman
+ :func:`py_fatigue.utils.compile_specialized_newton`
Raises
------
@@ -394,6 +394,8 @@ def goodman_haigh_mean_stress_correction( # pylint: disable=R0912 # noqa: C901,
for r_out_val in r_out:
# NOTE: Special cases for r_out = -1 and r_out = 0
# NOTE: These cases are solved analytically to improve performance
+ # NOTE: The case correction_exponent = 1 is also solved
+ # analytically, as the equation becomes linear in amp_out.
if r_out_val == -1:
# amp_out.append(
# amp_in / (1 - (mean_in / ult_s) ** correction_exponent)
@@ -404,6 +406,17 @@ def goodman_haigh_mean_stress_correction( # pylint: disable=R0912 # noqa: C901,
)
mean_out[r_out == -1, :] = np.zeros_like(amp_in)
continue
+ if correction_exponent == 1:
+ r_in_term = (1 + r_in) / (1 - r_in) * amp_in / ult_s
+ r_out_term = (1 + r_out_val) / (1 - r_out_val) / ult_s
+ amp_out_fsolve = ((1 - r_in_term) / amp_in + r_out_term) ** -1
+ amp_out[r_out == r_out_val, :] = np.clip(amp_out_fsolve, 0.0, ult_s)
+ mean_out[r_out == r_out_val, :] = (
+ amp_out[r_out == r_out_val, :]
+ * (1 + r_out_val)
+ / (1 - r_out_val)
+ )
+ continue
amp_out_fsolve = []
for i in range(len(initial_guess)):
# Solve the implicit equation using fsolve
@@ -425,7 +438,7 @@ def goodman_haigh_mean_stress_correction( # pylint: disable=R0912 # noqa: C901,
f"{r_out_val}")
amp_out_fsolve.append(np.nan)
else:
- amp_out_fsolve.append(sol if sol < ult_s else ult_s)
+ amp_out_fsolve.append(np.clip(sol, 0.0, ult_s))
# fmt: on
# amp_out.append(amp_out_fsolve)
# mean_out.append(
@@ -444,8 +457,8 @@ def goodman_haigh_mean_stress_correction( # pylint: disable=R0912 # noqa: C901,
else:
srt_idx = np.argsort(r_out)
srt_r_out = r_out[srt_idx]
- srt_amp_out = amp_out[srt_idx, :]
- srt_mean_out = mean_out[srt_idx, :]
+ srt_amp_out = cast(np.ndarray, amp_out[srt_idx, :])
+ srt_mean_out = cast(np.ndarray, mean_out[srt_idx, :])
if plot:
# Create figure with two subplots
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
diff --git a/py_fatigue/testing.py b/py_fatigue/testing.py
index 870eb6f..b996c72 100644
--- a/py_fatigue/testing.py
+++ b/py_fatigue/testing.py
@@ -1,5 +1,5 @@
from __future__ import annotations
-from typing import Any
+from typing import Any, cast
import warnings
import numpy as np
@@ -46,7 +46,10 @@ def get_sampled_time(
warnings.formatwarning = py_fatigue_formatwarning
warnings.warn(w_msg, UserWarning)
- return np.arange(start, start + duration, 1 / fs)
+ return cast(
+ np.ndarray[Any, np.dtype[np.float64]],
+ np.arange(start, start + duration, 1 / fs),
+ )
def get_random_data(
diff --git a/py_fatigue/utils.py b/py_fatigue/utils.py
index 510ee89..4ad16f6 100644
--- a/py_fatigue/utils.py
+++ b/py_fatigue/utils.py
@@ -6,8 +6,11 @@
# Packages from the Python Standard Library
from __future__ import annotations
+from contextlib import redirect_stdout
from dataclasses import dataclass
-from functools import wraps
+from functools import lru_cache, wraps
+from importlib import import_module
+from io import StringIO
from types import FunctionType
from typing import (
Any,
@@ -25,6 +28,7 @@
)
import copy
import logging
+import os
# Packages from non-standard libraries
# from pydantic.fields import ModelField
@@ -32,6 +36,9 @@
import matplotlib.pyplot as plt
import numba as nb
import numpy as np
+from numba.extending import overload
+
+NUMBA_SPECIALIZED_ATTR = "__py_fatigue_numba_specialized__"
# Decorator
@@ -613,7 +620,7 @@ def validate(cls, val, field):
def to_numba_dict(
data: dict, key_type: type = str, val_type: type = float
-) -> nb.types.DictType:
+) -> Any:
"""Converts a dictionary to a numba typed dict, provided the output
key and value types.
@@ -631,13 +638,22 @@ def to_numba_dict(
nb.types.DictType
The numba typed dictionary.
"""
+ filtered_items = {
+ key: value
+ for key, value in data.items()
+ if isinstance(value, val_type) and isinstance(key, key_type)
+ }
+ if (
+ getattr(nb.config, "DISABLE_JIT", False)
+ or os.environ.get("NUMBA_DISABLE_JIT") == "1"
+ ):
+ return filtered_items
+
dct = nb.typed.Dict.empty(
- key_type=nb.types.string,
- value_type=nb.float64,
+ key_type=nb.types.unicode_type,
+ value_type=nb.types.float64,
)
- for key, value in data.items():
- if not isinstance(value, val_type) or not isinstance(key, key_type):
- continue
+ for key, value in filtered_items.items():
dct[key] = value
return dct
@@ -735,18 +751,88 @@ def calc_slope_intercept(
# return nb.njit()(python_bisect)
+def scalarize_numba_result(value):
+ """Normalize scalar-like numba residual returns to a float."""
+
+ return value
+
+
+# pylint: disable-next=unreachable
+@nb.njit(nb.float64())
+def raise_non_scalar_numba_result():
+ """Raise the scalar residual error from nopython code."""
+
+ raise TypeError("Compiled function must return a scalar value")
+ return 0.0 # pylint: disable=unreachable
+
+
+@overload(scalarize_numba_result)
+# pylint: disable-next=too-many-return-statements
+def scalarize_numba_result_overload(value):
+ """Compile scalar-like residual return normalization."""
+
+ if isinstance(value, nb.types.Number):
+
+ def scalar_impl(value):
+ return float(value)
+
+ return scalar_impl
+ if isinstance(value, nb.types.UniTuple) and value.count == 1:
+
+ def unituple_scalar_impl(value):
+ return float(value[0])
+
+ return unituple_scalar_impl
+ if isinstance(value, nb.types.UniTuple):
+
+ def unituple_reject_impl(value):
+ if len(value) != 1:
+ return raise_non_scalar_numba_result()
+ return float(value[0])
+
+ return unituple_reject_impl
+ if isinstance(value, nb.types.Tuple) and len(value) == 1:
+
+ def tuple_scalar_impl(value):
+ return float(value[0])
+
+ return tuple_scalar_impl
+ if isinstance(value, nb.types.Tuple):
+
+ def tuple_reject_impl(value):
+ if len(value) != 1:
+ return raise_non_scalar_numba_result()
+ return float(value[0])
+
+ return tuple_reject_impl
+ if isinstance(value, nb.types.Array):
+
+ def array_impl(value):
+ if value.size != 1:
+ raise TypeError("Compiled function must return a scalar value")
+ return float(value.flat[0])
+
+ return array_impl
+ if isinstance(value, (nb.types.List, nb.types.ListType)):
+
+ def list_impl(value):
+ if len(value) != 1:
+ raise TypeError("Compiled function must return a scalar value")
+ return float(value[0])
+
+ return list_impl
+ return None
+
+
+@lru_cache(maxsize=None)
def compile_specialized_bisect(fun):
"""
Returns a compiled bisection implementation for `f`.
"""
compiled_f = nb.njit()(fun)
- def python_bisect(a, b, tol, mxiter, *args):
+ def python_bisect_compat(a, b, tol, mxiter, *args):
its = 0
- len_args = len(args)
-
- if len_args > 3:
- raise ValueError("Too many extra arguments for compiled function")
def ensure_scalar(value):
if isinstance(value, (tuple, list)):
@@ -766,15 +852,7 @@ def ensure_scalar(value):
return value
def evaluate(point):
- if len_args == 0:
- result = compiled_f(point)
- elif len_args == 1:
- result = compiled_f(point, args[0])
- elif len_args == 2:
- result = compiled_f(point, args[0], args[1])
- else:
- result = compiled_f(point, args[0], args[1], args[2])
- return ensure_scalar(result)
+ return ensure_scalar(compiled_f(point, *args))
fa = evaluate(a)
fb = evaluate(b)
@@ -799,9 +877,44 @@ def evaluate(point):
fc = evaluate(c)
return c
- return python_bisect
+ if (
+ getattr(nb.config, "DISABLE_JIT", False)
+ or os.environ.get("NUMBA_DISABLE_JIT") == "1"
+ ):
+ setattr(python_bisect_compat, NUMBA_SPECIALIZED_ATTR, True)
+ return python_bisect_compat
+
+ def python_bisect(a, b, tol, mxiter, *args):
+ its = 0
+ left = float(a)
+ right = float(b)
+ fa = scalarize_numba_result(compiled_f(left, *args))
+ fb = scalarize_numba_result(compiled_f(right, *args))
+
+ if abs(fa) < tol:
+ return left
+ if abs(fb) < tol:
+ return right
+
+ c = (left + right) / 2.0
+ fc = scalarize_numba_result(compiled_f(c, *args))
+
+ while abs(fc) > tol and its < mxiter:
+ its += 1
+ if fa * fc < 0:
+ right = c
+ fb = fc
+ else:
+ left = c
+ fa = fc
+ c = (left + right) / 2.0
+ fc = scalarize_numba_result(compiled_f(c, *args))
+ return c
+
+ return nb.njit()(python_bisect)
+@lru_cache(maxsize=None)
def compile_specialized_newton(fun):
"""
Returns a compiled Newton–Raphson implementation for f that accepts extra
@@ -882,7 +995,9 @@ def numba_bisect(fun, a, b, tol, mxiter, *args):
"""
A wrapper that compiles `f` if it is a regular Python function.
"""
- if isinstance(fun, FunctionType):
+ if isinstance(fun, FunctionType) and not getattr(
+ fun, NUMBA_SPECIALIZED_ATTR, False
+ ):
jit_bisect_func = compile_specialized_bisect(fun)
return jit_bisect_func(a, b, tol, mxiter, *args)
return fun(a, b, tol, mxiter, *args)
@@ -893,12 +1008,128 @@ def numba_newton(fun, x0, tol, mxiter, *args):
A wrapper that compiles f if it is a regular Python function and calls the
Newton–Raphson routine with extra arguments.
"""
- if isinstance(fun, FunctionType):
+ if isinstance(fun, FunctionType) and not getattr(
+ fun, NUMBA_SPECIALIZED_ATTR, False
+ ):
jit_newton_func = compile_specialized_newton(fun)
return jit_newton_func(x0, tol, mxiter, *args)
return fun(x0, tol, mxiter, *args)
+def warmup_numba() -> None:
+ """Compile the common numba dispatchers for the current Python process.
+
+ This function is intentionally opt-in. Calling it can take noticeable time,
+ but it moves first-call compilation cost to application startup. Functions
+ decorated with ``cache=True`` can also reuse numba's disk cache in later
+ processes; jitclasses and dynamically specialized root finders still need an
+ in-process warmup.
+ """
+
+ if (
+ getattr(nb.config, "DISABLE_JIT", False)
+ or os.environ.get("NUMBA_DISABLE_JIT") == "1"
+ ):
+ return
+
+ rainflow_module = import_module(
+ ".cycle_count.rainflow", package=__package__
+ )
+ crack_growth_module = import_module(
+ ".damage.crack_growth", package=__package__
+ )
+ stress_life_module = import_module(
+ ".damage.stress_life", package=__package__
+ )
+ cylinder_module = import_module(".geometry.cylinder", package=__package__)
+ crack_growth_curve_module = import_module(
+ ".material.crack_growth_curve", package=__package__
+ )
+ sn_curve_module = import_module(".material.sn_curve", package=__package__)
+ mean_stress_module = import_module(
+ ".mean_stress.corrections", package=__package__
+ )
+
+ findcross = rainflow_module.findcross
+ findtp = rainflow_module.findtp
+ calc_crack_growth = crack_growth_module.CalcCrackGrowth
+ get_sif = crack_growth_module.get_sif
+ calc_theil_sn_damage = stress_life_module.calc_theil_sn_damage
+ find_sn_curve_intersection = stress_life_module.find_sn_curve_intersection
+ f_hol_cyl_01 = cylinder_module.f_hol_cyl_01
+ paris_curve_cls = crack_growth_curve_module.ParisCurve
+ sn_curve_cls = sn_curve_module.SNCurve
+ goodman_correction = mean_stress_module.goodman_haigh_mean_stress_correction
+
+ sn_curve = sn_curve_cls([3.0, 5.0], [12.0, 15.0], endurance=1e9)
+ stress_range = np.ascontiguousarray([90.0, 120.0, 180.0], dtype=np.float64)
+ cycles = np.ascontiguousarray([1e5, 1e6, 1e7], dtype=np.float64)
+ sn_curve.get_cycles(stress_range)
+ sn_curve.get_stress(cycles)
+
+ paris_curve = paris_curve_cls(
+ slope=[2.88, 5.1],
+ intercept=[1e-16, 1e-20],
+ threshold=20.0,
+ critical=2_000.0,
+ )
+ sif_range = np.ascontiguousarray([20.0, 100.0, 1_000.0], dtype=np.float64)
+ growth_rate = np.ascontiguousarray([1e-10, 1e-8, 1e-6], dtype=np.float64)
+ paris_curve.get_growth_rate(sif_range)
+ paris_curve.get_sif(growth_rate)
+
+ signal = np.ascontiguousarray([0.0, 1.0, -1.0, 2.0, -0.5], dtype=np.float64)
+ findcross(signal)
+ findtp(signal)
+
+ geometry = to_numba_dict(
+ {
+ "initial_depth": 1.0,
+ "outer_diameter": 100.0,
+ "thickness": 10.0,
+ "height": 200.0,
+ "width_to_depth_ratio": 2.0,
+ }
+ )
+ f_hol_cyl_01(1.0, geometry)
+ get_sif(100.0, 1.0, "HOL_CYL_01", geometry)
+
+ inf_geometry = to_numba_dict({"initial_depth": 1.0})
+ with redirect_stdout(StringIO()):
+ calc_crack_growth(
+ stress_range,
+ np.ones_like(stress_range),
+ np.ascontiguousarray([3.0], dtype=np.float64),
+ np.ascontiguousarray([1e-12], dtype=np.float64),
+ 0.0,
+ 1e9,
+ "INF_SUR_00",
+ inf_geometry,
+ )
+
+ calc_theil_sn_damage(stress_range, np.ones_like(stress_range), sn_curve)
+ find_sn_curve_intersection(
+ sn_curve.slope,
+ sn_curve.intercept,
+ sn_curve.endurance,
+ 100.0,
+ 0.0,
+ 1e3,
+ 1e9,
+ )
+
+ amp_in = np.ascontiguousarray([50.0, 100.0], dtype=np.float64)
+ mean_in = np.ascontiguousarray([0.0, 20.0], dtype=np.float64)
+ r_out = np.ascontiguousarray([-1.0, 0.0], dtype=np.float64)
+ goodman_correction(
+ amp_in,
+ mean_in,
+ r_out,
+ 1_000.0,
+ 3.0,
+ )
+
+
def _plot_damage_accumulation( # pragma: no cover
cumsum_nl_dmg: np.ndarray,
cumsum_pm_dmg: np.ndarray,
@@ -952,14 +1183,41 @@ class CustomFormatter(logging.Formatter):
italic = "\033[3m"
reset = "\033[0m"
level = "\033[1m%(levelname)-8s → \033[22m"
+ lineofile = "\033[3m%(filename)s:%(lineno)d\033[0m - "
+ newline = "\n"
message = "%(message)s"
FORMATS = {
logging.DEBUG: grey + "🐞 " + level + italic + message + reset,
logging.INFO: blue + "ℹ️ " + level + italic + message + reset,
- logging.WARNING: yellow + "⚠️ " + level + italic + message + reset,
- logging.ERROR: red + "⛔ " + level + italic + message + reset,
- logging.CRITICAL: red + "🆘 " + level + bold + italic + message + reset,
+ logging.WARNING: yellow
+ + "⚠️ "
+ + level
+ + newline
+ + lineofile
+ + italic
+ + newline
+ + message
+ + reset,
+ logging.ERROR: red
+ + "⛔ "
+ + level
+ + newline
+ + lineofile
+ + italic
+ + newline
+ + message
+ + reset,
+ logging.CRITICAL: red
+ + "🆘 "
+ + level
+ + newline
+ + lineofile
+ + bold
+ + italic
+ + newline
+ + message
+ + reset,
}
def format(self, record):
diff --git a/py_fatigue/version.py b/py_fatigue/version.py
index d59e028..8fa452e 100644
--- a/py_fatigue/version.py
+++ b/py_fatigue/version.py
@@ -14,6 +14,7 @@
__version__ (str): Current version of py_fatigue package.
"""
+
import re
from typing import Optional, NamedTuple
diff --git a/pyproject.toml b/pyproject.toml
index b1b124b..34ee59a 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -23,16 +23,18 @@ classifiers=["Development Status :: 5 - Production/Stable",
"Programming Language :: Python :: 3.11",
"Programming Language :: Python :: 3.12",
"Programming Language :: Python :: 3.13",
+ "Programming Language :: Python :: 3.14",
"Topic :: Database"]
readme = "README.md"
license = {file = "LICENSE"}
-requires-python = ">=3.10,<3.14"
+requires-python = ">=3.10,<3.15"
dependencies = [
"numpy>=1.24",
"plotly",
"pandas>2.2",
- "numba>=0.61,<0.62",
+ "numba>=0.61,<0.66; python_version < '3.14'",
+ "numba>=0.65,<0.66; python_version >= '3.14'",
"matplotlib",
"pydantic<3.0.0",
]
@@ -52,6 +54,7 @@ dev = [
"hypothesis>=6.14.0",
"notebook>=6.4",
"rich>=13.9.4",
+ "rich-argparse>=1.8.0",
"jsbeautifier>=1.15.4",
]
docs = [
diff --git a/rst_docs/api/material/crack_growth_curve.rst b/rst_docs/api/material/crack_growth_curve.rst
index 1c29023..09c4538 100644
--- a/rst_docs/api/material/crack_growth_curve.rst
+++ b/rst_docs/api/material/crack_growth_curve.rst
@@ -1,18 +1,16 @@
crack_growth_curve
==================
-The crack_growth_curve module contains the ParisCurve class
-and all the funcitons strictly related to its correct functioning.
-The function relies on the AbstractCrackGrowthCurve class, which can be used
-as a abstract class for the definition of new growth models such
-as NASGRO law.
-
-The ParisCurve class
-********************
+The :mod:`py_fatigue.material.crack_growth_curve` module contains the
+:class:`py_fatigue.ParisCurve` and :class:`py_fatigue.WalkerCurve` classes,
+together with the helpers used to evaluate crack-growth curves.
.. autoclass:: py_fatigue.ParisCurve
:members:
:private-members:
:special-members:
- .. autoclasstoc::
+.. autoclass:: py_fatigue.WalkerCurve
+ :members:
+ :private-members:
+ :special-members:
diff --git a/rst_docs/api/material/sn_curve.rst b/rst_docs/api/material/sn_curve.rst
index 8ccf1e3..0a54c5b 100644
--- a/rst_docs/api/material/sn_curve.rst
+++ b/rst_docs/api/material/sn_curve.rst
@@ -1,18 +1,10 @@
sn_curve
===================
-The sn_curve module contains the SNCurve class
-and all the funcitons strictly related to its correct functioning.
-The function relies on the AbstractSNCurve class, which can be used
-as a abstract class for the definition of new SN curve models such
-as smart SN curves.
-
-The SNCurve class
-*****************
+The :mod:`py_fatigue.material.sn_curve` module contains the
+:class:`py_fatigue.SNCurve` class and the helpers used to evaluate S-N curves.
.. autoclass:: py_fatigue.SNCurve
:members:
:private-members:
:special-members:
-
- .. autoclasstoc::
diff --git a/rst_docs/user/01-absolute-noob.rst b/rst_docs/user/01-absolute-noob.rst
index e48582d..11337b9 100644
--- a/rst_docs/user/01-absolute-noob.rst
+++ b/rst_docs/user/01-absolute-noob.rst
@@ -241,11 +241,11 @@ Cycle-count matrix
import matplotlib.pyplot as plt
fig, axs = plt.subplots(1, 2, figsize=(12, 4.5))
cycle_count.plot_histogram(fig=fig, ax=axs[0], plot_type="mean-range",
- marker="d", s=2, cmap=mpl.cm.get_cmap("coolwarm"))
+ marker="d", s=2, cmap=plt.get_cmap("coolwarm"))
axs[0].set_title("Cycle-count from signal")
cycle_count_d.plot_histogram(fig=fig, ax=axs[1], plot_type="mean-range",
marker="s", s=10, edgecolors="#222",
- cmap=mpl.cm.get_cmap("coolwarm"), linewidth=0.25)
+ cmap=plt.get_cmap("coolwarm"), linewidth=0.25)
axs[1].set_title("Cycle-count from matrix")
plt.show()
diff --git a/rst_docs/user/additional/glossary.rst b/rst_docs/user/additional/glossary.rst
index 927b4d0..1d2e245 100644
--- a/rst_docs/user/additional/glossary.rst
+++ b/rst_docs/user/additional/glossary.rst
@@ -273,7 +273,7 @@ Glossary
See :term:`residuals`.
SHM
- `Stuctural health monitoring `_.
+ `Structural health monitoring `_.
Timestamp
A timestamp is a time information that is associated with a particular
@@ -310,8 +310,7 @@ Glossary
we also have to deal with gaps in the data. Even if the ambition were to
concatenate everything, we would still make errors as such.
- `Marsh *et al.* `_
- have introduced a brilliant approach that significantly
+ Marsh *et al.* have introduced a brilliant approach that significantly
reduces calculation time without losing accuracy in the final spectrum
histogram, as it can retrieve all the hysteresis cycles caused by LFFD without
needing the a-priori signal concatenation. We have modified the approach to
diff --git a/rst_docs/user/examples/02-sn_curve.rst b/rst_docs/user/examples/02-sn_curve.rst
index 0946370..ba3d6ee 100644
--- a/rst_docs/user/examples/02-sn_curve.rst
+++ b/rst_docs/user/examples/02-sn_curve.rst
@@ -59,10 +59,10 @@ a. Multiple SN curves
---------------------
.. note::
- In this example we define four SN curves for free corrosion as per
- `DNVGL-RP-C203 `_ and plot them using
- `matplotlib `_ and `plotly `_.
- We additionally define a random gaussian stress range-cycles history to
+ In this example we define four SN curves for free corrosion as per
+ DNVGL-RP-C203 and plot them using `matplotlib `_
+ and `plotly `_. We additionally define a random
+ Gaussian stress range-cycles history to
plot against the SN curves defined.
.. code-block:: python
diff --git a/rst_docs/user/examples/03-cg_curve.rst b/rst_docs/user/examples/03-cg_curve.rst
index 2b73eda..2c51169 100644
--- a/rst_docs/user/examples/03-cg_curve.rst
+++ b/rst_docs/user/examples/03-cg_curve.rst
@@ -114,11 +114,10 @@ a. Definition of multiple Paris' laws
-------------------------------------
.. note::
- In this example we define four SN curves for free corrosion as per
- `DNVGL-RP-C203 `_ and plot them using
- `matplotlib `_ and `plotly `_.
- We additionally define a random gaussian stress range-cycles history to
- plot against the SN curves defined.
+ In this example we define four Paris curves for free corrosion as per
+ DNVGL-RP-C203 and plot them using `matplotlib `_
+ and `plotly `_. We additionally define a random
+ Gaussian stress range-cycles history to plot against the curves defined.
.. code-block:: python
:linenos:
diff --git a/rst_docs/user/examples/04-cycle_counting.rst b/rst_docs/user/examples/04-cycle_counting.rst
index c9ab3eb..5065a63 100644
--- a/rst_docs/user/examples/04-cycle_counting.rst
+++ b/rst_docs/user/examples/04-cycle_counting.rst
@@ -261,7 +261,7 @@ the cumulative (or non-cumulative) rainflow matric both for `cc` and `cc_dct`:
cc.plot_histogram(fig=fig, ax=axs[1][0], plot_type="counts-range-cumsum", s=30)
cc_dct.plot_histogram(fig=fig, ax=axs[1][1], plot_type="counts-range",
marker='s', s=20,
- cmap=matplotlib.cm.get_cmap("gnuplot2")) # chg cmap
+ cmap=plt.get_cmap("gnuplot2")) # chg cmap
axs[1][1].set_xscale("log")
axs[1][0].set_xscale("log")
plt.show()
diff --git a/rst_docs/user/examples/06-mean_stress.rst b/rst_docs/user/examples/06-mean_stress.rst
index 4f84f8e..13e2991 100644
--- a/rst_docs/user/examples/06-mean_stress.rst
+++ b/rst_docs/user/examples/06-mean_stress.rst
@@ -37,7 +37,7 @@ DNVGL-RP-C203 correction
---------------------------
Calculates the mean stress correction according to par.
-2.5 of `DNVGL-RP-C203 `_ which includes
+2.5 of DNVGL-RP-C203, which includes
an attenuation factor :math:`p` for the stress ranges if the
following cases:
diff --git a/scripts/benchmark_numba_speedups.py b/scripts/benchmark_numba_speedups.py
new file mode 100644
index 0000000..db5976a
--- /dev/null
+++ b/scripts/benchmark_numba_speedups.py
@@ -0,0 +1,677 @@
+"""Benchmark numba call patterns used by py-fatigue.
+
+The script keeps package logic unchanged. It compares current public call
+paths with equivalent call patterns that hoist JIT compilation out of tight
+loops, and it inventories the numba-decorated functions found in the package.
+"""
+
+from __future__ import annotations
+
+import argparse
+import ast
+import gc
+import json
+import statistics
+import sys
+import time
+from contextlib import redirect_stdout
+from dataclasses import asdict, dataclass
+from io import StringIO
+from pathlib import Path
+from typing import Any, Callable, Iterable
+
+import numpy as np
+from rich.console import Console
+from rich.table import Table
+from rich_argparse import RichHelpFormatter
+
+from py_fatigue.cycle_count.rainflow import findcross, findtp
+from py_fatigue.damage.crack_growth import CalcCrackGrowth, get_sif
+from py_fatigue.damage.stress_life import (
+ calc_theil_sn_damage,
+ find_sn_curve_intersection,
+)
+from py_fatigue.geometry.cylinder import f_hol_cyl_01
+from py_fatigue.material.crack_growth_curve import ParisCurve
+from py_fatigue.material.sn_curve import SNCurve
+from py_fatigue.mean_stress.corrections import (
+ goodman_haigh_mean_stress_correction,
+)
+from py_fatigue.utils import (
+ compile_specialized_bisect,
+ compile_specialized_newton,
+ numba_bisect,
+ numba_newton,
+ py_bisect,
+ py_newton,
+ to_numba_dict,
+)
+
+PROJECT_ROOT = Path(__file__).resolve().parents[1]
+PACKAGE_ROOT = PROJECT_ROOT / "py_fatigue"
+
+
+@dataclass(frozen=True)
+class BenchmarkResult:
+ """A single benchmark row."""
+
+ group: str
+ name: str
+ calls: int
+ seconds: float | None
+ per_call_ms: float | None
+ speedup_vs_group_base: float | None = None
+ status: str = "ok"
+ note: str = ""
+
+
+@dataclass(frozen=True)
+class InventoryItem:
+ """A numba-related symbol discovered by AST inspection."""
+
+ file: str
+ line: int
+ kind: str
+ name: str
+ mechanism: str
+
+
+def quadratic_bisect(x_value: float) -> float:
+ """Simple scalar residual for bisection benchmarks."""
+
+ return x_value * x_value - 2.0
+
+
+def quadratic_newton(x_value: float, target: float) -> float:
+ """Simple scalar residual for Newton benchmarks."""
+
+ return x_value * x_value - target
+
+
+def time_repeated(
+ func: Callable[[], Any],
+ calls: int,
+ repeats: int,
+) -> float:
+ """Return the median wall time for repeated calls."""
+
+ timings: list[float] = []
+ was_enabled = gc.isenabled()
+ gc.disable()
+ try:
+ for _ in range(repeats):
+ start = time.perf_counter()
+ for _ in range(calls):
+ func()
+ timings.append(time.perf_counter() - start)
+ finally:
+ if was_enabled:
+ gc.enable()
+ return statistics.median(timings)
+
+
+def add_result(
+ results: list[BenchmarkResult],
+ group: str,
+ name: str,
+ calls: int,
+ seconds: float | None,
+ baseline: float | None = None,
+ status: str = "ok",
+ note: str = "",
+) -> None:
+ """Append a benchmark result with derived per-call and speedup values."""
+
+ per_call_ms = None if seconds is None else seconds / calls * 1_000
+ speedup = None
+ if seconds is not None and baseline is not None and seconds > 0:
+ speedup = baseline / seconds
+ results.append(
+ BenchmarkResult(
+ group=group,
+ name=name,
+ calls=calls,
+ seconds=seconds,
+ per_call_ms=per_call_ms,
+ speedup_vs_group_base=speedup,
+ status=status,
+ note=note,
+ )
+ )
+
+
+def measure_case(
+ results: list[BenchmarkResult],
+ group: str,
+ name: str,
+ func: Callable[[], Any],
+ calls: int,
+ repeats: int,
+ baseline: float | None = None,
+) -> float | None:
+ """Run a timed case and append the result."""
+
+ try:
+ elapsed = time_repeated(func, calls=calls, repeats=repeats)
+ except Exception as exc: # pylint: disable=broad-exception-caught
+ add_result(
+ results,
+ group,
+ name,
+ calls,
+ None,
+ baseline=baseline,
+ status="error",
+ note=f"{type(exc).__name__}: {exc}",
+ )
+ return None
+
+ add_result(results, group, name, calls, elapsed, baseline=baseline)
+ return elapsed
+
+
+def benchmark_root_finders(
+ calls: int,
+ repeats: int,
+) -> list[BenchmarkResult]:
+ """Benchmark current root-finder wrappers against hoisted compilation."""
+
+ results: list[BenchmarkResult] = []
+
+ bisect_base = measure_case(
+ results,
+ "bisect-wrapper",
+ "current numba_bisect(py_func)",
+ lambda: numba_bisect(quadratic_bisect, 0.0, 2.0, 1e-8, 100),
+ calls,
+ repeats,
+ )
+
+ solver = compile_specialized_bisect(quadratic_bisect)
+ solver(0.0, 2.0, 1e-8, 100)
+ measure_case(
+ results,
+ "bisect-wrapper",
+ "hoisted compile_specialized_bisect",
+ lambda: solver(0.0, 2.0, 1e-8, 100),
+ calls,
+ repeats,
+ baseline=bisect_base,
+ )
+ measure_case(
+ results,
+ "bisect-wrapper",
+ "pure Python py_bisect",
+ lambda: py_bisect(quadratic_bisect, 0.0, 2.0, 1e-8, 100),
+ calls,
+ repeats,
+ baseline=bisect_base,
+ )
+
+ newton_base = measure_case(
+ results,
+ "newton-wrapper",
+ "current numba_newton(py_func)",
+ lambda: numba_newton(quadratic_newton, 1.0, 1e-8, 100, 2.0),
+ calls,
+ repeats,
+ )
+
+ newton_solver = compile_specialized_newton(quadratic_newton)
+ newton_solver(1.0, 1e-8, 100, 2.0)
+ measure_case(
+ results,
+ "newton-wrapper",
+ "hoisted compile_specialized_newton",
+ lambda: newton_solver(1.0, 1e-8, 100, 2.0),
+ calls,
+ repeats,
+ baseline=newton_base,
+ )
+ measure_case(
+ results,
+ "newton-wrapper",
+ "pure Python py_newton",
+ lambda: py_newton(quadratic_newton, 1.0, 1e-8, 100, 2.0),
+ calls,
+ repeats,
+ baseline=newton_base,
+ )
+
+ return results
+
+
+def benchmark_public_paths(
+ calls: int,
+ repeats: int,
+) -> list[BenchmarkResult]:
+ """Benchmark public workflows that exercise each numba area."""
+
+ results: list[BenchmarkResult] = []
+
+ sn_curve = SNCurve([4, 5], [15.117, 17.146], endurance=1e9)
+ stress_range = np.linspace(80.0, 240.0, 10_000)
+ cycles = np.logspace(5, 9, 10_000)
+ sn_curve.get_cycles(stress_range)
+ sn_curve.get_stress(cycles)
+
+ sn_base = measure_case(
+ results,
+ "direct-dispatchers",
+ "SNCurve.get_cycles",
+ lambda: sn_curve.get_cycles(stress_range),
+ calls,
+ repeats,
+ )
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "SNCurve.get_stress",
+ lambda: sn_curve.get_stress(cycles),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+
+ paris_curve = ParisCurve(
+ slope=[2.88, 5.1, 8.16, 5.1, 2.88],
+ intercept=[1e-16, 1e-20, 1e-27, 1e-19, 1e-13],
+ threshold=20.0,
+ critical=2_000.0,
+ )
+ sif_range = np.linspace(20.0, 2_000.0, 10_000)
+ growth_rate = np.logspace(-10, -4, 10_000)
+ paris_curve.get_growth_rate(sif_range)
+ paris_curve.get_sif(growth_rate)
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "ParisCurve.get_growth_rate",
+ lambda: paris_curve.get_growth_rate(sif_range),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "ParisCurve.get_sif",
+ lambda: paris_curve.get_sif(growth_rate),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+
+ signal = np.sin(np.linspace(0.0, 2_000.0, 20_000))
+ signal += 0.15 * np.sin(np.linspace(0.0, 40_000.0, 20_000))
+ findcross(signal)
+ findtp(signal)
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "rainflow.findcross",
+ lambda: findcross(signal),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "rainflow.findtp",
+ lambda: findtp(signal),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+
+ geometry = to_numba_dict(
+ {
+ "initial_depth": 1.0,
+ "outer_diameter": 100.0,
+ "thickness": 10.0,
+ "height": 200.0,
+ "width_to_depth_ratio": 2.0,
+ }
+ )
+ f_hol_cyl_01(1.0, geometry)
+ get_sif(100.0, 1.0, "HOL_CYL_01", geometry)
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "f_hol_cyl_01",
+ lambda: f_hol_cyl_01(1.0, geometry),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "crack_growth.get_sif",
+ lambda: get_sif(100.0, 1.0, "HOL_CYL_01", geometry),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+
+ cg_stress = np.full(500, 100.0, dtype=np.float64)
+ count_cycle = np.ones(500, dtype=np.float64)
+ slope = np.array([3.0], dtype=np.float64)
+ intercept = np.array([1e-12], dtype=np.float64)
+ inf_geometry = to_numba_dict({"initial_depth": 1.0})
+
+ def build_crack_growth() -> CalcCrackGrowth:
+ with redirect_stdout(StringIO()):
+ return CalcCrackGrowth(
+ cg_stress,
+ count_cycle,
+ slope,
+ intercept,
+ 0.0,
+ 1e9,
+ "INF_SUR_00",
+ inf_geometry,
+ )
+
+ build_crack_growth()
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "CalcCrackGrowth construction",
+ build_crack_growth,
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+
+ small_stress = np.array([120.0, 180.0, 90.0], dtype=np.float64)
+ small_counts = np.array([10.0, 8.0, 12.0], dtype=np.float64)
+ calc_theil_sn_damage(small_stress, small_counts, sn_curve)
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "calc_theil_sn_damage",
+ lambda: calc_theil_sn_damage(small_stress, small_counts, sn_curve),
+ calls,
+ repeats,
+ baseline=sn_base,
+ )
+
+ amp_in = np.linspace(20.0, 300.0, 40)
+ mean_in = np.linspace(-10.0, 60.0, 40)
+ r_out = np.array([-1.0, -0.5, 0.0, 0.5])
+ goodman_haigh_mean_stress_correction(
+ amp_in,
+ mean_in,
+ r_out,
+ 1_000.0,
+ 3.0,
+ )
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "goodman_haigh_mean_stress_correction",
+ lambda: goodman_haigh_mean_stress_correction(
+ amp_in,
+ mean_in,
+ r_out,
+ 1_000.0,
+ 3.0,
+ ),
+ max(1, calls // 5),
+ repeats,
+ baseline=sn_base,
+ )
+
+ try:
+ find_sn_curve_intersection(
+ sn_curve.slope,
+ sn_curve.intercept,
+ sn_curve.endurance,
+ 0.1,
+ 10.0,
+ 1.0,
+ 1e15,
+ )
+ except Exception as exc: # pylint: disable=broad-exception-caught
+ add_result(
+ results,
+ "direct-dispatchers",
+ "find_sn_curve_intersection",
+ 1,
+ None,
+ baseline=sn_base,
+ status="error",
+ note=f"{type(exc).__name__}: {exc}",
+ )
+ else:
+ measure_case(
+ results,
+ "direct-dispatchers",
+ "find_sn_curve_intersection",
+ lambda: find_sn_curve_intersection(
+ sn_curve.slope,
+ sn_curve.intercept,
+ sn_curve.endurance,
+ 0.1,
+ 10.0,
+ 1.0,
+ 1e15,
+ ),
+ max(1, calls // 10),
+ repeats,
+ baseline=sn_base,
+ )
+
+ return results
+
+
+def decorator_text(node: ast.AST) -> str:
+ """Return best-effort source text for a decorator or call."""
+
+ try:
+ return ast.unparse(node)
+ except Exception: # pragma: no cover
+ return type(node).__name__
+
+
+def iter_python_files(paths: Iterable[Path]) -> Iterable[Path]:
+ """Yield Python files below the provided paths."""
+
+ for path in paths:
+ if path.is_file() and path.suffix == ".py":
+ yield path
+ elif path.is_dir():
+ yield from sorted(path.rglob("*.py"))
+
+
+def inspect_inventory() -> list[InventoryItem]:
+ """Inspect package sources for numba decorators and wrapper call sites."""
+
+ items: list[InventoryItem] = []
+ for py_file in iter_python_files([PACKAGE_ROOT]):
+ tree = ast.parse(py_file.read_text(encoding="utf-8"))
+ rel_path = py_file.relative_to(PROJECT_ROOT).as_posix()
+ for node in ast.walk(tree):
+ if isinstance(
+ node, (ast.FunctionDef, ast.AsyncFunctionDef, ast.ClassDef)
+ ):
+ for decorator in node.decorator_list:
+ text = decorator_text(decorator)
+ if any(
+ token in text
+ for token in ("njit", "jitclass", "vectorize")
+ ):
+ items.append(
+ InventoryItem(
+ file=rel_path,
+ line=node.lineno,
+ kind=type(node).__name__,
+ name=node.name,
+ mechanism=f"@{text}",
+ )
+ )
+ if isinstance(node, ast.Assign) and isinstance(
+ node.value, ast.Call
+ ):
+ text = decorator_text(node.value.func)
+ if text.endswith(
+ (
+ "compile_specialized_bisect",
+ "compile_specialized_newton",
+ )
+ ):
+ for target in node.targets:
+ if isinstance(target, ast.Name):
+ items.append(
+ InventoryItem(
+ file=rel_path,
+ line=node.lineno,
+ kind="Assign",
+ name=target.id,
+ mechanism=text,
+ )
+ )
+ if isinstance(node, ast.Call):
+ text = decorator_text(node.func)
+ if text in {"numba_bisect", "numba_newton"}:
+ items.append(
+ InventoryItem(
+ file=rel_path,
+ line=node.lineno,
+ kind="Call",
+ name=text,
+ mechanism="generic wrapper call",
+ )
+ )
+ return sorted(items, key=lambda item: (item.file, item.line, item.name))
+
+
+def render_inventory(console: Console, items: list[InventoryItem]) -> None:
+ """Print the numba inventory."""
+
+ table = Table(title="Numba inventory")
+ table.add_column("File")
+ table.add_column("Line", justify="right")
+ table.add_column("Kind")
+ table.add_column("Name")
+ table.add_column("Mechanism")
+ for item in items:
+ table.add_row(
+ item.file,
+ str(item.line),
+ item.kind,
+ item.name,
+ item.mechanism,
+ )
+ console.print(table)
+
+
+def render_results(console: Console, results: list[BenchmarkResult]) -> None:
+ """Print benchmark results."""
+
+ table = Table(title="Numba benchmark results")
+ table.add_column("Group")
+ table.add_column("Case")
+ table.add_column("Calls", justify="right")
+ table.add_column("Total (s)", justify="right")
+ table.add_column("ms/call", justify="right")
+ table.add_column("Speedup", justify="right")
+ table.add_column("Status")
+ table.add_column("Note")
+ for result in results:
+ total = "-" if result.seconds is None else f"{result.seconds:.6f}"
+ per_call = (
+ "-" if result.per_call_ms is None else f"{result.per_call_ms:.3f}"
+ )
+ speedup = (
+ "-"
+ if result.speedup_vs_group_base is None
+ else f"{result.speedup_vs_group_base:.1f}x"
+ )
+ table.add_row(
+ result.group,
+ result.name,
+ str(result.calls),
+ total,
+ per_call,
+ speedup,
+ result.status,
+ result.note,
+ )
+ console.print(table)
+
+
+def parse_args(argv: list[str]) -> argparse.Namespace:
+ """Parse command-line arguments."""
+
+ parser = argparse.ArgumentParser(
+ description=__doc__,
+ formatter_class=RichHelpFormatter,
+ )
+ parser.add_argument(
+ "--calls",
+ type=int,
+ default=10,
+ help="Calls per benchmark repeat. Higher values amplify differences.",
+ )
+ parser.add_argument(
+ "--repeats",
+ type=int,
+ default=3,
+ help="Number of repeats used to compute the median wall time.",
+ )
+ parser.add_argument(
+ "--skip-inventory",
+ action="store_true",
+ help="Do not print the static numba inventory.",
+ )
+ parser.add_argument(
+ "--skip-benchmarks",
+ action="store_true",
+ help="Do not run timing benchmarks.",
+ )
+ parser.add_argument(
+ "--json-output",
+ type=Path,
+ default=None,
+ help="Optional JSON file for machine-readable results.",
+ )
+ return parser.parse_args(argv)
+
+
+def main(argv: list[str] | None = None) -> int:
+ """Run the inventory and benchmarks."""
+
+ args = parse_args(sys.argv[1:] if argv is None else argv)
+ console = Console()
+
+ inventory: list[InventoryItem] = []
+ results: list[BenchmarkResult] = []
+
+ if not args.skip_inventory:
+ inventory = inspect_inventory()
+ render_inventory(console, inventory)
+
+ if not args.skip_benchmarks:
+ results.extend(benchmark_root_finders(args.calls, args.repeats))
+ results.extend(benchmark_public_paths(args.calls, args.repeats))
+ render_results(console, results)
+
+ if args.json_output is not None:
+ payload = {
+ "inventory": [asdict(item) for item in inventory],
+ "results": [asdict(result) for result in results],
+ }
+ args.json_output.write_text(
+ json.dumps(payload, indent=2),
+ encoding="utf-8",
+ )
+ console.print(f"Wrote JSON results to {args.json_output}")
+
+ return 0
+
+
+if __name__ == "__main__":
+ raise SystemExit(main())
diff --git a/scripts/numba_benchmark_results.json b/scripts/numba_benchmark_results.json
new file mode 100644
index 0000000..f85d780
--- /dev/null
+++ b/scripts/numba_benchmark_results.json
@@ -0,0 +1,326 @@
+{
+ "inventory": [
+ {
+ "file": "py_fatigue/cycle_count/rainflow.py",
+ "line": 209,
+ "kind": "FunctionDef",
+ "name": "_findcross",
+ "mechanism": "@njit(int64(int64[:], int8[:]), cache=True)"
+ },
+ {
+ "file": "py_fatigue/damage/crack_growth.py",
+ "line": 47,
+ "kind": "ClassDef",
+ "name": "CalcCrackGrowth",
+ "mechanism": "@jitclass(spec)"
+ },
+ {
+ "file": "py_fatigue/damage/crack_growth.py",
+ "line": 333,
+ "kind": "FunctionDef",
+ "name": "get_geometry_factor",
+ "mechanism": "@nb.njit(nb.float64(nb.float64, nb.types.unicode_type, nb.types.DictType(nb.types.unicode_type, nb.float64)), fastmath=True, cache=True)"
+ },
+ {
+ "file": "py_fatigue/damage/crack_growth.py",
+ "line": 360,
+ "kind": "FunctionDef",
+ "name": "get_sif",
+ "mechanism": "@nb.njit(nb.types.UniTuple(nb.float64, 2)(nb.float64, nb.float64, nb.types.unicode_type, nb.types.DictType(nb.types.unicode_type, nb.float64)), fastmath=True, cache=True)"
+ },
+ {
+ "file": "py_fatigue/damage/stress_life.py",
+ "line": 1269,
+ "kind": "FunctionDef",
+ "name": "_calc_theil_sn_damage",
+ "mechanism": "@nb.njit(fastmath=True, cache=True)"
+ },
+ {
+ "file": "py_fatigue/geometry/cylinder.py",
+ "line": 246,
+ "kind": "FunctionDef",
+ "name": "f_hol_cyl_01",
+ "mechanism": "@nb.njit(nb.float64(nb.float64, nb.types.DictType(nb.types.unicode_type, nb.float64)), cache=True, fastmath=True)"
+ },
+ {
+ "file": "py_fatigue/material/crack_growth_curve.py",
+ "line": 1077,
+ "kind": "FunctionDef",
+ "name": "get_array_min",
+ "mechanism": "@nb.njit(cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/crack_growth_curve.py",
+ "line": 1093,
+ "kind": "FunctionDef",
+ "name": "_calc_growth_rate",
+ "mechanism": "@nb.njit(fastmath=False, cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/crack_growth_curve.py",
+ "line": 1139,
+ "kind": "FunctionDef",
+ "name": "_calc_sif",
+ "mechanism": "@nb.njit(fastmath=False, cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/sn_curve.py",
+ "line": 956,
+ "kind": "FunctionDef",
+ "name": "get_sn_array_min",
+ "mechanism": "@nb.njit(cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/sn_curve.py",
+ "line": 972,
+ "kind": "FunctionDef",
+ "name": "_calc_cycles",
+ "mechanism": "@nb.njit(fastmath=False, cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/sn_curve.py",
+ "line": 996,
+ "kind": "FunctionDef",
+ "name": "_calc_cycles_2",
+ "mechanism": "@nb.njit(fastmath=False, cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/sn_curve.py",
+ "line": 1051,
+ "kind": "FunctionDef",
+ "name": "_calc_stress",
+ "mechanism": "@nb.njit(fastmath=False, cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/sn_curve.py",
+ "line": 1078,
+ "kind": "FunctionDef",
+ "name": "_calc_stress_2",
+ "mechanism": "@nb.njit(fastmath=False, cache=True)"
+ },
+ {
+ "file": "py_fatigue/material/sn_curve.py",
+ "line": 1183,
+ "kind": "Assign",
+ "name": "__jit_sn_curve_residuals",
+ "mechanism": "compile_specialized_bisect"
+ },
+ {
+ "file": "py_fatigue/mean_stress/corrections.py",
+ "line": 423,
+ "kind": "Call",
+ "name": "numba_newton",
+ "mechanism": "generic wrapper call"
+ },
+ {
+ "file": "py_fatigue/mean_stress/corrections.py",
+ "line": 582,
+ "kind": "Assign",
+ "name": "__jit_goodman_equation",
+ "mechanism": "compile_specialized_newton"
+ },
+ {
+ "file": "py_fatigue/utils.py",
+ "line": 762,
+ "kind": "FunctionDef",
+ "name": "raise_non_scalar_numba_result",
+ "mechanism": "@nb.njit(nb.float64())"
+ },
+ {
+ "file": "py_fatigue/utils.py",
+ "line": 1001,
+ "kind": "Assign",
+ "name": "jit_bisect_func",
+ "mechanism": "compile_specialized_bisect"
+ },
+ {
+ "file": "py_fatigue/utils.py",
+ "line": 1014,
+ "kind": "Assign",
+ "name": "jit_newton_func",
+ "mechanism": "compile_specialized_newton"
+ }
+ ],
+ "results": [
+ {
+ "group": "bisect-wrapper",
+ "name": "current numba_bisect(py_func)",
+ "calls": 8,
+ "seconds": 1.2343996786512434e-05,
+ "per_call_ms": 0.0015429995983140543,
+ "speedup_vs_group_base": null,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "bisect-wrapper",
+ "name": "hoisted compile_specialized_bisect",
+ "calls": 8,
+ "seconds": 6.1780010582879186e-06,
+ "per_call_ms": 0.0007722501322859898,
+ "speedup_vs_group_base": 1.9980567614103437,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "bisect-wrapper",
+ "name": "pure Python py_bisect",
+ "calls": 8,
+ "seconds": 9.328400483354926e-05,
+ "per_call_ms": 0.011660500604193658,
+ "speedup_vs_group_base": 0.13232704586963617,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "newton-wrapper",
+ "name": "current numba_newton(py_func)",
+ "calls": 8,
+ "seconds": 1.6310994396917522e-05,
+ "per_call_ms": 0.0020388742996146902,
+ "speedup_vs_group_base": null,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "newton-wrapper",
+ "name": "hoisted compile_specialized_newton",
+ "calls": 8,
+ "seconds": 8.942995918914676e-06,
+ "per_call_ms": 0.0011178744898643345,
+ "speedup_vs_group_base": 1.823884808268707,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "newton-wrapper",
+ "name": "pure Python py_newton",
+ "calls": 8,
+ "seconds": 4.595100472215563e-05,
+ "per_call_ms": 0.005743875590269454,
+ "speedup_vs_group_base": 0.3549649130751879,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "SNCurve.get_cycles",
+ "calls": 8,
+ "seconds": 0.0033126759954029694,
+ "per_call_ms": 0.4140844994253712,
+ "speedup_vs_group_base": null,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "SNCurve.get_stress",
+ "calls": 8,
+ "seconds": 0.0035862310032825917,
+ "per_call_ms": 0.44827887541032396,
+ "speedup_vs_group_base": 0.9237207509418025,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "ParisCurve.get_growth_rate",
+ "calls": 8,
+ "seconds": 0.010180105993640609,
+ "per_call_ms": 1.2725132492050761,
+ "speedup_vs_group_base": 0.32540682753915906,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "ParisCurve.get_sif",
+ "calls": 8,
+ "seconds": 0.017416592003428377,
+ "per_call_ms": 2.177074000428547,
+ "speedup_vs_group_base": 0.19020230793434698,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "rainflow.findcross",
+ "calls": 8,
+ "seconds": 0.0010041150089818984,
+ "per_call_ms": 0.1255143761227373,
+ "speedup_vs_group_base": 3.299100168577092,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "rainflow.findtp",
+ "calls": 8,
+ "seconds": 0.0013474489969667047,
+ "per_call_ms": 0.1684311246208381,
+ "speedup_vs_group_base": 2.4584796922631313,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "f_hol_cyl_01",
+ "calls": 8,
+ "seconds": 2.6398993213661015e-05,
+ "per_call_ms": 0.003299874151707627,
+ "speedup_vs_group_base": 125.48493681526907,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "crack_growth.get_sif",
+ "calls": 8,
+ "seconds": 4.6574990847148e-05,
+ "per_call_ms": 0.0058218738558935,
+ "speedup_vs_group_base": 71.12563921428703,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "CalcCrackGrowth construction",
+ "calls": 8,
+ "seconds": 0.005733437996241264,
+ "per_call_ms": 0.716679749530158,
+ "speedup_vs_group_base": 0.5777817772817458,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "calc_theil_sn_damage",
+ "calls": 8,
+ "seconds": 0.00019102801161352545,
+ "per_call_ms": 0.02387850145169068,
+ "speedup_vs_group_base": 17.341310143063964,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "goodman_haigh_mean_stress_correction",
+ "calls": 1,
+ "seconds": 0.0008232050022343174,
+ "per_call_ms": 0.8232050022343174,
+ "speedup_vs_group_base": 4.024120342334907,
+ "status": "ok",
+ "note": ""
+ },
+ {
+ "group": "direct-dispatchers",
+ "name": "find_sn_curve_intersection",
+ "calls": 1,
+ "seconds": 6.591700366698205e-05,
+ "per_call_ms": 0.06591700366698205,
+ "speedup_vs_group_base": 50.25525753777996,
+ "status": "ok",
+ "note": ""
+ }
+ ]
+}
\ No newline at end of file
diff --git a/tasks/test.py b/tasks/test.py
index ec834c6..fa2f0d8 100644
--- a/tasks/test.py
+++ b/tasks/test.py
@@ -1,5 +1,7 @@
"""Test tasks.""" # pylint: disable=R0801
+import shutil
+
from invoke import task
from .colors import Color, colorize
@@ -28,6 +30,7 @@ def run(c_r, test=None, pytest_args="-v -W ignore::UserWarning"):
test_command = f" {test}"
_command = (
+ "MPLBACKEND=Agg PY_FATIGUE_TEST_NO_PLOTS=1 "
f"pytest {pytest_args} "
f"--cov={c_r.project_slug} --cov-report=term:skip-covered "
f"--cov-report=html --cov-report=html:{COV_DOC_BUILD_DIR} {test_command}"
@@ -45,6 +48,16 @@ def coverage(c_r):
COV_PORT = c_r.start_port + 2 # pylint: disable=C0103
if SYSTEM in [OperatingSystem.LINUX, OperatingSystem.MAC]:
+ if shutil.which("screen") is None:
+ print(
+ colorize(
+ "screen is not installed; skipping coverage webserver. "
+ "Use `python -m http.server --directory "
+ f"{COV_DOC_BUILD_DIR} {COV_PORT}` to serve manually.",
+ color=Color.WARNING,
+ )
+ )
+ return
_command = (
f"screen -d -S {COV_SCREEN_NAME} "
"-m python -m http.server --bind localhost "
@@ -86,9 +99,12 @@ def stop(c_r):
"\nStopping coverage server...\n", color=Color.HEADER, bold=True
)
print(tmp_str)
- _command = f"kill $(lsof -ti:{COV_PORT})"
+ _command = (
+ f"pids=$(lsof -ti:{COV_PORT}); "
+ 'if [ -n "$pids" ]; then kill $pids; fi'
+ )
print(f"{colorize('>>> ' + _command, color=Color.OKBLUE)}\n")
- c_r.run(_command)
+ c_r.run(_command, warn=True)
elif SYSTEM == OperatingSystem.WINDOWS:
print(
@@ -113,7 +129,6 @@ def stop(c_r):
def all(
c_r, test=None, pytest_args="-v -W ignore::UserWarning"
): # pylint: disable=W0622
- """Run all tests and start coverage report webserver."""
+ """Run all tests."""
stop(c_r)
run(c_r, test, pytest_args)
- coverage(c_r)
diff --git a/tests/conftest.py b/tests/conftest.py
index eb901ff..e2cc7d1 100644
--- a/tests/conftest.py
+++ b/tests/conftest.py
@@ -57,3 +57,14 @@ def datadir(tmpdir, request):
dir_util.copy_tree(test_dir, str(tmpdir))
return tmpdir
+
+
+@pytest.fixture(autouse=True)
+def disable_plot_show(monkeypatch):
+ """Disable interactive plot rendering for task-based test runs."""
+ if os.getenv("PY_FATIGUE_TEST_NO_PLOTS") != "1":
+ return
+
+ import matplotlib.pyplot as plt
+
+ monkeypatch.setattr(plt, "show", lambda *args, **kwargs: None)
diff --git a/tests/damage/test_stress_life.py b/tests/damage/test_stress_life.py
index 2332ea1..9d54f5a 100644
--- a/tests/damage/test_stress_life.py
+++ b/tests/damage/test_stress_life.py
@@ -4,13 +4,13 @@
damage calculation methods in the stress-life approach.
"""
-
# Standard imports
import datetime as dt
import os
import sys
# Non-standard imports
+import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import pytest
@@ -18,6 +18,7 @@
from hypothesis import strategies as hy
import py_fatigue.damage as damage
+
# Local imports
from py_fatigue import CycleCount, SNCurve
from tests.cycle_count.test_cycle_count import CC_RF_1
@@ -173,7 +174,6 @@ def test_palmgren_miner_constant_load(
cc_obj.unit = "m"
damage.get_pm(cc_obj, sn_curve)
-
@pytest.mark.parametrize(
"sn_curve", [(DNV_B1A), (DNV_B1A_END), (DNV_B1W), (DNV_B1C)]
)
@@ -198,9 +198,7 @@ def test_palmgren_miner_variable_load(self, sn_curve: SNCurve):
)
assert damage_pm == pytest.approx(damage_ref, 1e-12)
- @pytest.mark.parametrize(
- "cc", [(CC_TS_1), (CC_TS_2), (CC_TS_3), (CC_RF_1)]
- )
+ @pytest.mark.parametrize("cc", [(CC_TS_1), (CC_TS_2), (CC_TS_3), (CC_RF_1)])
@pytest.mark.parametrize("exponent", [(3), (4), (5)])
@pytest.mark.parametrize("eq_cycles", [(1e6), (2e6), (1e7)])
# fmt: on
@@ -266,7 +264,7 @@ def test_damage_equivalent_moment(
assert "outer_radius must be greater" in ve.value.args[0]
a_i = np.pi / 4 * (r_o**4 - r_i**4)
dem = damage.get_dem(r_o, r_i, cc, exponent, eq_cycles)
- dem_ref = a_i * 1E6 * damage.get_des(cc, exponent, eq_cycles) / r_o
+ dem_ref = a_i * 1e6 * damage.get_des(cc, exponent, eq_cycles) / r_o
assert dem == pytest.approx(dem_ref, 1e-12)
# fmt: off
@@ -309,10 +307,14 @@ def test_miner_pandas_accessor_constant_load(
assert df_d.sn_curve == sn_curve
assert df_d._metadata["name"] == cc_obj.name
assert df_d._metadata["timestamp"] == cc_obj.timestamp
- assert df_d._metadata["mean_stress_corrected"] == \
- cc_obj.mean_stress_corrected
- assert df_d._metadata["stress_concentration_factor"] == \
- cc_obj.stress_concentration_factor
+ assert (
+ df_d._metadata["mean_stress_corrected"]
+ == cc_obj.mean_stress_corrected
+ )
+ assert (
+ df_d._metadata["stress_concentration_factor"]
+ == cc_obj.stress_concentration_factor
+ )
assert df_d._metadata["nr_small_cycles"] == cc_obj.nr_small_cycles
assert df_d._metadata["lffd_solved"] == cc_obj.lffd_solved
assert isinstance(df_d, pd.DataFrame)
@@ -327,9 +329,7 @@ def test_miner_pandas_accessor_constant_load(
@pytest.mark.parametrize(
"sn_curve", [DNV_B1A, DNV_B1A_END, DNV_B1W, DNV_B1C]
)
- @pytest.mark.parametrize(
- "cc", [CC_TS_1, CC_TS_3]
- )
+ @pytest.mark.parametrize("cc", [CC_TS_1, CC_TS_3])
# fmt: on
def test_miner_pandas_accessor_variable_load(
self, cc: CycleCount, sn_curve: SNCurve
@@ -353,23 +353,29 @@ def test_miner_pandas_accessor_variable_load(
assert df_d.sn_curve == sn_curve
assert df_d._metadata["name"] == cc.name
assert df_d._metadata["timestamp"] == cc.timestamp
- assert df_d._metadata["mean_stress_corrected"] == \
- cc.mean_stress_corrected
- assert df_d._metadata["stress_concentration_factor"] == \
- cc.stress_concentration_factor
+ assert (
+ df_d._metadata["mean_stress_corrected"] == cc.mean_stress_corrected
+ )
+ assert (
+ df_d._metadata["stress_concentration_factor"]
+ == cc.stress_concentration_factor
+ )
assert df_d._metadata["nr_small_cycles"] == cc.nr_small_cycles
assert df_d._metadata["lffd_solved"] == cc.lffd_solved
assert isinstance(df_d, pd.DataFrame)
assert df_d["pm_damage"].sum() == pytest.approx(
np.sum(damage.get_pm(cc, sn_curve)), 1e-12
)
+ fig, ax = plt.subplots()
+ df_d.miner.plot_histogram(fig=fig, ax=ax)
+ plt.close(fig)
+ with pytest.raises(AttributeError):
+ df_d.miner.damage(sn_curve)
- @pytest.mark.parametrize(
- "cc,", [CC_TS_1, CC_TS_3]
- )
+ @pytest.mark.parametrize("cc,", [CC_TS_1, CC_TS_3])
@given(
slope=hy.floats(min_value=3, max_value=20),
- n_eq=hy.floats(min_value=1E5, max_value=1e10),
+ n_eq=hy.floats(min_value=1e5, max_value=1e10),
)
# fmt: on
def test_des_pandas_accessor_variable_load(
@@ -394,28 +400,28 @@ def test_des_pandas_accessor_variable_load(
assert isinstance(df, pd.DataFrame)
assert df._metadata["name"] == cc.name
assert df._metadata["timestamp"] == cc.timestamp
- assert df._metadata["mean_stress_corrected"] == \
- cc.mean_stress_corrected
- assert df._metadata["stress_concentration_factor"] == \
- cc.stress_concentration_factor
+ assert df._metadata["mean_stress_corrected"] == cc.mean_stress_corrected
+ assert (
+ df._metadata["stress_concentration_factor"]
+ == cc.stress_concentration_factor
+ )
assert df._metadata["nr_small_cycles"] == cc.nr_small_cycles
assert df._metadata["lffd_solved"] == cc.lffd_solved
assert isinstance(df, pd.DataFrame)
- assert df.miner.des(slope=slope, equivalent_cycles=n_eq) == \
- pytest.approx(
- damage.get_des(cc, slope, equivalent_cycles=n_eq), 1e-12
- )
+ assert df.miner.des(
+ slope=slope, equivalent_cycles=n_eq
+ ) == pytest.approx(
+ damage.get_des(cc, slope, equivalent_cycles=n_eq), 1e-12
+ )
- @pytest.mark.parametrize(
- "cc,", [CC_TS_1, CC_TS_3]
- )
+ @pytest.mark.parametrize("cc,", [CC_TS_1, CC_TS_3])
@given(
slope=hy.floats(min_value=3, max_value=20),
- n_eq=hy.floats(min_value=1E5, max_value=1e10),
+ n_eq=hy.floats(min_value=1e5, max_value=1e10),
)
@pytest.mark.parametrize(
"r_i, r_o",
- [(3, 3.5), (4000, 4100), pytest.param(2, 1, marks=pytest.mark.xfail)]
+ [(3, 3.5), (4000, 4100), pytest.param(2, 1, marks=pytest.mark.xfail)],
)
# fmt: on
def test_dem_pandas_accessor_variable_load(
@@ -444,10 +450,11 @@ def test_dem_pandas_accessor_variable_load(
assert isinstance(df, pd.DataFrame)
assert df._metadata["name"] == cc.name
assert df._metadata["timestamp"] == cc.timestamp
- assert df._metadata["mean_stress_corrected"] == \
- cc.mean_stress_corrected
- assert df._metadata["stress_concentration_factor"] == \
- cc.stress_concentration_factor
+ assert df._metadata["mean_stress_corrected"] == cc.mean_stress_corrected
+ assert (
+ df._metadata["stress_concentration_factor"]
+ == cc.stress_concentration_factor
+ )
assert df._metadata["nr_small_cycles"] == cc.nr_small_cycles
assert df._metadata["lffd_solved"] == cc.lffd_solved
assert isinstance(df, pd.DataFrame)
@@ -455,23 +462,23 @@ def test_dem_pandas_accessor_variable_load(
outer_radius=r_o,
inner_radius=r_i,
slope=slope,
- equivalent_cycles=n_eq
+ equivalent_cycles=n_eq,
) == pytest.approx(
damage.get_dem(
outer_radius=r_o,
inner_radius=r_i,
cycle_count=cc,
slope=slope,
- equivalent_cycles=n_eq
+ equivalent_cycles=n_eq,
),
- 1e-12
+ 1e-12,
)
+
class TestGassner:
- """Test the shift factor calculation related with the Gassner curve
- """
+ """Test the shift factor calculation related with the Gassner curve"""
- # fmt: off
+ # fmt: off
@settings(deadline=None)
@pytest.mark.parametrize(
"sn_curve", [DNV_B1C, DNV_C_C, DNV_E_C,
@@ -499,6 +506,7 @@ def test_g_pandas_accessor_constant_load(
The length of the history to use.
"""
import py_fatigue as pf
+
time_reversals = [peak * (-1) ** _ for _ in range(len_hist)]
cc_obj = pf.CycleCount.from_timeseries(
time_reversals,
@@ -515,26 +523,34 @@ def test_g_pandas_accessor_constant_load(
assert df_g.sn_curve == sn_curve
assert df_g._metadata["name"] == cc_obj.name
assert df_g._metadata["timestamp"] == cc_obj.timestamp
- assert df_g._metadata["mean_stress_corrected"] == \
- cc_obj.mean_stress_corrected
- assert df_g._metadata["stress_concentration_factor"] == \
- cc_obj.stress_concentration_factor
+ assert (
+ df_g._metadata["mean_stress_corrected"]
+ == cc_obj.mean_stress_corrected
+ )
+ assert (
+ df_g._metadata["stress_concentration_factor"]
+ == cc_obj.stress_concentration_factor
+ )
assert df_g._metadata["nr_small_cycles"] == cc_obj.nr_small_cycles
assert df_g._metadata["lffd_solved"] == cc_obj.lffd_solved
assert isinstance(df_g, pd.DataFrame)
assert df_g["shift_factor"].sum() == pytest.approx(1, 1e-12)
+ @pytest.mark.parametrize("cc", [CC_TS_1, CC_TS_3])
@pytest.mark.parametrize(
- "cc", [CC_TS_1, CC_TS_3]
- )
- @pytest.mark.parametrize(
- "sn_curve", [DNV_B1C, DNV_C_C, DNV_E_C,
- pytest.param(DNV_B1A, marks=pytest.mark.xfail)
- ]
+ "sn_curve",
+ [
+ DNV_B1C,
+ DNV_C_C,
+ DNV_E_C,
+ pytest.param(DNV_B1A, marks=pytest.mark.xfail),
+ ],
)
# fmt: on
def test_g_pandas_accessor_variable_load(
- self, cc: CycleCount, sn_curve: SNCurve,
+ self,
+ cc: CycleCount,
+ sn_curve: SNCurve,
):
"""Test the Gassner shift factor that has to be less than one
for variable amplitude stress histories.
@@ -551,19 +567,24 @@ def test_g_pandas_accessor_variable_load(
assert isinstance(df, pd.DataFrame)
assert df._metadata["name"] == cc.name
assert df._metadata["timestamp"] == cc.timestamp
- assert df._metadata["mean_stress_corrected"] == \
- cc.mean_stress_corrected
- assert df._metadata["stress_concentration_factor"] == \
- cc.stress_concentration_factor
+ assert df._metadata["mean_stress_corrected"] == cc.mean_stress_corrected
+ assert (
+ df._metadata["stress_concentration_factor"]
+ == cc.stress_concentration_factor
+ )
assert df._metadata["nr_small_cycles"] == cc.nr_small_cycles
assert df._metadata["lffd_solved"] == cc.lffd_solved
assert df["shift_factor"].sum() < 1
+
@pytest.mark.parametrize("sn_curve", [DNV_B1C, DNV_C_C, DNV_E_C, DNV_B1A])
-@pytest.mark.parametrize("load", [
- [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 3, -3, 3, -3, 3, -3, 3, -3],
- [3, -3, 3, -3, 3, -3, 3, -3, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1]
-])
+@pytest.mark.parametrize(
+ "load",
+ [
+ [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 3, -3, 3, -3, 3, -3, 3, -3],
+ [3, -3, 3, -3, 3, -3, 3, -3, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1],
+ ],
+)
@pytest.mark.parametrize("base_exponent", [0.99, 1, 1.01])
def test_leve_damage_rule(load: list, sn_curve: SNCurve, base_exponent: float):
"""Test the nonlinear damage calculation.
@@ -585,7 +606,7 @@ def test_leve_damage_rule(load: list, sn_curve: SNCurve, base_exponent: float):
name="Test_CC",
)
d_nl = damage.get_nonlinear_damage(
- 'leve', cc, sn_curve, base_exponent=base_exponent
+ "leve", cc, sn_curve, base_exponent=base_exponent
)
d_l = np.sum(damage.get_pm(cc, sn_curve))
if base_exponent == 1:
@@ -595,6 +616,7 @@ def test_leve_damage_rule(load: list, sn_curve: SNCurve, base_exponent: float):
if base_exponent > 1:
assert d_nl[-1] < d_l
+
@pytest.mark.parametrize("sn_curve", [DNV_B1C, DNV_C_C, DNV_E_C, DNV_B1A])
@pytest.mark.parametrize("rule", ["pavlou", "manson", "si jian"])
@pytest.mark.filterwarnings("ignore::UserWarning")
@@ -650,10 +672,12 @@ class TestDamageExponents:
)
def test_calc_damage_exponents_no_kwargs(self, damage_rule):
"""Test the _calc_damage_exponents function with no kwargs."""
- stress_range = np.array([100., 200., 300.])
+ stress_range = np.array([100.0, 200.0, 300.0])
if "manson" in damage_rule:
with pytest.raises(ValueError, match="sn_curve must be provided"):
- damage.stress_life._calc_damage_exponents(damage_rule, stress_range)
+ damage.stress_life._calc_damage_exponents(
+ damage_rule, stress_range
+ )
else:
if damage_rule == "pavlou":
with pytest.warns(UserWarning, match="base_exponent"):
@@ -699,7 +723,10 @@ def test_calc_damage_exponents_manson(self):
base_exponent = 0.5
sn_curve = SNCurve([3, 5, 7], [10.970, 13.617, 16])
exponents = damage.stress_life._calc_damage_exponents(
- "manson", stress_range, sn_curve=sn_curve, base_exponent=base_exponent
+ "manson",
+ stress_range,
+ sn_curve=sn_curve,
+ base_exponent=base_exponent,
)
assert isinstance(exponents, np.ndarray)
assert np.allclose(
@@ -717,14 +744,20 @@ def test_calc_damage_exponents_leve(self):
"leve", stress_range, base_exponent=base_exponent
)
assert isinstance(exponents, np.ndarray)
- assert np.allclose(exponents, base_exponent * np.ones(len(stress_range)))
+ assert np.allclose(
+ exponents, base_exponent * np.ones(len(stress_range))
+ )
def test_calc_damage_exponents_si_jian(self):
"""Test the _calc_damage_exponents function with si jian rule."""
stress_range = np.array([100, 200, 300])
- exponents = damage.stress_life._calc_damage_exponents("si jian", stress_range)
+ exponents = damage.stress_life._calc_damage_exponents(
+ "si jian", stress_range
+ )
assert isinstance(exponents, np.ndarray)
- assert np.allclose(exponents, damage.stress_life.calc_si_jian_exponents(stress_range))
+ assert np.allclose(
+ exponents, damage.stress_life.calc_si_jian_exponents(stress_range)
+ )
def test_calc_damage_exponents_unknown_rule(self):
"""Test the _calc_damage_exponents function with unknown rule."""
@@ -732,15 +765,19 @@ def test_calc_damage_exponents_unknown_rule(self):
with pytest.raises(ValueError, match="Unknown damage rule: unknown"):
damage.stress_life._calc_damage_exponents("unknown", stress_range)
+
@pytest.mark.parametrize("sn_curve", [DNV_B1C, DNV_C_C, DNV_E_C, DNV_B1A])
-@pytest.mark.parametrize("load", [
- [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 3, -3, 3, -3, 3, -3, 3, -3],
- [3, -3, 3, -3, 3, -3, 3, -3, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1]
-])
-@pytest.mark.parametrize("damage_bands", [
- [0, 0.2, 0.4, 0.6, 0.8, 1],
- [0, 0.1, 0.3, 0.5, 0.7, 0.9, 1]
-])
+@pytest.mark.parametrize(
+ "load",
+ [
+ [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 3, -3, 3, -3, 3, -3, 3, -3],
+ [3, -3, 3, -3, 3, -3, 3, -3, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1],
+ ],
+)
+@pytest.mark.parametrize(
+ "damage_bands",
+ [[0, 0.2, 0.4, 0.6, 0.8, 1], [0, 0.1, 0.3, 0.5, 0.7, 0.9, 1]],
+)
def test_nonlinear_damage_dca(
load: list, sn_curve: SNCurve, damage_bands: list
):
@@ -765,15 +802,19 @@ def test_nonlinear_damage_dca(
name="Test_CC",
)
d_nl, _, _, _ = damage.stress_life.get_nonlinear_damage_with_dca(
- 'pavlou', cc, sn_curve, np.asarray(damage_bands)
+ "pavlou", cc, sn_curve, np.asarray(damage_bands)
)
assert isinstance(d_nl, np.ndarray)
+
@pytest.mark.parametrize("sn_curve", [DNV_B1C, DNV_C_C, DNV_E_C, DNV_B1A])
-@pytest.mark.parametrize("load", [
- [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 3, -3, 3, -3, 3, -3, 3, -3],
- [3, -3, 3, -3, 3, -3, 3, -3, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1]
-])
+@pytest.mark.parametrize(
+ "load",
+ [
+ [1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 3, -3, 3, -3, 3, -3, 3, -3],
+ [3, -3, 3, -3, 3, -3, 3, -3, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1],
+ ],
+)
def test_theil_damage_rule(load: list, sn_curve: SNCurve):
"""Test the nonlinear damage calculation with Theil's method.
@@ -794,7 +835,6 @@ def test_theil_damage_rule(load: list, sn_curve: SNCurve):
name="Test_CC",
)
d_nl, _, _, _ = damage.stress_life.get_nonlinear_damage_with_dca(
- 'theil', cc, sn_curve, damage_bands=np.array([0, 0.2, 0.4, 0.6, 0.8, 1])
+ "theil", cc, sn_curve, damage_bands=np.array([0, 0.2, 0.4, 0.6, 0.8, 1])
)
assert isinstance(d_nl, np.ndarray)
-
diff --git a/tests/material/test_crack_growth_curve.py b/tests/material/test_crack_growth_curve.py
index 20c14c0..9775084 100644
--- a/tests/material/test_crack_growth_curve.py
+++ b/tests/material/test_crack_growth_curve.py
@@ -10,6 +10,7 @@
nb.config.DISABLE_JIT = True
+
def test_pure_paris_curve():
paris_pure = ParisCurve(slope=3.1, intercept=1.2e-14)
assert paris_pure.linear
@@ -17,7 +18,9 @@ def test_pure_paris_curve():
assert paris_pure.intercept == 1.2e-14
assert paris_pure.threshold == 0
assert paris_pure.critical == np.inf
- assert paris_pure.threshold_growth_rate == paris_pure.get_growth_rate(paris_pure.threshold)
+ assert paris_pure.threshold_growth_rate == paris_pure.get_growth_rate(
+ paris_pure.threshold
+ )
assert np.isinf(paris_pure.critical_growth_rate)
assert np.isinf(paris_pure.get_growth_rate(paris_pure.critical))
assert len(paris_pure.get_knee_growth_rate()) == 0
@@ -38,10 +41,12 @@ def test_pure_paris_curve():
def test_paris_curve_initialization():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- pc = ParisCurve(slope=slope, intercept=intercept, threshold=threshold, critical=critical)
+ pc = ParisCurve(
+ slope=slope, intercept=intercept, threshold=threshold, critical=critical
+ )
assert not pc.linear
assert pc.__str__() == pc.name
@@ -80,54 +85,80 @@ def test_paris_curve_initialization():
_calc_growth_rate(
knee_sif, pc.slope, pc.intercept, pc.threshold, pc.critical
),
- rtol=1e-2
+ rtol=1e-2,
)
assert np.allclose(
pc.get_sif(knee_gr),
- _calc_sif(
- knee_gr, pc.slope, pc.intercept, pc.threshold, pc.critical
- ),
- rtol=1e-2
+ _calc_sif(knee_gr, pc.slope, pc.intercept, pc.threshold, pc.critical),
+ rtol=1e-2,
)
with pytest.raises(TypingError):
- assert _calc_growth_rate(
- threshold - 1e-2 , pc.slope, pc.intercept, pc.threshold, pc.critical
- ) == 0
- assert _calc_growth_rate(
- np.array([threshold - 1e-2]),
- pc.slope,
- pc.intercept,
- pc.threshold,
- pc.critical
- ) == 0
+ assert (
+ _calc_growth_rate(
+ threshold - 1e-2,
+ pc.slope,
+ pc.intercept,
+ pc.threshold,
+ pc.critical,
+ )
+ == 0
+ )
+ assert (
+ _calc_growth_rate(
+ np.array([threshold - 1e-2]),
+ pc.slope,
+ pc.intercept,
+ pc.threshold,
+ pc.critical,
+ )
+ == 0
+ )
with pytest.raises(TypingError):
- assert _calc_growth_rate(
- critical + 1e-2, pc.slope, pc.intercept, pc.threshold, pc.critical
- ) == np.inf
- assert _calc_growth_rate(
- np.array([critical + 1e-2]),
- pc.slope,
- pc.intercept,
- pc.threshold,
- pc.critical
- ) == np.inf
+ assert (
+ _calc_growth_rate(
+ critical + 1e-2,
+ pc.slope,
+ pc.intercept,
+ pc.threshold,
+ pc.critical,
+ )
+ == np.inf
+ )
+ assert (
+ _calc_growth_rate(
+ np.array([critical + 1e-2]),
+ pc.slope,
+ pc.intercept,
+ pc.threshold,
+ pc.critical,
+ )
+ == np.inf
+ )
+
def test_paris_curve_from_knee_points():
knee_sif = [20, 100, 500, 1000]
- knee_growth_rate = [1E-10, 1E-7, 1E-6, 1E-4]
- pc = ParisCurve.from_knee_points(knee_sif=knee_sif, knee_growth_rate=knee_growth_rate)
+ knee_growth_rate = [1e-10, 1e-7, 1e-6, 1e-4]
+ pc = ParisCurve.from_knee_points(
+ knee_sif=knee_sif, knee_growth_rate=knee_growth_rate
+ )
assert np.allclose(pc.get_knee_sif(), knee_sif[1:-1], rtol=1e-2)
- assert np.allclose(pc.get_knee_growth_rate(), knee_growth_rate[1:-1], rtol=1e-2)
+ assert np.allclose(
+ pc.get_knee_growth_rate(), knee_growth_rate[1:-1], rtol=1e-2
+ )
+
def test_paris_curve_get_growth_rate():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- pc = ParisCurve(slope=slope, intercept=intercept, threshold=threshold, critical=critical)
+ pc = ParisCurve(
+ slope=slope, intercept=intercept, threshold=threshold, critical=critical
+ )
sif_range = np.linspace(20, 2000, 10)
growth_rates = pc.get_growth_rate(sif_range)
@@ -135,12 +166,15 @@ def test_paris_curve_get_growth_rate():
assert len(growth_rates) == len(sif_range)
assert np.all(growth_rates >= 0)
+
def test_paris_curve_get_sif():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- pc = ParisCurve(slope=slope, intercept=intercept, threshold=threshold, critical=critical)
+ pc = ParisCurve(
+ slope=slope, intercept=intercept, threshold=threshold, critical=critical
+ )
growth_rate = np.logspace(-10, -4, 10)
sif = pc.get_sif(growth_rate)
@@ -149,20 +183,25 @@ def test_paris_curve_get_sif():
assert np.all(sif >= threshold)
assert np.all(sif <= critical)
+
def test_paris_curve_plot():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- pc = ParisCurve(slope=slope, intercept=intercept, threshold=threshold, critical=critical)
+ pc = ParisCurve(
+ slope=slope, intercept=intercept, threshold=threshold, critical=critical
+ )
fig, ax = pc.plot()
assert fig is not None
assert ax is not None
+
def test_pure_walker_curve():
- walker_pure = WalkerCurve(slope=3.1, intercept=1.2e-14, load_ratio=0.1,
- walker_exponent=0.5)
+ walker_pure = WalkerCurve(
+ slope=3.1, intercept=1.2e-14, load_ratio=0.1, walker_exponent=0.5
+ )
assert walker_pure.linear
assert walker_pure.slope == 3.1
assert walker_pure.intercept == 1.2e-14
@@ -170,7 +209,9 @@ def test_pure_walker_curve():
assert walker_pure.critical == np.inf
assert walker_pure.walker_exponent == 0.5
assert walker_pure.load_ratio == 0.1
- assert walker_pure.threshold_growth_rate == walker_pure.get_growth_rate(walker_pure.threshold)
+ assert walker_pure.threshold_growth_rate == walker_pure.get_growth_rate(
+ walker_pure.threshold
+ )
assert np.isinf(walker_pure.critical_growth_rate)
assert np.isinf(walker_pure.get_growth_rate(walker_pure.critical))
assert len(walker_pure.get_knee_growth_rate()) == 0
@@ -189,13 +230,49 @@ def test_pure_walker_curve():
walker_pure.get_knee_sif(check_knee=np.array([0]))
+def test_paris_curve_accepts_walker_parameters():
+ slope = np.array([3.1])
+ intercept = np.array([1.2e-14])
+
+ curve = ParisCurve(
+ slope=slope,
+ intercept=intercept,
+ load_ratio=0.1,
+ walker_exponent=0.5,
+ )
+
+ assert curve.load_ratio == 0
+ assert curve.walker_exponent == 0
+ np.testing.assert_allclose(curve.walker_intercept, intercept)
+
+
+def test_paris_curve_plot_filters_knees_above_critical():
+ curve = ParisCurve(
+ slope=[3.0, 4.0, 5.0],
+ intercept=[1e-12, 1e-15, 1e-20],
+ threshold=10.0,
+ critical=50.0,
+ )
+
+ fig, ax = curve.plot()
+
+ assert fig is not None
+ assert ax is not None
+
+
def test_walker_curve_initialization():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=0.5, load_ratio=0.1)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=0.5,
+ load_ratio=0.1,
+ )
assert not wc.linear
assert wc.__str__() == wc.name
@@ -234,55 +311,87 @@ def test_walker_curve_initialization():
_calc_growth_rate(
knee_sif, wc.slope, wc.walker_intercept, wc.threshold, wc.critical
),
- rtol=1e-2
+ rtol=1e-2,
)
assert np.allclose(
wc.get_sif(knee_gr),
_calc_sif(
knee_gr, wc.slope, wc.walker_intercept, wc.threshold, wc.critical
),
- rtol=1e-2
+ rtol=1e-2,
)
with pytest.raises(TypingError):
- assert _calc_growth_rate(
- wc.threshold - 1e-2 , wc.slope, wc.intercept, wc.threshold, wc.critical
- ) == 0
- assert _calc_growth_rate(
- np.array([wc.threshold - 1e-2]),
- wc.slope,
- wc.intercept,
- wc.threshold,
- wc.critical
- ) == 0
+ assert (
+ _calc_growth_rate(
+ wc.threshold - 1e-2,
+ wc.slope,
+ wc.intercept,
+ wc.threshold,
+ wc.critical,
+ )
+ == 0
+ )
+ assert (
+ _calc_growth_rate(
+ np.array([wc.threshold - 1e-2]),
+ wc.slope,
+ wc.intercept,
+ wc.threshold,
+ wc.critical,
+ )
+ == 0
+ )
with pytest.raises(TypingError):
- assert _calc_growth_rate(
- critical + 1e-2, wc.slope, wc.intercept, wc.threshold, wc.critical
- ) == np.inf
- assert _calc_growth_rate(
- np.array([critical + 1e-2]),
- wc.slope,
- wc.intercept,
- wc.threshold,
- wc.critical
- ) == np.inf
+ assert (
+ _calc_growth_rate(
+ critical + 1e-2,
+ wc.slope,
+ wc.intercept,
+ wc.threshold,
+ wc.critical,
+ )
+ == np.inf
+ )
+ assert (
+ _calc_growth_rate(
+ np.array([critical + 1e-2]),
+ wc.slope,
+ wc.intercept,
+ wc.threshold,
+ wc.critical,
+ )
+ == np.inf
+ )
+
def test_walker_curve_from_knee_points():
knee_sif = [20, 100, 500, 1000]
- knee_growth_rate = [1E-10, 1E-7, 1E-6, 1E-4]
- wc = WalkerCurve.from_knee_points(knee_sif=knee_sif, knee_growth_rate=knee_growth_rate)
+ knee_growth_rate = [1e-10, 1e-7, 1e-6, 1e-4]
+ wc = WalkerCurve.from_knee_points(
+ knee_sif=knee_sif, knee_growth_rate=knee_growth_rate
+ )
assert np.allclose(wc.get_knee_sif(), knee_sif[1:-1], rtol=1e-2)
- assert np.allclose(wc.get_knee_growth_rate(), knee_growth_rate[1:-1], rtol=1e-2)
+ assert np.allclose(
+ wc.get_knee_growth_rate(), knee_growth_rate[1:-1], rtol=1e-2
+ )
+
def test_walker_curve_get_growth_rate():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=0.5, load_ratio=0.2)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=0.5,
+ load_ratio=0.2,
+ )
sif_range = np.linspace(20, 2000, 10)
growth_rates = wc.get_growth_rate(sif_range)
@@ -290,13 +399,20 @@ def test_walker_curve_get_growth_rate():
assert len(growth_rates) == len(sif_range)
assert np.all(growth_rates >= 0)
+
def test_walker_curve_get_sif():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=0.4, load_ratio=0.3)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=0.4,
+ load_ratio=0.3,
+ )
growth_rate = np.logspace(-10, -4, 10)
sif = wc.get_sif(growth_rate)
@@ -305,66 +421,96 @@ def test_walker_curve_get_sif():
assert np.all(sif >= threshold)
assert np.all(sif <= critical)
+
def test_walker_curve_plot():
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=0.5, load_ratio=0.1)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=0.5,
+ load_ratio=0.1,
+ )
fig, ax = wc.plot()
assert fig is not None
assert ax is not None
+
def test_walker_curve_multiple_load_ratios():
"""Assert that increasing the load ratios reduces the threshold and critical
SIF"""
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
wc_old = None
w_exp = 0.5
for r in [0, 0.1, 0.3, 0.5, 0.7, 0.9]:
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=w_exp, load_ratio=r)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=w_exp,
+ load_ratio=r,
+ )
if wc_old is not None:
assert np.all(wc.walker_intercept >= wc_old.walker_intercept)
assert wc.critical <= wc_old.critical
wc_old = wc
+
def test_walker_curve_change_load_ratios():
"""Assert that increasing the load ratios reduces the threshold and critical
SIF"""
import copy
+
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
wc_old = None
w_exp = 0.5
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=w_exp, load_ratio=0)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=w_exp,
+ load_ratio=0,
+ )
for r in [0.1, 0.3, 0.5, 0.7, 0.9]:
wc_old = copy.copy(wc)
wc.load_ratio = r
assert np.all(wc.walker_intercept >= wc_old.walker_intercept)
assert wc.critical <= wc_old.critical
+
def test_walker_curve_change_exponent():
"""Assert that increasing the Walker exponent towards 1 moves the curve
closer to the Paris curve"""
import copy
+
slope = [2.88, 5.1, 8.16, 5.1, 2.88]
- intercept = [1E-16, 1E-20, 1E-27, 1E-19, 1E-13]
+ intercept = [1e-16, 1e-20, 1e-27, 1e-19, 1e-13]
threshold = 20
critical = 2000
wc_old = None
w_exp = 0.5
- wc = WalkerCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical, walker_exponent=w_exp, load_ratio=0)
+ wc = WalkerCurve(
+ slope=slope,
+ intercept=intercept,
+ threshold=threshold,
+ critical=critical,
+ walker_exponent=w_exp,
+ load_ratio=0,
+ )
for w_e in [0, 0.25, 0.5, 0.75, 1]:
wc_old = copy.copy(wc)
wc.walker_exponent = w_e
@@ -372,15 +518,17 @@ def test_walker_curve_change_exponent():
assert wc.critical >= wc_old.critical
assert wc.threshold >= wc_old.threshold
- pc = ParisCurve(slope=slope, intercept=intercept, threshold=threshold,
- critical=critical)
+ pc = ParisCurve(
+ slope=slope, intercept=intercept, threshold=threshold, critical=critical
+ )
# Assert that the two curves are the same when the exponent is 1
assert wc.threshold == pc.threshold
assert wc.critical == pc.critical
assert np.allclose(wc.get_knee_sif(), pc.get_knee_sif(), rtol=1e-2)
- assert np.allclose(wc.get_knee_growth_rate(), pc.get_knee_growth_rate(),
- rtol=1e-2)
+ assert np.allclose(
+ wc.get_knee_growth_rate(), pc.get_knee_growth_rate(), rtol=1e-2
+ )
nb.config.DISABLE_JIT = False
diff --git a/tests/material/test_sn_curve.py b/tests/material/test_sn_curve.py
index 57f2390..2cf3af4 100644
--- a/tests/material/test_sn_curve.py
+++ b/tests/material/test_sn_curve.py
@@ -5,6 +5,7 @@
import numba as nb
import numpy as np
import pytest
+from py_fatigue.material.sn_curve import _calc_cycles, _calc_stress
# os.environ["NUMBA_DISABLE_JIT"] = "1"
@@ -462,4 +463,76 @@ def test_from_knee_points():
assert sn_knee.name == sn.name
+def _kernel(fn):
+ """Return python-callable implementation for numba kernels."""
+ return getattr(fn, "py_func", fn)
+
+
+@pytest.mark.parametrize("sn", [DNV_B1C, DNV_B1A, EXOTIC])
+def test_cycles_kernel_boundary_regression(sn):
+ kernel = _kernel(_calc_cycles)
+
+ stress_values = [1e-12, sn.get_stress(1e4)[0], sn.get_stress(1e7)[0]]
+ stress_values.extend(list(sn.get_knee_stress()))
+ stress_values.extend([s * (1 - 1e-12) for s in sn.get_knee_stress()])
+ stress_values.extend([s * (1 + 1e-12) for s in sn.get_knee_stress()])
+ if sn.endurance < np.inf:
+ end_stress = sn.get_stress(sn.endurance)[0]
+ stress_values.extend(
+ [
+ end_stress,
+ end_stress * (1 - 1e-12),
+ end_stress * (1 + 1e-12),
+ ]
+ )
+ stress = np.asarray(stress_values, dtype=np.float64)
+
+ kernel_out = kernel(stress, sn.slope, sn.intercept, sn.endurance)
+ api_out = sn.get_cycles(stress)
+ np.testing.assert_allclose(kernel_out, api_out, rtol=1e-11, atol=0.0)
+
+ zero_stress = np.array([0.0], dtype=np.float64)
+ kernel_zero = kernel(zero_stress, sn.slope, sn.intercept, sn.endurance)
+ api_zero = sn.get_cycles(zero_stress)
+ np.testing.assert_allclose(kernel_zero, api_zero, rtol=0.0, atol=0.0)
+
+ with pytest.raises(AssertionError):
+ kernel(np.array([-1.0]), sn.slope, sn.intercept, sn.endurance)
+ with pytest.raises(AssertionError):
+ _ = sn.get_cycles(-1)
+
+
+@pytest.mark.parametrize("sn", [DNV_B1C, DNV_B1A, EXOTIC])
+def test_stress_kernel_boundary_regression(sn):
+ kernel = _kernel(_calc_stress)
+
+ cycle_values = [1e-12, 1e4, 1e7]
+ cycle_values.extend(list(sn.get_knee_cycles()))
+ cycle_values.extend([c * (1 - 1e-12) for c in sn.get_knee_cycles()])
+ cycle_values.extend([c * (1 + 1e-12) for c in sn.get_knee_cycles()])
+ if sn.endurance < np.inf:
+ cycle_values.extend(
+ [
+ sn.endurance,
+ sn.endurance * (1 - 1e-12),
+ sn.endurance * (1 + 1e-12),
+ ]
+ )
+
+ cycles = np.asarray(cycle_values, dtype=np.float64)
+ kernel_out = kernel(cycles, sn.slope, sn.intercept, sn.endurance)
+ api_out = sn.get_stress(cycles)
+ np.testing.assert_allclose(kernel_out, api_out, rtol=1e-11, atol=0.0)
+
+ with pytest.raises(AssertionError):
+ kernel(np.array([0.0]), sn.slope, sn.intercept, sn.endurance)
+ with pytest.raises(AssertionError):
+ _ = sn.get_stress(0)
+
+ with pytest.raises(AssertionError):
+ kernel(np.array([-1.0]), sn.slope, sn.intercept, sn.endurance)
+ with pytest.raises(AssertionError):
+ _ = sn.get_stress(-1)
+
+
nb.config.DISABLE_JIT = False
diff --git a/tests/mean_stress/test_corrections.py b/tests/mean_stress/test_corrections.py
index 3629fbd..9913e2e 100644
--- a/tests/mean_stress/test_corrections.py
+++ b/tests/mean_stress/test_corrections.py
@@ -86,19 +86,31 @@ def test_invalid_detail_factor(self):
@given(
mean_stress=hy.lists(
- hy.floats(min_value=-1000, max_value=1000, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=20
+ hy.floats(
+ min_value=-1000,
+ max_value=1000,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=20,
),
stress_amp=hy.lists(
- hy.floats(min_value=1, max_value=1000, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=20
- )
+ hy.floats(
+ min_value=1,
+ max_value=1000,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=20,
+ ),
)
def test_property_based_dnvgl_08(self, mean_stress, stress_amp):
"""Property-based test for DNVGL with detail_factor=0.8."""
if len(mean_stress) != len(stress_amp):
- mean_stress = mean_stress[:min(len(mean_stress), len(stress_amp))]
- stress_amp = stress_amp[:min(len(mean_stress), len(stress_amp))]
+ mean_stress = mean_stress[: min(len(mean_stress), len(stress_amp))]
+ stress_amp = stress_amp[: min(len(mean_stress), len(stress_amp))]
mean_stress = np.array(mean_stress)
stress_amp = np.array(stress_amp)
@@ -120,19 +132,31 @@ def test_property_based_dnvgl_08(self, mean_stress, stress_amp):
@given(
mean_stress=hy.lists(
- hy.floats(min_value=-1000, max_value=1000, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=20
+ hy.floats(
+ min_value=-1000,
+ max_value=1000,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=20,
),
stress_amp=hy.lists(
- hy.floats(min_value=1, max_value=1000, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=20
- )
+ hy.floats(
+ min_value=1,
+ max_value=1000,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=20,
+ ),
)
def test_property_based_dnvgl_06(self, mean_stress, stress_amp):
"""Property-based test for DNVGL with detail_factor=0.6."""
if len(mean_stress) != len(stress_amp):
- mean_stress = mean_stress[:min(len(mean_stress), len(stress_amp))]
- stress_amp = stress_amp[:min(len(mean_stress), len(stress_amp))]
+ mean_stress = mean_stress[: min(len(mean_stress), len(stress_amp))]
+ stress_amp = stress_amp[: min(len(mean_stress), len(stress_amp))]
mean_stress = np.array(mean_stress)
stress_amp = np.array(stress_amp)
@@ -166,7 +190,9 @@ def test_edge_cases(self):
result = dnvgl_mean_stress_correction(
np.array([0]), np.array([50]), detail_factor=0.8
)
- assert result[0] < 100 # Should be less than 2*amplitude for mixed loading
+ assert (
+ result[0] < 100
+ ) # Should be less than 2*amplitude for mixed loading
assert result[0] >= 80 # But at least detail_factor * 2 * amplitude
# Very small amplitude with high tensile mean
@@ -218,7 +244,9 @@ def test_gamma_variations(self):
np.testing.assert_allclose(result_1, stress_amplitude, rtol=1e-10)
@given(
- gamma=hy.floats(min_value=0.0, max_value=1.0, allow_nan=False, allow_infinity=False)
+ gamma=hy.floats(
+ min_value=0.0, max_value=1.0, allow_nan=False, allow_infinity=False
+ )
)
def test_gamma_range(self, gamma):
"""Test gamma parameter in valid range."""
@@ -258,19 +286,31 @@ def test_negative_max_stress_handling(self):
@given(
mean_stress=hy.lists(
- hy.floats(min_value=-500, max_value=500, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=10
+ hy.floats(
+ min_value=-500,
+ max_value=500,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=10,
),
stress_amp=hy.lists(
- hy.floats(min_value=1, max_value=500, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=10
- )
+ hy.floats(
+ min_value=1,
+ max_value=500,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=10,
+ ),
)
def test_property_based_walker(self, mean_stress, stress_amp):
"""Property-based test for Walker correction."""
if len(mean_stress) != len(stress_amp):
- mean_stress = mean_stress[:min(len(mean_stress), len(stress_amp))]
- stress_amp = stress_amp[:min(len(mean_stress), len(stress_amp))]
+ mean_stress = mean_stress[: min(len(mean_stress), len(stress_amp))]
+ stress_amp = stress_amp[: min(len(mean_stress), len(stress_amp))]
mean_stress = np.array(mean_stress)
stress_amp = np.array(stress_amp)
@@ -322,19 +362,31 @@ def test_plotting_functionality(self):
@given(
mean_stress=hy.lists(
- hy.floats(min_value=-500, max_value=500, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=10
+ hy.floats(
+ min_value=-500,
+ max_value=500,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=10,
),
stress_amp=hy.lists(
- hy.floats(min_value=1, max_value=500, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=10
- )
+ hy.floats(
+ min_value=1,
+ max_value=500,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=10,
+ ),
)
def test_property_based_swt(self, mean_stress, stress_amp):
"""Property-based test for SWT correction."""
if len(mean_stress) != len(stress_amp):
- mean_stress = mean_stress[:min(len(mean_stress), len(stress_amp))]
- stress_amp = stress_amp[:min(len(mean_stress), len(stress_amp))]
+ mean_stress = mean_stress[: min(len(mean_stress), len(stress_amp))]
+ stress_amp = stress_amp[: min(len(mean_stress), len(stress_amp))]
mean_stress = np.array(mean_stress)
stress_amp = np.array(stress_amp)
@@ -362,7 +414,10 @@ def test_basic_functionality(self):
amp_in, mean_in, r_out, ult_s, correction_exponent
)
- assert amp_out.shape == (1, 3) # r_out as scalar creates (1, len(amp_in))
+ assert amp_out.shape == (
+ 1,
+ 3,
+ ) # r_out as scalar creates (1, len(amp_in))
assert mean_out.shape == (1, 3)
assert np.all(amp_out >= 0)
assert np.all(np.isfinite(amp_out))
@@ -394,12 +449,38 @@ def test_r_out_minus_one_special_case(self):
)
# For r_out = -1, mean_out should be zeros
- np.testing.assert_allclose(mean_out[0], np.zeros_like(amp_in), rtol=1e-10)
+ np.testing.assert_allclose(
+ mean_out[0], np.zeros_like(amp_in), rtol=1e-10
+ )
# Analytical solution: amp_out = amp_in / (1 - (mean_in / ult_s)^n)
- expected_amp_out = amp_in / (1 - (mean_in / ult_s) ** correction_exponent)
+ expected_amp_out = amp_in / (
+ 1 - (mean_in / ult_s) ** correction_exponent
+ )
np.testing.assert_allclose(amp_out[0], expected_amp_out, rtol=1e-6)
+ def test_goodman_non_reversed_load_ratio_analytical_solution(self):
+ """Test Goodman correction with non-fully-reversed output ratios."""
+ amp_in = np.array([80.7335222, 189.5302955, 258.9631694])
+ mean_in = np.array([13.33359442, 122.13036772, 52.69749382])
+ r_out = np.array([-3.0, 0.0])
+ ult_s = 900.0
+
+ amp_out, mean_out = goodman_haigh_mean_stress_correction(
+ amp_in, mean_in, r_out, ult_s, correction_exponent=1.0
+ )
+
+ r_in = (mean_in - amp_in) / (mean_in + amp_in)
+ for idx, r_out_val in enumerate(r_out):
+ expected_amp = (
+ (1 - ((1 + r_in) / (1 - r_in) * amp_in / ult_s)) / amp_in
+ + ((1 + r_out_val) / (1 - r_out_val) / ult_s)
+ ) ** -1
+ expected_mean = expected_amp * (1 + r_out_val) / (1 - r_out_val)
+ np.testing.assert_allclose(amp_out[idx], expected_amp, rtol=1e-12)
+ np.testing.assert_allclose(mean_out[idx], expected_mean, rtol=1e-12)
+ assert np.all(amp_out[idx] >= 0)
+
def test_different_correction_exponents(self):
"""Test different correction exponents (Goodman=1, Gerber=2)."""
amp_in = np.array([100])
@@ -514,25 +595,47 @@ def test_initial_guess_parameter(self):
initial_guess = np.array([120, 180])
amp_out, mean_out = goodman_haigh_mean_stress_correction(
- amp_in, mean_in, r_out, ult_s, correction_exponent,
- initial_guess=initial_guess
+ amp_in,
+ mean_in,
+ r_out,
+ ult_s,
+ correction_exponent,
+ initial_guess=initial_guess,
)
assert amp_out.shape == (1, 2)
@given(
amp_in=hy.lists(
- hy.floats(min_value=10, max_value=200, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=3
+ hy.floats(
+ min_value=10,
+ max_value=200,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=3,
),
mean_in=hy.lists(
- hy.floats(min_value=-50, max_value=50, allow_nan=False, allow_infinity=False),
- min_size=1, max_size=3
+ hy.floats(
+ min_value=-50,
+ max_value=50,
+ allow_nan=False,
+ allow_infinity=False,
+ ),
+ min_size=1,
+ max_size=3,
+ ),
+ ult_s=hy.floats(
+ min_value=800, max_value=2000, allow_nan=False, allow_infinity=False
+ ),
+ correction_exponent=hy.floats(
+ min_value=1.0, max_value=2.0, allow_nan=False, allow_infinity=False
),
- ult_s=hy.floats(min_value=800, max_value=2000, allow_nan=False, allow_infinity=False),
- correction_exponent=hy.floats(min_value=1.0, max_value=2.0, allow_nan=False, allow_infinity=False)
)
- def test_property_based_goodman_haigh(self, amp_in, mean_in, ult_s, correction_exponent):
+ def test_property_based_goodman_haigh(
+ self, amp_in, mean_in, ult_s, correction_exponent
+ ):
"""Property-based test for Goodman-Haigh correction."""
if len(amp_in) != len(mean_in):
min_len = min(len(amp_in), len(mean_in))
@@ -605,7 +708,11 @@ def test_edge_case_stress_combinations(self):
zero_mean = np.array([0])
small_amp = np.array([1e-6])
- for correction in [dnvgl_mean_stress_correction, walker_mean_stress_correction, swt_mean_stress_correction]:
+ for correction in [
+ dnvgl_mean_stress_correction,
+ walker_mean_stress_correction,
+ swt_mean_stress_correction,
+ ]:
if correction == dnvgl_mean_stress_correction:
result = correction(zero_mean, small_amp, detail_factor=0.8)
else:
@@ -625,4 +732,4 @@ def test_large_stress_values(self):
for result in [result_dnvgl, result_walker, result_swt]:
assert np.all(np.isfinite(result))
- assert np.all(result >= 0)
\ No newline at end of file
+ assert np.all(result >= 0)
diff --git a/tests/test_utils.py b/tests/test_utils.py
index 35990d8..4ec4259 100644
--- a/tests/test_utils.py
+++ b/tests/test_utils.py
@@ -10,6 +10,7 @@
from unittest.mock import Mock
# Non-standard imports
+import numba as nb
import numpy as np
import matplotlib.pyplot as plt
import pytest
@@ -94,6 +95,8 @@ def test_split_full_cycles_and_residuals(
(0, 1, 99, *expected(99.0)),
(0, 1, 98.9, *expected(99.0)),
]
+
+
@pytest.mark.parametrize(bub_names, bub_data)
def test_bin_upper_bound(
bin_lower_bound: float,
@@ -172,6 +175,8 @@ def test_calc_half_cycles(
cfc_names = "count_cycle, expected, error"
+
+
# fmt: off
cfc_data = [
(COUNTS_1, *expected(np.array(
@@ -231,27 +236,19 @@ def test_faulty_fatigue_stress(self, lists: Tuple[List[int], List[int]]):
with self.assertRaises(ValueError) as ve:
pu.FatigueStress(_counts=cts, _values=vals[1:], bin_width=0.5)
self.assertEqual(
- "counts and values must have the same length",
- str(ve.exception)
+ "counts and values must have the same length", str(ve.exception)
)
@given(same_len_lists())
def test_empty_fatigue_stress(self, lists: Tuple[List[int], List[int]]):
- """Assert that an error is raised when the counts or values is empty
- """
+ """Assert that an error is raised when the counts or values is empty"""
vals, cts = lists
with self.assertRaises(ValueError) as ve:
pu.FatigueStress(_counts=cts, _values=[], bin_width=0.5)
- self.assertEqual(
- "No data provided",
- str(ve.exception)
- )
+ self.assertEqual("No data provided", str(ve.exception))
with self.assertRaises(ValueError) as ve:
pu.FatigueStress(_counts=[], _values=vals, bin_width=0.5)
- self.assertEqual(
- "No data provided",
- str(ve.exception)
- )
+ self.assertEqual("No data provided", str(ve.exception))
@pytest.mark.parametrize(
@@ -377,18 +374,17 @@ def test_make_axes():
plt.close(fig)
# Test with invalid figure type
- with pytest.raises(TypeError, match="fig must be a matplotlib.figure.Figure"):
+ with pytest.raises(
+ TypeError, match="fig must be a matplotlib.figure.Figure"
+ ):
pu.make_axes(fig="not_a_figure")
def test_to_numba_dict():
"""Test to_numba_dict function"""
import numba as nb
- data = {
- "key1": 1.0,
- "key2": 2.0,
- "key3": 3.0
- }
+
+ data = {"key1": 1.0, "key2": 2.0, "key3": 3.0}
result = pu.to_numba_dict(data)
# The next line fails since os.environ["NUMBA_DISABLE_JIT"] = "1"
@@ -559,6 +555,7 @@ def test_fatigue_stress_setters():
def test_py_bisect():
"""Test py_bisect function"""
+
def test_func(x):
return x**2 - 4
@@ -575,6 +572,7 @@ def simple_func(x):
def test_py_newton():
"""Test py_newton function"""
+
def test_func(x):
return x**2 - 4
@@ -584,6 +582,7 @@ def test_func(x):
def test_numba_bisect():
"""Test numba_bisect wrapper"""
+
def test_func(x):
return x**2 - 4
@@ -593,6 +592,7 @@ def test_func(x):
def test_numba_newton():
"""Test numba_newton wrapper"""
+
def test_func(x):
return x**2 - 4
@@ -602,6 +602,7 @@ def test_func(x):
def test_compile_specialized_bisect():
"""Test compile_specialized_bisect function"""
+
def test_func(x):
return (x**2 - 4,)
@@ -610,8 +611,105 @@ def test_func(x):
assert np.isclose(root, 2.0, atol=1e-5)
+def test_compile_specialized_bisect_many_args():
+ """Test specialized bisection with more than three extra arguments."""
+
+ def test_func(x, arg_1, arg_2, arg_3, arg_4, arg_5):
+ return x - (arg_1 + arg_2 + arg_3 + arg_4 + arg_5)
+
+ compiled_bisect = pu.compile_specialized_bisect(test_func)
+ root = compiled_bisect(
+ 0,
+ 20,
+ 1e-6,
+ 100,
+ 1.0,
+ 2.0,
+ 3.0,
+ 4.0,
+ 5.0,
+ )
+ assert np.isclose(root, 15.0, atol=1e-5)
+
+
+@pytest.mark.parametrize(
+ "test_func",
+ [
+ lambda x: x**2 - 4,
+ lambda x: (x**2 - 4,),
+ lambda x: np.array([x**2 - 4]),
+ lambda x: [x**2 - 4],
+ lambda x: np.float64(x**2 - 4),
+ ],
+)
+def test_compile_specialized_bisect_scalar_like_returns(test_func):
+ """Test bisection residuals that return scalar-like values."""
+
+ compiled_bisect = pu.compile_specialized_bisect(test_func)
+ root = compiled_bisect(0, 5, tol=1e-6, mxiter=100)
+ assert np.isclose(root, 2.0, atol=1e-5)
+
+
+@pytest.mark.parametrize(
+ "test_func",
+ [
+ lambda x: (x**2 - 4, x),
+ lambda x: np.array([x**2 - 4, x]),
+ lambda x: [x**2 - 4, x],
+ ],
+)
+def test_compile_specialized_bisect_rejects_non_scalar_returns(test_func):
+ """Test bisection rejects residuals that return multiple values."""
+
+ compiled_bisect = pu.compile_specialized_bisect(test_func)
+ with pytest.raises((TypeError, nb.core.errors.TypingError)):
+ compiled_bisect(0, 5, tol=1e-6, mxiter=100)
+
+
+def test_compile_specialized_bisect_uses_numba_dispatcher(monkeypatch):
+ """Test bisection uses a numba dispatcher when JIT is enabled."""
+
+ monkeypatch.delenv("NUMBA_DISABLE_JIT", raising=False)
+ monkeypatch.setattr(pu.nb.config, "DISABLE_JIT", False)
+ pu.compile_specialized_bisect.cache_clear()
+
+ def test_func(x):
+ return (x**2 - 4,)
+
+ compiled_bisect = pu.compile_specialized_bisect(test_func)
+ root = compiled_bisect(0, 5, tol=1e-6, mxiter=100)
+
+ assert compiled_bisect.__class__.__name__ == "CPUDispatcher"
+ assert np.isclose(root, 2.0, atol=1e-5)
+
+ pu.compile_specialized_bisect.cache_clear()
+
+
+def test_numba_bisect_accepts_specialized_function():
+ """Test numba_bisect does not recompile specialized bisection functions."""
+
+ def test_func(x, arg_1, arg_2, arg_3, arg_4, arg_5):
+ return x - (arg_1 + arg_2 + arg_3 + arg_4 + arg_5)
+
+ compiled_bisect = pu.compile_specialized_bisect(test_func)
+ root = pu.numba_bisect(
+ compiled_bisect,
+ 0,
+ 20,
+ 1e-6,
+ 100,
+ 1.0,
+ 2.0,
+ 3.0,
+ 4.0,
+ 5.0,
+ )
+ assert np.isclose(root, 15.0, atol=1e-5)
+
+
def test_compile_specialized_newton():
"""Test compile_specialized_newton function"""
+
def test_func(x):
return x**2 - 4
@@ -620,6 +718,12 @@ def test_func(x):
assert np.isclose(root, 2.0, atol=1e-5)
+def test_warmup_numba():
+ """Test warmup_numba compiles supported numba call paths."""
+
+ pu.warmup_numba()
+
+
def test_custom_formatter():
"""Test CustomFormatter class"""
formatter = pu.CustomFormatter()
@@ -629,8 +733,13 @@ def test_custom_formatter():
# Test DEBUG level
record_debug = logging.LogRecord(
- name="test", level=logging.DEBUG, pathname="test.py",
- lineno=1, msg="Debug message", args=(), exc_info=None
+ name="test",
+ level=logging.DEBUG,
+ pathname="test.py",
+ lineno=1,
+ msg="Debug message",
+ args=(),
+ exc_info=None,
)
formatted_debug = formatter.format(record_debug)
assert "Debug message" in formatted_debug
@@ -638,8 +747,13 @@ def test_custom_formatter():
# Test INFO level
record_info = logging.LogRecord(
- name="test", level=logging.INFO, pathname="test.py",
- lineno=1, msg="Info message", args=(), exc_info=None
+ name="test",
+ level=logging.INFO,
+ pathname="test.py",
+ lineno=1,
+ msg="Info message",
+ args=(),
+ exc_info=None,
)
formatted_info = formatter.format(record_info)
assert "Info message" in formatted_info
@@ -647,8 +761,13 @@ def test_custom_formatter():
# Test WARNING level
record_warning = logging.LogRecord(
- name="test", level=logging.WARNING, pathname="test.py",
- lineno=1, msg="Warning message", args=(), exc_info=None
+ name="test",
+ level=logging.WARNING,
+ pathname="test.py",
+ lineno=1,
+ msg="Warning message",
+ args=(),
+ exc_info=None,
)
formatted_warning = formatter.format(record_warning)
assert "Warning message" in formatted_warning
@@ -656,8 +775,13 @@ def test_custom_formatter():
# Test ERROR level
record_error = logging.LogRecord(
- name="test", level=logging.ERROR, pathname="test.py",
- lineno=1, msg="Error message", args=(), exc_info=None
+ name="test",
+ level=logging.ERROR,
+ pathname="test.py",
+ lineno=1,
+ msg="Error message",
+ args=(),
+ exc_info=None,
)
formatted_error = formatter.format(record_error)
assert "Error message" in formatted_error
@@ -665,8 +789,13 @@ def test_custom_formatter():
# Test CRITICAL level
record_critical = logging.LogRecord(
- name="test", level=logging.CRITICAL, pathname="test.py",
- lineno=1, msg="Critical message", args=(), exc_info=None
+ name="test",
+ level=logging.CRITICAL,
+ pathname="test.py",
+ lineno=1,
+ msg="Critical message",
+ args=(),
+ exc_info=None,
)
formatted_critical = formatter.format(record_critical)
assert "Critical message" in formatted_critical
@@ -789,7 +918,9 @@ def test_fatigue_stress_bin_properties():
"""Test FatigueStress bin-related properties"""
counts = np.array([1.0, 2.0, 3.0])
values = np.array([10.0, 20.0, 30.0])
- fs = pu.FatigueStress(_counts=counts, _values=values, bin_width=5.0, _bin_lb=5.0, _bin_ub=35.0)
+ fs = pu.FatigueStress(
+ _counts=counts, _values=values, bin_width=5.0, _bin_lb=5.0, _bin_ub=35.0
+ )
# Test bin_edges
edges = fs.bin_edges
@@ -820,7 +951,9 @@ def test_plot_damage_accumulation():
cumsum_pm_dmg = np.array([0.05, 0.2, 0.5, 0.8, 1.1])
limit_damage = 1.0
- fig, ax = pu._plot_damage_accumulation(cumsum_nl_dmg, cumsum_pm_dmg, limit_damage)
+ fig, ax = pu._plot_damage_accumulation(
+ cumsum_nl_dmg, cumsum_pm_dmg, limit_damage
+ )
assert isinstance(fig, plt.Figure)
assert isinstance(ax, plt.Axes)
plt.close(fig)
@@ -895,5 +1028,7 @@ def test_json_encoders():
# Test in actual JSON encoding
data = {"array": arr}
- json_str = json.dumps(data, default=lambda x: pu.JSON_ENCODERS.get(type(x), str)(x))
- assert "[1, 2, 3]" in json_str
\ No newline at end of file
+ json_str = json.dumps(
+ data, default=lambda x: pu.JSON_ENCODERS.get(type(x), str)(x)
+ )
+ assert "[1, 2, 3]" in json_str