-
Notifications
You must be signed in to change notification settings - Fork 0
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Showing
1 changed file
with
113 additions
and
0 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,113 @@ | ||
{ | ||
"nbformat": 4, | ||
"nbformat_minor": 0, | ||
"metadata": { | ||
"colab": { | ||
"private_outputs": true, | ||
"provenance": [], | ||
"authorship_tag": "ABX9TyNYWx+I2cFrt4q7k9DFe8h8", | ||
"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/1DElement.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": "EVc6sM8n_6WC" | ||
}, | ||
"outputs": [], | ||
"source": [ | ||
"import numpy as np\n", | ||
"import matplotlib.pyplot as plt\n", | ||
"\n", | ||
"# Define the function to assemble the global stiffness matrix\n", | ||
"def assemble_stiffness_matrix(nodal_coords, nodal_connectivity, material_props):\n", | ||
" num_nodes = len(nodal_coords)\n", | ||
" K_global = np.zeros((num_nodes, num_nodes))\n", | ||
"\n", | ||
" for element in nodal_connectivity:\n", | ||
" node1, node2 = element\n", | ||
" x1, x2 = nodal_coords[node1], nodal_coords[node2]\n", | ||
" length = x2 - x1\n", | ||
" k = material_props['E'] * material_props['A'] / length\n", | ||
"\n", | ||
" # Element stiffness matrix for 1D linear element\n", | ||
" K_local = (k / length) * np.array([[1, -1], [-1, 1]])\n", | ||
"\n", | ||
" # Assembly into the global stiffness matrix\n", | ||
" K_global[node1:node2+1, node1:node2+1] += K_local\n", | ||
"\n", | ||
" return K_global\n", | ||
"\n", | ||
"# Apply boundary conditions (Dirichlet)\n", | ||
"def apply_boundary_conditions(K, F, bc):\n", | ||
" for node, value in bc.items():\n", | ||
" # Modify stiffness matrix and force vector to enforce boundary condition\n", | ||
" K[node, :] = 0\n", | ||
" K[:, node] = 0\n", | ||
" K[node, node] = 1\n", | ||
" F[node] = value\n", | ||
" return K, F\n", | ||
"\n", | ||
"# Define function to solve FEM problem\n", | ||
"def fem_1d(nodal_coords, nodal_connectivity, material_props, boundary_conditions, load):\n", | ||
" num_nodes = len(nodal_coords)\n", | ||
"\n", | ||
" # Step 1: Assemble global stiffness matrix\n", | ||
" K_global = assemble_stiffness_matrix(nodal_coords, nodal_connectivity, material_props)\n", | ||
"\n", | ||
" # Step 2: Assemble global force vector\n", | ||
" F_global = np.zeros(num_nodes)\n", | ||
" for node, value in load.items():\n", | ||
" F_global[node] = value\n", | ||
"\n", | ||
" # Step 3: Apply boundary conditions\n", | ||
" K_global, F_global = apply_boundary_conditions(K_global, F_global, boundary_conditions)\n", | ||
"\n", | ||
" # Step 4: Solve the system of equations\n", | ||
" displacements = np.linalg.solve(K_global, F_global)\n", | ||
"\n", | ||
" return displacements\n", | ||
"\n", | ||
"# Example Input Data\n", | ||
"nodal_coords = np.array([0.0, 0.5, 1.0]) # Coordinates of nodes\n", | ||
"nodal_connectivity = np.array([[0, 1], [1, 2]]) # Elements defined by node numbers\n", | ||
"material_props = {'E': 210e9, 'A': 1e-4} # Material properties: E (Young's modulus) and A (Cross-sectional area)\n", | ||
"boundary_conditions = {0: 0.0} # Boundary condition: node 0 is fixed (Dirichlet condition)\n", | ||
"load = {2: 1000.0} # Load applied at node 2\n", | ||
"\n", | ||
"# Solve the FEM problem\n", | ||
"displacements = fem_1d(nodal_coords, nodal_connectivity, material_props, boundary_conditions, load)\n", | ||
"\n", | ||
"# Output displacements\n", | ||
"print(\"Nodal Displacements:\", displacements)\n", | ||
"\n", | ||
"# Plot the displacements\n", | ||
"plt.plot(nodal_coords, displacements, '-o')\n", | ||
"plt.xlabel('Position (m)')\n", | ||
"plt.ylabel('Displacement (m)')\n", | ||
"plt.title('1D FEM Nodal Displacements')\n", | ||
"plt.grid(True)\n", | ||
"plt.show()\n" | ||
] | ||
} | ||
] | ||
} |