Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- import math
- import operator
- import numpy as np
- from sklearn import linear_model
- import numpy as np
- from scipy.stats import norm
- import matplotlib.pyplot as plt
- import numpy as np
- from scipy.optimize import curve_fit
- from numpy import arange,array,ones,polyfit,polyval#,random,linalg
- from pylab import plot,show,figure,legend
- from scipy import stats
- import math
- #
- #
- # stress = [70, 70, 70, 70, 70, 45, 45, 45, 300, 300, 300, 300,\
- # 300, 300, 300, 300, 300, 300, 200, 200, 200, 200, 200]
- # temperature = [523, 523, 523, 523, 523, 573, 573, 573, 423,\
- # 423, 423, 423, 423, 423, 423, 423, 423, 423, 473, 473, 473, 473, 473]
- # time = [254, 253, 255, 256, 252, 160, 162, 158, 165, 155, 157,\
- # 163, 164, 156, 153, 167, 155, 165, 17.3, 17, 16.7, 17.5, 16.5]
- import pandas as pd
- parameters = ['Stress', 'Temperature', 'MinimumCreepRate', 'Priorplasticdeformation', 'Durationofprimarycreep',\
- 'Timetothirdcreepperiod', 'Timetorupture','Deformation', 'Time', 'MinimumCreepRate2', 'Time3', 'MinimumCreepRate3']
- path = "C:\\PostGraduate\\Experiments_new.csv"
- def read_csv():
- data = pd.read_csv(path, sep=',',
- names=parameters, skiprows=1).astype('float')
- return data
- def process_calc():
- data = read_csv()
- result = {}
- n, h, a0 = get_first_params(data, result)
- l = get_l(data, result)
- m, p, b0 = get_second_params(data, result)
- data['a'] = (data['MinimumCreepRate2'] / (data['Stress'].pow(n) * np.exp(-h/data['Temperature'])))
- data['b'] = (1/((l+1)*data['Timetorupture']*np.exp(-p/data['Temperature'])*data['Stress'].pow(m)))
- data['def_first'] = data['MinimumCreepRate2'] * data['Durationofprimarycreep']
- # data['teta'] = data['a'] * np.exp(-h/data['Temperature']) * data['Stress'].pow(n) * data['Durationofprimarycreep']
- temp = 1-(l+1) * data['b'] * np.exp(-p/data['Temperature']) * data['Stress'].pow(m) * data['Durationofprimarycreep']
- data['temp'] = temp
- data['teta'] = data['a'] * np.exp((p-h)/data['Temperature']) * data['Stress'].pow(n-m) * \
- ( data['temp'].pow((l - n +1) / (l + 1)) - 1) / (data['b'] * (n - l - 1))
- return result, data
- def get_first_params(data, result):
- x1 = np.log(data['Stress'])
- x2 = - data['Temperature'].pow(-1)
- y = np.log(data['MinimumCreepRate2'])
- # linregress doesn't do multi regression, so we use sklearn
- x = np.vstack([x1, x2]).T # don't need ones for
- regr = linear_model.LinearRegression()
- regr.fit( x, y )
- coef = regr.coef_
- intercept = regr.intercept_
- result['a0'] = np.exp(intercept)
- result['n'] = coef[0]
- result['h'] = coef[1]
- return coef[0], coef[1], np.exp(intercept)
- def get_l(data, result):
- n = result['n']
- a0 = result['a0']
- h = result['h']
- data['Minus_Temperature_reverse'] = - data['Temperature'].pow(-1)
- data['l'] = n - 1 + n/((data['Deformation']/\
- (a0*np.exp(h*data['Minus_Temperature_reverse'])*data['Stress'].pow(n)*data['Timetorupture'])) - 1)
- l = data['l'].mean()
- result['l'] = l
- return l
- def get_second_params(data, result):
- x1 = - np.log(data['Stress'])
- x2 = data['Temperature'].pow(-1)
- y = np.log(data['Timetorupture'])
- l = result['l']
- # linregress doesn't do multi regression, so we use sklearn
- x = np.vstack([x1, x2]).T # don't need ones for
- regr = linear_model.LinearRegression()
- regr.fit( x, y )
- coef = regr.coef_
- intercept = regr.intercept_
- result['b0'] = np.exp(-intercept)/(l+1)
- result['m'] = coef[0]
- result['p'] = coef[1]
- return coef[0], coef[1], np.exp(-intercept)/(l+1)
- def get_characteristics(data):
- means = []
- stds = []
- stresses = []
- temperatures = []
- for k, group in data.groupby(['Stress', 'Temperature']):
- stress = k[0]
- temperature = k[1]
- data = group['b']
- # mean = data.mean()
- mean, std = norm.fit(data)
- # var = data.var()
- print('stress = ', stress, 'temperature = ', temperature, \
- 'mean = ', mean, 'std = ', std)
- stresses.append(stress)
- temperatures.append(temperature)
- means.append(mean)
- stds.append(std)
- return stresses, temperatures, means, stds
- def get_params_multiple_regressions(depend_var1, depend_var2, independ_var):
- x = np.vstack([depend_var1, depend_var2]).T # don't need ones for
- regr = linear_model.LinearRegression()
- regr.fit(x, independ_var)
- return regr.coef_, regr.intercept_
- def get_params_multiple_regressions_many(depend_var1, depend_var2, independ_var):
- # depend_var1 - stress
- # depend_var2 - temperature
- depend_var1_sq = [dep*dep for dep in depend_var1]
- depend_var2_sq = [dep*dep for dep in depend_var2]
- depend_var1_var2_sq = [depend_var1[index]*depend_var2[index] for index in range(len(depend_var1))]
- x = np.vstack([depend_var1_sq, depend_var2_sq, depend_var1_var2_sq, depend_var1, depend_var2]).T # don't need ones for
- regr = linear_model.LinearRegression()
- regr.fit(x, independ_var)
- return regr.coef_, regr.intercept_
- def get_value_by_params(coefficient, intercept, params):
- summa = 0
- for i in range(len(coefficient)):
- summa+=params[i]*coefficient[i]
- summa+=intercept
- return summa
- def calculate_time(mean_b, std_b, result, stress, temperature):
- from numpy import random
- r = random.normal(mean_b, std_b, size = 1000)
- return 1/((result['l']+1)*r*np.exp(-result['p']/temperature)*math.pow(stress, result['m']))
- def plot_graphs(time):
- plot_pdf(time)
- plot_cdf(time)
- def plot_pdf(time):
- mu, std = norm.fit(time)
- # Plot the histogram.
- plt.hist(time, bins=25, normed=True, alpha=0.6, color='g')
- # Plot the PDF.
- xmin, xmax = plt.xlim()
- x = np.linspace(xmin, xmax, 100)
- p = norm.pdf(x, mu, std)
- plt.plot(x, p, 'k', linewidth=2)
- title = "Fit results: mu = %.2f, std = %.2f" % (mu, std)
- plt.title(title)
- plt.show()
- def plot_cdf(time):
- mu, std = norm.fit(time)
- # Plot the histogram.
- plt.hist(time, bins=25, normed=True, alpha=0.6, color='g')
- # Plot the CDF.
- xmin, xmax = plt.xlim()
- x = np.linspace(xmin, xmax, 100)
- c = norm.cdf(x, mu, std)
- plt.plot(x, c, 'k', linewidth=2)
- title = "Fit results: mu = %.2f, std = %.2f" % (mu, std)
- plt.title(title)
- plt.show()
- def get_selected_stresses(result, data):
- temperature_2_selected_stress = {}
- for temperature, data_temp in data.groupby(['Temperature']):
- # data_filtered = data[data['Temperature']==temperature]
- median_deform = data_temp['def_first'].median()
- selected_stress = get_selected_stress(data_temp, median_deform)
- temperature_2_selected_stress[temperature] = selected_stress
- return list(set(temperature_2_selected_stress.values()))
- def get_params_c_k(result, data):
- c_list = []
- k_list = []
- selected_stresses = get_selected_stresses(result, data)
- data_filtered = data[data['Stress'].isin(selected_stresses)]
- for stress, group in data_filtered.groupby(['Stress']):
- def_first = group['def_first'].values
- teta = group['teta'].values
- c, k = get_estimations_params(def_first, teta)
- c_list.append(c)
- k_list.append(k)
- filtered_c_list = [c for c in c_list if c != 0]
- filtered_k_list = [k for k in k_list if k != 0]
- return np.mean(filtered_c_list), np.mean(filtered_k_list)
- def get_charact_c_k_2(result, data):
- stresses = []
- temperatures = []
- c_list = []
- k_list = []
- for k, group in data.groupby(['Stress', 'Temperature']):
- stress = k[0]
- temperature = k[1]
- def_first = group['def_first'].values
- teta = group['teta'].values
- c, k = get_estimations_params(def_first, teta)
- if c != 0 or k != 0:
- stresses.append(stress)
- temperatures.append(temperature)
- c_list.append(c)
- k_list.append(k)
- return stresses, temperatures, c_list, k_list
- def get_params_c_k_2(stress, temperature, result, data):
- stresses, temperatures, c_list, k_list = get_charact_c_k_2(result, data)
- params = [stress * stress, temperature * temperature, stress * temperature, stress, temperature]
- # estimate c
- mean_coef_c, mean_intercept_c = get_params_multiple_regressions_many(stresses, temperatures, c_list)
- mean_c = get_value_by_params(mean_coef_c, mean_intercept_c, params)
- # estimate k
- mean_coef_k, mean_intercept_k = get_params_multiple_regressions_many(stresses, temperatures, k_list)
- mean_k = get_value_by_params(mean_coef_k, mean_intercept_k, params)
- return mean_c, mean_k
- def estimation_function(x, c, k):
- return k*np.log((1+c)*np.exp(x/k) - c)
- def func(x, p1,p2):
- return p1*np.cos(p2*x) + p2*np.sin(p1*x)
- def get_estimations_params(def_first, teta):
- try:
- popt, pcov = curve_fit(estimation_function, teta, def_first)
- c = popt[0]
- k = popt[1]
- data_checked = estimation_function(teta, c, k)
- residuals = def_first - data_checked
- error = np.sqrt(sum(residuals**2))
- return c, k
- except:
- return 0, 0
- def get_selected_stress(data_filtered, median_deform):
- stress_2_distance = {}
- for stress, data_stress in data_filtered.groupby(['Stress']):
- # deform_first = data_filtered[data_filtered['Stress'] == stress]
- median = data_stress['def_first'].median()
- stress_2_distance[stress] = math.fabs(median_deform - median)
- stress_res = sorted(stress_2_distance.items(), key=operator.itemgetter(1))
- return stress_res[0][0]
- def calculate_time_to_failure(data, result, stress, temperature):
- stresses, temperatures, means, stds = get_characteristics(data)
- # mean_coef, mean_intercept = get_params_multiple_regressions(stresses, temperatures, means)
- # std_coef, std_intercept = get_params_multiple_regressions(stresses, temperatures, stds)
- # params = [stress, temperature]
- mean_coef, mean_intercept = get_params_multiple_regressions_many(stresses, temperatures, means)
- std_coef, std_intercept = get_params_multiple_regressions_many(stresses, temperatures, stds)
- params = [stress * stress, temperature * temperature, stress * temperature, stress, temperature]
- mean_b = get_value_by_params(mean_coef, mean_intercept, params)
- std_b = get_value_by_params(std_coef, std_intercept, params)
- time = calculate_time(mean_b, std_b, result, stress, temperature)
- return time
- def mean_calculation(stress, temperature):
- result, data = process_calc()
- time = calculate_time_to_failure(data, result, stress, temperature)
- # plot_graphs(time)
- # print('mean_Time = ', np.mean(time), 'std_Time = ', np.std(time))
- # c, k = get_params_c_k(result, data)
- c, k = get_params_c_k_2(stress, temperature, result, data)
- result['c'] = c
- result['k'] = k
- print('res')
- if __name__ == '__main__':
- # stress = 45
- # temperature = 573
- stresses_test = [70, 45, 300, 200]
- temperatures_test = [523, 573, 423, 473]
- for i in range(len(stresses_test)):
- mean_calculation(stresses_test[i], temperatures_test[i])
- # mean_Time = 202.977268144 std_Time = 35.3658215247
- # mean_Time = 229.291468263 std_Time = 23.530169309
- # mean_Time = 36.3147799399 std_Time = 21.4642344076
- # mean_Time = 362.608885147 std_Time = 48.3941602839
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement