7-5　2D Diffusion Equation　[python, julia, fortran]

• English (United States)
• 日本語

This time, we did two things: “validation to confirm that the code is running properly” and “checking difference between python and julia fortran”.

1. validation

$$\frac{\partial}{\partial t}(\rho C T)= \frac{\partial}{\partial x} \left( \lambda \frac{\partial T}{\partial x} \right ) + \frac{\partial}{\partial y} \left ( \lambda \frac{\partial T}{\partial y} \right )$$

The density \(\rho\), specific heat \(C\), and thermal conductivity \(\lambda\)constant conditions are imposed as calculation conditions.

Then, using the diffusion coefficient (D), the equation can be written as:

$$\frac{\partial T}{\partial t}=\frac{\lambda}{\rho C} \left( \frac{\partial T}{\partial x^2}+\frac{\partial T}{\partial y^2} \right)  $$

$$\frac{\partial T}{\partial t}=D \left( \frac{\partial T}{\partial x^2}+\frac{\partial T}{\partial y^2} \right)  $$

Now, I calculated using python under the following conditions.

There are five things to do: “one-dimensional exact solution”, “one-dimensional explicit method”, “one-dimensional implicit method”, “two-dimensional explicit method”, and “two-dimensional implicit method”.

The number of time steps and time step size are set under each condition.

Initial Temperture：\( T(x) = \left\{ \begin{array}{c} 2x \ \ \ \ \ (0\le x \le 1/2) \\ 2(1-x) \ \ \ \ \ (1/2 \le x \le 1) \end{array} \right.\)

Boundary Condition：\( T(0,t)=T(L,t)=0\)

Diffusion Coefficient：\(D=1.0\)

Number of Cells in x Direction：\(nx_0=100\)

Number of Cells in y Direction：\(ny_0=100\)

Step Size in x Direction：\(dx=0.01\)

Step Size in y Direction：\(dy=0.01\)

⇒The width is determined as follows from the above conditions．
\(L_x=1.0, L_y=1.0\)

The one-dimensional exact solution under the above conditions is

$$T(x,t)=\sum ^{\inf}_{n=1} \frac{8D}{n^2 \pi ^2} sin(\frac{n \pi}{2}) sin(\frac{n \pi x}{L}) e^{(-D \frac{n^2 \pi ^2}{L^2}t)}$$

\(D=1.0, L=1.0\)

$$T(x,t)=\sum ^{\inf}_{n=1} \frac{8}{n^2 \pi ^2} sin(\frac{n \pi}{2}) sin(n \pi x) e^{(-n^2 \pi ^2t)}$$

$$ \simeq \frac{8}{\pi ^2} \left[ sin(\pi x) e^{-\pi ^2 t}- \frac{1}{3^2} sin(3 \pi x) e^{-3^2 \pi ^2 t} + \frac{1}{5^2} sin(5 \pi x) e^{-5^2 \pi ^2 t} \right]$$

Calculation of exact solution is troublesome, so I wrote only the result.

There are too many codes, so I put the code in [validation] on Github below.

https://github.com/hide-dog/2d-heat-equation

This is the result.

I guess the calculation is working well even in 2D.

In this case, \(\Delta x, \Delta y\) of the same size is used in all conditions, and \(\Delta t\) is changed for accuracy improvement.

I think that the two-dimensional calculation has a visible error because \(\Delta x, \Delta y, \Delta t\) is still coarse. If you are interested, please reduce \(\Delta t\). Notice! The calculation time will increase.

The error is easily reduced.

2. Comparison

Next, We will compare python, julia and fortran.

That’s because it took a lot of time when I increased the the number of time steps with Python.

The conditions were set as follows.

Number of Time Steps：(nt=1000)

Step Size：\(nt=1000\)

Initial Temperature：\(T(x,y,0)=293 (20℃)\)

Boundary Condition ： \( \left\{ \begin{array}{c} \frac{\partial T(x,0,t)}{\partial y}= \frac{\partial T(x,L_y,t)}{\partial y} =0\ \\ T(0,y,t)=573 (300℃) ,\ \ \ \ \ T(L_x,y,t)=293 (20℃) \end{array} \right.\)

Thermal Conductivity of Alumina：\(\lambda=32.0\)

Density of Alumina：\(\rho=3.94e3\)

Specific Heat of Alumina：\(c=0.78e3\)

Number of Cells in x Direction：\(nx_0=100\)

Number of Cells in y Direction：\(ny_0=100 \)

Step Size in x Direction：\(dx=0.001\)

Step Size in y Direction：\(dy=0.001\)

I use 3 languages, “python”, “julia” and “fortran”．

There are too many codes, so I put them in [comparison] on the following Github.

https://github.com/hide-dog/2d-heat-equation

“julia” is tuned so that global variables are reduced as much as possible.

Fortran is compiled with the option “-O4” so it is optimized at compile time.

The result is as follows. (A full plot is at the end.)

Although Python uses global variables, dynamically typed julia can be over 100 times faster with a little tweaking.

I haven’t tuned fortran because it’s statically typed, but it’s fast enough.

After all, python is for data processing and seems not suitable for numerical analysis.

Next time, the number of grids will increase and the calculation will be heavy, so I will use julia or fortran.

Thank you !

by hide

*おまけ

The results of increasing the number of steps under the condition of “2 compariosn” are attached.