diff --git a/.gitignore b/.gitignore new file mode 100644 index 00000000..68b5a9a3 --- /dev/null +++ b/.gitignore @@ -0,0 +1,6 @@ +.venv +__pycache__ +test_env +build +gemaakmd_diffusion2d.egg-info +dist diff --git a/README.md b/README.md index fe548354..30c5db86 100644 --- a/README.md +++ b/README.md @@ -8,12 +8,48 @@ The code used in this exercise is based on [Chapter 7 of the book "Learning Scie ## Project description +This code solves the diffusion equation in 2D over a square domain which is at a certain temperature and a circular disc at the center which is at a higher temperature. This code solves the diffusion equation using the Finite Difference Method. The thermal diffusivity and initial conditions of the system can be changed by the user. The code produces four plots at various timepoints of the simulation. The diffusion process can be clearly observed in these plots. + ## Installing the package +To install from TestPyPI (at the time this README is written, the latest version on TestPyPI is 0.0.14): + +```pip3 install -i https://test.pypi.org/simple/ gemaakmd-diffusion2d===0.0.14``` + +or you can install without specifcying the version: + +```pip3 install -i https://test.pypi.org/simple/ gemaakmd-diffusion2d``` + ### Using pip3 to install from PyPI +At the time this README is written, latest version on PyPI is 0.0.14. + +```pip3 install gemaakmd-diffusion2d===0.0.14``` + +or you can install without specifcying the version: + +```pip3 install gemaakmd-diffusion2d``` + ### Required dependencies +It is recommended to use the latest version of matplotlib. + +Matplotlib: + +```pip3 install matplotlib``` + +Numpy: + +```pip3 install numpy``` + ## Running this package +You can include this line to use the solve() method: + +```from gemaakmd_diffusion2d import diffusion2d``` + +``` diffusion2d.solve()``` + ## Citing + + diff --git a/__init__.py b/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/diffusion2d.py b/diffusion2d.py deleted file mode 100644 index c0c6083a..00000000 --- a/diffusion2d.py +++ /dev/null @@ -1,81 +0,0 @@ -""" -Solving the two-dimensional diffusion equation - -Example acquired from https://scipython.com/book/chapter-7-matplotlib/examples/the-two-dimensional-diffusion-equation/ -""" - -import numpy as np -import matplotlib.pyplot as plt - -# plate size, mm -w = h = 10. -# intervals in x-, y- directions, mm -dx = dy = 0.1 -# Thermal diffusivity of steel, mm^2/s -D = 4. - -# Initial cold temperature of square domain -T_cold = 300 - -# Initial hot temperature of circular disc at the center -T_hot = 700 - -# Number of discrete mesh points in X and Y directions -nx, ny = int(w / dx), int(h / dy) - -# Computing a stable time step -dx2, dy2 = dx * dx, dy * dy -dt = dx2 * dy2 / (2 * D * (dx2 + dy2)) - -print("dt = {}".format(dt)) - -u0 = T_cold * np.ones((nx, ny)) -u = u0.copy() - -# Initial conditions - circle of radius r centred at (cx,cy) (mm) -r = min(h, w) / 4.0 -cx = w / 2.0 -cy = h / 2.0 -r2 = r ** 2 -for i in range(nx): - for j in range(ny): - p2 = (i * dx - cx) ** 2 + (j * dy - cy) ** 2 - if p2 < r2: - u0[i, j] = T_hot - - -def do_timestep(u_nm1, u, D, dt, dx2, dy2): - # Propagate with forward-difference in time, central-difference in space - u[1:-1, 1:-1] = u_nm1[1:-1, 1:-1] + D * dt * ( - (u_nm1[2:, 1:-1] - 2 * u_nm1[1:-1, 1:-1] + u_nm1[:-2, 1:-1]) / dx2 - + (u_nm1[1:-1, 2:] - 2 * u_nm1[1:-1, 1:-1] + u_nm1[1:-1, :-2]) / dy2) - - u_nm1 = u.copy() - return u_nm1, u - - -# Number of timesteps -nsteps = 101 -# Output 4 figures at these timesteps -n_output = [0, 10, 50, 100] -fig_counter = 0 -fig = plt.figure() - -# Time loop -for n in range(nsteps): - u0, u = do_timestep(u0, u, D, dt, dx2, dy2) - - # Create figure - if n in n_output: - fig_counter += 1 - ax = fig.add_subplot(220 + fig_counter) - im = ax.imshow(u.copy(), cmap=plt.get_cmap('hot'), vmin=T_cold, vmax=T_hot) # image for color bar axes - ax.set_axis_off() - ax.set_title('{:.1f} ms'.format(n * dt * 1000)) - -# Plot output figures -fig.subplots_adjust(right=0.85) -cbar_ax = fig.add_axes([0.9, 0.15, 0.03, 0.7]) -cbar_ax.set_xlabel('$T$ / K', labelpad=20) -fig.colorbar(im, cax=cbar_ax) -plt.show() diff --git a/gemaakmd_diffusion2d/__init__.py b/gemaakmd_diffusion2d/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/gemaakmd_diffusion2d/diffusion2d.py b/gemaakmd_diffusion2d/diffusion2d.py new file mode 100644 index 00000000..c825e9de --- /dev/null +++ b/gemaakmd_diffusion2d/diffusion2d.py @@ -0,0 +1,74 @@ +""" +Solving the two-dimensional diffusion equation + +Example acquired from https://scipython.com/book/chapter-7-matplotlib/examples/the-two-dimensional-diffusion-equation/ +""" + +import numpy as np +import matplotlib.pyplot as plt +from gemaakmd_diffusion2d.output import create_plots, output_plots + +def do_timestep(u_nm1, u, D, dt, dx2, dy2): + # Propagate with forward-difference in time, central-difference in space + u[1:-1, 1:-1] = u_nm1[1:-1, 1:-1] + D * dt * ( + (u_nm1[2:, 1:-1] - 2 * u_nm1[1:-1, 1:-1] + u_nm1[:-2, 1:-1]) / dx2 + + (u_nm1[1:-1, 2:] - 2 * u_nm1[1:-1, 1:-1] + u_nm1[1:-1, :-2]) / dy2) + + u_nm1 = u.copy() + return u_nm1, u + +def solve(dx = 0.1, dy = 0.1, D = 4): + # plate size, mms + w = h = 10. + # intervals in x-, y- directions, mm + # dx = dy = 0.1 + # Thermal diffusivity of steel, mm^2/s + # D = 4. + + # Initial cold temperature of square domain + T_cold = 300 + + # Initial hot temperature of circular disc at the center + T_hot = 700 + + # Number of discrete mesh points in X and Y directions + nx, ny = int(w / dx), int(h / dy) + + # Computing a stable time step + dx2, dy2 = dx * dx, dy * dy + dt = dx2 * dy2 / (2 * D * (dx2 + dy2)) + + print("dt = {}".format(dt)) + + u0 = T_cold * np.ones((nx, ny)) + u = u0.copy() + + # Initial conditions - circle of radius r centred at (cx,cy) (mm) + r = min(h, w) / 4.0 + cx = w / 2.0 + cy = h / 2.0 + r2 = r ** 2 + for i in range(nx): + for j in range(ny): + p2 = (i * dx - cx) ** 2 + (j * dy - cy) ** 2 + if p2 < r2: + u0[i, j] = T_hot + + # Number of timesteps + nsteps = 101 + # Output 4 figures at these timesteps + n_output = [0, 10, 50, 100] + fig_counter = 0 + fig = plt.figure() + + # Time loop + for n in range(nsteps): + u0, u = do_timestep(u0, u, D, dt, dx2, dy2) + + # Create figure + if n in n_output: + fig_counter += 1 + im = create_plots(fig, u, n, dt, T_cold, T_hot, fig_counter) + + # Plot output figures + output_plots(fig, im) diff --git a/gemaakmd_diffusion2d/output.py b/gemaakmd_diffusion2d/output.py new file mode 100644 index 00000000..b8f9f9bd --- /dev/null +++ b/gemaakmd_diffusion2d/output.py @@ -0,0 +1,16 @@ +import matplotlib.pyplot as plt + +def output_plots(fig, im): + fig.subplots_adjust(right=0.85) + cbar_ax = fig.add_axes([0.9, 0.15, 0.03, 0.7]) + cbar_ax.set_xlabel('$T$ / K', labelpad=20) + fig.colorbar(im, cax=cbar_ax) + plt.show() + return + +def create_plots(fig, u, n, dt, T_cold, T_hot, fig_counter): + ax = fig.add_subplot(220 + fig_counter) + im = ax.imshow(u.copy(), cmap=plt.get_cmap('hot'), vmin=T_cold, vmax=T_hot) # image for color bar axes + ax.set_axis_off() + ax.set_title('{:.1f} ms'.format(n * dt * 1000)) + return im \ No newline at end of file diff --git a/setup.py b/setup.py new file mode 100644 index 00000000..b8a55c9d --- /dev/null +++ b/setup.py @@ -0,0 +1,25 @@ +from setuptools import setup +import setuptools + +# Read the README.md file. so this shows up as the Project Description in the PyPI/TestPyPI page. +with open("README.md", "r", encoding="utf-8") as fh: + long_description = fh.read() + +setup( + name="gemaakmd_diffusion2d", + version="0.0.14", + author="Muhammad Gema Akbar", + description="SSE Python Exercise", + long_description=long_description, + long_description_content_type="text/markdown", + url="https://github.com/mgemaakbar/diffusion2D", + package_dir={"": "."}, + packages=setuptools.find_packages(where="."), + install_requires=[ + "matplotlib>=3.9.0", + "numpy>=2.1.0" + ] + # entry_points={ + # 'console_scripts': ['package-import-name = '] + # } +) \ No newline at end of file