[Scientific / technical calculation by Python] Monte Carlo simulation of thermodynamics of 2D Ising spin system by Metropolis method

Introduction

As an application of the Monte Carlo method to statistical mechanics, simulation of spin magnetism using the Ising model is often performed. ** In this paper, Monte Carlo simulation of 2D Ising model is performed using Python. ** **

The outline of the theory and method is summarized in Addendum, so please refer to it if you are interested. ..

** Contents **

$ N_x \ times N_y $ The Hamiltonian of the two-dimensional Ising model on the lattice is given below [1,2].

H=-J \ \sum_{i=1,\ldots,N_x} \sum_{j=1,\ldots,N_y} s_{i,j}\ s_{i,i+1} - h \sum_i s_i \tag{1} The first term is ** spin-spin interaction ([exchange interaction](https://ja.wikipedia.org/wiki/%E4%BA%A4%E6%8F%9B%E7%9B%B8%E4] % BA% 92% E4% BD% 9C% E7% 94% A8)) **, the second term is the interaction with the magnetic field (Zeeman interaction E3% 82% BC% E3% 83% BC% E3% 83% 9E% E3% 83% B3% E5% 8A% B9% E6% 9E% 9C))) are shown respectively. ** $ J $ is a spin-spin interaction coupling constant, which is an important parameter in characterizing the magnetic properties of the system. $ J \ gt 0 $ describes the ferromagnetic interaction and $ J \ lt 0 $ describes the diamagnetic interaction. ** ** Here, $ (i, j) $ is an index of points on the two-dimensional lattice, and the sum considers the closest interaction. $ s_ {i, j} $ takes $ s_ {i, j} = 1 $ (upspin) or $ s_ {i, j} = -1 $ (downspin). $ h $ is the external magnetic field.

** The purpose of this paper is to determine the temperature dependence of the thermal energy, magnetization, and specific heat of this two-dimensional Ising spin system by Monte Carlo simulation by the Metropolis method. </ font> **

[Calculation condition] In the simulation, a grid of $ 20 \ times 20 $ was used, and $ J = 1 $ and $ h = 0 $. The total number of Monte Carlo steps was set to 40,000, and the number of samplings for thermal near value calculation was set to about 8000 (Monte Carlo steps averaged from about 32000).

In addition, the spin-spin interaction of spins at the boundary is [Periodic boundary condition](https://ja.wikipedia.org/wiki/%E5%91%A8%E6%9C%9F%E7%9A%84% E5% A2% 83% E7% 95% 8C% E6% 9D% A1% E4% BB% B6) was used for evaluation.


** Calculation code **

Uses Python 3.5.1. Implemented and confirmed operation on Jupyter notebook 5.0.0.

"""
Monte Carlo simulation by Metropolis method of 2D Ising model

"""

from random import random, randrange
import numpy as np
import matplotlib.pyplot as plt



#Define a function to calculate the energy for a given spin arrangement alis
def Ecalc2(alis2, Nx,Ny):
    dum=0
    for i in range(-1,Nx):
        for j in range(0,Ny):
            
            l=i
            m=j
            
            if l==-1:
                l=Nx
            if l == Nx:
                l=0
            if m ==-1:
                m=Ny
            if m== Ny:
                m=0
            
            ll=i+1
            if ll == Nx:
                ll=0
            dum+=alis2[l, m]*alis2[ll, m] 
            
    for i in range(0,Nx):
        for j in range(-1,Ny):
            
            l=i
            m=j
            
            if l==-1:
                l=Nx
            if l == Nx:
                l=0
            if m ==-1:
                m=Ny
            if m== Ny:
                m=0
            
            mm=j+1
            if mm == Ny:
                mm=0
            dum+=alis2[l, m]*alis2[l, mm] 
        
    return dum 



#Random initial state generation:The number of upspins and the number of downspins are the same. M=0
def Initial_rand(s, Nx,Ny):
    NN=int(Nx*Ny/2)
    for k in range(NN):
        i=randrange(Nx-1) #A random number is chosen from 1 to N
        j=randrange(Ny-1) #A random number is chosen from 1 to N
        s[i,j]=-1*s[i,j]
    
    return s


Nx= 20  #Divided into 100 in x direction
Ny =20  #Divided into 100 in the y direction
Ntot=Nx*Ny


KBT_lis=np.linspace(0.001,6, 201) #Temperature 0.From 001 to 6(In kBT units)It fluctuates in steps of 201.
ETplot=[]
MTplot=[]
CTplot=[]
for KBT in KBT_lis:
    
    J = 1 #Spin coupling constant
    B = 0.0 #External magnetic field
    steps =40000 #MC step
    avsteps=int(steps/5)
    
        
    #Initial state generation:Random spin placement
    s= np.ones([Nx,Ny], int)  #Quantum number for N spins(Sz = +1 or -1)Set as 1
    s=Initial_rand(s, Nx,Ny)
 
    E = -J* Ecalc2(s, Nx,Ny) -B*np.sum(s) #(initial)Calculate energy.
    E2 = E**2  # E^Store 2
    
    #setup
    eplot = []  #Define a list to store energy values under each temperature kBT
    Eavplot=[] #Define a list to store the thermal mean of energy
    
    e2plot = []#Define a list to store the square of the energy value under each temperature kBT. Used in the calculation of specific heat.
    
    eplot.append(E/Ntot)
    Eavplot.append(E/Ntot)
    
    e2plot.append((E/Ntot)**2)
    
    Cvplot = [] ##Define a list to store the specific heat under each temperature kBT
    
    
    M = np.sum(s)  #Magnetization calculation
    Mplot = []
    Mavplot=[]
    Mplot.append(M/Ntot) #Thermal mean value of magnetization
    Mavplot.append(M/Ntot) #Define a list to store the thermal mean of magnetization

    #Main
    for k in range(steps):
        i=randrange(Nx-1) #0 to N-A random number is chosen up to 1
        j=randrange(Ny-1) #0 to N-A random number is chosen up to 1
    
        s_trial=s.copy()
        s_trial[i,j]= -1*s[i,j]
        delta_E=2*s_trial[i,j]*-1*J*(s[i+1,j]+s[i-1,j]+s[i,j+1]+s[i,j-1])-B*(s_trial[i,j]-s[i,j])

        E_trial =E+  delta_E

        #Status update by metropolis method
        if E_trial < E : 
            s = s_trial
            E = E_trial
        else :
            if random() < np.exp(-(delta_E)/KBT):
                s = s_trial
                E = E_trial
   
        eplot.append(E/Ntot)
        e2plot.append((E/Ntot)**2)
        M = np.sum(s)
        Mplot.append(M/Ntot)
    
    ETplot.append(np.sum(eplot[-avsteps:])/avsteps)
    MTplot.append(np.sum(Mplot[-avsteps:])/avsteps)
    CTplot.append((np.sum(e2plot[-avsteps:])/avsteps-((np.sum(eplot[-avsteps:]))/avsteps)**2)/KBT**2)  #Calculation of specific heat



plt.plot(KBT_lis,ETplot)
plt.ylabel("Totla energy per spin")
plt.xlabel("kBT")
plt.show()

plt.plot(KBT_lis,MTplot)
plt.ylabel("Total magnetic moment per spin")
plt.hlines([0], 0, KBT_lis[-1], linestyles="-")  # y=Draw a line at 0.
plt.xlabel("kBT")
plt.show()


plt.plot(KBT_lis,CTplot)

plt.ylabel("Heat capacity per spin")
plt.xlabel("kBT")
plt.show()

** Result (1) Time evolution of energy and magnetization **

An example of how the amount of substance changes with the Monte Carlo step during a Monte Carlo simulation under isothermal conditions is shown below. This cannot be obtained directly from the calculation code. To display it, write the following plot statement without looping the temperature in the calculation code (for KBT in KBT_lis :).


plt.plot(eplot)
plt.ylabel("Totla energy per spin")
plt.xlabel("kBT")
plt.show()

plt.plot(Mplot)
plt.ylabel("Total magnetic moment per spin")
#plt.hlines([0], 0, steps, linestyles="-")  # y=Draw a line at 0.
plt.xlabel("kBT")
plt.show()

$ k_BT = 1 $. It can be seen that the physics converges as the number of steps increases (the state of reaching thermal equilibrium).

en.png Time evolution of energy per spin.

m.png Time evolution of magnetization per spin.


** Result (2) Temperature dependence **

The figure obtained by the simulation is shown below. If you want to obtain more accurate results, you can increase the number of grids in the simulation. The theoretical value of $ T_c $, the phase transition temperature from ferromagnetism to paramagnetism, is $ T_c = 2.2691 $ [see Addendum (4-1) (12)].

ene.png Total energy per spin. As the temperature rises, pairs with opposite spins occur, increasing energy.

m.png Temperature dependence of magnetization. Magnetization disappears with heating (phase transition from ferromagnet to paramagnetic). Although it is difficult to see in the figure, strictly speaking, it changes discontinuously at $ T = T_c $ (Secondary phase transition. % BB% A2% E7% A7% BB # .E7.AC.AC.E4.BA.8C.E7.A8.AE.E7.9B.B8.E8.BB.A2.E7.A7.BB)) is there.

cv.png Temperature dependence of specific heat. Magnetic transition is [secondary phase transition](https://ja.wikipedia.org/wiki/%E7%9B%B8%E8%BB%A2%E7%A7%BB#.E7.AC.AC.E4 .BA.8C.E7.A8.AE.E7.9B.B8.E8.BB.A2.E7.A7.BB), and the specific heat dissipates at $ T = T_c $. This tendency can be seen from this figure as well.


Addendum


** Addendum (0): Time evolution of spin state **

The change in spin state during the simulation is shown in animation below. The initial spin states were all upward (upspin, $ s_ {i, j} = 1 $), and $ k_BT = 30 $. ** At this temperature, the paramagnetic state becomes stable, and the spin state changes so that the upward and downward spins are about the same **.

output2.gif

Change in spin state of 100x100 grid. ** Pink represents upspin ($ s_ {i, j} = 1 ) and light blue represents downspin ( s_ {i, j} = -1 $). Downspin increases with time **.


** Addendum (1): Theory **

According to quantum statistical mechanics, the thermal mean value of the physical quantity (observable operator) $ \ hat A $ in the thermal equilibrium state of temperature T is $ \ <A > $, where the Hamilton operator of the system is $ \ hat H $. ,

\< A\> = Tr (\hat ρ \hat A) = e^{-\beta \hat H \hat A} /Z \tag{2}

Given in. Here, Tr is the diagonal sum, and $ \ beta = 1 / (k_B T) $ is the quantity called inverse temperature.

$ \ Hat ρ $ is the Density Matrix Operator And in thermal equilibrium \hat ρ = e^{-\beta \hat H} /Z \tag{3} Was used to be.

If you use the display that diagonalizes $ \ hat H $ and $ \ hat A $ at the same time, the diagonal sum is

\< A\> =\sum_{i=1}^{n} A_i \exp(\{ -\beta E_{i} \}), \ \ \ Z =\sum_{i=1}^{n} \exp(\{ -\beta E_{i} \})\tag{4}

Call. $ E_ {i} $ is the $ i $ th energy eigenvalue of $ \ hat H $, and $ Z $ is the partition function E9% 85% 8D% E9% 96% A2% E6% 95% B0).


** Addendum (2): Monte Carlo simulation by Metropolis method **

From Addendum (1), it was found that the thermal mean value of physical quantities under a finite temperature can be theoretically described. This theoretical system will be one of the monuments in theoretical physics. However, from the standpoint of ** actually calculating **, it is generally extremely difficult to calculate the thermal mean value. The reason is that the number of solutions (number of states) of the Schledinger-equation is often enormous and cannot be summed for the states of ** all **.

** Therefore, multiple appropriate states are sampled using random numbers, and $ \ <A > $ is approximately calculated as the average of them (Monte Carlo simulation). ** ** Assuming that the number of samples is N, the $ \ <A_N > _ {MC} $ obtained in this way becomes all pairs as N becomes infinite when the following conditions (1) and (2) are satisfied. It is known that the sum of angles gradually approaches $ \ <A > $, which is calculated exactly. Therefore, if a sufficiently large N is taken, it is possible to accurately obtain the thermal mean value of the physical quantity.

\lim_{N\to\infty} \< A_N \>_{MC} = \< A\> \tag{5}

conditions (1): Monte Carlo simulation is ergodic (2): ** Satisfy "Detailed Balance Conditions" **

(1) is that the state of the system changes with time evolution starting from an arbitrary initial state, but it is possible to go through all possible states (Ergodic condition. .org / wiki /% E3% 82% A8% E3% 83% AB% E3% 82% B4% E3% 83% BC% E3% 83% 89% E7% 90% 86% E8% AB% 96)).

** (2) sets the probability that the state $ i $ is realized in the thermal equilibrium state $ \ pi (i) $ and the probability of transition from the state $ i $ to the state $ j $ (transition probability) $ w (i, j). ) $ \pi(i)\ w(i,j)=\pi(j) \ w(j,i) To meet. ** **

The choice of $ w (i, j) $ is arbitrary. In the metropolis method, the transition probability $ w (i, j) $ from the state i to the state j is set as follows.

w(i,j) =$ 1 (E_j \le E_i)\ $, $ \ e^{-\beta \ (E_j-E_i)} \ (E_j > E_i) \tag{6}$

** This creates various states one after another, and it has been proved that after a sufficient time, the possible states of the system follow the Boltzmann distribution for the canonical ensemble in statistical mechanics [4]. The method of selecting $ w (i, j) $ by the metropolis method therefore has a high affinity with statistical mechanics and is indispensable in computational statistical physics. ** **


** Addendum (3): Procedure for Monte Carlo simulation by Metropolis method of Ising model **

  1. Make an arbitrary spin array $ A_k $

  2. Randomly pick up the i-th spin to generate a new spin placement $ A_ {k + 1} $

  3. Create $ A_ {trial} with the spin direction of spin i inverted (spin flip) to create a trial arrangement, and calculate $ energy $ E_ {tri} $.

  4. If $ E_i \ \ le E_ {trial} , accept the trial placement ( A_ {k + 1} $ = $ A_ {trial}) $

  5. If $ E_i \ \ gt E_ {trial} $, accept with probability $ P_ {tr} = exp [-(\ beta \ (E_ {trial} -E_i))] $ ($ A_ {k + 1} $ = $ A_ {trial} $)

Where 5 is 5-1. Generate a uniform random number $ r $$ (0 \ le r \ le 1) $ 5-2. Compare the magnitude relationship between $ r $ and $ P_ {tr} $, and if $ P_ {tr} \ ge r $, then $ A_ {k + 1} $ = $ A_ {trial} , otherwise If so, reject the state change (( A_ {k + 1} $ = $ A_k $))

And. 6. Repeat steps 1-5 above until the physical quantity is not so different (until the thermal equilibrium is reached). Energy E, magnetization M, and (equal volume) specific heat Cv are, respectively.

E= -J \ \sum_{i=1}^{n} s_i s_{i+1}\tag{7} M= -J \ \sum_{i=1}^{n} s_i \tag{8} $C_v =\frac{E^2-<E>^2}{k_B T^2}\tag{9} $ Can be given. $ n $ is the number of particles. 7. After reaching thermal equilibrium, calculate the average physical quantity in multiple steps. The number of steps depends on the temperature and the number of particles, but the higher the temperature, the more steps may be required (because the energy fluctuation becomes larger).


** Addendum (4): Other **

** Addendum (4-1) Exact solution of 2D Ising model **

There is an exact solution in the no magnetic field (h = 0) 2D Ising model. For the two-dimensional Hamiltonian given by Eq. (1), the asymptotic solutions of energy E, specific heat Cv, and magnetization M when the number of lattice points is increased are as follows [3]. \frac{E}{N}=-J\ coth(j)\ [1+\frac{2}{\pi}\kappa'\ K(\kappa)] \tag{10}

\frac{Cv}{Nk_B}=\frac{1}{\pi} \frac{j}{sinh(j)}[2K(\kappa)-\frac{2}{\pi}-cosh^2(j) \ E(\kappa)+\kappa' sinh^2(j)\ K(\kappa)]\tag{11}

Magnetization M is $ T \ lt T_c (= \ frac {J} {0.4407 \ k_B}) $

M=[\frac{cosh^2(j)}{sinh^4(j)}\ (sinh^2j-1)]^{\frac{1}{8}} \tag{12} And in $ T \ gt T_c $, $ M = 0 $. ** This $ T_c $ is the phase transition temperature from ferromagnetism to paramagnetism, which simply shows that the two-dimensional Ising model can describe the magnetic phase transition. ** In one dimension, $ M $ is a continuous function and there is no magnetic transition (same thing, but it is sometimes written as $ T_c $ = 0K $).

Where j, $ \ kappa $, $ \ kappa'$ are respectively

j=\frac{2J}{k_BT}, \kappa=2\frac{sinh(j)}{cosh^2(j)}, \kappa'=2\frac{sinh^2(j)-1}{cosh^2(j)} \tag{13} Given in. Also, $ K (\ kappa) $ and $ E (\ kappa') $ are first-class and second-class perfect elliptic integrals, respectively.

{K(\kappa)=\int_0^{\pi/2}{\frac{1}{\sqrt{1-k^2\sin^2\theta}}}d\theta}\tag{13}

{E(k)=\int_0^{\pi/2}{\sqrt{1-k^2\sin^2\theta}}d\theta}\tag{14}

(See [5] for calculations in Python).


** Addendum (4-2) Magnetism of matter **

It is in the behavior of "spin", which is one of the quantum mechanical degrees of freedom of electrons, and is the result of spin-spin competition (spin-spin interaction). Therefore, quantum mechanics is necessary to understand the essence of magnets [1]. The governing equation is the Schrodinger equation (in the non-relativistic limit). Therefore, (A) ** In order to simulate the magnetism of matter, it is necessary to solve the Schrodinger equation considering the spin-spin interaction. ** **

Statistical mechanics is a discipline that deals with the thermodynamic properties of matter at finite temperatures. According to it, quantitative information about magnetism at finite temperature is (B) ** knowing all the solutions of the (many-body) Schrodinger equation with many spins, and averaging with some weight ** You can get it.

Unfortunately, it is extremely difficult to do (A) and (B). The following three points are the main reasons.

** (C) The substance of spin-spin interaction is called exchange interaction derived from the statistical nature of electrons, and the functional form of this interaction is unknown **,

** (D) It is difficult to solve the many-body Schrodinger equation even if the functional form of the interaction is known **

** (E) Even if the Schrodinger equation can be solved, it is very difficult to calculate the thermal statistical mean (calculation of partition function, etc.) **

Because.

Therefore, in the study of magnetism, ** model interaction (model Hamiltonian) is set for (C) and (D), and the Schrodinger equation considering it is solved by numerical calculation. For (E), this paper It has often been done for a long time to reduce the amount of calculation by performing ** efficient sampling ** as dealt with in.

Phase by many-body interaction (electronic correlation) Descripting the transition is still a very challenging task. The Ising model and the Heisenberg model (two-dimensional XY model) that handles it quantum mechanically have played an important role in discussing the thermodynamics (phase transition, etc.) of magnetic materials together with the (quantum) Monte Carlo method [2]. ..


References

  1. Kei Yoshida, ["Magnetic"](https://www.amazon.co.jp/%E7%A3%81%E6%80%A7-%E8%8A%B3%E7%94%B0-% E5% A5% 8E / dp / 4000054422), Iwanami Publishing, 1991.

  2. Seiji Miyashita, ["Thermal and Statistical Dynamics"](https://www.amazon.co.jp/%E7%86%B1-%E7%B5%B1%E8%A8%88%E5%8A % 9B% E5% AD% A6-% E7% 89% A9% E7% 90% 86% E5% AD% A6% E5% 9F% BA% E7% A4% 8E% E3% 82% B7% E3% 83% AA% E3% 83% BC% E3% 82% BA-% E5% AE% AE% E4% B8% 8B-% E7% B2% BE% E4% BA% 8C / dp / 4563023841 / ref = sr_1_6? S = books & ie = UTF8 & qid = 1502402477 & sr = 1-6 & keywords =% E5% AE% AE% E4% B8% 8B% E7% B2% BE% E4% BA% 8C), Bakufukan, 1993.

  3. Akira Morita and Yutaka Okabe, ["Simulation Statistical Physics"](https://www.amazon.co.jp/Windows%E3%83%A6%E3%83%BC%E3%82%B6%E3% 83% BC% E3% 81% AE% E3% 81% 9F% E3% 82% 81% E3% 81% AE-% E3% 82% B7% E3% 83% 9F% E3% 83% A5% E3% 83 % AC% E3% 83% BC% E3% 82% B7% E3% 83% A7% E3% 83% B3% E7% B5% B1% E8% A8% 88% E7% 89% A9% E7% 90% 86 -% E6% A3% AE% E7% 94% B0-% E7% AB% A0 / dp / 4621042033 / ref = sr_1_1? S = books & ie = UTF8 & qid = 1502402499 & sr = 1-1 & keywords =% E3% 82% B7% E3% 83% 9F% E3% 83% A5% E3% 83% AC% E3% 83% BC% E3% 82% B7% E3% 83% A7% E3% 83% B3% E7% B5% B1% E8% A8% 88% E7% 89% A9% E7% 90% 86), Maruzen Co., Ltd., 1996.

  4. Qiita's article by kaityo256 "Detailed balance conditions in Monte Carlo method and how to determine transition probability" is very easy to understand and helpful. It was. recommendation.

  5. Qiita article "[Scientific / technical calculation by Python] List of usage of (special) functions used in physics by using scipy" by Python I am writing about the calculation of perfect elliptic integrals.

Postscript

August 12, 2017 We received the following valuable comments from Mr. kaityo256. It will be helpful. Thank you very much.

*** "As you may know, spin-based Monte Carlo simulations generally use cluster update methods such as the Swendsen-Wang and Wolff algorithms. It's a super-powerful method compared to single-spin flips. It's a method. I wrote a simple article about Swendsen-Wang, so I hope you find it helpful. ***

http://qiita.com/kaityo256/items/6539261993e282edc5aa "

Recommended Posts

[Scientific / technical calculation by Python] Monte Carlo simulation of thermodynamics of 2D Ising spin system by Metropolis method
[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] Drawing of 3D curved surface, surface, wireframe, visualization, matplotlib
[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] Basic operation of arrays, numpy
[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] 2D random walk (drunk walk problem), numerical calculation
Scientific / technical calculation by Python] Drawing and visualization of 3D isosurface and its cross section using mayavi
[Scientific / technical calculation by Python] Fitting by nonlinear function, equation of state, scipy
[Scientific / technical calculation by Python] Calculation of matrix product by @ operator, python3.5 or later, numpy
[Scientific / technical calculation by Python] Inverse matrix calculation, numpy
[Scientific / technical calculation by Python] Histogram, visualization, matplotlib
[Scientific / technical calculation by Python] Lagrange interpolation, numerical calculation
[Scientific / technical calculation by Python] Drawing animation of parabolic motion with locus, matplotlib
[Scientific / technical calculation by Python] Solving second-order ordinary differential equations by Numerov method, numerical calculation
[Scientific / technical calculation by Python] Numerical calculation to find the value of derivative (differential)
[Scientific / technical calculation by Python] Analytical solution to find the solution of equation sympy
[Scientific / technical calculation by Python] Plot, visualize, matplotlib 2D data with error bars
[Scientific / technical calculation by Python] Solving one-dimensional Newton equation by the 4th-order Runge-Kutta method
Calculation of the shortest path using the Monte Carlo method
Speed comparison of each language by Monte Carlo method
[Scientific / technical calculation by Python] Logarithmic graph, visualization, matplotlib
[Scientific / technical calculation by Python] Numerical solution of one-dimensional and two-dimensional wave equations by FTCS method (explicit method), hyperbolic partial differential equations
[Scientific / technical calculation by Python] Polar graph, visualization, matplotlib
[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 / technical calculation by Python] List of usage of (special) functions used in physics by using scipy
[Scientific and technical calculation by Python] Drawing of fractal figures [Sierpinski triangle, Bernsley fern, fractal tree]
[Scientific / technical calculation by Python] 3rd order spline interpolation, scipy
[Scientific / technical calculation by Python] Numerical integration, trapezoidal / Simpson law, numerical calculation, scipy
[Scientific / technical calculation by Python] Solving simultaneous linear equations, numerical calculation, numpy
[Scientific / technical calculation by Python] Vector field visualization example, electrostatic field, matplotlib
[Scientific / technical calculation by Python] 1D-3D discrete fast Fourier transform, scipy
Find the ratio of the area of Lake Biwa by the Monte Carlo method
Spin map generation of Ising model by Metropolis-Hastings method (M-H algorithm)
[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] Derivation of analytical solutions for quadratic and cubic equations, mathematical formulas, sympy
Simulate Monte Carlo method in Python
Estimating π by Monte Carlo method
[Scientific / technical calculation by Python] Solving ordinary differential equations, mathematical formulas, sympy
[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 one-dimensional Schrodinger equation in stationary state by shooting method (1), well-type potential, quantum mechanics
[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
Example of 3D skeleton analysis by Python
[Scientific / technical calculation by Python] Solving one-dimensional Schrodinger equation in stationary state by shooting method (2), harmonic oscillator potential, quantum mechanics