Kurvenanpassung einer nicht linearen Gleichung ⇐ Python

[Anonymous]

Ich versuche, meine thermische Leitfähigkeit in die Debye-Callaway-Gleichung einzubeziehen. Einer meiner Parameter ist jedoch negativ. Ich habe verschiedene anfängliche Vermutungen ausprobiert. Daher füge ich einen Code mit thermischen Leitfähigkeitsdaten aus der Literatur und den Werten bei, die sie verstehen mussten, wo mein Modell schief geht. Was mich betrifft, ist, dass es nicht sehr empfindlich gegenüber den anfänglichen Parametern ist. >
b=5.38e-18
Die Parameter, die ich erhielt: < /p>
Angepasste Parameter: 1.417070339751493e-44, B = 1.075485244260637E-15. Gleichungen " src = "https://i.static.net/nww2hypn.png"/>
< /p>
Bearbeiten: Ich habe einen zweiten Code hinzugefügt, der den zweiten Faktor noch zurückgibt negativ. Ich bin mir nicht sicher, ob es daran liegt, dass der Faktor irgendwo auf einem lokalen Minimum stecken bleibt?

Code: Select all

import numpy as np
from scipy.integrate import quad
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

# Constants
kb = 1.380649e-23  # Boltzmann constant (J/K)
hbar = 1.0545718e-34  # Reduced Planck's constant (J*s)
theta_D = 280  # Debye temperature (K)
v = 2300  # Average phonon group velocity (m/s) -

# Define the inverse relaxation times
def tau_PD_inv(omega, A):
return A * omega**4

def tau_U_inv(omega, B, T):
return B * omega**2 * np.exp(-theta_D / (3 * T))

def tau_GB_inv(omega):
D = 0.7e-6  # Grain boundary characteristic length (m)
return v / D

def tau_total_inv(omega, A, B, T):
return tau_PD_inv(omega, A) + tau_U_inv(omega, B, T) + tau_GB_inv(omega)

# Debye Callaway model for thermal conductivity
def kappa_lattice(T, A, B):
def integrand(x):
omega = x * kb * T / hbar
tau_total = 1 / tau_total_inv(omega, A, B, T)
return x**4 * np.exp(x) / (np.exp(x) - 1)**2 * tau_total

integral, _ = quad(integrand, 0, theta_D / T)
prefactor = kb**4 * T**3 / (2 * np.pi**2 * v * hbar**3)
return prefactor * integral

# Vectorize the model for curve fitting
def kappa_lattice_vectorized(T, A, B ):
return np.array([kappa_lattice(Ti, A, B ) for Ti in T])

# Example thermal conductivity and temperature data for Ge WQ
T_data = np.array([3.0748730372823156, 3.1962737114707203, 3.4441443065638015, 3.785528026732028,
4.043999986888003, 4.524033859711473, 4.902329259251676, 5.270E+00,
5.902534836759174, 6.721645553370004, 7.690683028568884, 9.031E+00,
1.103E+01, 12.811754075061538, 16.07373842258417, 19.67675447459597,
23.951202854256202, 29.710517518825778, 3.684E+01, 45.030707655028536,
56.92469139605285, 6.929E+01, 85.62796640159988, 105.41803279520694,
129.5096325800293, 159.77670578646953, 200.9314386671018, 240.34130109997898,
291.63127697234717])  # Temperature in K
kappa_data = np.array([ 0.5021456742832973, 0.5907837911587943, 0.7331309897165834, 0.9570228747931655,
1.1281866396926783, 1.4904675008427906, 1.8682223847710366, 2.1856825810901013,
2.8515011626042894, 4.057791636291224, 5.328111104734368, 7.645627510566391,
11.045409412114584, 14.277055520040301, 19.54065663430087, 24.738891215975602,
28.776422188514683, 32.562485211039274, 33.37572711709385, 32.455084248250884,
28.56802589740088, 24.937473590597776, 20.887192910896978, 17.37073097716312,
14.783541083880793, 12.26010669954394, 9.991973525044475, 8.531955488387288,
7.347906481962733])  # Thermal conductivity in W/mK

# Initial guess for parameters A, B, and C
initial_guess = [4.75E-42, 1.054E-17]

# Curve fitting
popt, pcov = curve_fit(kappa_lattice_vectorized, T_data, kappa_data, p0=initial_guess, maxfev=8_000)
#popt, pcov = curve_fit(kappa_lattice_vectorized, T_data, kappa_data, p0=initial_guess, bounds=(0, np.inf))

# Extract fitted parameters
A_fit, B_fit = popt

# Print the fitted parameters
print(f"Fitted parameters:  A = {A_fit}, B = {B_fit}")

# Generate curve fit data
T_fit = np.linspace(min(T_data), max(T_data), 100)
kappa_fit = kappa_lattice_vectorized(T_fit, *popt)

# Plot the original data, fitted data, and extrapolated data
plt.figure(figsize=(8, 6))
plt.plot(T_data, kappa_data, 'o', label='Experimental data')
plt.plot(T_fit, kappa_fit, '-', label='Fitted model')
#plt.plot(T_extrap, kappa_extrap, '--', label='Extrapolated (400-600K)')
plt.xlabel('Temperature (K)')
plt.ylabel('Thermal Conductivity (W/mK)')
plt.legend()
plt.title('Thermal Conductivity Fitting ')
plt.show()
< /code>
import numpy as np
from scipy.integrate import quad
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

# Constants
kb = 1.380649e-23  # Boltzmann constant (J/K)
hbar = 1.0545718e-34  # Reduced Planck's constant (J*s)
theta_D = 483.0  # Debye temperature (K)
v = 4618.0  # Average phonon group velocity (m/s) -

# Define the inverse relaxation times
def tau_PD_inv(omega, A):
return A * omega**4

def tau_U_inv(omega, B, T):
return B * omega**2 * T * np.exp(-theta_D / (3 * T))

def tau_GB_inv(omega):
D = 1e-6  # Grain boundary characteristic length (m)
return v / D

def tau_new_inv(omega, C):
return C * omega**2

def tau_total_inv(omega, A, B, C, T):
return tau_PD_inv(omega, A) + tau_U_inv(omega, B, T) + tau_GB_inv(omega) + tau_new_inv(omega, C)

# Debye Callaway model for thermal conductivity
def kappa_lattice(T, A, B, C):
def integrand(x):
omega = x * kb * T / hbar
tau_total = 1 / tau_total_inv(omega, A, B, C, T)
return x**4 * np.exp(x) / (np.exp(x) - 1)**2 * tau_total

integral, _ = quad(integrand, 0, theta_D / T)
prefactor = kb**4 * T**3 / (2 * np.pi**2 * v * hbar**3)
return prefactor * integral

# Vectorize the model for curve fitting
def kappa_lattice_vectorized(T, A, B, C):
return np.array([kappa_lattice(Ti, A, B, C) for Ti in T])

# Example thermal conductivity and temperature data for Ge WQ
T_data = np.array([31.21198032, 31.63863233, 33.88415583, 36.02698918, 38.39872854,
39.56437842, 43.09065509, 44.50006634, 46.62080745, 43.20705304,
47.77194375, 51.1374173, 53.33894737, 55.68342411, 57.7839973,
60.44084481, 63.31137245, 66.34735373, 69.16983511, 74.09615366,
77.63998236, 82.66616204, 87.95740211, 93.19748533, 99.59297955,
105.6903923, 113.1334233, 120.6776065, 127.0279836, 134.1693093,
141.3911575, 148.1783285, 155.5512459, 162.631427, 169.7085111,
177.378102, 183.5898785, 190.9367188, 197.6917958, 205.7338436,
212.8065032, 220.0698325, 226.7609477, 234.0571384, 241.1287763,
247.8786187, 255.2722521, 262.3425521, 269.4163916, 276.4872816,
283.5510925, 297.1613034, 304.7919793, 311.8668132, 318.9352276,
326.0927208, 333.0893085, 340.164332, 347.2375209, 354.9542919,
362.0246034, 369.0955138, 375.8247965, 383.5799572])  # Temperature in K
kappa_data = np.array([ 1.085444844, 1.13592882, 1.176429149, 1.217086899,
1.258678341, 1.305290878, 1.357081371, 1.398857871, 1.26084173,
1.436622858, 1.483042548, 1.528031021, 1.561551609, 1.598710119,
1.632418606, 1.665656346, 1.703762245, 1.742520105, 1.774186181,
1.820109169, 1.855169951, 1.893680154, 1.932568555, 1.962172074,
1.999232142, 2.029595098, 2.0490363, 2.067859052, 2.092450212,
2.119053495, 2.131413546, 2.145488003, 2.158842932, 2.169084281,
2.169487956, 2.184424479, 2.196747698, 2.20620033, 2.217744549,
2.223537925, 2.236483003, 2.258522251, 2.273748326, 2.278514432,
2.28270469, 2.288247862, 2.291668986, 2.298648488, 2.302662675,
2.299560106, 2.311410793, 2.339970141, 2.347949378, 2.349474819,
2.354437976, 2.369869851, 2.378039597, 2.384365082, 2.39183344,
2.395266095, 2.399300875, 2.40766296, 2.409888354, 2.442162261])  # Thermal conductivity in W/mK

# Initial guess for parameters A, B, and C
initial_guess = [4.75E-43, 3.6979E-19, 1e-16]

# Curve fitting
popt, pcov = curve_fit(kappa_lattice_vectorized, T_data, kappa_data, p0=initial_guess,maxfev=8_000)
#popt, pcov = curve_fit(kappa_lattice_vectorized, T_data, kappa_data, p0=initial_guess, bounds=(0, np.inf))

# Extract fitted parameters
A_fit, B_fit, C_fit = popt

# Print the fitted parameters
print(f"Fitted parameters:  A = {A_fit}, B = {B_fit}, C = {C_fit}")

# Generate curve fit data
T_fit = np.linspace(min(T_data), max(T_data), 100)
kappa_fit = kappa_lattice_vectorized(T_fit, *popt)

# Plot the original data, fitted data, and extrapolated data
plt.figure(figsize=(8, 6))
plt.plot(T_data, kappa_data, 'o', label='Experimental data')
plt.plot(T_fit, kappa_fit, '-', label='Fitted model')
#plt.plot(T_extrap, kappa_extrap, '--', label='Extrapolated (400-600K)')
plt.xlabel('Temperature (K)')
plt.ylabel('Thermal Conductivity (W/mK)')
plt.legend()
plt.title('Thermal Conductivity Fitting ')
plt.show()