Skip to content

Commit

Permalink
Created using Colab
Browse files Browse the repository at this point in the history
  • Loading branch information
GeoCalcul committed Nov 20, 2024
1 parent 7388fb8 commit b5d2185
Showing 1 changed file with 103 additions and 0 deletions.
103 changes: 103 additions & 0 deletions The coupled equation in the time domain.ipynb
Original file line number Diff line number Diff line change
@@ -0,0 +1,103 @@
{
"nbformat": 4,
"nbformat_minor": 0,
"metadata": {
"colab": {
"private_outputs": true,
"provenance": [],
"authorship_tag": "ABX9TyN0cY6Dqpal4FkyyfXXQdFk",
"include_colab_link": true
},
"kernelspec": {
"name": "python3",
"display_name": "Python 3"
},
"language_info": {
"name": "python"
}
},
"cells": [
{
"cell_type": "markdown",
"metadata": {
"id": "view-in-github",
"colab_type": "text"
},
"source": [
"<a href=\"https://colab.research.google.com/github/GEORMC/Nnumerical_Methods_Course/blob/main/The%20coupled%20equation%20in%20the%20time%20domain.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"id": "boju8AYKJmLI"
},
"outputs": [],
"source": [
"import numpy as np\n",
"\n",
"# Define parameters\n",
"theta = 0.5 # Time-stepping parameter\n",
"delta_t = 0.01 # Time step size\n",
"num_steps = 100 # Number of time steps\n",
"\n",
"# Example matrices (replace with your specific problem matrices)\n",
"n_dof = 3 # Degrees of freedom (example size)\n",
"K_e = np.array([[10, -2, 0], [-2, 10, -2], [0, -2, 10]]) # Example stiffness matrix\n",
"Q = np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]]) # Coupling matrix\n",
"S = np.eye(n_dof) # Storage matrix\n",
"H = 0.1 * np.eye(n_dof) # Hydraulic conductivity matrix\n",
"\n",
"# Initial conditions\n",
"u_n = np.zeros(n_dof) # Initial displacement\n",
"p_n = np.zeros(n_dof) # Initial pore pressure\n",
"f_u = np.zeros(n_dof) # External force\n",
"f_p = np.zeros(n_dof) # Fluid source\n",
"\n",
"# Time-stepping loop\n",
"u_history = [u_n] # Displacement history\n",
"p_history = [p_n] # Pressure history\n",
"\n",
"# Assemble constant matrices\n",
"LHS = np.block([\n",
" [theta * K_e, -theta * Q],\n",
" [Q.T, S + delta_t * theta * H]\n",
"])\n",
"\n",
"for step in range(num_steps):\n",
" # Calculate RHS\n",
" RHS_u = (1 - theta) * K_e @ u_n + f_u\n",
" RHS_p = S @ p_n - (1 - theta) * delta_t * H @ p_n + delta_t * f_p\n",
" RHS = np.hstack([RHS_u, RHS_p])\n",
"\n",
" # Solve for next step\n",
" solution = np.linalg.solve(LHS, RHS)\n",
" u_n1, p_n1 = solution[:n_dof], solution[n_dof:]\n",
"\n",
" # Update history\n",
" u_history.append(u_n1)\n",
" p_history.append(p_n1)\n",
"\n",
" # Update for next iteration\n",
" u_n, p_n = u_n1, p_n1\n",
"\n",
"# Post-process results\n",
"u_history = np.array(u_history)\n",
"p_history = np.array(p_history)\n",
"\n",
"# Example: plot results\n",
"import matplotlib.pyplot as plt\n",
"\n",
"plt.figure(figsize=(12, 6))\n",
"plt.plot(u_history, label=\"Displacement (u)\")\n",
"plt.plot(p_history, label=\"Pore Pressure (p)\", linestyle=\"--\")\n",
"plt.xlabel(\"Time Step\")\n",
"plt.ylabel(\"Value\")\n",
"plt.legend()\n",
"plt.title(\"Coupled System Solution\")\n",
"plt.show()\n"
]
}
]
}

0 comments on commit b5d2185

Please sign in to comment.