FTCS prototyping

This commit is contained in:
philippun 2023-06-22 11:57:46 +02:00
parent 72509832d4
commit be68985600

View File

@ -0,0 +1,249 @@
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"# imports\n",
"import numpy as np"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {},
"outputs": [],
"source": [
"# environment\n",
"dim = 2\n",
"grid_size = {'x': 20, 'y': 20} # to start out\n",
"domain_size = {'x': 20, 'y': 20} # to start out\n",
"\n",
"delta_x = domain_size['x'] / grid_size['x']\n",
"delta_y = domain_size['y'] / grid_size['y'] \n",
"\n",
"alpha_x = np.ones((grid_size['x'],grid_size['y']))\n",
"alpha_y = np.ones((grid_size['x'],grid_size['y']))\n",
"\n",
"max_stable_time_step_x = delta_x**2 / (2 * np.max(alpha_x))\n",
"max_stable_time_step_y = delta_y**2 / (2 * np.max(alpha_y))\n",
"max_stable_time_step = max_stable_time_step_x if max_stable_time_step_x <= max_stable_time_step_y else max_stable_time_step_y\n",
"\n",
"C_t = np.zeros((grid_size['x'],grid_size['y']))"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {},
"outputs": [],
"source": [
"# initialization\n",
"C_t[1,1] = 1"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"def alpha_interblock(alpha1, alpha2, harmonic=False):\n",
" if not harmonic:\n",
" return 0.5 * (alpha1 + alpha2)\n",
" else:\n",
" return 2 / ((1/alpha1) + (1/alpha2))\n"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {},
"outputs": [],
"source": [
"# simulating one time step\n",
"def simulate(C_t): \n",
" C_t1 = np.zeros((grid_size['x'],grid_size['y']))\n",
"\n",
" # inner cells\n",
" for i in range(1, grid_size['x']-1):\n",
" for j in range(1, grid_size['y']-1):\n",
" C_t1[i,j] = C_t[i,j] \\\n",
" + max_stable_time_step/delta_x**2 * (alpha_interblock(alpha_x[i+1,j], alpha_x[i,j]) * C_t[i+1,j]\n",
" - (alpha_interblock(alpha_x[i+1,j], alpha_x[i,j]) + alpha_interblock(alpha_x[i-1,j], alpha_x[i,j])) * C_t[i,j] \n",
" + alpha_interblock(alpha_x[i-1,j], alpha_x[i,j]) * C_t[i-1,j]) \\\n",
" + max_stable_time_step/delta_y**2 * (alpha_interblock(alpha_y[i,j+1], alpha_y[i,j]) * C_t[i,j+1]\n",
" - (alpha_interblock(alpha_y[i,j+1], alpha_y[i,j]) + alpha_interblock(alpha_y[i,j-1], alpha_y[i,j])) * C_t[i,j] \n",
" + alpha_interblock(alpha_y[i,j-1], alpha_y[i,j]) * C_t[i,j-1])\n",
"\n",
" return C_t1"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[[ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 4.93210938e+02 -6.85390625e+02 5.61357422e+02\n",
" -3.21421875e+02 1.35512695e+02 -4.28554688e+01 1.01855469e+01\n",
" -1.79687500e+00 2.28515625e-01 -1.95312500e-02 9.76562500e-04\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 -6.85390625e+02 9.41660156e+02 -7.56210938e+02\n",
" 4.20468750e+02 -1.70039062e+02 5.07128906e+01 -1.10742188e+01\n",
" 1.71875000e+00 -1.75781250e-01 9.76562500e-03 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 5.61357422e+02 -7.56210938e+02 5.86933594e+02\n",
" -3.10078125e+02 1.16542969e+02 -3.12890625e+01 5.84472656e+00\n",
" -7.03125000e-01 4.39453125e-02 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 -3.21421875e+02 4.20468750e+02 -3.10078125e+02\n",
" 1.51593750e+02 -5.08593750e+01 1.15312500e+01 -1.64062500e+00\n",
" 1.17187500e-01 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 1.35512695e+02 -1.70039062e+02 1.16542969e+02\n",
" -5.08593750e+01 1.43554688e+01 -2.46093750e+00 2.05078125e-01\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 -4.28554688e+01 5.07128906e+01 -3.12890625e+01\n",
" 1.15312500e+01 -2.46093750e+00 2.46093750e-01 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 1.01855469e+01 -1.10742188e+01 5.84472656e+00\n",
" -1.64062500e+00 2.05078125e-01 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 -1.79687500e+00 1.71875000e+00 -7.03125000e-01\n",
" 1.17187500e-01 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 2.28515625e-01 -1.75781250e-01 4.39453125e-02\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 -1.95312500e-02 9.76562500e-03 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 9.76562500e-04 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]\n",
" [ 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00\n",
" 0.00000000e+00 0.00000000e+00 0.00000000e+00 0.00000000e+00]]\n"
]
}
],
"source": [
"for i in range(10):\n",
" C_t = simulate(C_t)\n",
"print(C_t)"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"0.5\n"
]
}
],
"source": [
"print(max_stable_time_step)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.1"
},
"orig_nbformat": 4
},
"nbformat": 4,
"nbformat_minor": 2
}