In [1]:

    from IPython.display import display, Markdown as md, HTML as html, FileLink
    import math
    import time
    import json
    import zipfile

    try:
        import pandas as pd
        import ipywidgets as widgets
        from ipywidgets import interact, interactive, interactive_output
        import numpy as np
        import matplotlib.pylab as plt
        import scipy
        from sklearn.preprocessing import minmax_scale
        import plotly.express as px
        from plotly.subplots import make_subplots
        import plotly.graph_objects as go
        import statsmodels.api as sm
        from matplotlib import pyplot
        from scipy.optimize import curve_fit

    except Exception as e:
        !pip install pandas ipywidgets matplotlib plotly pandas-bokeh scipy scikit-learn statsmodels
        import pandas as pd
        import ipywidgets as widgets
        from ipywidgets import interact, interactive
        import numpy as np
        import matplotlib.pylab as plt
        import scipy
        from sklearn.preprocessing import minmax_scale
        import plotly.express as px
        from plotly.subplots import make_subplots
        import plotly.graph_objects as go
        import statsmodels.api as sm
        from matplotlib import pyplot
        from scipy.optimize import curve_fit

    #pd.options.plotting.backend = "plotly"

In [2]:

    display(md("# Wöhlerlinie"))
    df = pd.read_csv('woehlerlinie.csv')

    df["kN"] = df[["Prüflast"]].multiply(0.6, axis=1)
    df["Gesamt Lastwechsel"] = df[["Lastwechsel", "Vorbelastung"]].sum(axis=1)

# Wöhlerlinie¶

In [3]:

    display(df.info())
    display(df.head())
    display(df.describe())

    <class 'pandas.core.frame.DataFrame'>
    RangeIndex: 27 entries, 0 to 26
    Data columns (total 9 columns):
     #   Column              Non-Null Count  Dtype
    ---  ------              --------------  -----
     0   Serie               27 non-null     int64
     1   Prüfkörper          27 non-null     object
     2   Prüflast            27 non-null     int64
     3   Vorbelastung        27 non-null     int64
     4   Lastwechsel         27 non-null     int64
     5   Bemerkung           1 non-null      object
     6   Gebrochen           27 non-null     int64
     7   kN                  27 non-null     float64
     8   Gesamt Lastwechsel  27 non-null     int64
    dtypes: float64(1), int64(6), object(2)
    memory usage: 2.0+ KB

    None

| Serie | Prüfkörper | Prüflast | Vorbelastung | Lastwechsel | Bemerkung | Gebrochen | kN | Gesamt Lastwechsel
--- | --- | --- | --- | --- | --- | --- | --- | --- | ---
0 | 1 | Pk 1 | 125 | 0 | 3000000 | NaN | 0 | 75.0 | 3000000
1 | 1 | Pk 14 | 125 | 0 | 3000000 | NaN | 0 | 75.0 | 3000000
2 | 1 | Pk 19 | 125 | 0 | 3000000 | NaN | 0 | 75.0 | 3000000
3 | 1 | Pk 4 | 150 | 0 | 718000 | NaN | 1 | 90.0 | 718000
4 | 1 | Pk 5 | 150 | 0 | 3000000 | NaN | 0 | 90.0 | 3000000

| Serie | Prüflast | Vorbelastung | Lastwechsel | Gebrochen | kN | Gesamt Lastwechsel
--- | --- | --- | --- | --- | --- | --- | ---
count | 27.00000 | 27.000000 | 2.700000e+01 | 2.700000e+01 | 27.000000 | 27.000000 | 2.700000e+01
mean | 2.00000 | 151.851852 | 3.620000e+05 | 1.306313e+06 | 0.703704 | 91.111111 | 1.668313e+06
std | 0.83205 | 23.944975 | 9.533156e+05 | 1.985510e+06 | 0.465322 | 14.366985 | 2.354422e+06
min | 1.00000 | 125.000000 | 0.000000e+00 | 1.000000e+00 | 0.000000 | 75.000000 | 4.830000e+03
25% | 1.00000 | 125.000000 | 0.000000e+00 | 2.824700e+04 | 0.000000 | 75.000000 | 4.615000e+04
50% | 2.00000 | 150.000000 | 0.000000e+00 | 1.712000e+05 | 1.000000 | 90.000000 | 3.422500e+05
75% | 3.00000 | 175.000000 | 2.300000e+04 | 2.203000e+06 | 1.000000 | 105.000000 | 3.000000e+06
max | 3.00000 | 200.000000 | 3.000000e+06 | 6.000000e+06 | 1.000000 | 120.000000 | 9.000000e+06

In [4]:

    #df.set_index('Prüflast', inplace=True)
    df.reset_index()
    #df.groupby("Serie")["Lastwechsel"].plot()
    fig, ax = pyplot.subplots()
    df.plot.scatter(x="Lastwechsel", y="Prüflast",ax=ax)

    #model = sm.formula.ols(formula='Prüflast ~ Lastwechsel', data=df)
    #res = model.fit()
    #df.assign(fit=res.fittedvalues).plot(x='Lastwechsel', y='Prüflast', ax=ax)

Out[4]:

    <Axes: xlabel='Lastwechsel', ylabel='Prüflast'>

In [5]:

    #df.set_index('Prüflast', inplace=True)
    df.reset_index()

    # loc[(df["Prüflast"] != 175)]
    _exp = df.groupby("Prüflast")["Gesamt Lastwechsel"].mean().reset_index()
    print(_exp)
    fig, ax = pyplot.subplots()
    _exp.plot.scatter(x="Gesamt Lastwechsel", y="Prüflast",ax=ax)

    x = _exp["Gesamt Lastwechsel"].to_numpy()
    y = _exp["Prüflast"].to_numpy()

    fit = np.polyfit(x, np.log(y), 1)

    #plt. plot (x, p(x), color="purple")

    print("log(y) =", fit[1], fit[0], "x")
    #print("y =", np.exp(fit[1] + fit[0] * 10000000))

    print("x =", (np.log(100) - fit[1]) / fit[0])

    #def func(x, a, b, c):
    #    return a * np.exp(-b * x) + c
    #popt, pcov = curve_fit(func, x, y)
    #print(popt)
    #plt.plot(x, func(x, *popt), 'r-', label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(popt))

    def monoExp(x, m, t, b):
        return m * np.exp(-t * x) + b
    # perform the fit
    p0 = ( fit[1], -fit[0], 50) # start with values near those we expect
    params, cv = scipy.optimize.curve_fit(monoExp, x, y, p0)
    m, t, b = params
    sampleRate = 20_000 # Hz
    tauSec = (1 / t) / sampleRate

    # determine quality of the fit
    squaredDiffs = np.square(y - monoExp(x, m, t, b))
    squaredDiffsFromMean = np.square(y - np.mean(y))
    rSquared = 1 - np.sum(squaredDiffs) / np.sum(squaredDiffsFromMean)
    print(f"R² = {rSquared}")

    # plot the results
    plt.plot(x, monoExp(x, m, t, b), '--', label="fitted")

    # inspect the parameters
    print(f"Y = {m} * e^(-{t} * x) + {b}")
    print(f"Tau = {tauSec * 1e6} µs")

    #for i in [3000000,10000000,20000000]:
    #    print(monoExp(i, m, t, b))

    #def _exp_func(y):
    #    return (np.log(y) - fit[1]) / fit[0]
    def _exp_func(x):
        return np.exp(fit[1]) * np.exp(fit[0] * x)
    #print(_exp_func(6625238))
    plt. plot (x, _exp_func(x), color="purple")

    x = df["Gesamt Lastwechsel"].to_numpy()
    y = df["Prüflast"].to_numpy()

    fit = np.polyfit(x, np.log(y), 1)
    print("log(y) =", fit[1], fit[0], "x")
    #print("y =", np.exp(fit[1] + fit[0] * 10000000))

    print("x =", (np.log(100) - fit[1]) / fit[0])

    #model = sm.formula.ols(formula='Prüflast ~ Lastwechsel', data=df)
    #res = model.fit()
    #df.assign(fit=res.fittedvalues).plot(x='Lastwechsel', y='Prüflast', ax=ax)

       Prüflast  Gesamt Lastwechsel
    0       125        4.286822e+06
    1       150        6.858833e+05
    2       175        3.977914e+04
    3       200        5.823000e+03
    log(y) = 5.1853843533743715 -8.757634212738333e-08 x
    x = 6625238.657973803
    R² = 0.9229652299165378
    Y = 64.92298636146543 * e^(-1.5016737334552002e-06 * x) + 125.30034798807925
    Tau = 33296180.712274313 µs
    log(y) = 5.092285574940918 -4.864156719896677e-08 x
    x = 10014385.164859032

In [6]:

    import plotly.express as px
    #fig, ax = pyplot.subplots()

    fig = px.scatter(df, x="Lastwechsel", y="Prüflast", trendline="ols", trendline_options=dict(log_x=True), trendline_color_override="green")
    #ax.fill_between(x, y1, 1)
    fig.show()

In [7]:

    _df = df.sort_values(by=["Lastwechsel"])
    x = _df["Lastwechsel"].to_numpy()
    y = _df["Prüflast"].to_numpy()
    z = np. polyfit (x, y, 1)
    p = np. poly1d (z)
    print(z)
    print(p)
    fig, ax = pyplot.subplots()

    ax.fill_between(x, p(x), 100, where=(p(x)>100), color="green", interpolate=True)#, alpha=0.3)
    ax.fill_between(x, p(x), 130, where=(p(x)>130), color="yellow", interpolate=True)#, alpha=0.3)
    ax.fill_between(x, p(x), 155, where=(p(x)>155), color="red", interpolate=True)#, alpha=0.3)

    plt.scatter (x,y, color="purple")
    #calculate equation for trendline
    #add trendline to plot
    plt. plot (x, p(x), color="purple")

    print("Wert fuer Prueflast100:", (100-z[1])/z[0])

    [-7.79554891e-06  1.62035279e+02]

    -7.796e-06 x + 162
    Wert fuer Prueflast100: 7957781.993149513

In [8]:

    import statsmodels.api as sm

    fig, ax = pyplot.subplots()

    #fig = sm.graphics.abline_plot(model_results=results, color='red')
    #fig.axes[0].set(ylim=(-1,0), xlim=(1,2))

    _df = df.sort_values(by=["Lastwechsel"])
    x = _df["Lastwechsel"].to_numpy()
    y = _df["kN"].to_numpy()

    A = np.vstack([np.log(x), np.ones(len(x))]).T
    m, c = np.linalg.lstsq(A, y, rcond=None)[0]

    def p(x):
        return m * np.log(x) + c

    _ = plt.plot(x, y, 'o', label='Meßwerte', markersize=10)
    _ = plt.plot(x, p(x), 'r', label='Wöhlerlinie')
    _ = plt.legend()

    #ax.fill_between(x, p(x), 100 * 0.6, where=(p(x)>110*0.6), color="green", interpolate=True)#, alpha=0.3)
    #ax.fill_between(x, p(x), 130 * 0.6, where=(p(x)>130*0.6), color="yellow", interpolate=True)#, alpha=0.3)
    #ax.fill_between(x, p(x), 155 * 0.6, where=(p(x)>155*0.6), color="red", interpolate=True)#, alpha=0.3)

    ax.set_xlabel("Lastzyklen in Mio.")
    ax.set_ylabel("Last in kN")
    #plt.scatter (x,y, color="purple")
    #calculate equation for trendline
    #add trendline to plot
    #plt. plot (x, p(x), color="purple")
    print(m, c)
    def p_rev(y):
        return math.exp((y-c)/m)
    #print(p(703474306382))
    print("y = 60 :", "%.0f" % p_rev(60))

    print(703474306382 / 20 / 19 / 365)

    -2.0258903908369246 115.26486604656881
    y = 60 : 703474306382
    5071912.807368421

In [9]:

    _df["Berechnete Lastwechsel"] = _df["kN"].apply(p_rev)
    _df["Berechnete Lastwechsel"] = _df["Berechnete Lastwechsel"].map('{:.0f}'.format)
    _df2 = _df.sort_values(by=["Serie", "Prüflast"])
    _df2.to_csv('woehlerlinie_berechnete_lw.csv', index=False)
    _df2

Out[9]:

| Serie | Prüfkörper | Prüflast | Vorbelastung | Lastwechsel | Bemerkung | Gebrochen | kN | Gesamt Lastwechsel | Berechnete Lastwechsel
--- | --- | --- | --- | --- | --- | --- | --- | --- | --- | ---
2 | 1 | Pk 19 | 125 | 0 | 3000000 | NaN | 0 | 75.0 | 3000000 | 428219033
1 | 1 | Pk 14 | 125 | 0 | 3000000 | NaN | 0 | 75.0 | 3000000 | 428219033
0 | 1 | Pk 1 | 125 | 0 | 3000000 | NaN | 0 | 75.0 | 3000000 | 428219033
5 | 1 | Pk 9 | 150 | 0 | 342250 | NaN | 1 | 90.0 | 342250 | 260666
3 | 1 | Pk 4 | 150 | 0 | 718000 | NaN | 1 | 90.0 | 718000 | 260666
4 | 1 | Pk 5 | 150 | 0 | 3000000 | NaN | 0 | 90.0 | 3000000 | 260666
6 | 1 | Pk 8 | 175 | 0 | 31340 | NaN | 1 | 105.0 | 31340 | 159
8 | 1 | Pk 7 | 200 | 0 | 4830 | NaN | 1 | 120.0 | 4830 | 0
7 | 1 | Pk 6 | 200 | 0 | 6816 | NaN | 1 | 120.0 | 6816 | 0
9 | 2 | Pk 18 | 125 | 0 | 1406000 | NaN | 1 | 75.0 | 1406000 | 428219033
10 | 2 | Pk 2 | 125 | 0 | 6000000 | NaN | 0 | 75.0 | 6000000 | 428219033
11 | 2 | Pk 11 | 125 | 0 | 6000000 | NaN | 0 | 75.0 | 6000000 | 428219033
13 | 2 | Pk 17 | 150 | 0 | 107600 | NaN | 1 | 90.0 | 107600 | 260666
12 | 2 | Pk 12 | 150 | 0 | 114700 | NaN | 1 | 90.0 | 114700 | 260666
14 | 2 | Pk 15 | 150 | 0 | 171200 | NaN | 1 | 90.0 | 171200 | 260666
16 | 2 | Pk 13 | 175 | 0 | 11220 | NaN | 1 | 105.0 | 11220 | 159
15 | 2 | Pk 10 | 175 | 0 | 26230 | NaN | 1 | 105.0 | 26230 | 159
17 | 2 | Pk 16 | 175 | 0 | 30264 | NaN | 1 | 105.0 | 30264 | 159
18 | 3 | Pk 1 | 125 | 3000000 | 1 | * Nach dem FTW wurde an der Spannklemme einen ... | 0 | 75.0 | 3000001 | 428219033
20 | 3 | Pk 19 | 125 | 3000000 | 1175400 | NaN | 1 | 75.0 | 4175400 | 428219033
19 | 3 | Pk 14 | 125 | 3000000 | 6000000 | NaN | 0 | 75.0 | 9000000 | 428219033
22 | 3 | Pk 21 | 150 | 235000 | 72900 | NaN | 1 | 90.0 | 307900 | 260666
21 | 3 | Pk 20 | 150 | 235000 | 281800 | NaN | 1 | 90.0 | 516800 | 260666
23 | 3 | Pk 22 | 150 | 235000 | 659500 | NaN | 1 | 90.0 | 894500 | 260666
26 | 3 | Pk 3 | 175 | 23000 | 20800 | NaN | 1 | 105.0 | 43800 | 159
24 | 3 | Pk 23 | 175 | 23000 | 25500 | NaN | 1 | 105.0 | 48500 | 159
25 | 3 | Pk 24 | 175 | 23000 | 64100 | NaN | 1 | 105.0 | 87100 | 159

In [10]:

    import statsmodels.api as sm

    for i in range(1,4):

        fig, ax = pyplot.subplots()

        #fig = sm.graphics.abline_plot(model_results=results, color='red')
        #fig.axes[0].set(ylim=(-1,0), xlim=(1,2))

        _df = df.loc[(df["Serie"] == i)].sort_values(by=["Lastwechsel"])
        x = _df["Lastwechsel"].to_numpy()
        y = _df["Prüflast"].to_numpy()

        #X = sm.add_constant(x)
        #model = sm.OLS(y,X)
        #results = model.fit()
        #print(results.summary())

        A = np.vstack([np.log(x), np.ones(len(x))]).T
        m, c = np.linalg.lstsq(A, y, rcond=None)[0]

        def p(x):
            return m*np.log(x) + c

        _ = plt.plot(x, y, 'o', label='Original data', markersize=10)
        _ = plt.plot(x, p(x), 'r', label='Fitted line')
        _ = plt.legend()

        ax.fill_between(x, p(x), 100, where=(p(x)>100), color="green", interpolate=True)#, alpha=0.3)
        ax.fill_between(x, p(x), 130, where=(p(x)>130), color="yellow", interpolate=True)#, alpha=0.3)
        ax.fill_between(x, p(x), 155, where=(p(x)>155), color="red", interpolate=True)#, alpha=0.3)

        ax.set_xlabel("Lastzyklen in Mio.")
        ax.set_ylabel("Last in %")
    #plt.scatter (x,y, color="purple")
    #calculate equation for trendline
    #add trendline to plot
    #plt. plot (x, p(x), color="purple")

In [11]:

    S_M = 0
    sum_kN = 0

    for index, row in df.iterrows(): #loc[(df["Serie"] == 1)].iterrows():
        #sum_kN += row['kN']
        #if row['Gesamt Lastwechsel'] < 3000000:
        #    S_M += row['kN'] / row['Gesamt Lastwechsel']
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

        #sum_kN += row['Lastwechsel']
        #S_M += row['Lastwechsel'] / row['kN']
        #print(S_M, row['Lastwechsel'], row['kN'])

        # Berechnung für N_OM
        sum_kN += row['Gesamt Lastwechsel']
        S_M += row['Gesamt Lastwechsel'] / row['kN']
        print(sum_kN, S_M, row['Lastwechsel'],row['Vorbelastung'], row['kN'])

        #sum_kN += float(row['Gesamt Lastwechsel'])
        #if row['Gesamt Lastwechsel'] < 3000000:
        #    S_M += float(row['Gesamt Lastwechsel'] / 3000000.0)
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

        #sum_kN += 3000000.0
        #S_M += float(3000000.0 / row['Gesamt Lastwechsel'])
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

    N_M = float(sum_kN / S_M)
    print(sum_kN, S_M, N_M)

    3000000 40000.0 3000000 0 75.0
    6000000 80000.0 3000000 0 75.0
    9000000 120000.0 3000000 0 75.0
    9718000 127977.77777777778 718000 0 90.0
    12718000 161311.11111111112 3000000 0 90.0
    13060250 165113.8888888889 342250 0 90.0
    13091590 165412.3650793651 31340 0 105.0
    13098406 165469.16507936508 6816 0 120.0
    13103236 165509.41507936508 4830 0 120.0
    14509236 184256.08174603173 1406000 0 75.0
    20509236 264256.08174603176 6000000 0 75.0
    26509236 344256.08174603176 6000000 0 75.0
    26623936 345530.5261904762 114700 0 90.0
    26731536 346726.08174603176 107600 0 90.0
    26902736 348628.303968254 171200 0 90.0
    26928966 348878.11349206354 26230 0 105.0
    26940186 348984.9706349207 11220 0 105.0
    26970450 349273.19920634927 30264 0 105.0
    29970451 389273.21253968263 1 3000000 75.0
    38970451 509273.21253968263 6000000 3000000 75.0
    43145851 564945.2125396826 1175400 3000000 75.0
    43662651 570687.4347619049 281800 235000 90.0
    43970551 574108.545873016 72900 235000 90.0
    44865051 584047.4347619049 659500 235000 90.0
    44913551 584509.3395238096 25500 23000 105.0
    45000651 585338.8633333334 64100 23000 105.0
    45044451 585756.0061904762 20800 23000 105.0
    45044451 585756.0061904762 76.89968267325361

0 75.0 3000000 3000000 3000000 0 75.0 3000000 3000000 3000000 0 75.0 3000000 3000000 3000000 0.0001253481894150418 90.0 718000 718000 3000000 0.0001253481894150418 90.0 3000000 3000000 3000000 0.0003883138577861155 90.0 342250 342250 3000000 0.00373866484693736 105.0 31340 31340 3000000 0.021344298649754262 120.0 6816 6816 3000000 0.04618901914664867 120.0 4830 4830 3000000 840.0 0.04618901914664867 18186.140678437587

In [12]:

    display(md("# Test Miner Original"))

    d = {'oai': [200,160,120,80,40], 'ni': [10,90,1000,9000,10000], 'Ni': [31250,95367,401878,1000000,1000000]}

    _df = pd.DataFrame(data=d)

    S_M = 0
    sum_kN = 0

    for index, row in _df.iterrows():
        #sum_kN += row['kN']
        #S_M += row['kN'] / row['Gesamt Lastwechsel']
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'])

        #sum_kN += row['Lastwechsel']
        #S_M += row['Lastwechsel'] / row['kN']
        #print(S_M, row['Lastwechsel'], row['kN'])

        #sum_kN += row['Gesamt Lastwechsel']
        #S_M += row['Gesamt Lastwechsel'] / row['kN']
        #print(sum_kN, S_M, row['Lastwechsel'],row['Vorbelastung'], row['kN'])

        # Gesamtlebensdauer N_OM
        sum_kN += float(row['ni'])
        if row['Ni'] < 1000000:
            S_M += float(row['ni'] / row['Ni'])
        print(S_M, row['oai'], row['ni'],row['Ni'], 1000000)

        #sum_kN += 3000000.0
        #S_M += float(3000000.0 / row['Gesamt Lastwechsel'])
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

    N_M = float(sum_kN / S_M)
    print(sum_kN, S_M, N_M)

# Test Miner Original¶

    0.00032 200 10 31250 1000000
    0.0012637226713643084 160 90 95367 1000000
    0.0037520400214058635 120 1000 401878 1000000
    0.0037520400214058635 80 9000 1000000 1000000
    0.0037520400214058635 40 10000 1000000 1000000
    20100.0 0.0037520400214058635 5357085.714791674

In [13]:

    S_M = 0
    sum_kN = 0

    for index, row in _df2.iterrows():
        #sum_kN += row['kN']
        #if row['Gesamt Lastwechsel'] < 3000000:
        #    S_M += row['kN'] / row['Gesamt Lastwechsel']
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

        #sum_kN += row['Lastwechsel']
        #S_M += row['Lastwechsel'] / row['kN']
        #print(S_M, row['Lastwechsel'], row['kN'])

        # Berechnung für N_OM
        sum_kN += int(row['Berechnete Lastwechsel'])
        if row['Gebrochen']:
            S_M += int(row['Berechnete Lastwechsel']) / row['Gesamt Lastwechsel']

        #print(sum_kN, S_M, row['Gesamt Lastwechsel'], row['Berechnete Lastwechsel'])

        #sum_kN += float(row['Gesamt Lastwechsel'])
        #if row['Gesamt Lastwechsel'] < 3000000:
        #    S_M += float(row['Gesamt Lastwechsel'] / 3000000.0)
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

        #sum_kN += 3000000.0
        #S_M += float(3000000.0 / row['Gesamt Lastwechsel'])
        #print(S_M, row['kN'], row['Gesamt Lastwechsel'],row['Lastwechsel'], 3000000)

    N_M = float(sum_kN / S_M)
    print(sum_kN, S_M, N_M)

    3856318404 416.14713268646807 9266718.670161815

In [ ]: