From 2a97d9fde1bd42178d647acd35284a1b8c9b6431 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Wed, 11 Dec 2024 15:36:29 +0100 Subject: [PATCH 01/25] Introduce draft of NotchApproxBinner Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 71 ++++++ .../test_notch_approximation_law.py | 209 +++++++++++++++++- 2 files changed, 274 insertions(+), 6 deletions(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index 2978e293..8069f70f 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -198,6 +198,19 @@ def load(self, stress, *, rtol=1e-4, tol=1e-4): ) return load + def primary(self, load): + load = np.asarray(load) + stress = self.stress(load) + strain = self.strain(stress, None) + return np.stack([stress, strain], axis=len(load.shape)) + + def secondary(self, delta_load): + delta_load = np.asarray(delta_load) + delta_stress = self.stress_secondary_branch(delta_load) + delta_strain = self.strain_secondary_branch(delta_stress, None) + return np.stack([delta_stress, delta_strain], axis=len(delta_load.shape)) + + def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): """Calculate the stress on secondary branches in the stress-strain diagram at a given elastic-plastic stress (load), from a FE computation. @@ -906,3 +919,61 @@ def _create_bins_multiple_assessment_points(self): self._lut_secondary_branch.delta_strain \ = self._notch_approximation_law.strain_secondary_branch( self._lut_secondary_branch.delta_stress, self._lut_secondary_branch.delta_load) + + + +class NotchApproxBinner: + + def __init__(self, notch_approximation_law, number_of_bins=100): + self._n_bins = number_of_bins + self._notch_approximation_law = notch_approximation_law + self.ramberg_osgood_relation = notch_approximation_law.ramberg_osgood_relation + self._max_load_rep = None + + def initialize(self, max_load): + max_load = np.asarray(max_load) + self._max_load_rep, _ = self._rep_abs_and_sign(max_load) + + load = self._param_for_lut(self._n_bins, max_load) + self._lut_primary = self._notch_approximation_law.primary(load) + + delta_load = self._param_for_lut(2 * self._n_bins, 2.0*max_load) + self._lut_secondary = self._notch_approximation_law.secondary(delta_load) + + return self + + def primary(self, load): + self._raise_if_uninitialized() + load_rep, sign = self._rep_abs_and_sign(load) + + if load_rep > self._max_load_rep: + msg = f"Requested load `{load_rep}`, higher than initialized maximum load `{self._max_load_rep}`" + raise ValueError(msg) + + idx = int(np.ceil(load_rep / self._max_load_rep * self._n_bins)) - 1 + return sign * self._lut_primary[idx, :] + + def secondary(self, delta_load): + self._raise_if_uninitialized() + delta_load_rep, sign = self._rep_abs_and_sign(delta_load) + + if delta_load_rep > 2.0 * self._max_load_rep: + msg = f"Requested load `{delta_load_rep}`, higher than initialized maximum delta load `{2.0*self._max_load_rep}`" + raise ValueError(msg) + + idx = int(np.ceil(delta_load_rep / (2.0*self._max_load_rep) * 2*self._n_bins)) - 1 + return sign * self._lut_secondary[idx, :] + + def _raise_if_uninitialized(self): + if self._max_load_rep is None: + raise RuntimeError("NotchApproxBinner not initialized.") + + def _param_for_lut(self, number_of_bins, max_val): + scale = np.linspace(0.0, 1.0, number_of_bins + 1)[1:] + max_val, scale_m = np.meshgrid(max_val, scale) + return (max_val * scale_m) + + def _rep_abs_and_sign(self, value): + value = np.asarray(value) + value_rep = value if len(value.shape) == 0 else value[0] + return np.abs(value_rep), np.sign(value_rep) diff --git a/tests/materiallaws/test_notch_approximation_law.py b/tests/materiallaws/test_notch_approximation_law.py index 982b0ac2..2b1ede8c 100644 --- a/tests/materiallaws/test_notch_approximation_law.py +++ b/tests/materiallaws/test_notch_approximation_law.py @@ -22,7 +22,7 @@ import pandas as pd import pylife.strength.damage_parameter -from pylife.materiallaws.notch_approximation_law import ExtendedNeuber +from pylife.materiallaws.notch_approximation_law import ExtendedNeuber, NotchApproxBinner from .data import * def test_extended_neuber_example_1(): @@ -32,7 +32,7 @@ def test_extended_neuber_example_1(): K = 1184 # [MPa] n = 0.187 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - + L = pd.Series([100, -200, 100, -250, 200, 0, 200, -200]) c = 1.4 gamma_L = (250+6.6)/250 @@ -171,10 +171,11 @@ def test_derivatives(stress, load): assert np.isclose(numeric_derivative, derivative) + @pytest.mark.parametrize('E, K, n, L', [ - (260e3, 1184, 0.187, pd.Series([100, -200, 100, -250, 200, 100, 200, -200])), - (100e3, 1500, 0.4, pd.Series([-100, 100, -200])), - (200e3, 1000, 0.2, pd.Series([100, 10])), + (260e3, 1184, 0.187, np.array([100, -200, 100, -250, 200, 100, 200, -200])), + (100e3, 1500, 0.4, np.array([-100, 100, -200])), + (200e3, 1000, 0.2, np.array([100, 10])), ]) def test_load(E, K, n, L): c = 1.4 @@ -184,7 +185,7 @@ def test_load(E, K, n, L): # initialize notch approximation law and damage parameter notch_approximation_law = ExtendedNeuber(E, K, n, K_p=3.5) - # The "load" method is the inverse operation of "stress", + # The "load" method is the inverse operation of "stress", # i.e., ``L = load(stress(L))`` and ``S = stress(load(stress))``. stress = notch_approximation_law.stress(L) load = notch_approximation_law.load(stress) @@ -193,3 +194,199 @@ def test_load(E, K, n, L): np.testing.assert_allclose(L, load, rtol=1e-3) np.testing.assert_allclose(stress, stress2, rtol=1e-3) + +@pytest.mark.parametrize("L, expected", [ + (150.0, [148.5, 7.36e-4]), + (175.0, [171.7, 8.67e-4]), + (200.0, [193.7, 1.00e-3]) +]) +def test_load_primary_scalar(L, expected): + notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + result = notch_approximation_law.primary(L) + + assert result.shape == (2, ) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_load_primary_vectorized(): + notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + + result = notch_approximation_law.primary([150.0, 175.0, 200.0]) + expected = [[148.4, 7.36e-4], [171.7, 8.67e-4], [193.7, 1.00e-3]] + + assert result.shape == (3, 2) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (100.0, [100.0, 4.88e-4]), + (400.0, [386.0, 1.99e-3]), + (600.0, [533.0, 3.28e-3]) +]) +def test_load_secondary_scalar(L, expected): + notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + result = notch_approximation_law.secondary(L) + + assert result.shape == (2, ) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_load_secondary_vectorized(): + notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + + result = notch_approximation_law.secondary([100.0, 400.0, 600.0]) + expected = [[100.0, 4.88e-4], [386.0, 1.99e-3], [533.0, 3.28e-3]] + + assert result.shape == (3, 2) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_binner_uninitialized(): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned) + + with pytest.raises(RuntimeError, match="NotchApproxBinner not initialized."): + binned.primary(100.0) + + with pytest.raises(RuntimeError, match="NotchApproxBinner not initialized."): + binned.secondary(100.0) + + +@pytest.mark.parametrize("L, expected", [ + (120.0, [148.0, 7.36e-4]), + (160.0, [193.7, 1.00e-3]), + (200.0, [193.7, 1.00e-3]), +]) +def test_binner_initialized_five_points_primary_scalar(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=4).initialize(max_load=200.0) + + result = binned.primary(L) + + np.testing.assert_allclose(result, np.array([expected]), rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (-120.0, [-148.0, -7.36e-4]), + (-160.0, [-193.7, -1.00e-3]), + (-200.0, [-193.7, -1.00e-3]), +]) +def test_binner_initialized_five_points_primary_scalar_symmetry(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=4).initialize(max_load=200.0) + + result = binned.primary(L) + + np.testing.assert_allclose(result, np.array([expected]), rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (120.0, [124.0, 6.10e-4]), + (160.0, [171.7, 8.67e-4]), + (200.0, [193.7, 1.00e-3]) +]) +def test_binner_initialized_nine_points_primary_scalar(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=8).initialize(max_load=200.0) + + result = binned.primary(L) + + np.testing.assert_allclose(result, np.array([expected]), rtol=1e-2) + + +def test_binner_initialized_nine_points_primary_out_of_scale(): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=8).initialize(max_load=200.0) + + with pytest.raises( + ValueError, + match="Requested load `400.0`, higher than initialized maximum load `200.0`", + ): + binned.primary(400.0) + + +@pytest.mark.parametrize("L, expected", [ + (120.0, [[148.0, 7.36e-4], [182.9, 9.34e-4], [214.3, 1.15e-3]]), + (160.0, [[193.7, 1.00e-3], [233.3, 1.30e-3], [266.7, 1.64e-3]]), + (200.0, [[193.7, 1.00e-3], [233.3, 1.30e-3], [266.7, 1.64e-3]]) +]) +def test_binner_initialized_five_points_primary_vectorized(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=4).initialize( + max_load=[200.0, 250.0, 300.0] + ) + + result = binned.primary([L, 1.25*L, 1.5*L]) + + assert result.shape == (3, 2) + + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (20.0, [297.0, 1.48e-3]), + (340.0, [533.0, 3.28e-3]), + (600.0, [533.0, 3.28e-3]) +]) +def test_binner_initialized_one_bin_secondary_scalar(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=1).initialize(max_load=300.0) + + result = binned.secondary(L) + + np.testing.assert_allclose(result, np.array([expected]), rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (-20.0, [-297.0, -1.48e-3]), + (-340.0, [-533.0, -3.28e-3]), + (-600.0, [-533.0, -3.28e-3]) +]) +def test_binner_initialized_one_bin_secondary_scalar_symmetry(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=1).initialize(max_load=300.0) + + result = binned.secondary(L) + + np.testing.assert_allclose(result, np.array([expected]), rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (20.0, [100.0, 4.88e-4]), + (400.0, [386.0, 1.99e-3]), + (600.0, [533.0, 3.28e-3]) +]) +def test_binner_initialized_three_points_secondary_scalar(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=3).initialize(max_load=300.0) + + result = binned.secondary(L) + + np.testing.assert_allclose(result, np.array([expected]), rtol=1e-2) + + +def test_binner_initialized_three_points_secondary_out_of_scale(): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=3).initialize(max_load=300.0) + + with pytest.raises( + ValueError, + match="Requested load `700.0`, higher than initialized maximum delta load `600.0`", + ): + binned.secondary(700.0) + + +@pytest.mark.parametrize("L, expected", [ + (20.0, [[100.0, 4.88e-4], [133.3, 6.47e-4], [200.0, 9.74e-4]]), + (400.0, [[386.0, 1.99e-3], [490.0, 2.81e-3], [638.2, 4.90e-3]]), + (600.0, [[533.0, 3.28e-3], [638.2, 4.90e-03], [784.8, 9.26e-3]]) +]) +def test_binner_initialized_three_points_secondary_vectorized(L, expected): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=3).initialize(max_load=[300.0, 400.0, 600.0]) + + result = binned.secondary(([L, 1.25*L, 1.5*L])) + + assert result.shape == (3, 2) + + np.testing.assert_allclose(result, expected, rtol=1e-2) From 7ccf6cab2eed3ec5cd77bd6fb9377aafa615b6c7 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Wed, 11 Dec 2024 16:40:59 +0100 Subject: [PATCH 02/25] Adjust Rainflow HCM to NotchApproxBinner Signed-off-by: Johannes Mueller --- src/pylife/stress/rainflow/fkm_nonlinear.py | 16 +-- tests/stress/rainflow/test_fkm_nonlinear.py | 137 ++++++++++---------- 2 files changed, 75 insertions(+), 78 deletions(-) diff --git a/src/pylife/stress/rainflow/fkm_nonlinear.py b/src/pylife/stress/rainflow/fkm_nonlinear.py index 7c45fde6..22deabbf 100644 --- a/src/pylife/stress/rainflow/fkm_nonlinear.py +++ b/src/pylife/stress/rainflow/fkm_nonlinear.py @@ -318,21 +318,15 @@ def _update_residuals(self, record_vals, turning_point, load_turning_points_rep) def _collect_record(self, load_turning_points, num_turning_points, record): def primary(_prev, load): - sigma = self._notch_approximation_law.stress(load) - epsilon = self._notch_approximation_law.strain(sigma, load) - return np.array([load, sigma, epsilon]) + return self._notch_approximation_law.primary(load) def secondary(prev, load): prev_load = prev[LOAD] delta_L = load - prev_load - delta_sigma = self._notch_approximation_law.stress_secondary_branch(delta_L) - delta_epsilon = self._notch_approximation_law.strain_secondary_branch(delta_sigma, delta_L) + delta = self._notch_approximation_law.secondary(delta_L) - sigma = prev[STRESS] + delta_sigma - epsilon = prev[STRAIN] + delta_epsilon - - return np.array([load, sigma, epsilon]) + return prev[STRESS:STRAIN+1].T + delta def determine_prev_record(prev_idx): if prev_idx < 0: @@ -372,7 +366,9 @@ def determine_prev_record(prev_idx): return record_vals.set_index("turning_point", drop=True, append=True) def _process_deformation(self, deformation_func, result_buf, load, prev_record): - result_buf[:3] = deformation_func(prev_record, load) + result_buf[0] = load + res = deformation_func(prev_record, load).T + result_buf[1:3] = res old_load = self._last_record[LOAD, 0] diff --git a/tests/stress/rainflow/test_fkm_nonlinear.py b/tests/stress/rainflow/test_fkm_nonlinear.py index 325f023a..ab32aac9 100644 --- a/tests/stress/rainflow/test_fkm_nonlinear.py +++ b/tests/stress/rainflow/test_fkm_nonlinear.py @@ -25,8 +25,8 @@ from pylife.stress.rainflow.fkm_nonlinear import FKMNonlinearDetector import pylife.stress.rainflow.recorders as RFR -import pylife.materiallaws.notch_approximation_law -import pylife.materiallaws.notch_approximation_law_seegerbeste +import pylife.materiallaws.notch_approximation_law as NAL +from pylife.materiallaws.notch_approximation_law_seegerbeste import SeegerBeste @pytest.fixture(autouse=True) @@ -54,12 +54,13 @@ def setUp(self): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) @@ -117,12 +118,13 @@ def setUp(self): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) @@ -194,12 +196,13 @@ def setUp(self): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) @@ -269,12 +272,13 @@ def setUp(self): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(self.signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) @@ -338,13 +342,14 @@ def setUp(self): np.testing.assert_allclose(max(abs(signal)), 1013) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) @@ -416,12 +421,12 @@ def test_edge_case_value_in_sample_tail_simple_signal(vals, expected_loads_min, n = 0.187 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100 + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load ) detector = FKMNonlinearDetector( @@ -466,12 +471,12 @@ def test_edge_case_value_in_sample_tail(vals, expected_loads_min, expected_loads n = 0.187 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100 + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load ) detector = FKMNonlinearDetector( @@ -512,12 +517,14 @@ def test_flush_edge_case_load(): n = 0.07 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = max(abs(pd.concat([signal_1, signal_2]))) + maximum_absolute_load = pd.concat([signal_1, signal_2]).abs().groupby("node_id").max() - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100 + print(maximum_absolute_load) + + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load ) detector = FKMNonlinearDetector( @@ -570,11 +577,7 @@ def test_flush_edge_case_load(): def test_flush_edge_case_load_simple_signal(): - signal_1 = np.array([0.0, 143.0, -287.0, 143.0, -359.0, 287.0, 0.0, 287.0, -287.0]) - - mi_2 = pd.MultiIndex.from_product([range(9, 17), range(3)], names=["load_step", "node_id"]) - signal_2 = np.array([143.0, -287.0, 143.0, -359.0, 287.0, 0.0, 287.0, -287.0]) E = 206e3 # [MPa] Young's modulus @@ -583,12 +586,12 @@ def test_flush_edge_case_load_simple_signal(): n = 0.07 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) maximum_absolute_load = max(abs(np.concatenate([signal_1, signal_2]))) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100 + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load ) detector = FKMNonlinearDetector( @@ -613,9 +616,7 @@ def test_flush_edge_case_load_simple_signal(): def test_flush_edge_case_S_simple_signal(): - signal_1 = np.array([0.0, 143.0, -287.0, 143.0, -359.0, 287.0, 0.0, 287.0, -287.0]) - signal_2 = np.array([143.0, -287.0, 143.0, -359.0, 287.0, 0.0, 287.0, -287.0]) E = 206e3 # [MPa] Young's modulus @@ -624,12 +625,12 @@ def test_flush_edge_case_S_simple_signal(): n = 0.07 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) maximum_absolute_load = max(abs(np.concatenate([signal_1, signal_2]))) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100 + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load ) detector = FKMNonlinearDetector( @@ -672,12 +673,15 @@ def test_flush_edge_case_S(): n = 0.07 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) + + maximum_absolute_load = pd.concat([signal_1, signal_2]).abs().groupby("node_id").max() - maximum_absolute_load = max(abs(pd.concat([signal_1, signal_2]))) + print(maximum_absolute_load) + print(max(abs(pd.concat([signal_1, signal_2])))) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100 + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load ) detector = FKMNonlinearDetector( @@ -714,11 +718,8 @@ def test_flush_edge_case_S(): S_min = detector.recorder.S_min S_max = detector.recorder.S_max - pd.testing.assert_series_equal(S_min, expected_S_min, check_index=False) - pd.testing.assert_series_equal(S_max, expected_S_max, check_index=False) - - pd.testing.assert_series_equal(S_min, expected_S_min) - pd.testing.assert_series_equal(S_max, expected_S_max) + pd.testing.assert_series_equal(S_min, expected_S_min, rtol=1e-1) + pd.testing.assert_series_equal(S_max, expected_S_max, rtol=1e-1) @pytest.mark.parametrize('vals, num', [ @@ -768,15 +769,15 @@ def test_edge_case_value_in_sample_tail_compare_simple(vals, num): n = 0.187 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) maximum_absolute_load_simple = max(abs(vals)) maximum_absolute_load_multiple = signal.abs().groupby('node_id').max() print("single") - extended_neuber_binned_simple = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load_simple, 100 + extended_neuber_binned_simple = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load_simple ) detector_simple = FKMNonlinearDetector( recorder=RFR.FKMNonlinearRecorder(), @@ -785,8 +786,8 @@ def test_edge_case_value_in_sample_tail_compare_simple(vals, num): detector_simple.process(vals).process(vals) print("multiple") - extended_neuber_binned_multiple = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load_multiple, 100 + extended_neuber_binned_multiple = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load_multiple ) detector_multiindex = FKMNonlinearDetector( recorder=RFR.FKMNonlinearRecorder(), @@ -852,16 +853,16 @@ def test_hcm_first_second(vals, num): n = 0.187 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) maximum_absolute_load_simple = max(abs(vals)) maximum_absolute_load_multiple = signal.abs().groupby('node_id').max() - extended_neuber_binned_simple = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load_simple, 100 + extended_neuber_binned_simple = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load_simple ) - extended_neuber_binned_multiple = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load_multiple, 100 + extended_neuber_binned_multiple = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load_multiple ) detector_simple = FKMNonlinearDetector( recorder=RFR.FKMNonlinearRecorder(), @@ -897,10 +898,8 @@ def detector_seeger_beste(): n = 0.187 # [-] K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - seeger_beste = pylife.materiallaws.notch_approximation_law_seegerbeste.SeegerBeste(E, K, n, K_p) - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.Binned( - seeger_beste, 800, 100 - ) + seeger_beste = SeegerBeste(E, K, n, K_p) + seeger_beste_binned = NAL.NotchApproxBinner(seeger_beste).initialize(800) return FKMNonlinearDetector( recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=seeger_beste_binned @@ -1097,12 +1096,13 @@ def test_history_guideline_at_once(): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber_binned) @@ -1144,12 +1144,13 @@ def test_history_guideline_at_split(split_point): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) # initialize notch approximation law - extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K, n, K_p) + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( + maximum_absolute_load + ) # first run detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber_binned) From 026cf3565474d7682e2e0fb8056464c73b19d623 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 12 Dec 2024 11:04:17 +0100 Subject: [PATCH 03/25] Apply NotchApproxBinner to Seeger Beste Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 25 +++++----- .../notch_approximation_law_seegerbeste.py | 18 ++++--- .../test_notch_approximation_law.py | 8 +-- ...est_notch_approximation_law_seegerbeste.py | 50 ++++++++++++++++++- 4 files changed, 74 insertions(+), 27 deletions(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index 8069f70f..9e015b89 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -79,6 +79,18 @@ def K(self, value): """Set the strain hardening coefficient""" self.K_prime = value + def primary(self, load): + load = np.asarray(load) + stress = self.stress(load) + strain = self.strain(stress, None) + return np.stack([stress, strain], axis=len(load.shape)) + + def secondary(self, delta_load): + delta_load = np.asarray(delta_load) + delta_stress = self.stress_secondary_branch(delta_load) + delta_strain = self.strain_secondary_branch(delta_stress, None) + return np.stack([delta_stress, delta_strain], axis=len(delta_load.shape)) + class ExtendedNeuber(NotchApproximationLawBase): r"""Implementation of the extended Neuber notch approximation material relation. @@ -198,19 +210,6 @@ def load(self, stress, *, rtol=1e-4, tol=1e-4): ) return load - def primary(self, load): - load = np.asarray(load) - stress = self.stress(load) - strain = self.strain(stress, None) - return np.stack([stress, strain], axis=len(load.shape)) - - def secondary(self, delta_load): - delta_load = np.asarray(delta_load) - delta_stress = self.stress_secondary_branch(delta_load) - delta_strain = self.strain_secondary_branch(delta_stress, None) - return np.stack([delta_stress, delta_strain], axis=len(delta_load.shape)) - - def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): """Calculate the stress on secondary branches in the stress-strain diagram at a given elastic-plastic stress (load), from a FE computation. diff --git a/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py b/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py index c8ab1370..8c78c18e 100644 --- a/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py +++ b/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py @@ -96,7 +96,7 @@ def stress(self, load, *, rtol=1e-4, tol=1e-4): # or (value, converged, zero_der) for vector-valued invocation # only for multiple points at once, if some points diverged - if len(stress) == 3 and sum(stress[1]) < len(stress[1]): + if len(stress) == 3 and not stress[1].all(): stress = self._stress_fix_not_converged_values(stress, load, x0, rtol, tol) return stress[0] @@ -189,6 +189,7 @@ def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): ''' # initial value as given by correction document to FKM nonlinear + delta_load = np.asarray(delta_load) x0 = delta_load * (1 - (1 - 1/self._K_p)/1000) # suppress the divergence warnings @@ -197,7 +198,7 @@ def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): delta_stress = optimize.newton( func=self._stress_secondary_implicit, - x0=np.asarray(x0), + x0=x0, args=([delta_load]), full_output=True, rtol=rtol, tol=tol, maxiter=50 @@ -208,7 +209,9 @@ def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): # or (value, converged, zero_der) for vector-valued invocation # only for multiple points at once, if some points diverged - if len(delta_stress) == 3 and sum(delta_stress[1]) < len(delta_stress[1]): + + multidim = len(x0.shape) > 1 and x0.shape[1] > 1 + if multidim and x0.shape[1] > 1 and not delta_stress[1].all(): delta_stress = self._stress_secondary_fix_not_converged_values(delta_stress, delta_load, x0, rtol, tol) return delta_stress[0] @@ -482,22 +485,21 @@ def _stress_secondary_fix_not_converged_values(self, delta_stress, delta_load, x '''For the values that did not converge in the previous vectorized call to optimize.newton, call optimize.newton again on the scalar value. This usually finds the correct solution.''' - indices_diverged = [index for index, is_converged in enumerate(delta_stress[1]) if not is_converged] + indices_diverged = np.where(~delta_stress[1].all(axis=1))[0] x0_array = np.asarray(x0) delta_load_array = np.asarray(delta_load) # recompute previously failed points individually for index_diverged in indices_diverged: - x0_diverged = x0_array[index_diverged] - delta_load_diverged = delta_load_array[index_diverged] + x0_diverged = x0_array[index_diverged, 0] + delta_load_diverged = delta_load_array[index_diverged, 0] result = optimize.newton( func=self._stress_secondary_implicit, - x0=x0_diverged, + x0=np.asarray(x0_diverged), args=([delta_load_diverged]), full_output=True, rtol=rtol, tol=tol, maxiter=50 ) - if result[1].converged: delta_stress[0][index_diverged] = result[0] return delta_stress diff --git a/tests/materiallaws/test_notch_approximation_law.py b/tests/materiallaws/test_notch_approximation_law.py index 2b1ede8c..e7eae318 100644 --- a/tests/materiallaws/test_notch_approximation_law.py +++ b/tests/materiallaws/test_notch_approximation_law.py @@ -200,7 +200,7 @@ def test_load(E, K, n, L): (175.0, [171.7, 8.67e-4]), (200.0, [193.7, 1.00e-3]) ]) -def test_load_primary_scalar(L, expected): +def test_primary_scalar(L, expected): notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) result = notch_approximation_law.primary(L) @@ -208,7 +208,7 @@ def test_load_primary_scalar(L, expected): np.testing.assert_allclose(result, expected, rtol=1e-2) -def test_load_primary_vectorized(): +def test_primary_vectorized(): notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) result = notch_approximation_law.primary([150.0, 175.0, 200.0]) @@ -223,7 +223,7 @@ def test_load_primary_vectorized(): (400.0, [386.0, 1.99e-3]), (600.0, [533.0, 3.28e-3]) ]) -def test_load_secondary_scalar(L, expected): +def test_secondary_scalar(L, expected): notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) result = notch_approximation_law.secondary(L) @@ -231,7 +231,7 @@ def test_load_secondary_scalar(L, expected): np.testing.assert_allclose(result, expected, rtol=1e-2) -def test_load_secondary_vectorized(): +def test_secondary_vectorized(): notch_approximation_law = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) result = notch_approximation_law.secondary([100.0, 400.0, 600.0]) diff --git a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py index f94d0fcc..4e76f65a 100644 --- a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py +++ b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py @@ -176,7 +176,7 @@ def test_seeger_beste_load(E, K, n, L): # initialize notch approximation law and damage parameter notch_approximation_law = SeegerBeste(E, K, n, K_p=3.5) - # The "load" method is the inverse operation of "stress", + # The "load" method is the inverse operation of "stress", # i.e., ``L = load(stress(L))`` and ``S = stress(load(stress))``. stress = notch_approximation_law.stress(L) load = notch_approximation_law.load(stress) @@ -200,7 +200,7 @@ def test_seeger_beste_load_secondary_branch(E, K, n, L): # initialize notch approximation law and damage parameter notch_approximation_law = SeegerBeste(E, K, n, K_p=3.5) - # The "load" method is the inverse operation of "stress", + # The "load" method is the inverse operation of "stress", # i.e., ``L = load(stress(L))`` and ``S = stress(load(stress))``. stress = notch_approximation_law.stress_secondary_branch(L) load = notch_approximation_law.load_secondary_branch(stress) @@ -209,3 +209,49 @@ def test_seeger_beste_load_secondary_branch(E, K, n, L): np.testing.assert_allclose(L, load, rtol=1e-3) np.testing.assert_allclose(stress, stress2, rtol=1e-3) + + +@pytest.mark.parametrize("L, expected", [ + (7.18, [7.18, 3.50e-5]), + (179.6, [172.9, 8.73e-4]), + (359.2, [283.7, 1.86e-3]) +]) # Values from FKM NL Guideline p. 135 +def test_primary_scalar(L, expected): + notch_approximation_law = SeegerBeste(E=206e3, K=1184., n=0.187, K_p=3.5) + result = notch_approximation_law.primary(L) + + assert result.shape == (2, ) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_primary_vectorized(): + notch_approximation_law = SeegerBeste(E=206e3, K=1184., n=0.187, K_p=3.5) + + result = notch_approximation_law.primary([7.18, 179.6, 359.2]) + expected = [[7.18, 3.50e-5], [172.9, 8.73e-4], [283.7, 1.86e-3]] + + assert result.shape == (3, 2) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +@pytest.mark.parametrize("L, expected", [ + (359.2, [345.9, 1.75e-3]), + (362.8, [348.9, 1.77e-3]), + (718.5, [567.7, 3.72e-3]) +]) +def test_secondary_scalar(L, expected): + notch_approximation_law = SeegerBeste(E=206e3, K=1184., n=0.187, K_p=3.5) + result = notch_approximation_law.secondary(L) + + assert result.shape == (2, ) + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_secondary_vectorized(): + notch_approximation_law = SeegerBeste(E=206e3, K=1184., n=0.187, K_p=3.5) + + result = notch_approximation_law.secondary([359.2, 362.8, 718.5]) + expected = [[345.9, 1.75e-3], [348.9, 1.77e-3], [567.7, 3.72e-3]] + + assert result.shape == (3, 2) + np.testing.assert_allclose(result, expected, rtol=1e-2) From ad4131b96900437976b0345d1f634b57ba715fd5 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 12 Dec 2024 11:04:46 +0100 Subject: [PATCH 04/25] Adjust other modules to NotchApproxBinner Signed-off-by: Johannes Mueller --- .../assessment_nonlinear_standard.py | 12 ++++++---- tests/strength/test_damage_calculator.py | 24 ++++++++++++------- 2 files changed, 24 insertions(+), 12 deletions(-) diff --git a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py index 060a8249..9c27851e 100644 --- a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py +++ b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py @@ -398,8 +398,11 @@ def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolu extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K_prime, n_prime, K_p) # wrap the notch approximation law by a binning class, which precomputes the values - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = ( + pylife.materiallaws.notch_approximation_law.NotchApproxBinner( + extended_neuber + ).initialize(maximum_absolute_load) + ) # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() @@ -491,8 +494,9 @@ def _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolu seeger_beste = pylife.materiallaws.notch_approximation_law_seegerbeste.SeegerBeste(E, K_prime, n_prime, K_p) # wrap the notch approximation law by a binning class, which precomputes the values - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.Binned( - seeger_beste, maximum_absolute_load, 100) + seeger_beste_binned = pylife.materiallaws.notch_approximation_law.NotchApproxBinner( + seeger_beste + ).initialize(maximum_absolute_load) # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() diff --git a/tests/strength/test_damage_calculator.py b/tests/strength/test_damage_calculator.py index 02c99624..96910856 100644 --- a/tests/strength/test_damage_calculator.py +++ b/tests/strength/test_damage_calculator.py @@ -94,8 +94,11 @@ def test_woehler_curve_P_RAM_collective_has_no_index(): # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(np.abs(load_sequence_list)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = ( + pylife.materiallaws.notch_approximation_law.NotchApproxBinner( + extended_neuber + ).initialize(maximum_absolute_load) + ) # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() @@ -200,8 +203,11 @@ def test_woehler_curve_P_RAM_collective_has_MultiIndex(): # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(np.abs(load_sequence_list)) - extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned( - extended_neuber, maximum_absolute_load, 100) + extended_neuber_binned = ( + pylife.materiallaws.notch_approximation_law.NotchApproxBinner( + extended_neuber + ).initialize(maximum_absolute_load) + ) # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() @@ -295,8 +301,9 @@ def test_woehler_curve_P_RAJ_has_no_index(): # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(np.abs(load_sequence_list)) - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.Binned( - seeger_beste, maximum_absolute_load, 100) + seeger_beste_binned = pylife.materiallaws.notch_approximation_law.NotchApproxBinner( + seeger_beste + ).initialize(maximum_absolute_load) # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() @@ -402,8 +409,9 @@ def test_woehler_curve_P_RAJ_has_MultiIndex(): # wrap the notch approximation law by a binning class, which precomputes the values maximum_absolute_load = max(np.abs(load_sequence_list)) - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.Binned( - seeger_beste, maximum_absolute_load, 100) + seeger_beste_binned = pylife.materiallaws.notch_approximation_law.NotchApproxBinner( + seeger_beste + ).initialize(maximum_absolute_load) # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() From eabdbe00d798cbb04a3ecef4a956ccd56d07f416 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 12 Dec 2024 13:28:05 +0100 Subject: [PATCH 05/25] Adjust HCM rainflow counting unit tests to NotchApproxBinner Signed-off-by: Johannes Mueller --- .../reference-fkm-nonlinear/reference_first_second-1.json | 2 +- .../reference_process-process-1.json | 2 +- tests/stress/rainflow/test_fkm_nonlinear.py | 7 ++----- 3 files changed, 4 insertions(+), 7 deletions(-) diff --git a/tests/stress/rainflow/reference-fkm-nonlinear/reference_first_second-1.json b/tests/stress/rainflow/reference-fkm-nonlinear/reference_first_second-1.json index c9b4b5d0..af849f5b 100644 --- a/tests/stress/rainflow/reference-fkm-nonlinear/reference_first_second-1.json +++ b/tests/stress/rainflow/reference-fkm-nonlinear/reference_first_second-1.json @@ -468,7 +468,7 @@ "(13, 2)":5.754102661, "(14, 0)":4.9721593167, "(14, 1)":3.8025932597, - "(14, 2)":2.5366373199, + "(14, 2)":2.53663732, "(15, 0)":0.0, "(15, 1)":0.0, "(15, 2)":0.0, diff --git a/tests/stress/rainflow/reference-fkm-nonlinear/reference_process-process-1.json b/tests/stress/rainflow/reference-fkm-nonlinear/reference_process-process-1.json index ef3af193..f8a9be9f 100644 --- a/tests/stress/rainflow/reference-fkm-nonlinear/reference_process-process-1.json +++ b/tests/stress/rainflow/reference-fkm-nonlinear/reference_process-process-1.json @@ -468,7 +468,7 @@ "(13, 2)":5.754102661, "(14, 0)":4.9721593167, "(14, 1)":3.8025932597, - "(14, 2)":2.5366373199, + "(14, 2)":2.53663732, "(15, 0)":0.0, "(15, 1)":0.0, "(15, 2)":0.0, diff --git a/tests/stress/rainflow/test_fkm_nonlinear.py b/tests/stress/rainflow/test_fkm_nonlinear.py index ab32aac9..aba64c5c 100644 --- a/tests/stress/rainflow/test_fkm_nonlinear.py +++ b/tests/stress/rainflow/test_fkm_nonlinear.py @@ -677,9 +677,6 @@ def test_flush_edge_case_S(): maximum_absolute_load = pd.concat([signal_1, signal_2]).abs().groupby("node_id").max() - print(maximum_absolute_load) - print(max(abs(pd.concat([signal_1, signal_2])))) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( maximum_absolute_load ) @@ -694,9 +691,9 @@ def test_flush_edge_case_S(): expected_S_min = pd.Series( [ -49.096964, -57.761134, -11.552227, -96.749900, -115.522268, -21.660425, -9.610801e-11, - -1.143832e-10, -1.241333e-07, -96.749900, -115.522268, -21.660425, + -1.143832e-10, -1.549e-07, -96.749900, -115.522268, -21.660425, -96.749900, -115.522268, -21.660425, -121.298382, -144.402835, -27.436539, - -9.610801e-11, -1.143832e-10, -1.241333e-07, + -9.610801e-11, -1.143832e-10, -1.549e-07, ], index=pd.MultiIndex.from_product( [[1, 2, 6, 8, 10, 4, 14], range(3)], names=["load_step", "node_id"] From a5e22bb34a8b9a667a1810f390b161c697ae22a6 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 08:34:20 +0100 Subject: [PATCH 06/25] Add docstrings to NotchApproxBinner Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 79 +++++++++++++++++++ 1 file changed, 79 insertions(+) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index 9e015b89..6a5d97d9 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -922,6 +922,31 @@ def _create_bins_multiple_assessment_points(self): class NotchApproxBinner: + """Binning for notch approximation laws, as described in FKM nonlinear 2.5.8.2, p.55. + The implicitly defined stress function of the notch approximation law is precomputed + for various loads at a fixed number of equispaced `bins`. The values are stored in two + look-up tables for the primary and secondary branches of the stress-strain hysteresis + curves. When stress and strain values are needed for a given load, the nearest value + of the corresponding bin is retrived. This is faster than invoking the nonlinear + root finding algorithm for every new load. + + There are two variants of the data structure. + + * First, for a single assessment point, the lookup-table contains one load, + strain and stress value in every bin. + * Second, for vectorized assessment of multiple nodes at once, the lookup-table + contains at every load bin an array with stress and strain values for every node. + The representative load, stored in the lookup table and used for the lookup + is the first element of the given load array. + + Parameters + ---------- + notch_approximation_law : NotchApproximationLawBase + The law for the notch approximation to be used. + + number_of_bins : int, optional + The number of bins in the lookup table, default 100 + """ def __init__(self, notch_approximation_law, number_of_bins=100): self._n_bins = number_of_bins @@ -930,6 +955,18 @@ def __init__(self, notch_approximation_law, number_of_bins=100): self._max_load_rep = None def initialize(self, max_load): + """Initialize with a maximum expected load. + + Parameters + ---------- + max_load : array_like + The state of the maximum nominal load that is expected. The first + element is chosen as representative to caclulate the lookup table. + + Returns + ------- + self + """ max_load = np.asarray(max_load) self._max_load_rep, _ = self._rep_abs_and_sign(max_load) @@ -942,6 +979,27 @@ def initialize(self, max_load): return self def primary(self, load): + """Lookup the stress strain of the primary branch. + + Parameters + ---------- + load : array-like + The load as argument for the stress strain laws. + If non-scalar, the first element will be used to look it up in the + lookup table. + + Returns + ------- + stress strain : ndarray + The resulting stress strain data. + + If the argument is scalar, the resulting array is of the strucuture + ``[<σ>, <ε>]`` + + If the argument is an 1D-array with length `n`the resulting array is of the + structure ``[[<σ1>, <σ2>, <σ3>, ... <σn>], [<ε1>, <ε2>, <ε3>, ... <εn>]]`` + + """ self._raise_if_uninitialized() load_rep, sign = self._rep_abs_and_sign(load) @@ -953,6 +1011,27 @@ def primary(self, load): return sign * self._lut_primary[idx, :] def secondary(self, delta_load): + """Lookup the stress strain of the secondary branch. + + Parameters + ---------- + load : array-like + The load as argument for the stress strain laws. + If non-scalar, the first element will be used to look it up in the + lookup table. + + Returns + ------- + stress strain : ndarray + The resulting stress strain data. + + If the argument is scalar, the resulting array is of the strucuture + ``[<σ>, <ε>]`` + + If the argument is an 1D-array with length `n`the resulting array is of the + structure ``[[<σ1>, <σ2>, <σ3>, ..., <σn>], [<ε1>, <ε2>, <ε3>, ..., <εn>]]`` + + """ self._raise_if_uninitialized() delta_load_rep, sign = self._rep_abs_and_sign(delta_load) From 4c648a0d55e204f7431bd036269356ba046429e1 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 10:33:53 +0100 Subject: [PATCH 07/25] Drop Binned class for notch approximation and adjust tests accordingly Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 446 ------------------ .../test_notch_approximation_law.py | 112 +---- ...est_notch_approximation_law_seegerbeste.py | 100 ++-- 3 files changed, 53 insertions(+), 605 deletions(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index 6a5d97d9..37e67abd 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -474,452 +474,6 @@ def _d_load_secondary_implicit(self, delta_load, delta_stress): - delta_load/delta_stress * self.K_p * self._d_delta_e_star(delta_load) -class Binned: - """Binning for notch approximation laws, as described in FKM nonlinear 2.5.8.2, p.55. - The implicitly defined stress function of the notch approximation law is precomputed - for various loads at a fixed number of equispaced `bins`. The values are stored in two - look-up tables for the primary and secondary branches of the stress-strain hysteresis - curves. When stress and strain values are needed for a given load, the nearest value - of the corresponding bin is retrived. This is faster than invoking the nonlinear - root finding algorithm for every new load. - - There are two variants of the data structure. - - * First, for a single assessment point, the lookup-table contains one load, - strain and stress value in every bin. - * Second, for vectorized assessment of multiple nodes at once, the lookup-table - contains at every load bin a list with stress and strain values for every node. - The DataFrame has a multi-index over class_index and node_id. - - """ - - def __init__(self, notch_approximation_law, maximum_absolute_load, number_of_bins=100): - self._notch_approximation_law = notch_approximation_law - self._maximum_absolute_load = maximum_absolute_load - self._number_of_bins = number_of_bins - self._create_bins() - - @property - def ramberg_osgood_relation(self): - """The ramberg osgood relation object - """ - return self._notch_approximation_law.ramberg_osgood_relation - - def stress(self, load, *, rtol=1e-5, tol=1e-6): - """The stress of the primary path in the stress-strain diagram at a given load - by using the value of the look-up table. - - .. note:: - - The exact value would be computed by ``self._notch_approximation_law.stress(load)``. - - Parameters - ---------- - load : array-like float - The load, either a scalar value or a pandas DataFrame with RangeIndex (no named index) - rtol : float, optional - The relative tolerance to which the implicit formulation of the stress gets solved. - In this case for the `Binning` class, the parameter is not used. - tol : float, optional - The absolute tolerance to which the implicit formulation of the stress gets solved. - In this case for the `Binning` class, the parameter is not used. - - Returns - ------- - stress : array-like float - The resulting stress - """ - - - - # FIXME consolidate theese methods (duplicated code) - - load = np.asarray(load) - sign = np.sign(load) - - # if the assessment is performed for multiple points at once, i.e. load is a DataFrame with values for every node - if isinstance(self._lut_primary_branch.index, pd.MultiIndex): - - # the lut is a DataFrame with MultiIndex with levels class_index and node_id - - # find the corresponding class only for the first node, use the result for all nodes - first_node_id = self._lut_primary_branch.index.get_level_values("node_id")[0] - lut_for_first_node = self._lut_primary_branch.load[self._lut_primary_branch.index.get_level_values("node_id")==first_node_id] - first_abs_load = abs(load[0]) - - # get the class index of the corresponding bin/class - class_index = lut_for_first_node.searchsorted(first_abs_load) - - max_class_index = max(self._lut_primary_branch.index.get_level_values("class_index")) - - # raise error if requested load is higher than initialized maximum absolute load - if class_index+1 > max_class_index: - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of {first_abs_load} is requested (in stress()).") - - # get stress from matching class, "+1", because the next higher class is used - stress = self._lut_primary_branch[self._lut_primary_branch.index.get_level_values("class_index") == class_index+1].stress - - # multiply with sign - return sign * stress.to_numpy() - - # if the assessment is performed for multiple points at once, but only one lookup-table is used - elif isinstance(load, pd.Series): - - load = load.fillna(0) - sign = sign.fillna(0) - - index = self._lut_primary_branch.load.searchsorted(np.abs(load.values))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_primary_branch)): - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of |{load.max()}| is requested (in stress()).") - - return sign.values.flatten() * self._lut_primary_branch.iloc[(index+1).flatten()].stress.reset_index(drop=True) # "+1", because the next higher class is used - - # if the assessment is done only for one value, i.e. load is a scalar - else: - - index = self._lut_primary_branch.load.searchsorted(np.abs(load))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_primary_branch)): - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of |{load}| is requested (in stress()).") - - return sign * self._lut_primary_branch.iloc[index+1].stress # "+1", because the next higher class is used - - def strain(self, stress, load): - """Get the strain of the primary path in the stress-strain diagram at a given stress and load - by using the value of the look-up table. - - This method performs the task for for multiple points at once, - i.e. delta_load is a DataFrame with values for every node. - - Parameters - ---------- - load : array-like float - The load - - Returns - ------- - strain : array-like float - The resulting strain - """ - load = np.asarray(load) - sign = np.sign(load) - - # if the assessment is performed for multiple points at once, i.e. load is a DataFrame with values for every node - if isinstance(self._lut_primary_branch.index, pd.MultiIndex): - - # the lut is a DataFrame with MultiIndex with levels class_index and node_id - - # find the corresponding class only for the first node, use the result for all nodes - first_node_id = self._lut_primary_branch.index.get_level_values("node_id")[0] - lut_for_first_node = self._lut_primary_branch.load[self._lut_primary_branch.index.get_level_values("node_id")==first_node_id] - first_abs_load = abs(load[0]) - - # get the class index of the corresponding bin/class - class_index = lut_for_first_node.searchsorted(first_abs_load) - - max_class_index = max(self._lut_primary_branch.index.get_level_values("class_index")) - - # raise error if requested load is higher than initialized maximum absolute load - if class_index+1 > max_class_index: - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of {first_abs_load} is requested (in strain()).") - - # get strain from matching class, "+1", because the next higher class is used - strain = self._lut_primary_branch[self._lut_primary_branch.index.get_level_values("class_index") == class_index+1].strain - - # multiply with sign - return sign * strain.to_numpy() - - # if the assessment is performed for multiple points at once, but only one lookup-table is used - elif isinstance(load, pd.Series): - - load = load.fillna(0) - sign = sign.fillna(0) - - index = self._lut_primary_branch.load.searchsorted(np.abs(load.values))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_primary_branch)): - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of |{load.max()}| is requested (in strain()).") - - return sign.values.flatten() * self._lut_primary_branch.iloc[(index+1).flatten()].strain.reset_index(drop=True) # "+1", because the next higher class is used - - # if the assessment is done only for one value, i.e. load is a scalar - else: - - index = self._lut_primary_branch.load.searchsorted(np.abs(load))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_primary_branch)): - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of |{load}| is requested (in strain()).") - - return sign * self._lut_primary_branch.iloc[index+1].strain # "+1", because the next higher class is used - - def stress_secondary_branch(self, delta_load, *, rtol=1e-5, tol=1e-6): - """Get the stress on secondary branches in the stress-strain diagram at a given load - by using the value of the look-up table (lut). - - This method performs the task for for multiple points at once, - i.e. delta_load is a DataFrame with values for every node. - - Parameters - ---------- - delta_load : array-like float - The load increment of the hysteresis - rtol : float, optional - The relative tolerance to which the implicit formulation of the stress gets solved. - In this case for the `Binning` class, the parameter is not used. - tol : float, optional - The absolute tolerance to which the implicit formulation of the stress gets solved. - In this case for the `Binning` class, the parameter is not used. - - Returns - ------- - delta_stress : array-like float - The resulting stress increment within the hysteresis - """ - - delta_load = np.asarray(delta_load) - sign = np.sign(delta_load) - - # if the assessment is performed for multiple points at once, i.e. load is a DataFrame with values for every node - if isinstance(self._lut_primary_branch.index, pd.MultiIndex): - - # the lut is a DataFrame with MultiIndex with levels class_index and node_id - - # find the corresponding class only for the first node, use the result for all nodes - first_node_id = self._lut_primary_branch.index.get_level_values("node_id")[0] - lut_for_first_node = self._lut_secondary_branch.delta_load[self._lut_secondary_branch.index.get_level_values("node_id")==first_node_id] - first_abs_load = abs(delta_load[0]) - - # get the class index of the corresponding bin/class - class_index = lut_for_first_node.searchsorted(first_abs_load) - - max_class_index = max(self._lut_secondary_branch.index.get_level_values("class_index")) - - # raise error if requested load is higher than initialized maximum absolute load - if class_index+1 > max_class_index: - raise ValueError(f"Binned class is initialized with a maximum absolute delta_load load of {2*self._maximum_absolute_load}, "\ - f" but a higher absolute delta_load value of {first_abs_load} is requested (in stress_secondary_branch()).") - - # get stress from matching class, "+1", because the next higher class is used - delta_stress = self._lut_secondary_branch[self._lut_secondary_branch.index.get_level_values("class_index") == class_index+1].delta_stress - - # multiply with sign - return sign * delta_stress.to_numpy() - - # if the assessment is performed for multiple points at once, but only one lookup-table is used - elif isinstance(delta_load, pd.Series): - - delta_load = delta_load.fillna(0) - sign = sign.fillna(0) - - index = self._lut_secondary_branch.delta_load.searchsorted(np.abs(delta_load.values))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_secondary_branch)): - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of |{delta_load.max()}| is requested (in stress_secondary_branch()).") - - return sign.values.flatten() * self._lut_secondary_branch.iloc[(index+1).flatten()].delta_stress.reset_index(drop=True) # "+1", because the next higher class is used - - # if the assessment is done only for one value, i.e. load is a scalar - else: - - index = self._lut_secondary_branch.delta_load.searchsorted(np.abs(delta_load))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_secondary_branch)): - raise ValueError(f"Binned class is initialized for a maximum absolute delta_load of {2*self._maximum_absolute_load}, "\ - f" but a higher absolute delta_load value of |{delta_load}| is requested (in stress_secondary_branch()).") - - return sign * np.asarray(self._lut_secondary_branch.iloc[index+1].delta_stress) - - def strain_secondary_branch(self, delta_stress, delta_load): - """Get the strain on secondary branches in the stress-strain diagram at a given stress and load - by using the value of the look-up table (lut). - The lut is a DataFrame with MultiIndex with levels class_index and node_id. - - This method performs the task for for multiple points at once, - i.e. delta_load is a DataFrame with values for every node. - - Parameters - ---------- - delta_load : array-like float - The load increment - - Returns - ------- - strain : array-like float - The resulting strain - """ - #return self._notch_approximation_law.strain_secondary_branch(delta_stress, delta_load) - - delta_load = np.asarray(delta_load) - sign = np.sign(delta_load) - - # if the assessment is performed for multiple points at once, i.e. load is a DataFrame with values for every node - if isinstance(self._lut_primary_branch.index, pd.MultiIndex): - - # the lut is a DataFrame with MultiIndex with levels class_index and node_id - - # find the corresponding class only for the first node, use the result for all nodes - first_node_id = self._lut_secondary_branch.index.get_level_values("node_id")[0] - lut_for_first_node = self._lut_secondary_branch.delta_load[self._lut_secondary_branch.index.get_level_values("node_id")==first_node_id] - first_abs_load = abs(delta_load[0]) - - # get the class index of the corresponding bin/class - class_index = lut_for_first_node.searchsorted(first_abs_load) - - max_class_index = max(self._lut_secondary_branch.index.get_level_values("class_index")) - - # raise error if requested load is higher than initialized maximum absolute load - if class_index+1 > max_class_index: - raise ValueError(f"Binned class is initialized with a maximum absolute delta_load of {2*self._maximum_absolute_load}, "\ - f" but a higher absolute delta_load value of {first_abs_load} is requested (in strain_secondary_branch()).") - - # get strain from matching class, "+1", because the next higher class is used - delta_strain = self._lut_secondary_branch[self._lut_secondary_branch.index.get_level_values("class_index") == class_index+1].delta_strain - - # multiply with sign - return sign * delta_strain.to_numpy() - - # if the assessment is performed for multiple points at once, but only one lookup-table is used - elif isinstance(delta_load, pd.Series): - - delta_load = delta_load.fillna(0) - sign = sign.fillna(0) - - index = self._lut_secondary_branch.delta_load.searchsorted(np.abs(delta_load.values))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_secondary_branch)): - raise ValueError(f"Binned class is initialized with a maximum absolute load of {self._maximum_absolute_load}, "\ - f" but a higher absolute load value of |{delta_load.max()}| is requested (in strain_secondary_branch()).") - - return sign.values.flatten() * self._lut_secondary_branch.iloc[(index+1).flatten()].delta_strain.reset_index(drop=True) # "+1", because the next higher class is used - - # if the assessment is done only for one value, i.e. load is a scalar - else: - - index = self._lut_secondary_branch.delta_load.searchsorted(np.abs(delta_load))-1 # "-1", transform to zero-based indices - - # raise error if requested load is higher than initialized maximum absolute load - if np.any(index+1 >= len(self._lut_secondary_branch)): - raise ValueError(f"Binned class is initialized for a maximum absolute delta_load of {2*self._maximum_absolute_load}, "\ - f" but a higher absolute delta_load value of |{delta_load}| is requested (in strain_secondary_branch()).") - - return sign * np.asarray(self._lut_secondary_branch.iloc[index+1].delta_strain) # "-1", transform to zero-based indices - - def _create_bins(self): - """Initialize the lookup tables by precomputing the notch approximation law values. - """ - # for multiple assessment points at once use a Series with MultiIndex - if isinstance(self._maximum_absolute_load, pd.Series): - assert self._maximum_absolute_load.index.name == "node_id" - - self._create_bins_multiple_assessment_points() - - # for a single assessment point use the standard data structure - else: - self._create_bins_single_assessment_point() - - def _create_bins_single_assessment_point(self): - """Initialize the lookup tables by precomputing the notch approximation law values, - for the case of scalar variables, i.e., only a single assessment point.""" - - # create look-up table (lut) for the primary branch values, named PFAD in FKM nonlinear - self._lut_primary_branch = pd.DataFrame(0, - index=pd.Index(np.arange(1, self._number_of_bins+1), name="class_index"), - columns=["load", "strain", "stress"]) - - self._lut_primary_branch.load \ - = self._lut_primary_branch.index/self._number_of_bins * self._maximum_absolute_load - - self._lut_primary_branch.stress \ - = self._notch_approximation_law.stress(self._lut_primary_branch.load) - - self._lut_primary_branch.strain \ - = self._notch_approximation_law.strain( - self._lut_primary_branch.stress, self._lut_primary_branch.load) - - # create look-up table (lut) for the secondary branch values, named AST in FKM nonlinear - # Note that this time, we used twice the number of entries with the same bin width. - self._lut_secondary_branch = pd.DataFrame(0, - index=pd.Index(np.arange(1, 2*self._number_of_bins+1), name="class_index"), - columns=["delta_load", "delta_strain", "delta_stress"]) - - self._lut_secondary_branch.delta_load \ - = self._lut_secondary_branch.index/self._number_of_bins * self._maximum_absolute_load - - self._lut_secondary_branch.delta_stress \ - = self._notch_approximation_law.stress_secondary_branch(self._lut_secondary_branch.delta_load) - - self._lut_secondary_branch.delta_strain \ - = self._notch_approximation_law.strain_secondary_branch( - self._lut_secondary_branch.delta_stress, self._lut_secondary_branch.delta_load) - - def _create_bins_multiple_assessment_points(self): - """Initialize the lookup tables by precomputing the notch approximation law values, - for the case of vector-valued variables caused by an assessment on multiple points at once.""" - - self._maximum_absolute_load.name = "max_abs_load" - - # create look-up table (lut) for the primary branch values, named PFAD in FKM nonlinear - index = pd.MultiIndex.from_product([np.arange(1, self._number_of_bins+1), self._maximum_absolute_load.index], - names = ["class_index", "node_id"]) - - self._lut_primary_branch = pd.DataFrame(0, index=index, columns=["load", "strain", "stress"]) - - # create cartesian product of class index and max load - a = pd.DataFrame({"class_index": np.arange(1, self._number_of_bins+1)}) - class_index_with_max_load = a.merge(self._maximum_absolute_load, how='cross') - - # calculate load of each bin - load = class_index_with_max_load.class_index.astype(float)/self._number_of_bins * class_index_with_max_load.max_abs_load - load.index = index - self._lut_primary_branch.load = load - - self._lut_primary_branch.stress \ - = self._notch_approximation_law.stress(self._lut_primary_branch.load) - - self._lut_primary_branch.strain \ - = self._notch_approximation_law.strain( - self._lut_primary_branch.stress, self._lut_primary_branch.load) - - # ---------------- - # create look-up table (lut) for the secondary branch values, named AST in FKM nonlinear - # Note that this time, we used twice the number of entries with the same bin width. - index = pd.MultiIndex.from_product([np.arange(1, 2*self._number_of_bins+1), self._maximum_absolute_load.index], - names = ["class_index", "node_id"]) - - self._lut_secondary_branch = pd.DataFrame(0, index=index, columns=["delta_load", "delta_strain", "delta_stress"]) - - # create cartesian product of class index and max load - a = pd.DataFrame({"class_index": np.arange(1, 2*self._number_of_bins+1)}) - class_index_with_max_load = a.merge(self._maximum_absolute_load, how='cross') - - # calculate load of each bin - delta_load = class_index_with_max_load.class_index.astype(float)/self._number_of_bins * class_index_with_max_load.max_abs_load - delta_load.index = index - self._lut_secondary_branch.delta_load = delta_load - - self._lut_secondary_branch.delta_stress \ - = self._notch_approximation_law.stress_secondary_branch(self._lut_secondary_branch.delta_load) - - self._lut_secondary_branch.delta_strain \ - = self._notch_approximation_law.strain_secondary_branch( - self._lut_secondary_branch.delta_stress, self._lut_secondary_branch.delta_load) - - class NotchApproxBinner: """Binning for notch approximation laws, as described in FKM nonlinear 2.5.8.2, p.55. diff --git a/tests/materiallaws/test_notch_approximation_law.py b/tests/materiallaws/test_notch_approximation_law.py index e7eae318..cafffe0e 100644 --- a/tests/materiallaws/test_notch_approximation_law.py +++ b/tests/materiallaws/test_notch_approximation_law.py @@ -19,105 +19,45 @@ import pytest import numpy as np -import pandas as pd - -import pylife.strength.damage_parameter from pylife.materiallaws.notch_approximation_law import ExtendedNeuber, NotchApproxBinner from .data import * -def test_extended_neuber_example_1(): - """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ - - E = 206e3 # [MPa] Young's modulus - K = 1184 # [MPa] - n = 0.187 # [-] - K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - - L = pd.Series([100, -200, 100, -250, 200, 0, 200, -200]) - c = 1.4 - gamma_L = (250+6.6)/250 - L = c * gamma_L * L - - # initialize notch approximation law and damage parameter - notch_approximation_law = ExtendedNeuber(E, K, n, K_p) - - assert notch_approximation_law.E == E - assert notch_approximation_law.K == K - assert notch_approximation_law.n == n - assert notch_approximation_law.K_p == K_p - - maximum_absolute_load = max(abs(L)) - - # the FKM example seems to round here, real value is 359.24 - assert np.isclose(maximum_absolute_load, 359.3, rtol=1e-3) - - binned_notch_approximation_law = pylife.materiallaws.notch_approximation_law.Binned( - notch_approximation_law, maximum_absolute_load, 100) - - # some rows of PFAD are given in the FKM nonlinear example on p.76 - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[0], \ - pd.Series([3.592, 0.0017e-2, 3.592]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[1], \ - pd.Series([7.185, 0.0035e-2, 7.185]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - # Note that the original FKM document has an error at this row (it is row 49, not 50) - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[48], \ - pd.Series([176.057, 8.71e-4, 172.639]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[99], \ - pd.Series([359.3, 0.0021, 299.78]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - # matrix AST on page 162, chapter 3.4.1 - pd.testing.assert_frame_equal( - binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_162, rtol=1e-3, atol=1e-5) +@pytest.mark.parametrize('load, strain, stress', [ + (3.592, 0.0017e-2, 3.592), + (7.185, 0.0035e-2, 7.185), + (176.057, 8.71e-4, 172.639), + (359.3, 0.0021, 299.78) +]) +def test_extended_neuber_example_1(load, strain, stress): + """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ + notch_approximation_law = ExtendedNeuber(E=206e3, K=1184, n=0.187, K_p=3.5) - assert np.isclose(binned_notch_approximation_law._lut_primary_branch.load.max(), maximum_absolute_load) + assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) -def test_extended_neuber_example_2(): +@pytest.mark.parametrize('load, strain, stress', [ + (10.13, 0.0049e-2, 10.130), + (20.26, 0.0098e-2, 20.260), + (1013.00, 0.6035e-2, 829.681), +]) +def test_extended_neuber_example_2(load, strain, stress): """ example under 2.7.2, p.78 of FKM nonlinear, "Welle mit V-Kerbe" """ + notch_approximation_law = ExtendedNeuber(E=206e3, K=2650.5, n=0.187, K_p=3.5) - E = 206e3 # [MPa] Young's modulus - K = 2650.5 # [MPa] - n = 0.187 # [-] - K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - - L = 1266.25 * pd.Series([0.3, -0.3, 0.5, -0.5, 0.6, -0.6, 0.3, -0.3, 0.7, -0.7, 0.2, -0.2, 0.6, -0.6, 0.8, -0.8, 0.8, -0.8]) - - # initialize notch approximation law and damage parameter - notch_approximation_law = ExtendedNeuber(E, K, n, K_p) - - maximum_absolute_load = max(abs(L)) - - # the FKM example seems to round here, real value is 359.24 - assert np.isclose(maximum_absolute_load, 1013, rtol=1e-3) - - binned_notch_approximation_law = pylife.materiallaws.notch_approximation_law.Binned( - notch_approximation_law, maximum_absolute_load) - - # some rows of PFAD are given in the FKM nonlinear example on p.79 - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[0], \ - pd.Series([10.13, 0.0049e-2, 10.130]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[1], \ - pd.Series([20.26, 0.0098e-2, 20.260]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - # this row seems off in the FKM nonlinear example, error is as high as 5% - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[99], \ - pd.Series([1013.00, 0.6035e-2, 829.681]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - # test binning directly - assert np.isclose(binned_notch_approximation_law.stress(10.13), 10.130) - assert np.isclose(binned_notch_approximation_law.strain(10.30, 10.13), 0.0049174e-2) + assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) - assert np.isclose(binned_notch_approximation_law.stress_secondary_branch(10.13), 10.130) - assert np.isclose(binned_notch_approximation_law.strain_secondary_branch(10.30, 10.13), 0.0049174e-2) + assert np.isclose(notch_approximation_law.stress_secondary_branch(10.13), 10.130) + assert np.isclose(notch_approximation_law.strain_secondary_branch(10.30, 10.13), 0.0049174e-2, rtol=1e-1) # matrix AST on page 171, chapter 3.4.2 - pd.testing.assert_frame_equal( - binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_171, rtol=1e-3, atol=1e-5) + # pd.testing.assert_frame_equal( + # binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_171, rtol=1e-3, atol=1e-5) def test_extended_neuber_example_no_binning_vectorized(): diff --git a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py index 4e76f65a..1b065e4a 100644 --- a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py +++ b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py @@ -26,91 +26,45 @@ from pylife.materiallaws.notch_approximation_law_seegerbeste import SeegerBeste from .data import * - -def test_seeger_beste_example_1(): +@pytest.mark.parametrize('load, strain, stress', [ + (3.592, 0.0017e-2, 3.592), + (7.185, 0.0035e-2, 7.185), + (176.057, 8.71e-4, 172.639), + (359.3, 0.0019, 283.86) +]) +def test_seeger_beste_example_1(load, strain, stress): """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ + notch_approximation_law = SeegerBeste(E=206e3, K=1184, n=0.187, K_p=3.5) - E = 206e3 # [MPa] Young's modulus - K = 1184 # [MPa] - n = 0.187 # [-] - K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - - L = pd.Series([100, -200, 100, -250, 200, 0, 200, -200]) - c = 1.4 - gamma_L = (250+6.6)/250 - L = c * gamma_L * L - - # initialize notch approximation law and damage parameter - #notch_approximation_law = ExtendedNeuber(E, K, n, K_p) - notch_approximation_law = SeegerBeste(E, K, n, K_p) -# damage_parameter_type = pylife.strength.damage_parameter.P_RAM() - - maximum_absolute_load = max(abs(L)) -# print(maximum_absolute_load) - - # the FKM example seems to round here, real value is 359.24 - assert np.isclose(maximum_absolute_load, 359.3, rtol=1e-3) - - binned_notch_approximation_law = pylife.materiallaws.notch_approximation_law.Binned(notch_approximation_law, maximum_absolute_load) - - # some rows of PFAD are given in the FKM nonlinear example on p.76 - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[0], \ - pd.Series([3.59, 0.0017e-2, 3.59]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) + assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[1], \ - pd.Series([7.18, 0.0035e-2, 7.18]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - # this row seems off in the FKM nonlinear example, error is as high as 5% - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[49], \ - #pd.Series([179.62, 8.73e-4, 172.93]), check_names=False, check_index=False, rtol=10e-2, atol=1e-5) - pd.Series([179.62, 8.73e-4, 172.93]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[99], \ - pd.Series([359.24, 0.001859, 283.86]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) # matrix AST on page 162, chapter 3.4.1 - pd.testing.assert_frame_equal( - binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_162_seeger_beste, rtol=1e-3, atol=1e-5) + # pd.testing.assert_frame_equal( + # binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_162_seeger_beste, rtol=1e-3, atol=1e-5) - assert np.isclose(binned_notch_approximation_law._lut_primary_branch.load.max(), maximum_absolute_load) + # assert np.isclose(binned_notch_approximation_law._lut_primary_branch.load.max(), maximum_absolute_load) -def test_seeger_beste_example_2(): +@pytest.mark.parametrize('load, strain, stress', [ + (10.13, 0.0049e-2, 10.130), + (20.26, 0.0098e-2, 20.260), + (1013.00, 0.53e-2, 784.89), +]) +def test_seeger_beste_example_2(load, strain, stress): """ example under 2.7.2, p.78 of FKM nonlinear, "Welle mit V-Kerbe" """ + notch_approximation_law = SeegerBeste(E=206e3, K=2650.5, n=0.187, K_p=3.5) - E = 206e3 # [MPa] Young's modulus - K = 2650.5 # [MPa] - n = 0.187 # [-] - K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - - L = 1266.25 * pd.Series([0.3, -0.3, 0.5, -0.5, 0.6, -0.6, 0.3, -0.3, 0.7, -0.7, 0.2, -0.2, 0.6, -0.6, 0.8, -0.8, 0.8, -0.8]) - - # initialize notch approximation law and damage parameter - notch_approximation_law = SeegerBeste(E, K, n, K_p) - - maximum_absolute_load = max(abs(L)) - - # the FKM example seems to round here, real value is 359.24 - assert np.isclose(maximum_absolute_load, 1013, rtol=1e-3) - - binned_notch_approximation_law = pylife.materiallaws.notch_approximation_law.Binned( - notch_approximation_law, maximum_absolute_load) - - # some rows of PFAD are given in the FKM nonlinear example on p.79 - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[0], \ - pd.Series([10.13, 0.0049e-2, 10.130]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[1], \ - pd.Series([20.26, 0.01e-2, 20.260]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) - - # this row seems off in the FKM nonlinear example, error is as high as 5% - pd.testing.assert_series_equal(binned_notch_approximation_law._lut_primary_branch.iloc[99], \ - pd.Series([1013.00, 0.53e-2, 784.890]), check_names=False, check_index=False, rtol=1e-3, atol=1e-5) + assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) - # matrix AST on page 171, chapter 3.4.2 - pd.testing.assert_frame_equal( - binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_171_seeger_beste, rtol=2e-3, atol=1e-5) + # # matrix AST on page 171, chapter 3.4.2 + # pd.testing.assert_frame_equal( + # binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_171_seeger_beste, rtol=2e-3, atol=1e-5) def test_seeger_beste_example_no_binning(): From 7eb2d95ae7dbee793f17e9dded65129c4f6b6f91 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 10:43:05 +0100 Subject: [PATCH 08/25] Drop assertions that checked for the LUTs of former Binned The LUTs are implementation details and should not be accessed directly by tests. Signed-off-by: Johannes Mueller --- .../assessment_nonlinear_standard.py | 2 - .../fkm_nonlinear/test_fkm_nonlinear.py | 52 ------------------- 2 files changed, 54 deletions(-) diff --git a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py index 9c27851e..6f9592f7 100644 --- a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py +++ b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py @@ -630,7 +630,6 @@ def _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence """ detector_1st, detector, seeger_beste_binned, recorder = _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolute_load) - result["seeger_beste_binned"] = seeger_beste_binned result, damage_calculator = _compute_damage_and_lifetimes_RAJ(assessment_parameters, recorder, component_woehler_curve_P_RAJ, result) @@ -650,7 +649,6 @@ def _compute_lifetimes_P_RAM(assessment_parameters, result, scaled_load_sequence * Do statistical assessment and store all results in a dict. """ detector_1st, detector, extended_neuber_binned, recorder = _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolute_load) - result["extended_neuber_binned"] = extended_neuber_binned result, damage_calculator = _compute_damage_and_lifetimes_RAM(assessment_parameters, recorder, component_woehler_curve_P_RAM, result) diff --git a/tests/strength/fkm_nonlinear/test_fkm_nonlinear.py b/tests/strength/fkm_nonlinear/test_fkm_nonlinear.py index 1fe70d04..9edfe911 100644 --- a/tests/strength/fkm_nonlinear/test_fkm_nonlinear.py +++ b/tests/strength/fkm_nonlinear/test_fkm_nonlinear.py @@ -511,19 +511,6 @@ def test_assessment_multiple_points_p_ram(): # assert that the results are equal for i in range(3): - # assert that the lookup table in the notch approximation law is the same - # lut for primary branch - assert np.allclose( - result_multiple["extended_neuber_binned"]._lut_primary_branch[ - result_multiple["extended_neuber_binned"]._lut_primary_branch.index.get_level_values("node_id")==i].values, - result[i]["extended_neuber_binned"]._lut_primary_branch.values) - - # lut for secondary branch - assert np.allclose( - result_multiple["extended_neuber_binned"]._lut_secondary_branch[ - result_multiple["extended_neuber_binned"]._lut_secondary_branch.index.get_level_values("node_id")==i].values, - result[i]["extended_neuber_binned"]._lut_secondary_branch.values) - assert result[i]["P_RAM_is_life_infinite"] == result_multiple["P_RAM_is_life_infinite"][i] assert np.isclose(result[i]["P_RAM_lifetime_n_cycles"], result_multiple["P_RAM_lifetime_n_cycles"][i]) assert np.isclose(result[i]["P_RAM_lifetime_n_times_load_sequence"], result_multiple["P_RAM_lifetime_n_times_load_sequence"][i]) @@ -585,19 +572,6 @@ def test_assessment_multiple_points_p_raj(): # assert that the results are equal for i in range(3): - # assert that the lookup tabe in the notch approximation law is the same - # lut for primary branch - assert np.allclose( - result_multiple["seeger_beste_binned"]._lut_primary_branch[ - result_multiple["seeger_beste_binned"]._lut_primary_branch.index.get_level_values("node_id")==i].values, - result[i]["seeger_beste_binned"]._lut_primary_branch.values) - - # lut for secondary branch - assert np.allclose( - result_multiple["seeger_beste_binned"]._lut_secondary_branch[ - result_multiple["seeger_beste_binned"]._lut_secondary_branch.index.get_level_values("node_id")==i].values, - result[i]["seeger_beste_binned"]._lut_secondary_branch.values, rtol=1e-3) - assert result[i]["P_RAJ_is_life_infinite"] == result_multiple["P_RAJ_is_life_infinite"][i] assert np.isclose(result[i]["P_RAJ_lifetime_n_cycles"], result_multiple["P_RAJ_lifetime_n_cycles"][i]) assert np.isclose(result[i]["P_RAJ_lifetime_n_times_load_sequence"], result_multiple["P_RAJ_lifetime_n_times_load_sequence"][i]) @@ -666,19 +640,6 @@ def test_assessment_multiple_points_p_ram_with_spatially_varying_gradient(): # assert that the results are equal for i in range(3): - # assert that the lookup table in the notch approximation law is the same - # lut for primary branch - assert np.allclose( - result_multiple["extended_neuber_binned"]._lut_primary_branch[ - result_multiple["extended_neuber_binned"]._lut_primary_branch.index.get_level_values("node_id")==i].values, - result[i]["extended_neuber_binned"]._lut_primary_branch.values) - - # lut for secondary branch - assert np.allclose( - result_multiple["extended_neuber_binned"]._lut_secondary_branch[ - result_multiple["extended_neuber_binned"]._lut_secondary_branch.index.get_level_values("node_id")==i].values, - result[i]["extended_neuber_binned"]._lut_secondary_branch.values) - assert result[i]["P_RAM_is_life_infinite"] == result_multiple["P_RAM_is_life_infinite"][i] assert np.isclose(result[i]["P_RAM_lifetime_n_cycles"], result_multiple["P_RAM_lifetime_n_cycles"][i]) assert np.isclose(result[i]["P_RAM_lifetime_n_times_load_sequence"], result_multiple["P_RAM_lifetime_n_times_load_sequence"][i]) @@ -746,19 +707,6 @@ def test_assessment_multiple_points_p_raj_with_spatially_varying_gradient(): # assert that the results are equal for i in range(3): - # assert that the lookup tabe in the notch approximation law is the same - # lut for primary branch - assert np.allclose( - result_multiple["seeger_beste_binned"]._lut_primary_branch[ - result_multiple["seeger_beste_binned"]._lut_primary_branch.index.get_level_values("node_id")==i].values, - result[i]["seeger_beste_binned"]._lut_primary_branch.values) - - # lut for secondary branch - assert np.allclose( - result_multiple["seeger_beste_binned"]._lut_secondary_branch[ - result_multiple["seeger_beste_binned"]._lut_secondary_branch.index.get_level_values("node_id")==i].values, - result[i]["seeger_beste_binned"]._lut_secondary_branch.values, rtol=1e-3) - assert result[i]["P_RAJ_is_life_infinite"] == result_multiple["P_RAJ_is_life_infinite"][i] assert np.isclose(result[i]["P_RAJ_lifetime_n_cycles"], result_multiple["P_RAJ_lifetime_n_cycles"][i]) assert np.isclose(result[i]["P_RAJ_lifetime_n_times_load_sequence"], result_multiple["P_RAJ_lifetime_n_times_load_sequence"][i]) From 2c8ee98731adfa42943a090c56e060172e14cc8a Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 11:15:55 +0100 Subject: [PATCH 09/25] Add tests for secondary branch notch approximation Signed-off-by: Johannes Mueller --- tests/materiallaws/data.py | 10 ++--- .../test_notch_approximation_law.py | 40 +++++++++++++++---- ...est_notch_approximation_law_seegerbeste.py | 35 ++++++++++++---- 3 files changed, 65 insertions(+), 20 deletions(-) diff --git a/tests/materiallaws/data.py b/tests/materiallaws/data.py index 9b7af37a..ffe39c47 100644 --- a/tests/materiallaws/data.py +++ b/tests/materiallaws/data.py @@ -417,8 +417,8 @@ ] ) -# matrix AST on page 162, chapter 3.4.1 -expected_matrix_AST_162_seeger_beste = pd.DataFrame( +# matrix AST on page 180, chapter 3.5.1 +expected_matrix_AST_180_seeger_beste = pd.DataFrame( index=pd.Index(np.arange(1, 200+1), name="class_index"), columns=["delta_load", "delta_strain", "delta_stress"], data=[ @@ -625,8 +625,8 @@ ] ) -# matrix AST on page 171, chapter 3.4.2 -expected_matrix_AST_171_seeger_beste = pd.DataFrame( +# matrix AST on page 189, chapter 3.5.2 +expected_matrix_AST_189_seeger_beste = pd.DataFrame( index=pd.Index(np.arange(1, 200+1), name="class_index"), columns=["delta_load", "delta_strain", "delta_stress"], data=[ @@ -831,4 +831,4 @@ [2015.87, 0.010537, 1565.324], [2026.0, 0.010604, 1569.781] ] -) \ No newline at end of file +) diff --git a/tests/materiallaws/test_notch_approximation_law.py b/tests/materiallaws/test_notch_approximation_law.py index cafffe0e..793f448b 100644 --- a/tests/materiallaws/test_notch_approximation_law.py +++ b/tests/materiallaws/test_notch_approximation_law.py @@ -29,13 +29,28 @@ (176.057, 8.71e-4, 172.639), (359.3, 0.0021, 299.78) ]) -def test_extended_neuber_example_1(load, strain, stress): +def test_extended_neuber_primary_example_1(load, strain, stress): """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ notch_approximation_law = ExtendedNeuber(E=206e3, K=1184, n=0.187, K_p=3.5) - assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) + assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-3) assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) - assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) + assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-3) + + +@pytest.mark.parametrize('delta_load, delta_strain, delta_stress', [ + (100.587, 0.0488e-2, 100.578), + (201.174, 0.0978e-2, 200.794), + (402.349, 0.2019e-2, 389.430), + (700.518, 0.4051e-2, 590.289) +]) # values from matrix AST in FKM guideline nonlinear on page 162, chapter 3.4.1 +def test_extended_neuber_secondary_example_1(delta_load, delta_strain, delta_stress): + """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ + notch_approximation_law = ExtendedNeuber(E=206e3, K=1184, n=0.187, K_p=3.5) + + assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) @pytest.mark.parametrize('load, strain, stress', [ @@ -43,7 +58,7 @@ def test_extended_neuber_example_1(load, strain, stress): (20.26, 0.0098e-2, 20.260), (1013.00, 0.6035e-2, 829.681), ]) -def test_extended_neuber_example_2(load, strain, stress): +def test_extended_neuber_primary_example_2(load, strain, stress): """ example under 2.7.2, p.78 of FKM nonlinear, "Welle mit V-Kerbe" """ notch_approximation_law = ExtendedNeuber(E=206e3, K=2650.5, n=0.187, K_p=3.5) @@ -54,10 +69,21 @@ def test_extended_neuber_example_2(load, strain, stress): assert np.isclose(notch_approximation_law.stress_secondary_branch(10.13), 10.130) assert np.isclose(notch_approximation_law.strain_secondary_branch(10.30, 10.13), 0.0049174e-2, rtol=1e-1) - # matrix AST on page 171, chapter 3.4.2 - # pd.testing.assert_frame_equal( - # binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_171, rtol=1e-3, atol=1e-5) +@pytest.mark.parametrize('delta_load, delta_strain, delta_stress', [ + (101.30, 0.0492e-2, 101.300), + (506.50, 0.2462e-2, 505.784), + (1002.87, 0.4990e-2, 978.725), + (1499.24, 0.8017e-2, 1362.884), + (2005.74, 1.1897e-2, 1649.573) +]) # values from matrix AST in FKM guideline nonlinear on page 171, chapter 3.4.2 +def test_extended_neuber_secondary_example_2(delta_load, delta_strain, delta_stress): + """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ + notch_approximation_law = ExtendedNeuber(E=206e3, K=2650.5, n=0.187, K_p=3.5) + + assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) def test_extended_neuber_example_no_binning_vectorized(): diff --git a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py index 1b065e4a..b253f725 100644 --- a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py +++ b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py @@ -32,7 +32,7 @@ (176.057, 8.71e-4, 172.639), (359.3, 0.0019, 283.86) ]) -def test_seeger_beste_example_1(load, strain, stress): +def test_seeger_beste_primary_example_1(load, strain, stress): """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ notch_approximation_law = SeegerBeste(E=206e3, K=1184, n=0.187, K_p=3.5) @@ -41,11 +41,19 @@ def test_seeger_beste_example_1(load, strain, stress): assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) - # matrix AST on page 162, chapter 3.4.1 - # pd.testing.assert_frame_equal( - # binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_162_seeger_beste, rtol=1e-3, atol=1e-5) +@pytest.mark.parametrize('delta_load, delta_strain, delta_stress', [ + (100.59, 0.000488, 100.568), + (301.76, 0.001466, 295.844), + (499.34, 0.002456, 449.183), + (700.52, 0.003605, 559.407), +]) # values from matrix AST in FKM guideline nonlinear on page 180, chapter 3.5.1 +def test_seeger_beste_secondary_example_1(delta_load, delta_strain, delta_stress): + """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ + notch_approximation_law = SeegerBeste(E=206e3, K=1184, n=0.187, K_p=3.5) - # assert np.isclose(binned_notch_approximation_law._lut_primary_branch.load.max(), maximum_absolute_load) + assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) @pytest.mark.parametrize('load, strain, stress', [ @@ -62,9 +70,20 @@ def test_seeger_beste_example_2(load, strain, stress): assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) - # # matrix AST on page 171, chapter 3.4.2 - # pd.testing.assert_frame_equal( - # binned_notch_approximation_law._lut_secondary_branch, expected_matrix_AST_171_seeger_beste, rtol=2e-3, atol=1e-5) +@pytest.mark.parametrize('delta_load, delta_strain, delta_stress', [ + (101.3, 0.000492, 101.3), + (496.37, 0.00241, 495.095), + (1002.87, 0.004879, 960.634), + (1499.24, 0.007441, 1304.533), + (2005.74, 0.010471, 1560.847), +]) # values from matrix AST in FKM guideline nonlinear on page 189, chapter 3.5.2 +def test_seeger_beste_secondary_example_2(delta_load, delta_strain, delta_stress): + """ example under 2.7.1, p.74 of FKM nonlinear "Akademisches Beispiel" """ + notch_approximation_law = SeegerBeste(E=206e3, K=2650.5, n=0.187, K_p=3.5) + + assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) def test_seeger_beste_example_no_binning(): From 737102f1ce7817288a09b9f61369d4fbaf6c5d25 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 11:22:05 +0100 Subject: [PATCH 10/25] Drop `load` argument of notch approximation strain methods Signed-off-by: Johannes Mueller --- src/pylife/materiallaws/notch_approximation_law.py | 8 ++++---- .../notch_approximation_law_seegerbeste.py | 4 ++-- tests/materiallaws/test_notch_approximation_law.py | 11 ++++------- .../test_notch_approximation_law_seegerbeste.py | 8 ++++---- 4 files changed, 14 insertions(+), 17 deletions(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index 37e67abd..be9d5627 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -82,13 +82,13 @@ def K(self, value): def primary(self, load): load = np.asarray(load) stress = self.stress(load) - strain = self.strain(stress, None) + strain = self.strain(stress) return np.stack([stress, strain], axis=len(load.shape)) def secondary(self, delta_load): delta_load = np.asarray(delta_load) delta_stress = self.stress_secondary_branch(delta_load) - delta_strain = self.strain_secondary_branch(delta_stress, None) + delta_strain = self.strain_secondary_branch(delta_stress) return np.stack([delta_stress, delta_strain], axis=len(delta_load.shape)) @@ -151,7 +151,7 @@ def stress(self, load, *, rtol=1e-4, tol=1e-4): ) return stress - def strain(self, stress, load): + def strain(self, stress): """Calculate the strain of the primary path in the stress-strain diagram at a given stress and load. The formula is given by eq. 2.5-42 of FKM nonlinear. load / stress * self._K_p * e_star @@ -240,7 +240,7 @@ def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): ) return delta_stress - def strain_secondary_branch(self, delta_stress, delta_load): + def strain_secondary_branch(self, delta_stress): """Calculate the strain on secondary branches in the stress-strain diagram at a given stress and load. The formula is given by eq. 2.5-46 of FKM nonlinear. diff --git a/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py b/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py index 8c78c18e..cc2f0985 100644 --- a/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py +++ b/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py @@ -101,7 +101,7 @@ def stress(self, load, *, rtol=1e-4, tol=1e-4): return stress[0] - def strain(self, stress, load): + def strain(self, stress): '''Calculate the strain of the primary path in the stress-strain diagram at a given stress and load. The formula is given by eq. 2.8-39 of FKM nonlinear. load / stress * self._K_p * e_star @@ -216,7 +216,7 @@ def stress_secondary_branch(self, delta_load, *, rtol=1e-4, tol=1e-4): return delta_stress[0] - def strain_secondary_branch(self, delta_stress, delta_load): + def strain_secondary_branch(self, delta_stress): '''Calculate the strain on secondary branches in the stress-strain diagram at a given stress and load. The formula is given by eq. 2.8-43 of FKM nonlinear. diff --git a/tests/materiallaws/test_notch_approximation_law.py b/tests/materiallaws/test_notch_approximation_law.py index 793f448b..4d8f3e08 100644 --- a/tests/materiallaws/test_notch_approximation_law.py +++ b/tests/materiallaws/test_notch_approximation_law.py @@ -34,7 +34,7 @@ def test_extended_neuber_primary_example_1(load, strain, stress): notch_approximation_law = ExtendedNeuber(E=206e3, K=1184, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-3) - assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress), strain, rtol=1e-1) assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-3) @@ -49,7 +49,7 @@ def test_extended_neuber_secondary_example_1(delta_load, delta_strain, delta_str notch_approximation_law = ExtendedNeuber(E=206e3, K=1184, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) - assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress), delta_strain, rtol=1e-3) assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) @@ -63,12 +63,9 @@ def test_extended_neuber_primary_example_2(load, strain, stress): notch_approximation_law = ExtendedNeuber(E=206e3, K=2650.5, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) - assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress), strain, rtol=1e-1) assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) - assert np.isclose(notch_approximation_law.stress_secondary_branch(10.13), 10.130) - assert np.isclose(notch_approximation_law.strain_secondary_branch(10.30, 10.13), 0.0049174e-2, rtol=1e-1) - @pytest.mark.parametrize('delta_load, delta_strain, delta_stress', [ (101.30, 0.0492e-2, 101.300), @@ -82,7 +79,7 @@ def test_extended_neuber_secondary_example_2(delta_load, delta_strain, delta_str notch_approximation_law = ExtendedNeuber(E=206e3, K=2650.5, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) - assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress), delta_strain, rtol=1e-3) assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) diff --git a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py index b253f725..b6b2795b 100644 --- a/tests/materiallaws/test_notch_approximation_law_seegerbeste.py +++ b/tests/materiallaws/test_notch_approximation_law_seegerbeste.py @@ -37,7 +37,7 @@ def test_seeger_beste_primary_example_1(load, strain, stress): notch_approximation_law = SeegerBeste(E=206e3, K=1184, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) - assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress), strain, rtol=1e-1) assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) @@ -52,7 +52,7 @@ def test_seeger_beste_secondary_example_1(delta_load, delta_strain, delta_stress notch_approximation_law = SeegerBeste(E=206e3, K=1184, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) - assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress), delta_strain, rtol=1e-3) assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) @@ -66,7 +66,7 @@ def test_seeger_beste_example_2(load, strain, stress): notch_approximation_law = SeegerBeste(E=206e3, K=2650.5, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress(load), stress, rtol=1e-1) - assert np.isclose(notch_approximation_law.strain(stress, None), strain, rtol=1e-1) + assert np.isclose(notch_approximation_law.strain(stress), strain, rtol=1e-1) assert np.isclose(notch_approximation_law.load(stress), load, rtol=1e-1) @@ -82,7 +82,7 @@ def test_seeger_beste_secondary_example_2(delta_load, delta_strain, delta_stress notch_approximation_law = SeegerBeste(E=206e3, K=2650.5, n=0.187, K_p=3.5) assert np.isclose(notch_approximation_law.stress_secondary_branch(delta_load), delta_stress, rtol=1e-3) - assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress, None), delta_strain, rtol=1e-3) + assert np.isclose(notch_approximation_law.strain_secondary_branch(delta_stress), delta_strain, rtol=1e-3) assert np.isclose(notch_approximation_law.load_secondary_branch(delta_stress), delta_load, rtol=1e-3) From dd6c18f20e8d0250db5fce552cf0df98038692ca Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 13:28:15 +0100 Subject: [PATCH 11/25] Let FKMNonLinearDetector initialize NotchApproxBinner Signed-off-by: Johannes Mueller --- .../assessment_nonlinear_standard.py | 4 +- src/pylife/stress/rainflow/fkm_nonlinear.py | 19 +- tests/strength/test_damage_calculator.py | 12 +- tests/stress/rainflow/test_fkm_nonlinear.py | 172 +++++------------- 4 files changed, 70 insertions(+), 137 deletions(-) diff --git a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py index 6f9592f7..a59f4c03 100644 --- a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py +++ b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py @@ -409,7 +409,7 @@ def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolu # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=extended_neuber_binned) + recorder=recorder, notch_approximation_law=extended_neuber_binned, binner=None) # perform HCM algorithm, first run detector.process_hcm_first(scaled_load_sequence) @@ -503,7 +503,7 @@ def _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolu # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=seeger_beste_binned) + recorder=recorder, notch_approximation_law=seeger_beste_binned, binner=None) detector_1st = copy.deepcopy(detector) # perform HCM algorithm, first run diff --git a/src/pylife/stress/rainflow/fkm_nonlinear.py b/src/pylife/stress/rainflow/fkm_nonlinear.py index 22deabbf..db61d1dc 100644 --- a/src/pylife/stress/rainflow/fkm_nonlinear.py +++ b/src/pylife/stress/rainflow/fkm_nonlinear.py @@ -21,6 +21,7 @@ import pandas as pd import pylife.stress.rainflow.general as RFG +import pylife.materiallaws.notch_approximation_law as NAL INDEX = 0 LOAD_TYPE = 1 @@ -86,9 +87,15 @@ class FKMNonlinearDetector(RFG.AbstractDetector): """ - def __init__(self, recorder, notch_approximation_law): + def __init__(self, recorder, notch_approximation_law, binner=NAL.NotchApproxBinner): super().__init__(recorder) - self._notch_approximation_law = notch_approximation_law + + if binner is not None: + self._binner = binner(notch_approximation_law) + self._notch_approximation_law = self._binner + else: + self._binner = None + self._notch_approximation_law = notch_approximation_law if notch_approximation_law is not None: self._ramberg_osgood_relation = self._notch_approximation_law.ramberg_osgood_relation @@ -252,6 +259,14 @@ def _determine_load_turning_points(self, samples, flush): else: rep_samples = np.asarray(samples) + if self._binner is not None: + if have_multi_index: + load_max_idx = samples.groupby("load_step").first().abs().idxmax() + load_max = samples.xs(load_max_idx, level="load_step").abs() + else: + load_max = np.abs(samples).max() + self._binner.initialize(load_max) + loads_indices, load_turning_points = self._new_turns(rep_samples, flush) self._group_size = len(samples) // len(rep_samples) diff --git a/tests/strength/test_damage_calculator.py b/tests/strength/test_damage_calculator.py index 96910856..c7c620c5 100644 --- a/tests/strength/test_damage_calculator.py +++ b/tests/strength/test_damage_calculator.py @@ -105,7 +105,8 @@ def test_woehler_curve_P_RAM_collective_has_no_index(): # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=extended_neuber_binned) + recorder=recorder, notch_approximation_law=extended_neuber_binned, binner=None + ) # perform HCM algorithm, first run detector.process_hcm_first(load_sequence_list) @@ -214,7 +215,8 @@ def test_woehler_curve_P_RAM_collective_has_MultiIndex(): # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=extended_neuber_binned) + recorder=recorder, notch_approximation_law=extended_neuber_binned, binner=None + ) # perform HCM algorithm, first run detector.process_hcm_first(load_sequence_list) @@ -310,7 +312,8 @@ def test_woehler_curve_P_RAJ_has_no_index(): # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=seeger_beste_binned) + recorder=recorder, notch_approximation_law=seeger_beste_binned, binner=None + ) # perform HCM algorithm, first run detector.process_hcm_first(load_sequence_list) @@ -418,7 +421,8 @@ def test_woehler_curve_P_RAJ_has_MultiIndex(): # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=seeger_beste_binned) + recorder=recorder, notch_approximation_law=seeger_beste_binned, binner=None + ) # perform HCM algorithm, first run detector.process_hcm_first(load_sequence_list) diff --git a/tests/stress/rainflow/test_fkm_nonlinear.py b/tests/stress/rainflow/test_fkm_nonlinear.py index aba64c5c..c8b4fcb3 100644 --- a/tests/stress/rainflow/test_fkm_nonlinear.py +++ b/tests/stress/rainflow/test_fkm_nonlinear.py @@ -56,14 +56,8 @@ def setUp(self): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber) detector.process(signal) # second run @@ -103,6 +97,31 @@ def test_interpolation(self): pd.testing.assert_frame_equal(df, expected, rtol=1e-1) +def test_no_binning_class(): + signal = np.array([100, 0, 80, 20, 60, 40]) + + E = 206e3 # [MPa] Young's modulus + K = 2650 # 1184 [MPa] + n = 0.187 # [-] + K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) + + # initialize notch approximation law + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) + + detector = FKMNonlinearDetector( + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber, binner=None + ) + detector.process_hcm_first(signal).process_hcm_second(signal) + + recorder = detector.recorder + np.testing.assert_array_equal(recorder.loads_min, np.array([40., 20., 0.])) + np.testing.assert_array_equal(recorder.loads_max, np.array([60., 80., 100.])) + np.testing.assert_allclose(recorder.S_min, np.array([39.997581, 19.997582, -0.002388]), rtol=1e-3, atol=1e-6) + np.testing.assert_allclose(recorder.S_max, np.array([59.997581, 79.997574, 99.997488]), rtol=1e-3, atol=1e-6) + np.testing.assert_allclose(recorder.epsilon_min, np.array([1.941866e-04, 9.709922e-05, 1.169416e-08]), rtol=1e-3, atol=1e-6) + np.testing.assert_allclose(recorder.epsilon_max, np.array([0.000291, 0.000388, 0.000485]), rtol=1e-3, atol=1e-6) + + class TestFKMMemory1Inner(unittest.TestCase): """Example given in FKM nonlinear 3.2.1, p.147 """ @@ -120,14 +139,8 @@ def setUp(self): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber) detector.process_hcm_first(signal) # second run @@ -198,14 +211,8 @@ def setUp(self): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber) detector.process_hcm_first(signal) # second run @@ -274,14 +281,8 @@ def setUp(self): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(self.signal)) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber) detector.process_hcm_first(self.signal) self._detector_1st = copy.deepcopy(detector) @@ -344,15 +345,8 @@ def setUp(self): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(signal)) - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=self._recorder, notch_approximation_law=extended_neuber) detector.process_hcm_first(signal) # second run @@ -423,15 +417,8 @@ def test_edge_case_value_in_sample_tail_simple_signal(vals, expected_loads_min, extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = max(abs(signal)) - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - detector = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector.process(signal).process(signal) @@ -473,15 +460,8 @@ def test_edge_case_value_in_sample_tail(vals, expected_loads_min, expected_loads extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = max(abs(signal)) - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - detector = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector.process(signal).process(signal) @@ -519,17 +499,8 @@ def test_flush_edge_case_load(): extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = pd.concat([signal_1, signal_2]).abs().groupby("node_id").max() - - print(maximum_absolute_load) - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - detector = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector.process(signal_1, flush=True).process(signal_2, flush=True) @@ -588,15 +559,8 @@ def test_flush_edge_case_load_simple_signal(): extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = max(abs(np.concatenate([signal_1, signal_2]))) - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - detector = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector.process(signal_1, flush=True).process(signal_2, flush=True) @@ -627,15 +591,8 @@ def test_flush_edge_case_S_simple_signal(): extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = max(abs(np.concatenate([signal_1, signal_2]))) - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - detector = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector.process(signal_1, flush=True).process(signal_2, flush=True) @@ -675,15 +632,8 @@ def test_flush_edge_case_S(): extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load = pd.concat([signal_1, signal_2]).abs().groupby("node_id").max() - - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - detector = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector.process(signal_1, flush=True).process(signal_2, flush=True) @@ -768,27 +718,16 @@ def test_edge_case_value_in_sample_tail_compare_simple(vals, num): extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load_simple = max(abs(vals)) - maximum_absolute_load_multiple = signal.abs().groupby('node_id').max() - print("single") - extended_neuber_binned_simple = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load_simple - ) detector_simple = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned_simple + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector_simple.process(vals).process(vals) print("multiple") - extended_neuber_binned_multiple = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load_multiple - ) detector_multiindex = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned_multiple + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector_multiindex.process(signal).process(signal) @@ -852,24 +791,13 @@ def test_hcm_first_second(vals, num): extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - maximum_absolute_load_simple = max(abs(vals)) - maximum_absolute_load_multiple = signal.abs().groupby('node_id').max() - - extended_neuber_binned_simple = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load_simple - ) - extended_neuber_binned_multiple = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load_multiple - ) detector_simple = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned_simple + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector_simple.process_hcm_first(vals).process_hcm_second(vals) detector_multiindex = FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), - notch_approximation_law=extended_neuber_binned_multiple + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=extended_neuber ) detector_multiindex.process_hcm_first(signal).process_hcm_second(signal) @@ -896,10 +824,8 @@ def detector_seeger_beste(): K_p = 3.5 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) seeger_beste = SeegerBeste(E, K, n, K_p) - seeger_beste_binned = NAL.NotchApproxBinner(seeger_beste).initialize(800) - return FKMNonlinearDetector( - recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=seeger_beste_binned + recorder=RFR.FKMNonlinearRecorder(), notch_approximation_law=seeger_beste ) @@ -1095,14 +1021,8 @@ def test_history_guideline_at_once(): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) detector.process(signal) @@ -1143,14 +1063,8 @@ def test_history_guideline_at_split(split_point): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(abs(signal)) - extended_neuber_binned = NAL.NotchApproxBinner(extended_neuber).initialize( - maximum_absolute_load - ) - # first run - detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber_binned) + detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) detector.process(signal[:split_point]).process(signal[split_point:]) From bd04e545d45768399d809552ee2fa400edba816b Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Fri, 13 Dec 2024 13:59:27 +0100 Subject: [PATCH 12/25] Use automatic binning of FKMNonLinearDetector Signed-off-by: Johannes Mueller --- .../notch_approximation_law_seegerbeste.py | 11 +++--- .../assessment_nonlinear_standard.py | 24 ++++-------- tests/strength/test_damage_calculator.py | 37 ++----------------- 3 files changed, 17 insertions(+), 55 deletions(-) diff --git a/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py b/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py index cc2f0985..6551947e 100644 --- a/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py +++ b/src/pylife/materiallaws/notch_approximation_law_seegerbeste.py @@ -77,7 +77,7 @@ def stress(self, load, *, rtol=1e-4, tol=1e-4): The resulting stress ''' # initial value as given by correction document to FKM nonlinear - x0 = load * (1 - (1 - 1/self._K_p)/1000) + x0 = np.asarray(load * (1 - (1 - 1/self._K_p)/1000)) # suppress the divergence warnings with warnings.catch_warnings(): @@ -96,7 +96,8 @@ def stress(self, load, *, rtol=1e-4, tol=1e-4): # or (value, converged, zero_der) for vector-valued invocation # only for multiple points at once, if some points diverged - if len(stress) == 3 and not stress[1].all(): + multidim = len(x0.shape) > 1 and x0.shape[1] > 1 + if multidim and not stress[1].all(): stress = self._stress_fix_not_converged_values(stress, load, x0, rtol, tol) return stress[0] @@ -461,7 +462,7 @@ def _stress_fix_not_converged_values(self, stress, load, x0, rtol, tol): '''For the values that did not converge in the previous vectorized call to optimize.newton, call optimize.newton again on the scalar value. This usually finds the correct solution.''' - indices_diverged = [index for index, is_converged in enumerate(stress[1]) if not is_converged] + indices_diverged = np.where(~stress[1].all(axis=1))[0] x0_array = np.asarray(x0) load_array = np.asarray(load) @@ -471,13 +472,13 @@ def _stress_fix_not_converged_values(self, stress, load, x0, rtol, tol): load_diverged = load_array[index_diverged] result = optimize.newton( func=self._stress_implicit, - x0=x0_diverged, + x0=np.asarray(x0_diverged), args=([load_diverged]), full_output=True, rtol=rtol, tol=tol, maxiter=50 ) - if result[1].converged: + if result.converged.all(): stress[0][index_diverged] = result[0] return stress diff --git a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py index a59f4c03..393f0638 100644 --- a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py +++ b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py @@ -389,7 +389,7 @@ def _compute_component_woehler_curves(assessment_parameters): return assessment_parameters, component_woehler_curve_P_RAM, component_woehler_curve_P_RAJ -def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolute_load): +def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolute_): """Perform the HCM rainflow counting with the extended Neuber notch approximation. The HCM algorithm is executed twice, as described in the FKM nonlinear guideline.""" @@ -397,19 +397,13 @@ def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolu E, K_prime, n_prime, K_p = assessment_parameters[["E", "K_prime", "n_prime", "K_p"]] extended_neuber = pylife.materiallaws.notch_approximation_law.ExtendedNeuber(E, K_prime, n_prime, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - extended_neuber_binned = ( - pylife.materiallaws.notch_approximation_law.NotchApproxBinner( - extended_neuber - ).initialize(maximum_absolute_load) - ) - # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=extended_neuber_binned, binner=None) + recorder=recorder, notch_approximation_law=extended_neuber + ) # perform HCM algorithm, first run detector.process_hcm_first(scaled_load_sequence) @@ -418,7 +412,7 @@ def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolu # perform HCM algorithm, second run detector.process_hcm_second(scaled_load_sequence) - return detector_1st, detector, extended_neuber_binned, recorder + return detector_1st, detector, extended_neuber, recorder def _compute_damage_and_lifetimes_RAM(assessment_parameters, recorder, component_woehler_curve_P_RAM, result): @@ -493,17 +487,13 @@ def _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolu E, K_prime, n_prime, K_p = assessment_parameters[["E", "K_prime", "n_prime", "K_p"]] seeger_beste = pylife.materiallaws.notch_approximation_law_seegerbeste.SeegerBeste(E, K_prime, n_prime, K_p) - # wrap the notch approximation law by a binning class, which precomputes the values - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.NotchApproxBinner( - seeger_beste - ).initialize(maximum_absolute_load) - # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=seeger_beste_binned, binner=None) + recorder=recorder, notch_approximation_law=seeger_beste + ) detector_1st = copy.deepcopy(detector) # perform HCM algorithm, first run @@ -512,7 +502,7 @@ def _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolu # perform HCM algorithm, second run detector.process_hcm_second(scaled_load_sequence) - return detector_1st, detector, seeger_beste_binned, recorder + return detector_1st, detector, seeger_beste, recorder def _compute_damage_and_lifetimes_RAJ(assessment_parameters, recorder, component_woehler_curve_P_RAJ, result): diff --git a/tests/strength/test_damage_calculator.py b/tests/strength/test_damage_calculator.py index c7c620c5..8fcb29e8 100644 --- a/tests/strength/test_damage_calculator.py +++ b/tests/strength/test_damage_calculator.py @@ -92,20 +92,12 @@ def test_woehler_curve_P_RAM_collective_has_no_index(): load_sequence_list = pd.Series([143.696, -287.392, 143.696, -359.240, 287.392, 0.000, 287.392, -287.392]) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(np.abs(load_sequence_list)) - extended_neuber_binned = ( - pylife.materiallaws.notch_approximation_law.NotchApproxBinner( - extended_neuber - ).initialize(maximum_absolute_load) - ) - # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=extended_neuber_binned, binner=None + recorder=recorder, notch_approximation_law=extended_neuber ) # perform HCM algorithm, first run @@ -202,20 +194,12 @@ def test_woehler_curve_P_RAM_collective_has_MultiIndex(): load_sequence_list = pd.Series([143.696, -287.392, 143.696, -359.240, 287.392, 0.000, 287.392, -287.392]) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(np.abs(load_sequence_list)) - extended_neuber_binned = ( - pylife.materiallaws.notch_approximation_law.NotchApproxBinner( - extended_neuber - ).initialize(maximum_absolute_load) - ) - # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=extended_neuber_binned, binner=None + recorder=recorder, notch_approximation_law=extended_neuber ) # perform HCM algorithm, first run @@ -300,19 +284,12 @@ def test_woehler_curve_P_RAJ_has_no_index(): seeger_beste = pylife.materiallaws.notch_approximation_law_seegerbeste.SeegerBeste(E, K_prime, n_prime, K_p) load_sequence_list = pd.Series([143.696, -287.392, 143.696, -359.240, 287.392, 0.000, 287.392, -287.392]) - - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(np.abs(load_sequence_list)) - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.NotchApproxBinner( - seeger_beste - ).initialize(maximum_absolute_load) - # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=seeger_beste_binned, binner=None + recorder=recorder, notch_approximation_law=seeger_beste ) # perform HCM algorithm, first run @@ -410,18 +387,12 @@ def test_woehler_curve_P_RAJ_has_MultiIndex(): load_sequence_list = pd.Series([143.696, -287.392, 143.696, -359.240, 287.392, 0.000, 287.392, -287.392]) - # wrap the notch approximation law by a binning class, which precomputes the values - maximum_absolute_load = max(np.abs(load_sequence_list)) - seeger_beste_binned = pylife.materiallaws.notch_approximation_law.NotchApproxBinner( - seeger_beste - ).initialize(maximum_absolute_load) - # create recorder object recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder() # create detector object detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector( - recorder=recorder, notch_approximation_law=seeger_beste_binned, binner=None + recorder=recorder, notch_approximation_law=seeger_beste ) # perform HCM algorithm, first run From fc08bc0635931f502e9b29392aaa6d595572fe88 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Mon, 16 Dec 2024 09:52:25 +0100 Subject: [PATCH 13/25] Adjust jupyter notebook Signed-off-by: Johannes Mueller --- demos/fkm_nonlinear/fkm_nonlinear_full.ipynb | 20 +++++--------------- 1 file changed, 5 insertions(+), 15 deletions(-) diff --git a/demos/fkm_nonlinear/fkm_nonlinear_full.ipynb b/demos/fkm_nonlinear/fkm_nonlinear_full.ipynb index 5f0d80c2..04144c7e 100644 --- a/demos/fkm_nonlinear/fkm_nonlinear_full.ipynb +++ b/demos/fkm_nonlinear/fkm_nonlinear_full.ipynb @@ -215,7 +215,7 @@ "plt.rcParams.update({'font.size': 22})\n", "plt.plot(range(len(scaled_load_sequence)), load_sequence, \"o-\", c=\"g\", lw=3)\n", "plt.grid()\n", - "plt.ylabel(\"$\\sigma_V$ [MPa]\")" + "plt.ylabel(r\"$\\sigma_V$ [MPa]\")" ] }, { @@ -559,18 +559,13 @@ "load_sequence_list = scaled_load_sequence\n", "print(load_sequence_list)\n", "\n", - "# wrap the notch approximation law by a binning class, which precomputes the values\n", - "maximum_absolute_load = max(np.abs(load_sequence_list))\n", - "print(f\"maximum_absolute_load: {maximum_absolute_load}\")\n", - "extended_neuber_binned = pylife.materiallaws.notch_approximation_law.Binned(\n", - " extended_neuber, maximum_absolute_load, 100)\n", - "\n", "# create recorder object\n", "recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder()\n", "\n", "# create detector object\n", "detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector(\n", - " recorder=recorder, notch_approximation_law=extended_neuber_binned)\n", + " recorder=recorder, notch_approximation_law=extended_neuber\n", + ")\n", "\n", "# perform HCM algorithm, first run\n", "detector.process(load_sequence_list, flush=False)\n", @@ -759,18 +754,13 @@ "load_sequence_list = scaled_load_sequence\n", "print(load_sequence_list)\n", "\n", - "# wrap the notch approximation law by a binning class, which precomputes the values\n", - "maximum_absolute_load = max(np.abs(load_sequence_list))\n", - "print(f\"maximum_absolute_load: {maximum_absolute_load}\")\n", - "seeger_beste_binned = pylife.materiallaws.notch_approximation_law.Binned(\n", - " seeger_beste, maximum_absolute_load, 100)\n", - "\n", "# create recorder object\n", "recorder = pylife.stress.rainflow.recorders.FKMNonlinearRecorder()\n", "\n", "# create detector object\n", "detector = pylife.stress.rainflow.fkm_nonlinear.FKMNonlinearDetector(\n", - " recorder=recorder, notch_approximation_law=seeger_beste_binned)\n", + " recorder=recorder, notch_approximation_law=seeger_beste\n", + ")\n", "\n", "# perform HCM algorithm, first run\n", "detector.process(load_sequence_list, flush=False)\n", From 8983ff31864597858ea87d88a8e984645567460b Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Mon, 16 Dec 2024 10:14:28 +0100 Subject: [PATCH 14/25] Don't regenerate the test signal while collecting benchmark tests Signed-off-by: Johannes Mueller --- benchmarks/generate_time_signal.py | 34 ++++++++++++++++-------------- 1 file changed, 18 insertions(+), 16 deletions(-) diff --git a/benchmarks/generate_time_signal.py b/benchmarks/generate_time_signal.py index 1986ce28..05777e56 100644 --- a/benchmarks/generate_time_signal.py +++ b/benchmarks/generate_time_signal.py @@ -2,21 +2,23 @@ import numpy as np import pylife.stress.timesignal as TS -np.random.seed(23424711) -load_signal = TS.TimeSignalGenerator( - 10, - { - 'number': 50_000, - 'amplitude_median': 50.0, - 'amplitude_std_dev': 0.5, - 'frequency_median': 4, - 'frequency_std_dev': 3, - 'offset_median': 0, - 'offset_std_dev': 0.4, - }, - None, - None, -).query(1_000_000) +if __name__ == "__main__": + np.random.seed(23424711) -np.savetxt('load_signal.csv', load_signal) + load_signal = TS.TimeSignalGenerator( + 10, + { + 'number': 50_000, + 'amplitude_median': 50.0, + 'amplitude_std_dev': 0.5, + 'frequency_median': 4, + 'frequency_std_dev': 3, + 'offset_median': 0, + 'offset_std_dev': 0.4, + }, + None, + None, + ).query(1_000_000) + + np.savetxt('load_signal.csv', load_signal) From 8bd876382941836546ce4a2d75c671cf784595e0 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Mon, 16 Dec 2024 12:55:28 +0100 Subject: [PATCH 15/25] Add docstrings for .primary() and .secondary() Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 38 +++++++++++++++++++ 1 file changed, 38 insertions(+) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index be9d5627..e1ba8295 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -80,12 +80,50 @@ def K(self, value): self.K_prime = value def primary(self, load): + """Calculate stress and strain for primary branch. + + Parameters + ---------- + load : array-like + The load for which the stress and strain are to be calculated + + Returns + ------- + stress strain : ndarray + The resulting stress strain data. + + If the argument is scalar, the resulting array is of the strucuture + ``[<σ>, <ε>]`` + + If the argument is an 1D-array with length `n`the resulting array is of the + structure ``[[<σ1>, <σ2>, <σ3>, ... <σn>], [<ε1>, <ε2>, <ε3>, ... <εn>]]`` + + """ load = np.asarray(load) stress = self.stress(load) strain = self.strain(stress) return np.stack([stress, strain], axis=len(load.shape)) def secondary(self, delta_load): + """Calculate stress and strain for secondary branch. + + Parameters + ---------- + load : array-like + The load for which the stress and strain are to be calculated + + Returns + ------- + stress strain : ndarray + The resulting stress strain data. + + If the argument is scalar, the resulting array is of the strucuture + ``[<σ>, <ε>]`` + + If the argument is an 1D-array with length `n`the resulting array is of the + structure ``[[<σ1>, <σ2>, <σ3>, ... <σn>], [<ε1>, <ε2>, <ε3>, ... <εn>]]`` + + """ delta_load = np.asarray(delta_load) delta_stress = self.stress_secondary_branch(delta_load) delta_strain = self.strain_secondary_branch(delta_stress) From bd97d9646e392edad2dcf6f5299dd1a64c9e8083 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Tue, 17 Dec 2024 09:27:35 +0100 Subject: [PATCH 16/25] Add abstract base class methods to NotchApproximationLawBase Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 28 ++++++++++++++++++- 1 file changed, 27 insertions(+), 1 deletion(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index e1ba8295..d067cd23 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -17,13 +17,15 @@ __author__ = ["Benjamin Maier"] __maintainer__ = __author__ +from abc import ABC, abstractmethod + import numpy as np from scipy import optimize import pandas as pd import pylife.materiallaws.rambgood -class NotchApproximationLawBase: +class NotchApproximationLawBase(ABC): """This is a base class for any notch approximation law, e.g., the extended Neuber and the Seeger-Beste laws. It initializes the internal variables used by the derived classes and provides getters and setters. @@ -79,6 +81,30 @@ def K(self, value): """Set the strain hardening coefficient""" self.K_prime = value + @abstractmethod + def load(self, load, *, rtol=1e-4, tol=1e-4): + ... + + @abstractmethod + def stress(self, load, *, rtol=1e-4, tol=1e-4): + ... + + @abstractmethod + def strain(self, load): + ... + + @abstractmethod + def load_secondary_branch(self, load, *, rtol=1e-4, tol=1e-4): + ... + + @abstractmethod + def stress_secondary_branch(self, load, *, rtol=1e-4, tol=1e-4): + ... + + @abstractmethod + def strain_secondary_branch(self, load): + ... + def primary(self, load): """Calculate stress and strain for primary branch. From fff5b23bf6e1044b926857ba6864f5f6d0a8e74f Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Tue, 17 Dec 2024 09:30:39 +0100 Subject: [PATCH 17/25] Drop obsolete file Signed-off-by: Johannes Mueller --- .../notch_approximation_law_base.py | 77 ------------------- 1 file changed, 77 deletions(-) delete mode 100644 src/pylife/materiallaws/notch_approximation_law_base.py diff --git a/src/pylife/materiallaws/notch_approximation_law_base.py b/src/pylife/materiallaws/notch_approximation_law_base.py deleted file mode 100644 index c7f35054..00000000 --- a/src/pylife/materiallaws/notch_approximation_law_base.py +++ /dev/null @@ -1,77 +0,0 @@ -# Copyright (c) 2019-2023 - for information on the respective copyright owner -# see the NOTICE file and/or the repository -# https://github.com/boschresearch/pylife -# -# Licensed under the Apache License, Version 2.0 (the "License"); -# you may not use this file except in compliance with the License. -# You may obtain a copy of the License at -# -# http://www.apache.org/licenses/LICENSE-2.0 -# -# Unless required by applicable law or agreed to in writing, software -# distributed under the License is distributed on an "AS IS" BASIS, -# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. -# See the License for the specific language governing permissions and -# limitations under the License. - -__author__ = ["Benjamin Maier"] -__maintainer__ = __author__ - -import pylife.materiallaws.rambgood - -class NotchApproximationLawBase: - """This is a base class for any notch approximation law, e.g., the extended Neuber and the Seeger-Beste laws. - - It initializes the internal variables used by the derived classes and provides getters and setters. - """ - - def __init__(self, E, K, n, K_p=None): - self._E = E - self._K = K - self._n = n - self._K_p = K_p - - self._ramberg_osgood_relation = pylife.materiallaws.rambgood.RambergOsgood(E, K, n) - - @property - def E(self): - '''Get Young's Modulus''' - return self._E - - @property - def K(self): - '''Get the strength coefficient''' - return self._K - - @property - def n(self): - '''Get the strain hardening coefficient''' - return self._n - - @property - def K_p(self): - '''Get the shape factor (de: Traglastformzahl)''' - return self._K_p - - @property - def ramberg_osgood_relation(self): - '''Get the ramberg osgood relation object - ''' - return self._ramberg_osgood_relation - - @K_p.setter - def K_p(self, value): - """Set the shape factor value K_p (de: Traglastformzahl)""" - self._K_p = value - - @K.setter - def K_prime(self, value): - """Set the strain hardening coefficient""" - self._K = value - self._ramberg_osgood_relation = pylife.materiallaws.rambgood.RambergOsgood(self._E, self._K, self._n) - - @K.setter - def K(self, value): - """Set the strain hardening coefficient""" - self.K_prime = value - From a7105fdf57820b857edfc795ca1249ba9b5bd73b Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Mon, 13 Jan 2025 13:14:30 +0100 Subject: [PATCH 18/25] Simplifications Signed-off-by: Johannes Mueller --- src/pylife/stress/rainflow/fkm_nonlinear.py | 32 ++-- src/pylife/stress/rainflow/recorders.py | 176 +++++++------------- tests/stress/rainflow/test_recorders.py | 39 ++++- 3 files changed, 104 insertions(+), 143 deletions(-) diff --git a/src/pylife/stress/rainflow/fkm_nonlinear.py b/src/pylife/stress/rainflow/fkm_nonlinear.py index db61d1dc..8bbe686a 100644 --- a/src/pylife/stress/rainflow/fkm_nonlinear.py +++ b/src/pylife/stress/rainflow/fkm_nonlinear.py @@ -223,7 +223,7 @@ def process(self, samples, flush=False): self._store_recordings_for_history(record, record_vals, turning_point_idx, hysts) results = self._process_recording(load_turning_points_rep, record_vals, hysts) - results_min, results_max, epsilon_min_LF, epsilon_max_LF = results + results_min, results_max = results self._update_residuals(record_vals, turning_point_idx, load_turning_points_rep) @@ -234,14 +234,7 @@ def process(self, samples, flush=False): is_zero_mean_stress_and_strain = (hysts[:, 0] == MEMORY_3).tolist() self._recorder.record_values_fkm_nonlinear( - loads_min=results_min["loads_min"], - loads_max=results_max["loads_max"], - S_min=results_min["S_min"], - S_max=results_max["S_max"], - epsilon_min=results_min["epsilon_min"], - epsilon_max=results_max["epsilon_max"], - epsilon_min_LF=epsilon_min_LF, - epsilon_max_LF=epsilon_max_LF, + results_min, results_max, is_closed_hysteresis=is_closed_hysteresis, is_zero_mean_stress_and_strain=is_zero_mean_stress_and_strain, run_index=self._run_index @@ -433,15 +426,12 @@ def turn_memory_3(values, index): result_len = len(hysts) * self._group_size - results_min = np.zeros((result_len, 3)) + results_min = np.zeros((result_len, 4)) results_min_idx = np.zeros((result_len, signal_index_num), dtype=np.int64) - results_max = np.zeros((result_len, 3)) + results_max = np.zeros((result_len, 4)) results_max_idx = np.zeros((result_len, signal_index_num), dtype=np.int64) - epsilon_min_LF = np.zeros(result_len) - epsilon_max_LF = np.zeros(result_len) - for i, hyst in enumerate(hysts): idx = (hyst[FROM:CLOSE] + start) * self._group_size @@ -457,27 +447,27 @@ def turn_memory_3(values, index): beg = i * self._group_size end = beg + self._group_size - results_min[beg:end] = min_val[:, :3] - results_max[beg:end] = max_val[:, :3] + results_min[beg:end, :3] = min_val[:, :3] + results_max[beg:end, :3] = max_val[:, :3] results_min_idx[beg:end] = min_idx results_max_idx[beg:end] = max_idx - epsilon_min_LF[beg:end] = min_val[:, EPS_MIN_LF] - epsilon_max_LF[beg:end] = max_val[:, EPS_MAX_LF] + results_min[beg:end, -1] = min_val[:, EPS_MIN_LF] + results_max[beg:end, -1] = max_val[:, EPS_MAX_LF] results_min = pd.DataFrame( results_min, - columns=["loads_min", "S_min", "epsilon_min"], + columns=["loads_min", "S_min", "epsilon_min", "epsilon_min_LF"], index=pd.MultiIndex.from_arrays(results_min_idx.T, names=signal_index_names) ) results_max = pd.DataFrame( results_max, - columns=["loads_max", "S_max", "epsilon_max"], + columns=["loads_max", "S_max", "epsilon_max", "epsilon_max_LF"], index=pd.MultiIndex.from_arrays(results_max_idx.T, names=signal_index_names) ) - return results_min, results_max, pd.Series(epsilon_min_LF), pd.Series(epsilon_max_LF) + return results_min, results_max def _adjust_samples_and_flush_for_hcm_first_run(self, samples): diff --git a/src/pylife/stress/rainflow/recorders.py b/src/pylife/stress/rainflow/recorders.py index dac66416..e459ba8f 100644 --- a/src/pylife/stress/rainflow/recorders.py +++ b/src/pylife/stress/rainflow/recorders.py @@ -151,14 +151,14 @@ class FKMNonlinearRecorder(AbstractRecorder): def __init__(self): """Instantiate a FKMNonlinearRecorder.""" super().__init__() - self._loads_min = pd.Series(dtype=np.float64) - self._loads_max = pd.Series(dtype=np.float64) - self._S_min = pd.Series(dtype=np.float64) - self._S_max = pd.Series(dtype=np.float64) - self._epsilon_min = pd.Series(dtype=np.float64) - self._epsilon_max = pd.Series(dtype=np.float64) - self._epsilon_min_LF = pd.Series(dtype=np.float64) - self._epsilon_max_LF = pd.Series(dtype=np.float64) + self._results_min = pd.DataFrame( + columns=["loads_min", "S_min", "epsilon_min", "epsilon_min_LF"], + dtype=np.float64, + ) + self._results_max = pd.DataFrame( + columns=["loads_max", "S_max", "epsilon_max", "epsilon_max_LF"], + dtype=np.float64, + ) self._is_closed_hysteresis = [] self._is_zero_mean_stress_and_strain = [] self._run_index = [] @@ -166,38 +166,38 @@ def __init__(self): @property def loads_min(self): """1-D float array containing the start load values of the recorded hystereses.""" - return self._loads_min + return self._results_min["loads_min"] @property def loads_max(self): """1-D float array containing the end load values of the recorded hystereses, i.e., the values the loops go to before turning back.""" - return self._loads_max + return self._results_max["loads_max"] @property def S_min(self): """1-D float array containing the minimum stress values of the recorded hystereses.""" - return self._S_min + return self._results_min["S_min"] @property def S_max(self): """1-D float array containing the maximum stress values of the recorded hystereses.""" - return self._S_max + return self._results_max["S_max"] @property def epsilon_min(self): """1-D float array containing the minimum strain values of the recorded hystereses.""" - return self._epsilon_min + return self._results_min["epsilon_min"] @property def epsilon_max(self): """1-D float array containing the maximum strain values of the recorded hystereses.""" - return self._epsilon_max + return self._results_max["epsilon_max"] @property def S_a(self): """1-D numpy array containing the stress amplitudes of the recorded hystereses.""" - return 0.5 * (np.array(self._S_max) - np.array(self._S_min)) + return 0.5 * (np.array(self.S_max) - np.array(self.S_min)) @property def S_m(self): @@ -206,13 +206,13 @@ def S_m(self): Only for hystereses resulting from Memory 3, the FKM nonlinear document defines ``S_m`` to be zero (eq. 2.9-52). This is indicated by ``_is_zero_mean_stress_and_strain=True`̀ . For these hystereses, this function returns 0 instead of ``(S_min + S_max) / 2``. """ - median = 0.5 * (np.array(self._S_min) + np.array(self._S_max)) + median = 0.5 * (np.array(self.S_min) + np.array(self.S_max)) return np.where(self.is_zero_mean_stress_and_strain, 0, median) @property def epsilon_a(self): """1-D float array containing the strain amplitudes of the recorded hystereses.""" - return 0.5 * (np.array(self._epsilon_max) - np.array(self._epsilon_min)) + return 0.5 * (np.array(self.epsilon_max) - np.array(self.epsilon_min)) @property def epsilon_m(self): @@ -222,25 +222,13 @@ def epsilon_m(self): to be zero (eq. 2.9-53). This is indicated by ``_is_zero_mean_stress_and_strain=True`̀ . For these hystereses, this function returns 0 instead of ``(epsilon_min + epsilon_max) / 2``. """ return np.where(self.is_zero_mean_stress_and_strain, \ - 0, 0.5 * (np.array(self._epsilon_min) + np.array(self._epsilon_max))) - - def _get_for_every_node(self, boolean_array): - - # number of points, i.e., number of values for every load step - m = len(self._S_min[0]) - - # bring the array of boolean values to the right shape - # numeric_array contains only 0s and 1s for False and True - numeric_array = np.array(boolean_array).reshape(-1,1).dot(np.ones((1,m))) - - # transform the array to boolean type - return np.where(numeric_array == 1, True, False) + 0, 0.5 * (np.array(self.epsilon_min) + np.array(self.epsilon_max))) @property def is_zero_mean_stress_and_strain(self): # if the assessment is performed for multiple points at once - if len(self._S_min) > 0 and len(self._S_min.index.names) > 1: + if len(self.S_min) > 0 and len(self.S_min.index.names) > 1: return self._get_for_every_node(self._is_zero_mean_stress_and_strain) else: return self._is_zero_mean_stress_and_strain @@ -253,7 +241,7 @@ def R(self): (eq. 2.9-54). This is indicated by ``_is_zero_mean_stress_and_strain=True`̀ . For these hystereses, this function returns -1 instead of ``S_min / S_max``, which may be different. """ with np.errstate(all="ignore"): - R = np.array(self._S_min) / np.array(self._S_max) + R = np.array(self.S_min) / np.array(self.S_max) return np.where(self.is_zero_mean_stress_and_strain, -1, R) @property @@ -262,7 +250,7 @@ def is_closed_hysteresis(self): was recorded as a memory 3 hysteresis, which counts only half the damage in the FKM nonlinear procedure.""" # if the assessment is performed for multiple points at once - if len(self._S_min) > 0 and len(self._S_min.index.names) > 1: + if len(self.S_min) > 0 and len(self.S_min.index.names) > 1: return self._get_for_every_node(self._is_closed_hysteresis) else: return self._is_closed_hysteresis @@ -291,99 +279,61 @@ def collective(self): node_id starts, e.g., with 1. """ - # if the assessment is performed for multiple points at once - if len(self._S_min) > 0 and len(self._S_min.index.names) > 1: - R = self.R - n_nodes = self._S_min.groupby('node_id').first().count() - n_hystereses = int(len(self._S_min) / n_nodes) - - index = pd.MultiIndex.from_product([range(n_hystereses), range(n_nodes)], names=["hysteresis_index", "assessment_point_index"]) - - return pd.DataFrame( - index=index, - data={ - "loads_min": self._loads_min.to_numpy(), - "loads_max": self._loads_max.to_numpy(), - "S_min": self._S_min.to_numpy(), - "S_max": self._S_max.to_numpy(), - "R": self.R, - "epsilon_min": self._epsilon_min.to_numpy(), - "epsilon_max": self._epsilon_max.to_numpy(), - "S_a": self.S_a, - "S_m": self.S_m, - "epsilon_a": self.epsilon_a, - "epsilon_m": self.epsilon_m, - "epsilon_min_LF": self._epsilon_min_LF.to_numpy(), - "epsilon_max_LF": self._epsilon_max_LF.to_numpy(), - "is_closed_hysteresis": self.is_closed_hysteresis, - "is_zero_mean_stress_and_strain": self.is_zero_mean_stress_and_strain, - "run_index": np.array(self._run_index, dtype=np.int64) - }) + if len(self.S_min) > 0 and len(self.S_min.index.names) > 1: + n_nodes = self.S_min.groupby('node_id').first().count() + n_hystereses = int(len(self.S_min) / n_nodes) + index = pd.MultiIndex.from_product( + [range(n_hystereses), range(n_nodes)], + names=["hysteresis_index", "assessment_point_index"], + ) else: - # if the assessment is performed by a single point - n_hystereses = len(self._S_min) + n_hystereses = len(self.S_min) index = pd.MultiIndex.from_product([range(n_hystereses), [0]], names=["hysteresis_index", "assessment_point_index"]) - # determine load - loads_min = self._loads_min - loads_max = self._loads_max - - if len(self._loads_min) == 0 or len(self._loads_max) == 0: - loads_min = pd.Series(np.nan, index=index) - loads_max = pd.Series(np.nan, index=index) - - return pd.DataFrame( - index=index, - data={ - "loads_min": loads_min.values, - "loads_max": loads_max.values, - "S_min": self._S_min.values, - "S_max": self._S_max.values, - "R": self.R, # FIXME .values, - "epsilon_min": self._epsilon_min.values, - "epsilon_max": self._epsilon_max.values, - "S_a": self.S_a, # FIXME .values, - "S_m": self.S_m, # FIXME .values, - "epsilon_a": self.epsilon_a, # FIXME .values, - "epsilon_m": self.epsilon_m, # FIXME .values, - "epsilon_min_LF": self._epsilon_min_LF.values, - "epsilon_max_LF": self._epsilon_max_LF.values, - "is_closed_hysteresis": self._is_closed_hysteresis, # FIXME .values - "is_zero_mean_stress_and_strain": self._is_zero_mean_stress_and_strain, # FIXME .values, - "run_index": np.array(self._run_index, dtype=np.int64), - }, - ) - - def record_values_fkm_nonlinear(self, loads_min, loads_max, S_min, S_max, epsilon_min, epsilon_max, - epsilon_min_LF, epsilon_max_LF, - is_closed_hysteresis, is_zero_mean_stress_and_strain, run_index): + return pd.DataFrame( + index=index, + data={ + "loads_min": self.loads_min.to_numpy(), + "loads_max": self.loads_max.to_numpy(), + "S_min": self.S_min.to_numpy(), + "S_max": self.S_max.to_numpy(), + "R": self.R, + "epsilon_min": self.epsilon_min.to_numpy(), + "epsilon_max": self.epsilon_max.to_numpy(), + "S_a": self.S_a, + "S_m": self.S_m, + "epsilon_a": self.epsilon_a, + "epsilon_m": self.epsilon_m, + "epsilon_min_LF": self._results_min["epsilon_min_LF"].to_numpy(), + "epsilon_max_LF": self._results_max["epsilon_max_LF"].to_numpy(), + "is_closed_hysteresis": self.is_closed_hysteresis, + "is_zero_mean_stress_and_strain": self.is_zero_mean_stress_and_strain, + "run_index": np.array(self._run_index, dtype=np.int64) + }) + + def record_values_fkm_nonlinear( + self, + results_min, + results_max, + is_closed_hysteresis, + is_zero_mean_stress_and_strain, + run_index, + ): """Record the loop values.""" - if len(loads_min) > 0: - self._loads_min = loads_min if len(self._loads_min) == 0 else pd.concat([self._loads_min, loads_min]) - elif len(self._loads_min) == 0: - self._loads_min.index = loads_min.index - if len(loads_max) > 0: - self._loads_max = loads_max if len(self._loads_max) == 0 else pd.concat([self._loads_max, loads_max]) - elif len(self._loads_max) == 0: - self._loads_max.index = loads_max.index - self._S_min = S_min if len(self._S_min) == 0 else pd.concat([self._S_min, S_min]) - self._S_max = S_max if len(self._S_max) == 0 else pd.concat([self._S_max, S_max]) - self._epsilon_min = pd.concat([self._epsilon_min, epsilon_min]) - self._epsilon_max = pd.concat([self._epsilon_max, epsilon_max]) - self._epsilon_min_LF = pd.concat([self._epsilon_min_LF, epsilon_min_LF]) - self._epsilon_max_LF = pd.concat([self._epsilon_max_LF, epsilon_max_LF]) + self._results_min = results_min if len(self._results_min) == 0 else pd.concat([self._results_min, results_min]) + self._results_max = results_max if len(self._results_max) == 0 else pd.concat([self._results_max, results_max]) self._is_closed_hysteresis += is_closed_hysteresis self._is_zero_mean_stress_and_strain += is_zero_mean_stress_and_strain - self._run_index += [run_index] * len(S_min) + self._run_index += [run_index] * len(results_min) def _get_for_every_node(self, boolean_array): # number of points, i.e., number of values for every load step - li = self._S_min.index.to_frame()['load_step'] - m = self._S_min.groupby((li!=li.shift()).cumsum(), sort=False).count().iloc[0] + li = self.S_min.index.to_frame()['load_step'] + m = self.S_min.groupby((li!=li.shift()).cumsum(), sort=False).count().iloc[0] # bring the array of boolean values to the right shape # numeric_array contains only 0s and 1s for False and True numeric_array = np.array(boolean_array).reshape(-1,1).dot(np.ones((1,m))).flatten() diff --git a/tests/stress/rainflow/test_recorders.py b/tests/stress/rainflow/test_recorders.py index 1a07044d..5c4dc2a1 100644 --- a/tests/stress/rainflow/test_recorders.py +++ b/tests/stress/rainflow/test_recorders.py @@ -262,12 +262,21 @@ def test_fkm_nonlinear_recorder_record_two_values(): h1, h2 = 4., 5. fr = RFR.FKMNonlinearRecorder() - # arguments: loads_min, loads_max, S_min, S_max, epsilon_min, epsilon_max, epsilon_min_LF, epsilon_max_LF, - # is_closed_hysteresis, is_zero_mean_stress_and_strain, run_index - - args_1 = [pd.Series([v]) for v in [a1, b1, c1, d1, e1, f1, g1, h1]] + [[False], [False], 1] - args_2 = [pd.Series([v]) for v in [a2, b2, c2, d2, e2, f2, g2, h2]] + [[True], [False], 2] - args_3 = [pd.Series([v]) for v in [a2, b2, c2, d2, e2, f2, g2, h2]] + [[True], [True], 2] + results_min_1 = pd.DataFrame( + {"loads_min": [a1], "S_min": [c1], "epsilon_min": [e1], "epsilon_min_LF": [g1]}, + ) + results_max_1 = pd.DataFrame( + {"loads_max": [b1], "S_max": [d1], "epsilon_max": [f1], "epsilon_max_LF": [h1]}, + ) + results_min_2 = pd.DataFrame( + {"loads_min": [a2], "S_min": [c2], "epsilon_min": [e2], "epsilon_min_LF": [g2]}, + ) + results_max_2 = pd.DataFrame( + {"loads_max": [b2], "S_max": [d2], "epsilon_max": [f2], "epsilon_max_LF": [h2]}, + ) + args_1 = [results_min_1, results_max_1] + [[False], [False], 1] + args_2 = [results_min_2, results_max_2] + [[True], [False], 2] + args_3 = [results_min_2, results_max_2] + [[True], [True], 2] fr.record_values_fkm_nonlinear(*args_1) fr.record_values_fkm_nonlinear(*args_2) @@ -329,9 +338,21 @@ def test_fkm_nonlinear_recorder_two_non_zero_collective(): h1, h2 = 4., 5. fr = RFR.FKMNonlinearRecorder() - args_1 = [pd.Series([v]) for v in [a1, b1, c1, d1, e1, f1, g1, h1]] + [[False], [False], 1] - args_2 = [pd.Series([v]) for v in [a2, b2, c2, d2, e2, f2, g2, h2]] + [[True], [False], 2] - args_3 = [pd.Series([v]) for v in [a2, b2, c2, d2, e2, f2, g2, h2]] + [[True], [True], 2] + results_min_1 = pd.DataFrame( + {"loads_min": [a1], "S_min": [c1], "epsilon_min": [e1], "epsilon_min_LF": [g1]}, + ) + results_max_1 = pd.DataFrame( + {"loads_max": [b1], "S_max": [d1], "epsilon_max": [f1], "epsilon_max_LF": [h1]}, + ) + results_min_2 = pd.DataFrame( + {"loads_min": [a2], "S_min": [c2], "epsilon_min": [e2], "epsilon_min_LF": [g2]}, + ) + results_max_2 = pd.DataFrame( + {"loads_max": [b2], "S_max": [d2], "epsilon_max": [f2], "epsilon_max_LF": [h2]}, + ) + args_1 = [results_min_1, results_max_1] + [[False], [False], 1] + args_2 = [results_min_2, results_max_2] + [[True], [False], 2] + args_3 = [results_min_2, results_max_2] + [[True], [True], 2] fr.record_values_fkm_nonlinear(*args_1) fr.record_values_fkm_nonlinear(*args_2) From 02dadf30589d629c813c7b8eb4086e42e38786e0 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 16 Jan 2025 20:44:24 +0100 Subject: [PATCH 19/25] Handle an edge case in recording epslilon LF on open hysteresis The epsilon_min_LF and epsilon_max_LF values where wrongly picked for open hysteresis loops. Now the epsilon_(min|max)_LF calculation is identical to the state of bf10295a. Signed-off-by: Johannes Mueller --- src/pylife/stress/rainflow/fkm_nonlinear.py | 14 +++- src/pylife/stress/rainflow/recorders.py | 8 ++ tests/stress/rainflow/test_fkm_nonlinear.py | 93 +++++++++++++++++++++ 3 files changed, 112 insertions(+), 3 deletions(-) diff --git a/src/pylife/stress/rainflow/fkm_nonlinear.py b/src/pylife/stress/rainflow/fkm_nonlinear.py index 8bbe686a..d14a35f3 100644 --- a/src/pylife/stress/rainflow/fkm_nonlinear.py +++ b/src/pylife/stress/rainflow/fkm_nonlinear.py @@ -405,7 +405,11 @@ def turn_memory_1_2(values, index): def turn_memory_3(values, index): abs_point = np.abs(values[0]) - return (-abs_point, abs_point, index[0], index[0]) + point_min = -abs_point + point_max = abs_point + point_min[:, EPS_MIN_LF:] = values[1][:, EPS_MIN_LF:] + point_max[:, EPS_MIN_LF:] = values[1][:, EPS_MIN_LF:] + return (point_min, point_max, index[0], index[0]) memory_functions = [turn_memory_3, turn_memory_1_2] @@ -435,13 +439,17 @@ def turn_memory_3(values, index): for i, hyst in enumerate(hysts): idx = (hyst[FROM:CLOSE] + start) * self._group_size + hyst_type = hyst[IS_CLOSED] + beg0, beg1 = idx[0], idx[1] + vbeg1 = beg1 - self._group_size if hyst_type == MEMORY_3 else beg1 + end0, end1 = beg0 + self._group_size, beg1 + self._group_size + vend1 = vbeg1 + self._group_size - values = value_array[beg0:end0], value_array[beg1:end1] + values = value_array[beg0:end0], value_array[vbeg1:vend1] index = index_array[beg0:end0], index_array[beg1:end1] - hyst_type = hyst[IS_CLOSED] min_val, max_val, min_idx, max_idx = memory_functions[hyst_type](values, index) beg = i * self._group_size diff --git a/src/pylife/stress/rainflow/recorders.py b/src/pylife/stress/rainflow/recorders.py index e459ba8f..d77f0270 100644 --- a/src/pylife/stress/rainflow/recorders.py +++ b/src/pylife/stress/rainflow/recorders.py @@ -194,6 +194,14 @@ def epsilon_max(self): """1-D float array containing the maximum strain values of the recorded hystereses.""" return self._results_max["epsilon_max"] + @property + def epsilon_min_LF(self): + return self._results_min["epsilon_min_LF"] + + @property + def epsilon_max_LF(self): + return self._results_max["epsilon_max_LF"] + @property def S_a(self): """1-D numpy array containing the stress amplitudes of the recorded hystereses.""" diff --git a/tests/stress/rainflow/test_fkm_nonlinear.py b/tests/stress/rainflow/test_fkm_nonlinear.py index c8b4fcb3..3a60578d 100644 --- a/tests/stress/rainflow/test_fkm_nonlinear.py +++ b/tests/stress/rainflow/test_fkm_nonlinear.py @@ -181,6 +181,7 @@ def test_epsilon_LF(self): rtol=1e-2, ) + def test_interpolation(self): df = self._detector.interpolated_stress_strain_data(load_segment=5, n_points_per_branch=5) expected = pd.DataFrame( @@ -1087,3 +1088,95 @@ def test_history_guideline_at_split(split_point): ).set_index(["load_segment", "load_step", "run_index", "turning_point", "hyst_from", "hyst_to", "hyst_close"]) pd.testing.assert_frame_equal(df, expected, rtol=1e-1) + + +def test_epsilon_max_LF_unclosed_hysteresis(): + # Fig 2.3 FKM NL Guideline + + signal = pd.Series([0.0, -100.0, 100.0, 0.0, 200.0, 100.0]) + signal.index.name = "load_step" + + recorder = RFR.FKMNonlinearRecorder() + E = 1.0 # [MPa] Young's modulus + K = 1e6 # 1184 [MPa] + n = 1.0 # [-] + K_p = 1.0 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) + + # initialize notch approximation law + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) + + detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) + + detector.process(signal) + + print(detector.recorder.loads_min) + print(detector.recorder.loads_max) + print(detector.recorder.epsilon_min_LF) + print(detector.recorder.epsilon_max_LF) + + print(detector.recorder.collective) + print(detector.history()) + + print([np.diff(recorder.epsilon_min_LF) <= 0.0]) + assert np.all([np.diff(recorder.epsilon_min_LF) <= 0.0]) + + print([np.diff(recorder.epsilon_max_LF) >= 0.0]) + assert np.all([np.diff(recorder.epsilon_max_LF) >= 0.0]) + + +def test_epsilon_min_LF_flip(): + signal = pd.Series([10.0, -1.0, 260.0, -250.0, 60.0, -280.0, -100.0]) + signal.index.name = "load_step" + + recorder = RFR.FKMNonlinearRecorder() + E = 1.0 # [MPa] Young's modulus + K = 1e6 # 1184 [MPa] + n = 1.0 # [-] + K_p = 1.0 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) + + # initialize notch approximation law + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) + + detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) + + detector.process(signal) + + print(detector.recorder.loads_min) + print(detector.recorder.loads_max) + print(detector.recorder.epsilon_min_LF) + print(detector.recorder.epsilon_max_LF) + + print(detector.recorder.collective) + print(detector.history()) + + print([np.diff(recorder.epsilon_min_LF) <= 0.0]) + assert np.all([np.diff(recorder.epsilon_min_LF) <= 0.0]) + + +def test_epsilon_max_LF_flip(): + signal = -1.0 * pd.Series([10.0, -1.0, 260.0, -250.0, 60.0, -280.0, -100.0]) + signal.index.name = "load_step" + + recorder = RFR.FKMNonlinearRecorder() + E = 1.0 # [MPa] Young's modulus + K = 1e6 # 1184 [MPa] + n = 1.0 # [-] + K_p = 1.0 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) + + # initialize notch approximation law + extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) + + detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) + + detector.process(signal) + + print(detector.recorder.loads_min) + print(detector.recorder.loads_max) + print(detector.recorder.epsilon_min_LF) + print(detector.recorder.epsilon_max_LF) + + print(detector.recorder.collective) + print(detector.history()) + + print([np.diff(recorder.epsilon_max_LF) >= 0.0]) + assert np.all([np.diff(recorder.epsilon_max_LF) >= 0.0]) From 0e7c15e1e1320bf951708fcb896ffbcb499a6ad0 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Wed, 29 Jan 2025 09:45:59 +0100 Subject: [PATCH 20/25] Cleanups in tests Signed-off-by: Johannes Mueller --- tests/stress/rainflow/test_fkm_nonlinear.py | 80 ++++----------------- 1 file changed, 15 insertions(+), 65 deletions(-) diff --git a/tests/stress/rainflow/test_fkm_nonlinear.py b/tests/stress/rainflow/test_fkm_nonlinear.py index 3a60578d..5ca10722 100644 --- a/tests/stress/rainflow/test_fkm_nonlinear.py +++ b/tests/stress/rainflow/test_fkm_nonlinear.py @@ -1090,12 +1090,8 @@ def test_history_guideline_at_split(split_point): pd.testing.assert_frame_equal(df, expected, rtol=1e-1) -def test_epsilon_max_LF_unclosed_hysteresis(): - # Fig 2.3 FKM NL Guideline - - signal = pd.Series([0.0, -100.0, 100.0, 0.0, 200.0, 100.0]) - signal.index.name = "load_step" - +@pytest.fixture +def unitdetector(): recorder = RFR.FKMNonlinearRecorder() E = 1.0 # [MPa] Young's modulus K = 1e6 # 1184 [MPa] @@ -1105,78 +1101,32 @@ def test_epsilon_max_LF_unclosed_hysteresis(): # initialize notch approximation law extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) + return FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) - detector.process(signal) - print(detector.recorder.loads_min) - print(detector.recorder.loads_max) - print(detector.recorder.epsilon_min_LF) - print(detector.recorder.epsilon_max_LF) - - print(detector.recorder.collective) - print(detector.history()) +def test_epsilon_max_LF_unclosed_hysteresis(unitdetector): + signal = pd.Series([0.0, -100.0, 100.0, 0.0, 200.0, 100.0]) + signal.index.name = "load_step" - print([np.diff(recorder.epsilon_min_LF) <= 0.0]) - assert np.all([np.diff(recorder.epsilon_min_LF) <= 0.0]) + unitdetector.process(signal) - print([np.diff(recorder.epsilon_max_LF) >= 0.0]) - assert np.all([np.diff(recorder.epsilon_max_LF) >= 0.0]) + assert np.all([np.diff(unitdetector.recorder.epsilon_min_LF) <= 0.0]) + assert np.all([np.diff(unitdetector.recorder.epsilon_max_LF) >= 0.0]) -def test_epsilon_min_LF_flip(): +def test_epsilon_min_LF_flip(unitdetector): signal = pd.Series([10.0, -1.0, 260.0, -250.0, 60.0, -280.0, -100.0]) signal.index.name = "load_step" - recorder = RFR.FKMNonlinearRecorder() - E = 1.0 # [MPa] Young's modulus - K = 1e6 # 1184 [MPa] - n = 1.0 # [-] - K_p = 1.0 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) + unitdetector.process(signal) - # initialize notch approximation law - extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) + assert np.all([np.diff(unitdetector.recorder.epsilon_min_LF) <= 0.0]) - detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) - detector.process(signal) - - print(detector.recorder.loads_min) - print(detector.recorder.loads_max) - print(detector.recorder.epsilon_min_LF) - print(detector.recorder.epsilon_max_LF) - - print(detector.recorder.collective) - print(detector.history()) - - print([np.diff(recorder.epsilon_min_LF) <= 0.0]) - assert np.all([np.diff(recorder.epsilon_min_LF) <= 0.0]) - - -def test_epsilon_max_LF_flip(): +def test_epsilon_max_LF_flip(unitdetector): signal = -1.0 * pd.Series([10.0, -1.0, 260.0, -250.0, 60.0, -280.0, -100.0]) signal.index.name = "load_step" - recorder = RFR.FKMNonlinearRecorder() - E = 1.0 # [MPa] Young's modulus - K = 1e6 # 1184 [MPa] - n = 1.0 # [-] - K_p = 1.0 # [-] (de: Traglastformzahl) K_p = F_plastic / F_yield (3.1.1) - - # initialize notch approximation law - extended_neuber = NAL.ExtendedNeuber(E, K, n, K_p) - - detector = FKMNonlinearDetector(recorder=recorder, notch_approximation_law=extended_neuber) - - detector.process(signal) - - print(detector.recorder.loads_min) - print(detector.recorder.loads_max) - print(detector.recorder.epsilon_min_LF) - print(detector.recorder.epsilon_max_LF) - - print(detector.recorder.collective) - print(detector.history()) + unitdetector.process(signal) - print([np.diff(recorder.epsilon_max_LF) >= 0.0]) - assert np.all([np.diff(recorder.epsilon_max_LF) >= 0.0]) + assert np.all([np.diff(unitdetector.recorder.epsilon_max_LF) >= 0.0]) From 18675840b36ec635590f0d50011ad9e10f121ad8 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Wed, 5 Feb 2025 11:05:17 +0100 Subject: [PATCH 21/25] Improvements from review and some refactorings Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 19 +- src/pylife/stress/rainflow/fkm_nonlinear.py | 251 ++++++++++++------ 2 files changed, 188 insertions(+), 82 deletions(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index d067cd23..7a0286ac 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -569,7 +569,7 @@ class NotchApproxBinner: def __init__(self, notch_approximation_law, number_of_bins=100): self._n_bins = number_of_bins self._notch_approximation_law = notch_approximation_law - self.ramberg_osgood_relation = notch_approximation_law.ramberg_osgood_relation + self._ramberg_osgood_relation = notch_approximation_law.ramberg_osgood_relation self._max_load_rep = None def initialize(self, max_load): @@ -579,14 +579,14 @@ def initialize(self, max_load): ---------- max_load : array_like The state of the maximum nominal load that is expected. The first - element is chosen as representative to caclulate the lookup table. + element is chosen as representative to calculate the lookup table. Returns ------- self """ max_load = np.asarray(max_load) - self._max_load_rep, _ = self._rep_abs_and_sign(max_load) + self._max_load_rep, _ = self._representative_value_and_sign(max_load) load = self._param_for_lut(self._n_bins, max_load) self._lut_primary = self._notch_approximation_law.primary(load) @@ -596,6 +596,13 @@ def initialize(self, max_load): return self + @property + def ramberg_osgood_relation(self): + """the Ramberg-Osgood relation object, i.e., an object of type RambergOsgood + provided by the notch approximation law. + """ + return self._ramberg_osgood_relation + def primary(self, load): """Lookup the stress strain of the primary branch. @@ -619,7 +626,7 @@ def primary(self, load): """ self._raise_if_uninitialized() - load_rep, sign = self._rep_abs_and_sign(load) + load_rep, sign = self._representative_value_and_sign(load) if load_rep > self._max_load_rep: msg = f"Requested load `{load_rep}`, higher than initialized maximum load `{self._max_load_rep}`" @@ -651,7 +658,7 @@ def secondary(self, delta_load): """ self._raise_if_uninitialized() - delta_load_rep, sign = self._rep_abs_and_sign(delta_load) + delta_load_rep, sign = self._representative_value_and_sign(delta_load) if delta_load_rep > 2.0 * self._max_load_rep: msg = f"Requested load `{delta_load_rep}`, higher than initialized maximum delta load `{2.0*self._max_load_rep}`" @@ -669,7 +676,7 @@ def _param_for_lut(self, number_of_bins, max_val): max_val, scale_m = np.meshgrid(max_val, scale) return (max_val * scale_m) - def _rep_abs_and_sign(self, value): + def _representative_value_and_sign(self, value): value = np.asarray(value) value_rep = value if len(value.shape) == 0 else value[0] return np.abs(value_rep), np.sign(value_rep) diff --git a/src/pylife/stress/rainflow/fkm_nonlinear.py b/src/pylife/stress/rainflow/fkm_nonlinear.py index d14a35f3..4f000fb8 100644 --- a/src/pylife/stress/rainflow/fkm_nonlinear.py +++ b/src/pylife/stress/rainflow/fkm_nonlinear.py @@ -36,6 +36,7 @@ STRAIN = 2 EPS_MIN_LF = 3 EPS_MAX_LF = 4 +STRESS_AND_STRAIN = slice(STRESS, STRAIN+1) PRIMARY = 0 SECONDARY = 1 @@ -43,7 +44,7 @@ MEMORY_1_2 = 1 MEMORY_3 = 0 -HYSTORY_COLUMNS = ["load", "stress", "strain", "secondary_branch"] +HISTORY_COLUMNS = ["load", "stress", "strain", "secondary_branch"] HISTORY_INDEX_LEVELS = [ "load_segment", "load_step", "run_index", "turning_point", "hyst_from", "hyst_to", "hyst_close" ] @@ -106,7 +107,7 @@ def __init__(self, recorder, notch_approximation_law, binner=NAL.NotchApproxBinn self._load_max_seen = 0.0 # maximum seen load value self._run_index = 0 # which run through the load sequence is currently performed - self._last_record = None + self._last_deformation_record = None self._residuals_record = _ResidualsRecord() self._residuals = np.array([]) self._record_vals_residuals = pd.DataFrame() @@ -212,17 +213,17 @@ def process(self, samples, flush=False): load_turning_points.groupby(turning_point_idx, sort=False).first() ) - record, hysts = self._perform_hcm_algorithm(load_turning_points_rep) + deform_type_record, hysts = self._perform_hcm_algorithm(load_turning_points_rep) - if self._last_record is None: - self._last_record = np.zeros((5, self._group_size)) + if self._last_deformation_record is None: + self._last_deformation_record = np.zeros((5, self._group_size)) num_turning_points = len(load_turning_points_rep) - record_vals = self._collect_record(load_turning_points, num_turning_points, record) + record_vals = self._process_deformation(load_turning_points, num_turning_points, deform_type_record) - self._store_recordings_for_history(record, record_vals, turning_point_idx, hysts) + self._store_recordings_for_history(deform_type_record, record_vals, turning_point_idx, hysts) - results = self._process_recording(load_turning_points_rep, record_vals, hysts) + results = self._process_hysteresis(record_vals, hysts) results_min, results_max = results self._update_residuals(record_vals, turning_point_idx, load_turning_points_rep) @@ -324,7 +325,33 @@ def _update_residuals(self, record_vals, turning_point, load_turning_points_rep) self._residuals_record.reindex() - def _collect_record(self, load_turning_points, num_turning_points, record): + def _process_deformation(self, load_turning_points, num_turning_points, deform_type_record): + """Calculate the local stress and strain of all turning points + + In ._perform_hcm_algorithm we recorded which turning point is in + PRIMARY deformation regime and which in SECONDARY + + Now we use this information to calculate the local stress and strain + according to the notch approximation law for each turining point for + each point in the mesh. + + Parameters + ---------- + load_turning_points : pd.Series (N * self._group_size) float + The load distribution of all the turning points. It carries all th + index levels of the initial load signal. + + num_turning_points : int + The number of tunring_points + + deform_type_record : np.ndarray(N, 2) int + The inndex and deformation type of each turning point + (see _perform_hcm_algorithm) + + Returns pd.DataFrame + columns: ["load", "stress", "strain", "epsilon_min_LF", "epsilon_max_LF"] + index: the same like load_turning_points + """ def primary(_prev, load): return self._notch_approximation_law.primary(load) @@ -332,34 +359,44 @@ def secondary(prev, load): prev_load = prev[LOAD] delta_L = load - prev_load - delta = self._notch_approximation_law.secondary(delta_L) + delta_stress_strain = self._notch_approximation_law.secondary(delta_L) + + return prev[STRESS_AND_STRAIN].T + delta_stress_strain + + def prev_record_from_residuals(prev_idx): + idx = len(self._record_vals_residuals) + prev_idx*self._group_size + return self._record_vals_residuals.iloc[idx:idx+self._group_size].to_numpy().T - return prev[STRESS:STRAIN+1].T + delta + def prev_record_from_this_run(prev_idx): + idx = prev_idx * self._group_size + return record_vals[:, idx:idx+self._group_size] def determine_prev_record(prev_idx): if prev_idx < 0: - idx = len(self._record_vals_residuals) + prev_idx*self._group_size - return self._record_vals_residuals.iloc[idx:idx+self._group_size].to_numpy().T - if prev_idx < i: - idx = prev_idx * self._group_size - return record_vals[:, idx:idx+self._group_size] - return self._last_record + return prev_record_from_residuals(prev_idx) + if prev_idx == curr_idx: + return self._last_deformation_record + return prev_record_from_this_run(prev_idx) record_vals = np.empty((5, num_turning_points*self._group_size)) turning_points = load_turning_points.to_numpy() - for i in range(num_turning_points): - prev_record = determine_prev_record(record[i, INDEX]) + for curr_idx in range(num_turning_points): + prev_record = determine_prev_record(deform_type_record[curr_idx, INDEX]) - idx = i * self._group_size - load_turning_point = turning_points[idx:idx+self._group_size] + idx = curr_idx * self._group_size + load = turning_points[idx:idx+self._group_size] + + deformation_function = secondary if deform_type_record[curr_idx, SECONDARY] else primary - deformation_function = secondary if record[i, SECONDARY] else primary result_buf = record_vals[:, idx:idx+self._group_size] - self._process_deformation( - deformation_function, result_buf, load_turning_point, prev_record - ) + + result_buf[LOAD] = load + result_buf[STRESS_AND_STRAIN] = deformation_function(prev_record, load).T + + self._calculate_epsilon_LF(result_buf) + self._last_deformation_record = result_buf record_vals = pd.DataFrame( record_vals.T, @@ -373,31 +410,55 @@ def determine_prev_record(prev_idx): record_vals["turning_point"] = np.stack(tp_index).T.flatten() return record_vals.set_index("turning_point", drop=True, append=True) - def _process_deformation(self, deformation_func, result_buf, load, prev_record): - result_buf[0] = load - res = deformation_func(prev_record, load).T - result_buf[1:3] = res + def _calculate_epsilon_LF(self, deformation_record): + """Calculate epsilon_LF values for the current deformation record - old_load = self._last_record[LOAD, 0] + Parameters + ---------- + deformation_record : np.ndarray (5) + [:3] load, stress, strain + [3:] reserved for epsilon_min_LF and epsilon_max_LF + """ + + old_load = self._last_deformation_record[LOAD, 0] + current_load = deformation_record[LOAD, 0] - if old_load < load[0]: - result_buf[EPS_MAX_LF] = ( - self._last_record[EPS_MAX_LF] - if self._last_record[EPS_MAX_LF, 0] > result_buf[STRAIN, 0] - else result_buf[STRAIN, :] + if old_load < current_load: + deformation_record[EPS_MAX_LF] = ( + self._last_deformation_record[EPS_MAX_LF] + if self._last_deformation_record[EPS_MAX_LF, 0] > deformation_record[STRAIN, 0] + else deformation_record[STRAIN, :] ) - result_buf[EPS_MIN_LF] = self._last_record[EPS_MIN_LF] + deformation_record[EPS_MIN_LF] = self._last_deformation_record[EPS_MIN_LF] else: - result_buf[EPS_MIN_LF] = ( - self._last_record[EPS_MIN_LF] - if self._last_record[EPS_MIN_LF, 0] < result_buf[STRAIN, 0] - else result_buf[STRAIN, :] + deformation_record[EPS_MIN_LF] = ( + self._last_deformation_record[EPS_MIN_LF] + if self._last_deformation_record[EPS_MIN_LF, 0] < deformation_record[STRAIN, 0] + else deformation_record[STRAIN, :] ) - result_buf[EPS_MAX_LF] = self._last_record[EPS_MAX_LF] + deformation_record[EPS_MAX_LF] = self._last_deformation_record[EPS_MAX_LF] + - self._last_record = result_buf + def _process_hysteresis(self, record_vals, hysts): + """Calcuclate all the recorded hysteresis values - def _process_recording(self, turning_points, record_vals, hysts): + For each hysteresis we calculate the two records consisting of + load, stress, strain, epsilon_min_LF, epsilon_max_LF + + Parameters + ---------- + record_vals: pd.DataFrame + colimns: load, stress, strain, epsilon_min_LF, epsilon_max_LF + index of load_turning_points + + hysts: np.ndarray(N, 2) + the recorded hysteresis information + (see ._perform_hcm_algorithm) + + Returns + ------- + result_min, result_max: pd.DataFrame + """ def turn_memory_1_2(values, index): if values[0][0, 0] < values[1][0, 0]: return (values[0], values[1], index[0], index[1]) @@ -414,8 +475,6 @@ def turn_memory_3(values, index): memory_functions = [turn_memory_3, turn_memory_1_2] start = len(self._residuals) - if start: - turning_points = np.concatenate((self._residuals, turning_points)) record_vals_with_residuals = pd.concat([self._record_vals_residuals, record_vals]) @@ -516,64 +575,104 @@ def _adjust_samples_and_flush_for_hcm_first_run(self, samples): return samples, flush def _perform_hcm_algorithm(self, load_turning_points): - """Perform the entire HCM algorithm for all load samples""" + """Perform the entire HCM algorithm for all load samples - # iz: number of not yet closed branches - # ir: number of residual loads corresponding to hystereses that cannot be closed, - # because they contain parts of the primary branch + Parameters + ---------- + load_turning_points : np.ndarray + The representative tunring points of the load signal - # iterate over loads from the given list of samples + Returns + ------- + deform_type_record : np.ndarray (N,2) integer + first column: index in the turning point array + second column: indicate whether the deformation is in PRIMARY or SECONDARY regime + + hysts : np.ndarray (N,4) integer + first column: Type of hysteresis memort (MEMORY_1_2 or MEMORY_3) + second column: index in turning point array of the hysteresis origin + third column: index in turning point array of the hysteresis front + fourth column: index in turning point array of the hysteresis + closing point (-1 if hysteresis not closed) + """ hysts = np.zeros((len(load_turning_points), 4), dtype=np.int64) - hyst_index = 0 + hyst_ptr = 0 - record = -np.ones((len(load_turning_points), 2), dtype=np.int64) - rec_index = 0 + deform_type_record = -np.ones((len(load_turning_points), 2), dtype=np.int64) + deform_ptr = 0 for index, current_load in enumerate(load_turning_points): - hyst_index = self._hcm_process_sample(current_load, index, hysts, hyst_index, record, rec_index) + hyst_ptr = self._hcm_process_sample( + current_load, index, hysts, hyst_ptr, deform_type_record[deform_ptr, :] + ) if np.abs(current_load) > self._load_max_seen: self._load_max_seen = np.abs(current_load) self._iz += 1 - self._residuals_record.append(rec_index, current_load) + self._residuals_record.append(deform_ptr, current_load) + + deform_ptr += 1 + + hysts = hysts[:hyst_ptr, :] + return deform_type_record, hysts + + def _hcm_process_sample( + self, + current_load, + current_idx, + hysts, + hyst_ptr, + deform_type_record, + ): + """Process one sample in the HCM algorithm, i.e., one load value + + Parameters + ---------- + current_load: float + The current representative load - rec_index += 1 + current_index: int + The index of the current load in the turning point list - hysts = hysts[:hyst_index, :] - return record, hysts + hysts: np.ndarray (N,4) integer + The hysteresis record (see ._perform_hcm_algorithm) - def _hcm_process_sample(self, current_load, current_index, hysts, hyst_index, record, rec_index): - """ Process one sample in the HCM algorithm, i.e., one load value """ + hyst_ptr: int + The pointer to the next hysteresis record + + deform_type_record : np.ndarray (N,2) integer + The deformation type record (see ._perform_hcm_algorithm) + """ - record_index = current_index + record_idx = current_idx while True: if self._iz == self._ir: if np.abs(current_load) > self._load_max_seen: # case a) i, "Memory 3" - record[rec_index, :] = [record_index, PRIMARY] + deform_type_record[:] = [record_idx, PRIMARY] residuals_idx = self._residuals_record.current_index - hysts[hyst_index, :] = [MEMORY_3, residuals_idx, current_index, -1] - hyst_index += 1 + hysts[hyst_ptr, :] = [MEMORY_3, residuals_idx, current_idx, -1] + hyst_ptr += 1 self._ir += 1 else: - record[rec_index, :] = [record_index, SECONDARY] + deform_type_record[:] = [record_idx, SECONDARY] break if self._iz < self._ir: - record[rec_index, :] = [record_index, PRIMARY] + deform_type_record[:] = [record_idx, PRIMARY] break # here we have iz > ir: if self._residuals_record.will_remain_open_by(current_load): - record[rec_index, :] = [record_index, SECONDARY] + deform_type_record[:] = [record_idx, SECONDARY] break # no -> we have a new hysteresis @@ -582,7 +681,7 @@ def _hcm_process_sample(self, current_load, current_index, hysts, hyst_index, re prev_idx_0, prev_load_0 = self._residuals_record.pop() if len(self._residuals_record): - record_index = self._residuals_record.current_index + record_idx = self._residuals_record.current_index self._iz -= 2 @@ -594,8 +693,8 @@ def _hcm_process_sample(self, current_load, current_index, hysts, hyst_index, re # the primary branch is not yet reached, continue processing residual loads, potentially # closing even more hysteresis - hysts[hyst_index, :] = [MEMORY_1_2, prev_idx_0, prev_idx_1, current_index] - hyst_index += 1 + hysts[hyst_ptr, :] = [MEMORY_1_2, prev_idx_0, prev_idx_1, current_idx] + hyst_ptr += 1 continue @@ -609,11 +708,11 @@ def _hcm_process_sample(self, current_load, current_index, hysts, hyst_index, re # Proceed on primary path for the rest, which was not part of the closed hysteresis - record[rec_index, :] = [record_index, PRIMARY] - hysts[hyst_index, :] = [MEMORY_1_2, prev_idx_0, prev_idx_1, current_index] - hyst_index += 1 + deform_type_record[:] = [record_idx, PRIMARY] + hysts[hyst_ptr, :] = [MEMORY_1_2, prev_idx_0, prev_idx_1, current_idx] + hyst_ptr += 1 - return hyst_index + return hyst_ptr @property def strain_values(self): @@ -745,8 +844,8 @@ def history(self): ["turning_point", "load_step", "hyst_from", "hyst_to"], ] = -1 history.loc[turning_point_drop_idx, "secondary_branch"] = True - history.loc[negate, HYSTORY_COLUMNS] = -history.loc[ - negate, HYSTORY_COLUMNS + history.loc[negate, HISTORY_COLUMNS] = -history.loc[ + negate, HISTORY_COLUMNS ] history.loc[negate, "hyst_to"] = history.loc[negate + 1, "hyst_to"].to_numpy() history.loc[negate + 1, "hyst_to"] = -1 From 435217fc101be5cd8ca8fa2a246b9d151d29babb Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 6 Feb 2025 10:45:27 +0100 Subject: [PATCH 22/25] Drop maximum_absolute_load parameter for FKM NK rainflow counter Signed-off-by: Johannes Mueller --- .../assessment_nonlinear_standard.py | 34 +++++-------------- 1 file changed, 8 insertions(+), 26 deletions(-) diff --git a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py index 393f0638..612e6c74 100644 --- a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py +++ b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py @@ -213,17 +213,15 @@ def perform_fkm_nonlinear_assessment(assessment_parameters, load_sequence, calcu assessment_parameters, component_woehler_curve_P_RAM, component_woehler_curve_P_RAJ \ = _compute_component_woehler_curves(assessment_parameters) - maximum_absolute_load = _get_maximum_absolute_load(assessment_parameters, scaled_load_sequence) - result = {} # HCM rainflow counting and damage computation for P_RAM if calculate_P_RAM: - result = _compute_lifetimes_P_RAM(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAM, maximum_absolute_load) + result = _compute_lifetimes_P_RAM(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAM) # HCM rainflow counting and damage computation for P_RAJ if calculate_P_RAJ: - result = _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAJ, maximum_absolute_load) + result = _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAJ) # additional quantities result["assessment_parameters"] = assessment_parameters @@ -389,7 +387,7 @@ def _compute_component_woehler_curves(assessment_parameters): return assessment_parameters, component_woehler_curve_P_RAM, component_woehler_curve_P_RAJ -def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolute_): +def _compute_hcm_RAM(assessment_parameters, scaled_load_sequence): """Perform the HCM rainflow counting with the extended Neuber notch approximation. The HCM algorithm is executed twice, as described in the FKM nonlinear guideline.""" @@ -479,7 +477,7 @@ def _store_additional_objects_in_result_RAM(result, recorder, damage_calculator, return result -def _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolute_load): +def _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence): """Perform the HCM rainflow counting with the Seeger-Beste notch approximation. The HCM algorithm is executed twice, as described in the FKM nonlinear guideline.""" @@ -595,23 +593,7 @@ def _store_additional_objects_in_result_RAJ(result, recorder, damage_calculator, return result -def _get_maximum_absolute_load(assessment_parameters, scaled_load_sequence): - """Compute the maximum absolute load, separately for every node (max_load_independently_for_nodes=True). - # The value is used to initialize the binning in the notch approximation schemes. Using 'True' yields smaller bin sizes - # and higher accuracy but takes longer runtime. Using 'False' generates only one lookup-table for the whole assessment, which is - # most accurate for the nodes with the highest loads. - """ - - if "max_load_independently_for_nodes" not in assessment_parameters: - assessment_parameters["max_load_independently_for_nodes"] = False - - maximum_absolute_load = scaled_load_sequence.fkm_load_sequence.maximum_absolute_load( - max_load_independently_for_nodes=assessment_parameters.max_load_independently_for_nodes) - - return maximum_absolute_load - - -def _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAJ, maximum_absolute_load): +def _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAJ): """Compute the lifetimes using the given parameters and woehler curve, with P_RAJ. * Execute the HCM algorithm to detect closed hysteresis. @@ -619,7 +601,7 @@ def _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence * Do statistical assessment and store all results in a dict. """ - detector_1st, detector, seeger_beste_binned, recorder = _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence, maximum_absolute_load) + detector_1st, detector, seeger_beste_binned, recorder = _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence) result, damage_calculator = _compute_damage_and_lifetimes_RAJ(assessment_parameters, recorder, component_woehler_curve_P_RAJ, result) @@ -631,14 +613,14 @@ def _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence return result -def _compute_lifetimes_P_RAM(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAM, maximum_absolute_load): +def _compute_lifetimes_P_RAM(assessment_parameters, result, scaled_load_sequence, component_woehler_curve_P_RAM): """Compute the lifetimes using the given parameters and woehler curve, with P_RAM. * Execute the HCM algorithm to detect closed hysteresis. * Use the woehler curve and the damage parameter to predict lifetimes. * Do statistical assessment and store all results in a dict. """ - detector_1st, detector, extended_neuber_binned, recorder = _compute_hcm_RAM(assessment_parameters, scaled_load_sequence, maximum_absolute_load) + detector_1st, detector, extended_neuber_binned, recorder = _compute_hcm_RAM(assessment_parameters, scaled_load_sequence) result, damage_calculator = _compute_damage_and_lifetimes_RAM(assessment_parameters, recorder, component_woehler_curve_P_RAM, result) From b8e70068ec24959fda75f31bf9869d0cbbe5674b Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 6 Feb 2025 11:27:09 +0100 Subject: [PATCH 23/25] Handle case when first mesh load point is zero Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 22 ++++++++++++++++++- .../test_notch_approximation_law.py | 19 ++++++++++++++++ 2 files changed, 40 insertions(+), 1 deletion(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index 7a0286ac..eb136233 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -571,6 +571,7 @@ def __init__(self, notch_approximation_law, number_of_bins=100): self._notch_approximation_law = notch_approximation_law self._ramberg_osgood_relation = notch_approximation_law.ramberg_osgood_relation self._max_load_rep = None + self._max_load_index = None def initialize(self, max_load): """Initialize with a maximum expected load. @@ -678,5 +679,24 @@ def _param_for_lut(self, number_of_bins, max_val): def _representative_value_and_sign(self, value): value = np.asarray(value) - value_rep = value if len(value.shape) == 0 else value[0] + + single_point = len(value.shape) == 0 + + if self._max_load_index is None and single_point: + self._max_load_index = 0 + + value_rep = value if single_point else self._first_or_maximum_load_of_mesh(value) + return np.abs(value_rep), np.sign(value_rep) + + def _first_or_maximum_load_of_mesh(self, mesh_values): + if self._max_load_index is None: + if mesh_values[0] != 0.0: + self._max_load_index = 0 + else: + self._max_load_index = np.argmax(np.abs(mesh_values)) + if mesh_values[self._max_load_index] == 0.0: + raise ValueError( + "NotchApproxBinner must have at least one non zero point in max_load." + ) + return mesh_values[self._max_load_index] diff --git a/tests/materiallaws/test_notch_approximation_law.py b/tests/materiallaws/test_notch_approximation_law.py index 4d8f3e08..3a2ebe53 100644 --- a/tests/materiallaws/test_notch_approximation_law.py +++ b/tests/materiallaws/test_notch_approximation_law.py @@ -353,3 +353,22 @@ def test_binner_initialized_three_points_secondary_vectorized(L, expected): assert result.shape == (3, 2) np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_binner_zero_first_load(): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + binned = NotchApproxBinner(unbinned, number_of_bins=3).initialize(max_load=[0.0, 300.0]) + + result = binned.secondary(([0.0, 400.0])) + + expected = [[0.0, 0.0], [386.0, 1.99e-3]] + np.testing.assert_allclose(result, expected, rtol=1e-2) + + +def test_binner_total_zero_load(): + unbinned = ExtendedNeuber(E=206e3, K=1184., n=0.187, K_p=3.5) + with pytest.raises( + ValueError, + match="NotchApproxBinner must have at least one non zero point in max_load", + ): + NotchApproxBinner(unbinned, number_of_bins=3).initialize(max_load=[0.0, 0.0]) From b1361904618d17af838ad36ae4b380842b565fc1 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 6 Feb 2025 11:40:18 +0100 Subject: [PATCH 24/25] Add a couple of docstrings to abstract methods Signed-off-by: Johannes Mueller --- .../materiallaws/notch_approximation_law.py | 122 +++++++++++++++++- 1 file changed, 121 insertions(+), 1 deletion(-) diff --git a/src/pylife/materiallaws/notch_approximation_law.py b/src/pylife/materiallaws/notch_approximation_law.py index eb136233..310fd208 100644 --- a/src/pylife/materiallaws/notch_approximation_law.py +++ b/src/pylife/materiallaws/notch_approximation_law.py @@ -82,27 +82,147 @@ def K(self, value): self.K_prime = value @abstractmethod - def load(self, load, *, rtol=1e-4, tol=1e-4): + def load(self, stress, *, rtol=1e-4, tol=1e-4): + """Apply the notch-approximation law "backwards", i.e., compute the linear-elastic stress (called "load" or "L" in FKM nonlinear) + + This is to be reimplemented by derived classes + + Parameters + ---------- + stress : array-like float + The elastic-plastic stress as computed by the notch approximation + rtol : float, optional + The relative tolerance to which the implicit formulation of the load gets solved, + by default 1e-4 + tol : float, optional + The absolute tolerance to which the implicit formulation of the load gets solved, + by default 1e-4 + + Returns + ------- + load : array-like float + The resulting load or lienar-elastic stress. + + """ ... @abstractmethod def stress(self, load, *, rtol=1e-4, tol=1e-4): + r"""Calculate the stress of the primary path in the stress-strain diagram at a given + + This is to be reimplemented by derived classes + + Parameters + ---------- + load : array-like float + The elastic von Mises stress from a linear elastic FEA. + In the FKM nonlinear document, this is also called load "L", because it is derived + from a load-time series. Note that this value is scaled to match the actual loading + in the assessment, it equals the FEM solution times the transfer factor. + rtol : float, optional + The relative tolerance to which the implicit formulation of the stress gets solved, + by default 1e-4 + tol : float, optional + The absolute tolerance to which the implicit formulation of the stress gets solved, + by default 1e-4 + + Returns + ------- + stress : array-like float + The resulting elastic-plastic stress according to the notch-approximation law. + """ "Compute the local notch stress from the local nominal load." ... @abstractmethod def strain(self, load): + """Calculate the strain of the primary path in the stress-strain diagram at a given stress and load. + + This is to be reimplemented by derived classes + + Parameters + ---------- + stress : array-like float + The stress + load : array-like float + The load + + Returns + ------- + strain : array-like float + The resulting strain + """ ... @abstractmethod def load_secondary_branch(self, load, *, rtol=1e-4, tol=1e-4): + """Apply the notch-approximation law "backwards", i.e., compute the linear-elastic stress (called "load" or "L" in FKM nonlinear) from the elastic-plastic stress as from the notch approximation. + + This backward step is needed for the pfp FKM nonlinear surface layer & roughness. + + This is to be reimplemented by derived classes + + Parameters + ---------- + delta_stress : array-like float + The increment of the elastic-plastic stress as computed by the notch approximation + rtol : float, optional + The relative tolerance to which the implicit formulation of the stress gets solved, + by default 1e-4 + tol : float, optional + The absolute tolerance to which the implicit formulation of the stress gets solved, + by default 1e-4 + + Returns + ------- + delta_load : array-like float + The resulting load or lienar-elastic stress. + + """ ... @abstractmethod def stress_secondary_branch(self, load, *, rtol=1e-4, tol=1e-4): + """Calculate the stress on secondary branches in the stress-strain diagram at a given + elastic-plastic stress (load), from a FE computation. + + This is to be reimplemented by derived classes + + Parameters + ---------- + delta_load : array-like float + The load increment of the hysteresis + rtol : float, optional + The relative tolerance to which the implicit formulation of the stress gets solved, + by default 1e-4 + tol : float, optional + The absolute tolerance to which the implicit formulation of the stress gets solved, + by default 1e-4 + + Returns + ------- + delta_stress : array-like float + The resulting stress increment within the hysteresis + """ ... @abstractmethod def strain_secondary_branch(self, load): + """Calculate the strain on secondary branches in the stress-strain diagram at a given stress and load. + + This is to be reimplemented by derived classes + + Parameters + ---------- + delta_sigma : array-like float + The stress increment + delta_load : array-like float + The load increment + + Returns + ------- + strain : array-like float + The resulting strain + """ ... def primary(self, load): From a3d74527e8d3fb3585f2b38937d21eaafc5269c1 Mon Sep 17 00:00:00 2001 From: Johannes Mueller Date: Thu, 6 Feb 2025 12:48:46 +0100 Subject: [PATCH 25/25] Restore notch approximation law object in result output Signed-off-by: Johannes Mueller --- .../strength/fkm_nonlinear/assessment_nonlinear_standard.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py index 612e6c74..34e37890 100644 --- a/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py +++ b/src/pylife/strength/fkm_nonlinear/assessment_nonlinear_standard.py @@ -602,6 +602,7 @@ def _compute_lifetimes_P_RAJ(assessment_parameters, result, scaled_load_sequence """ detector_1st, detector, seeger_beste_binned, recorder = _compute_hcm_RAJ(assessment_parameters, scaled_load_sequence) + result["seeger_beste_binned"] = seeger_beste_binned result, damage_calculator = _compute_damage_and_lifetimes_RAJ(assessment_parameters, recorder, component_woehler_curve_P_RAJ, result) @@ -621,6 +622,7 @@ def _compute_lifetimes_P_RAM(assessment_parameters, result, scaled_load_sequence * Do statistical assessment and store all results in a dict. """ detector_1st, detector, extended_neuber_binned, recorder = _compute_hcm_RAM(assessment_parameters, scaled_load_sequence) + result["extended_neuber_binned"] = extended_neuber_binned result, damage_calculator = _compute_damage_and_lifetimes_RAM(assessment_parameters, recorder, component_woehler_curve_P_RAM, result)