[Scientific / technical calculation by Python] Fitting by nonlinear function, equation of state, scipy

Using the curvefit method of scipy.optimize, the XY data is fitted with a non-linear function to obtain the optimum fitting parameter.

Purpose: Prepare the volume (V) and pressure (P) data of the substance as XY data. It expresses the relationship between the volume and pressure of a solid [Equation of state](https://ja.wikipedia.org/wiki/%E7%8A%B6%E6%85%8B%E6%96%B9%E7%A8 Widely used as one of the% 8B% E5% BC% 8F_ (% E7% 86% B1% E5% 8A% 9B% E5% AD% A6)) 3rd order Birch-Marnagan equation of state Find the fitting parameters by fitting to (: //en.wikipedia.org/wiki/Birch%E2%80%93Murnaghan_equation_of_state). This equation has the following form:


P(V)=\frac{3}{2}K_0[(\frac{V}{V_0})^{-7/3}-(\frac{V}{V_0})^{-5/3}][1+\frac{3}{4}(K_d-4)((\frac{V}{V_0})^{-2/3}-1)]

$ K_0, V_0, K_d $ are the fitting parameters. The X-axis is V and the y-axis is P.

** Depending on the problem, the optimal fit parameter can be obtained for any functional form by changing the specification of the fit function form in the code. Of course, a linear fit is also possible. It is a versatile code except in special situations. ** **

import numpy as np
import scipy.optimize
import matplotlib.pylab as plt

"""
Curve fitting with non-linear function
Example:Fitting the pressure vs. volume relationship with a cubic Birch-Managan equation
"""



x_list=[] # x_define list(Create an empty list)
y_list=[] # y_define list(Create an empty list)

f=open('EoS.dat','rt') #R the file that contains the data you want to plot(Read) t(text)Read in mode

##Read the data, x_list and y_Store the value in list
for line in f:
    data = line[:-1].split(' ')
    x_list.append(float(data[0]))
    y_list.append(float(data[1]))
##
plt.plot(x_list, y_list, 'o',label='Raw data')

#Specifying fitting function
def fitfunc(V, V0, K0, K0d ):
    return (3.0/2.0)*K0*((V/V0)**(-7.0/3.0)-(V/V0)**(-5.0/3.0))*(1.0+(3.0/4.0)*(K0d-4.0)*((V/V0)**(-2.0/3.0)-1.0))


para_ini =[200.0, 100.0, 4.0] #Initial conditions for fitting parameters


para_opt, cov = scipy.optimize.curve_fit(fitfunc, x_list, y_list, para_ini)#Optimizing fitting parameters

#Display of the obtained fitting parameters
print ("Fitted V0 =", str(para_opt[0]), "(Å^3)")
print ("Fitted K0 =", str(para_opt[1]), "(GPa)")
print ("Fitted K0' =", str(para_opt[2]))

##
# for plot
V0=para_opt[0]
K0=para_opt[1]
K0d=para_opt[2]

V1 = np.arange(min(x_list) - 1, max(x_list) + 1, 0.1)
y_fit = (3.0/2.0)*K0*((V1/V0)**(-7/3)-(V1/V0)**(-5.0/3.0))*(1.0+(3.0/4.0)*(K0d-4.0)*((V1/V0)**(-2.0/3.0)-1.0))

plt.plot(V1, y_fit, color='RED',linewidth=2.0, label='Fitted curve')

plt.xticks(fontsize=18)
plt.yticks(fontsize=18)
plt.grid(True)
plt.title("Curve fitting for eqaution of states")
plt.legend(loc='upper right')
plt.xlabel('Volume (Å^3)', fontsize=18)
plt.ylabel('Pressure (GPa)', fontsize=18)

plt.show()

result

Fitted V0 = 163.627068001 (Å^3) Fitted K0 = 225.52189108 (GPa) Fitted K0' = 4.06160382625

t.png


[Appendix] Contents of the used EoS.dat file

164.26 0.0 
156.05 11.9 
147.83 27.9 
146.19 31.7 
144.55 35.7 
142.91 39.9 
141.26 44.4 
139.62 49.2 
137.98 54.2 
136.34 59.5 
134.69 65.1 
133.05 71.0 
131.41 77.3 
129.77 83.9 
128.12 90.9 
126.48 98.3 
124.84 106.2
123.20 114.4
121.55 123.2
119.91 132.5
118.27 142.3
116.62 152.6
114.98 163.6
113.34 175.2
111.70 187.5
110.05 200.6

Recommended Posts

[Scientific / technical calculation by Python] Fitting by nonlinear function, equation of state, scipy
[Scientific / technical calculation by Python] Analytical solution to find the solution of equation sympy
[Scientific / technical calculation by Python] 3rd order spline interpolation, scipy
[Scientific / technical calculation by Python] Basic operation of arrays, numpy
[Scientific / technical calculation by Python] Sum calculation, numerical calculation
[Scientific / technical calculation by Python] Numerical integration, trapezoidal / Simpson law, numerical calculation, scipy
[Scientific / technical calculation by Python] 1D-3D discrete fast Fourier transform, scipy
[Scientific / technical calculation by Python] List of usage of (special) functions used in physics by using scipy
[Scientific / technical calculation by Python] Histogram, visualization, matplotlib
[Scientific / technical calculation by Python] Lagrange interpolation, numerical calculation
[Scientific / technical calculation by Python] Calculation of matrix product by @ operator, python3.5 or later, numpy
[Scientific / technical calculation by Python] Solving (generalized) eigenvalue problem using numpy / scipy, using library
[Scientific / technical calculation by Python] Drawing animation of parabolic motion with locus, matplotlib
[Scientific / technical calculation by Python] Numerical calculation to find the value of derivative (differential)
[Scientific / technical calculation by Python] Drawing of 3D curved surface, surface, wireframe, visualization, matplotlib
[Scientific / technical calculation by Python] Solving one-dimensional Newton equation by the 4th-order Runge-Kutta method
[Scientific / technical calculation by Python] Logarithmic graph, visualization, matplotlib
[Scientific / technical calculation by Python] Polar graph, visualization, matplotlib
[Scientific / technical calculation by Python] Solving one-dimensional Schrodinger equation in stationary state by shooting method (1), well-type potential, quantum mechanics
[Scientific / technical calculation by Python] Plot, visualization, matplotlib of 2D data read from file
[Scientific / technical calculation by Python] Drawing, visualization, matplotlib of 2D (color) contour lines, etc.
[Scientific / technical calculation by Python] Monte Carlo integration, numerical calculation, numpy
[Scientific / technical calculation by Python] Generation of non-uniform random numbers giving a given probability density function, Monte Carlo simulation
[Scientific / technical calculation by Python] Solving one-dimensional Schrodinger equation in stationary state by shooting method (2), harmonic oscillator potential, quantum mechanics
[Scientific / technical calculation by Python] Numerical solution of second-order ordinary differential equations, initial value problem, numerical calculation
[Scientific / technical calculation by Python] List of matrices that appear in Hinpan in numerical linear algebra
[Scientific and technical calculation by Python] Drawing of fractal figures [Sierpinski triangle, Bernsley fern, fractal tree]
[Scientific / technical calculation by Python] Numerical solution of one-dimensional harmonic oscillator problem by velocity Verlet method
[Scientific / technical calculation by Python] Numerical solution of eigenvalue problem of matrix by power method, numerical linear algebra
[Scientific / technical calculation by Python] Vector field visualization example, electrostatic field, matplotlib
[Scientific / technical calculation by Python] 2D random walk (drunk walk problem), numerical calculation
[Scientific / technical calculation by Python] Monte Carlo simulation of thermodynamics of 2D Ising spin system by Metropolis method
[Scientific / technical calculation by Python] Solving ordinary differential equations, mathematical formulas, sympy
[Control engineering] Calculation of transfer functions and state space models by Python
[Scientific / technical calculation by Python] Derivation of analytical solutions for quadratic and cubic equations, mathematical formulas, sympy
[Scientific / technical calculation by Python] Numerical solution of two-dimensional Laplace-Poisson equation by Jacobi method for electrostatic potential, elliptic partial differential equation, boundary value problem
[Scientific / technical calculation by Python] Numerical solution of one-dimensional unsteady heat conduction equation by Crank-Nicholson method (implicit method) and FTCS method (positive solution method), parabolic partial differential equation
[Scientific / technical calculation by Python] Solving the boundary value problem of ordinary differential equations in matrix format, numerical calculation
[Scientific / technical calculation by Python] Solving second-order ordinary differential equations by Numerov method, numerical calculation
Calculation of technical indicators by TA-Lib and pandas
[Scientific / technical calculation by Python] Plot, visualize, matplotlib 2D data with error bars
[Scientific / technical calculation by Python] Solving the Schrodinger equation in the steady state in the 3D isotropic harmonic oscillator potential by the matrix method, boundary value problem, quantum mechanics
[Scientific / technical calculation by Python] Numerical solution of one-dimensional and two-dimensional wave equations by FTCS method (explicit method), hyperbolic partial differential equations
Transfer function / state space model of spring / mass / damper system and simulation by Python
[Python] Numerical calculation of two-dimensional thermal diffusion equation by ADI method (alternate direction implicit method)
[Science / technical calculation by Python] Numerical solution of first-order ordinary differential equations, initial value problem, numerical calculation
[Scientific / technical calculation by Python] Wave "beat" and group velocity, wave superposition, visualization, high school physics
[Python3] Call by dynamically specifying the keyword argument of the function
[python] Value of function object (?)
[Python] Etymology of python function names
Calculation of similarity by MinHash
Nonlinear function modeling in Python
[Science / technical calculation by Python] Comparison of convergence speeds of SOR method, Gauss-Seidel method, and Jacobi method for Laplace equation, partial differential equation, boundary value problem
Comparison of calculation speed by implementation of python mpmath (arbitrary precision calculation) (Note)