Skip to main content Skills Marktplatz Entdecken und erkunden Sie KI-Skills, die von der Community erstellt wurden.
Mit Codex oder Claude installieren Kopieren Sie diesen Prompt, fügen Sie ihn in Codex, Claude oder einen anderen Assistant ein und lassen Sie die Skill-Seite prüfen und installieren.
Prompt-Details anzeigen Ein direkter Befehl überspringt den Prüf-Prompt. Prüfen Sie die Quelle, bevor Sie ihn ausführen.
npx skills add https://github.com/xjtulyc/awesome-rosetta-skills --skill finite-element-analysisDer Befehl bleibt in einer Zeile. Scrollen Sie horizontal, um ihn vor dem Kopieren vollständig zu prüfen.
Sie bevorzugen eine lokale Kopie? Laden Sie die Dateien herunter, die SkillsMP derzeit vorliegen.
ZIP herunterladen Herunterladen... name finite-element-analysis description Use this Skill for FEA with FEniCSx or scikit-fem: mesh generation, boundary conditions, linear elasticity, heat conduction, and result visualization.
tags ["engineering","finite-element","fenics","pde","structural-analysis"] version 1.0.0 authors [{"name":"Rosetta Skills Contributors","github":"@xjtulyc"}] license MIT platforms ["claude-code","codex","gemini-cli","cursor"] dependencies {"python":["scikit-fem>=9.0","numpy>=1.24","scipy>=1.11","matplotlib>=3.7"]} last_updated 2026-03-17 status stable
Finite Element Analysis
One-line summary : Solve PDEs numerically using the Finite Element Method with scikit-fem: structural mechanics, heat conduction, and eigenvalue problems on 2D/3D meshes.
When to Use This Skill
When solving structural mechanics problems (stress, displacement, strain)
When computing temperature distributions in heat conduction problems
When solving Poisson, Laplace, or Helmholtz equations on complex geometries
When performing modal analysis (eigenfrequencies and mode shapes)
When working with boundary conditions: Dirichlet, Neumann, Robin
When post-processing FEA results (von Mises stress, heat flux)
Trigger keywords : finite element, FEA, FEM, FEniCS, scikit-fem, elasticity, heat conduction, PDE, mesh, stiffness matrix, boundary condition, Galerkin, stress analysis
Background & Key Concepts
Weak Formulation
For a PDE $\mathcal{L}u = f$ on domain $\Omega$, the weak form seeks $u \in V$:
$$
a(u, v) = L(v) \quad \forall v \in V
$$
where $a(\cdot,\cdot)$ is the bilinear form and $L(\cdot)$ is the linear functional.
Linear Elasticity
Equilibrium: $\nabla \cdot \boldsymbol{\sigma} + \mathbf{f} = \mathbf{0}$
Constitutive law (plane stress):
$$
\boldsymbol{\sigma} = \frac{E}{1-\nu^2}\begin{pmatrix} 1 & \nu & 0 \ \nu & 1 & 0 \ 0 & 0 & \frac{1-\nu}{2} \end{pmatrix} \boldsymbol{\varepsilon}
$$
Heat Equation (Steady-State)
$$
-\nabla \cdot (k \nabla T) = q \quad \text{in } \Omega
$$
Boundary conditions: Dirichlet $T=T_0$ on $\Gamma_D$, Neumann $k\partial T/\partial n = h(T-T_\infty)$ on $\Gamma_N$.
Environment Setup
Install Dependencies
pip install scikit-fem>=9.0 numpy>=1.24 scipy>=1.11 matplotlib>=3.7
pip install meshio
Verify Installation
import skfem
import numpy as np
from skfem import MeshTri, Basis, ElementTriP1
from skfem.models.poisson import laplace, mass
mesh = MeshTri.init_sqsymmetric(3 )
print
f"scikit-fem version: {skfem.__version__} "
print
f"Mesh: {mesh.nvertices} nodes, {mesh.nelements} elements"
Core Workflow
Step 1: Poisson Equation on Unit Square import numpy as np
import matplotlib.pyplot as plt
from skfem import MeshTri, Basis, ElementTriP1, BilinearForm, LinearForm
from skfem import enforce, solve
from skfem.helpers import grad, dot
mesh = MeshTri.init_sqsymmetric(6 )
mesh = mesh.refined(3 )
basis = Basis(mesh, ElementTriP1())
@BilinearForm
def stiffness (u, v, w ):
return dot(grad(u), grad(v))
@LinearForm
def load (v, w ):
x, y = w.x
f = 2 * np.pi**2 * np.sin(np.pi * x) * np.sin(np.pi * y)
return f * v
K = stiffness.assemble(basis)
F = load.assemble(basis)
boundary_dofs = basis.get_dofs().all ()
K, F = enforce(K, F, D=boundary_dofs)
u = solve(K, F)
x_nodes, y_nodes = mesh.p[0 ], mesh.p[1 ]
u_exact = np.sin(np.pi * x_nodes) * np.sin(np.pi * y_nodes)
L2_error = np.sqrt(np.mean((u - u_exact)**2 ))
print (f"L² error: {L2_error:.2 e} " )
fig, axes = plt.subplots(1 , 3 , figsize=(15 , 4 ))
im0 = axes[0 ].tricontourf(x_nodes, y_nodes, mesh.t.T, u, levels=20 , cmap='viridis' )
plt.colorbar(im0, ax=axes[0 ])
axes[0 ].set_title("FEM Solution u" )
axes[0 ].set_aspect('equal' )
im1 = axes[1 ].tricontourf(x_nodes, y_nodes, mesh.t.T, u_exact, levels=20 , cmap='viridis' )
plt.colorbar(im1, ax=axes[1 ])
axes[1 ].set_title("Exact Solution" )
axes[1 ].set_aspect('equal' )
error = np.abs (u - u_exact)
im2 = axes[2 ].tricontourf(x_nodes, y_nodes, mesh.t.T, error, levels=20 , cmap='hot_r' )
plt.colorbar(im2, ax=axes[2 ])
axes[2 ].set_title(f"|Error| (L²={L2_error:.2 e} )" )
axes[2 ].set_aspect('equal' )
plt.suptitle("Poisson Equation FEM Solution on Unit Square" )
plt.tight_layout()
plt.savefig("poisson_fem.png" , dpi=150 )
plt.show()
Step 2: Steady-State Heat Conduction import numpy as np
import matplotlib.pyplot as plt
from skfem import MeshTri, Basis, ElementTriP1, BilinearForm, LinearForm
from skfem import enforce, solve
from skfem.helpers import grad, dot
mesh = MeshTri.init_lshaped()
mesh = mesh.refined(3 )
basis = Basis(mesh, ElementTriP1())
k_cond = 50.0
h_conv = 25.0
T_inf = 20.0
q_vol = 1000.0
@BilinearForm
def conduction (u, v, w ):
return k_cond * dot(grad(u), grad(v))
@BilinearForm
def convection_bc (u, v, w ):
return h_conv * u * v
@LinearForm
def heat_source (v, w ):
return q_vol * v
@LinearForm
def convection_load (v, w ):
return h_conv * T_inf * v
K = conduction.assemble(basis)
F = heat_source.assemble(basis)
facets_top = mesh.facets_satisfying(lambda x: np.isclose(x[1 ], 1.0 ))
facets_right = mesh.facets_satisfying(lambda x: np.isclose(x[0 ], 1.0 ))
conv_facets = np.union1d(facets_top, facets_right)
if len (conv_facets) > 0 :
facet_basis = basis.boundary(facets=conv_facets)
K += convection_bc.assemble(facet_basis)
F += convection_load.assemble(facet_basis)
dofs_left = basis.get_dofs(mesh.facets_satisfying(lambda x: np.isclose(x[0 ], 0.0 ))).all ()
dofs_bottom = basis.get_dofs(mesh.facets_satisfying(lambda x: np.isclose(x[1 ], 0.0 ))).all ()
dofs_fixed = np.union1d(dofs_left, dofs_bottom)
K, F = enforce(K, F, D=dofs_fixed, overwrite=True )
T = solve(K, F)
print (f"Temperature range: {T.min ():.1 f} – {T.max ():.1 f} °C" )
print (f"Max temperature: {T.max ():.1 f} °C at node {T.argmax()} " )
fig, ax = plt.subplots(figsize=(7 , 6 ))
x_n, y_n = mesh.p[0 ], mesh.p[1 ]
im = ax.tricontourf(x_n, y_n, mesh.t.T, T, levels=30 , cmap='hot' )
plt.colorbar(im, ax=ax, label='Temperature (°C)' )
ax.triplot(x_n, y_n, mesh.t.T, color='gray' , linewidth=0.2 , alpha=0.3 )
ax.set_title("Steady-State Heat Conduction — L-shaped Domain" )
ax.set_aspect('equal' )
plt.tight_layout()
plt.savefig("heat_conduction.png" , dpi=150 )
plt.show()
Step 3: Plane-Stress Linear Elasticity import numpy as np
import matplotlib.pyplot as plt
from skfem import MeshTri, Basis, ElementTriP2, BilinearForm, LinearForm
from skfem import enforce, solve
from skfem.helpers import sym_grad, eye, transpose, dd, ddot, trace
mesh = MeshTri.init_symmetric()
mesh = mesh.refined(3 )
mesh = mesh.with_defaults()
E = 200e9
nu = 0.30
P = 1e6
lam = E * nu / ((1 + nu) * (1 - 2 *nu))
mu = E / (2 * (1 + nu))
@BilinearForm
def elasticity (u, v, w ):
eps_u = sym_grad(u)
eps_v = sym_grad(v)
return (lam * trace(eps_u) * trace(eps_v)
+ 2 * mu * ddot(eps_u, eps_v))
basis = Basis(mesh, ElementTriP2())
K = elasticity.assemble(basis)
F = np.zeros(K.shape[0 ])
facets_right = mesh.facets_satisfying(lambda x: np.isclose(x[0 ], 1.0 ))
if len (facets_right) > 0 :
facet_basis = basis.boundary(facets=facets_right)
@LinearForm
def traction_x (v, w ):
return P * v[0 ]
F += traction_x.assemble(facet_basis)
dofs_left = basis.get_dofs(mesh.facets_satisfying(lambda x: np.isclose(x[0 ], 0.0 )))
K, F = enforce(K, F, D=dofs_left.all (), overwrite=True )
u_vec = solve(K, F)
uh = basis.interpolate(u_vec)
u_x = u_vec[basis.nodal_dofs[0 ]]
u_y = u_vec[basis.nodal_dofs[1 ]]
print (f"Max displacement x: {u_x.max ()*1e3 :.4 f} mm" )
print (f"Max displacement y: {u_y.max ()*1e3 :.4 f} mm" )
u_analytical = P * 1.0 / E
print (f"Analytical u_x max: {u_analytical*1e3 :.4 f} mm" )
fig, axes = plt.subplots(1 , 2 , figsize=(12 , 4 ))
x_n, y_n = mesh.p[0 ], mesh.p[1 ]
for ax, data, title, unit in zip (
axes,
[u_x * 1e6 , u_y * 1e6 ],
["Displacement u_x" , "Displacement u_y" ],
["μm" , "μm" ]
):
im = ax.tricontourf(x_n, y_n, mesh.t.T, data, levels=20 , cmap='RdBu_r' )
plt.colorbar(im, ax=ax, label=f"Displacement ({unit} )" )
ax.set_title(title); ax.set_aspect('equal' )
ax.set_xlabel("x (m)" ); ax.set_ylabel("y (m)" )
plt.suptitle("Plane-Stress Linear Elasticity — Uniaxial Tension" )
plt.tight_layout()
plt.savefig("elasticity_fem.png" , dpi=150 )
plt.show()
Advanced Usage
Modal Analysis (Eigenvalue Problem) import numpy as np
import matplotlib.pyplot as plt
from skfem import MeshTri, Basis, ElementTriP2, BilinearForm
from skfem import enforce
from skfem.helpers import sym_grad, trace, ddot
from scipy.sparse.linalg import eigsh
mesh = MeshTri.init_sqsymmetric(4 )
mesh = mesh.refined(2 )
basis = Basis(mesh, ElementTriP2())
E, nu, rho = 70e9 , 0.33 , 2700.0
lam = E * nu / ((1 + nu) * (1 - 2 *nu))
mu = E / (2 * (1 + nu))
@BilinearForm
def stiffness (u, v, w ):
eps_u = sym_grad(u)
eps_v = sym_grad(v)
return lam * trace(eps_u) * trace(eps_v) + 2 *mu * ddot(eps_u, eps_v)
@BilinearForm
def mass_matrix (u, v, w ):
return rho * (u[0 ]*v[0 ] + u[1 ]*v[1 ])
K = stiffness.assemble(basis)
M = mass_matrix.assemble(basis)
boundary_dofs = basis.get_dofs().all ()
from skfem import enforce
K_free, _ = enforce(K.copy(), np.zeros(K.shape[0 ]), D=boundary_dofs, overwrite=True )
n_modes = 6
eigenvalues, eigenvectors = eigsh(K_free, M=M, k=n_modes, sigma=0 , which='LM' )
eigenvalues = np.abs (eigenvalues)
frequencies_hz = np.sqrt(eigenvalues) / (2 * np.pi)
print ("Natural frequencies (clamped plate):" )
for i, f in enumerate (frequencies_hz):
print (f" Mode {i+1 } : {f:.2 f} Hz" )
fig, axes = plt.subplots(2 , 2 , figsize=(12 , 8 ))
x_n, y_n = mesh.p[0 ], mesh.p[1 ]
for idx, ax in enumerate (axes.flat):
mode_x = eigenvectors[basis.nodal_dofs[0 ], idx]
im = ax.tricontourf(x_n, y_n, mesh.t.T, mode_x, levels=20 , cmap='RdBu' )
plt.colorbar(im, ax=ax)
ax.set_title(f"Mode {idx+1 } : f={frequencies_hz[idx]:.1 f} Hz" )
ax.set_aspect('equal' )
plt.suptitle("Modal Analysis — Clamped Square Plate" )
plt.tight_layout()
plt.savefig("modal_analysis.png" , dpi=150 )
plt.show()
Troubleshooting
Error: solve returns NaN or very large values Cause : Singular stiffness matrix — missing or insufficient Dirichlet BCs (rigid body modes).
from scipy.sparse.linalg import norm as spnorm
import numpy as np
boundary_dofs = basis.get_dofs().all ()
print (f"Fixed DOFs: {len (boundary_dofs)} " )
Error: enforce changes matrix shape unexpectedly Fix : Use overwrite=True to modify in-place, or capture the returned (K, F):
K, F = enforce(K, F, D=dofs, overwrite=True )
Poor convergence / large errors Cause : Coarse mesh or low-order elements.
mesh = mesh.refined(4 )
basis = Basis(mesh, ElementTriP2())
Version Compatibility Package Tested versions Notes scikit-fem 9.x API stable; 9.0 renamed some helpers scipy 1.11, 1.12 eigsh requires CSR sparse formatnumpy 1.24, 1.26 No known issues
External Resources
Official Documentation
Key Textbooks
Zienkiewicz, O.C. et al. (2005). The Finite Element Method (3 volumes). Elsevier.
Langtangen, H.P. & Logg, A. (2016). Solving PDEs in Python . Springer Open.
Examples
Example 1: Convergence Study — Mesh Refinement import numpy as np
import matplotlib.pyplot as plt
from skfem import MeshTri, Basis, ElementTriP1, BilinearForm, LinearForm
from skfem import enforce, solve
from skfem.helpers import grad, dot
def solve_poisson_refined (n_refinements ):
"""Return L² error for n levels of uniform mesh refinement."""
mesh = MeshTri.init_sqsymmetric(2 )
for _ in range (n_refinements):
mesh = mesh.refined()
basis = Basis(mesh, ElementTriP1())
@BilinearForm
def a (u, v, w ): return dot(grad(u), grad(v))
@LinearForm
def L (v, w ):
x, y = w.x
return 2 * np.pi**2 * np.sin(np.pi*x) * np.sin(np.pi*y) * v
K = a.assemble(basis)
F = L.assemble(basis)
K, F = enforce(K, F, D=basis.get_dofs().all ())
u = solve(K, F)
x_n, y_n = mesh.p[0 ], mesh.p[1 ]
u_ex = np.sin(np.pi * x_n) * np.sin(np.pi * y_n)
h = 1.0 / (2 **n_refinements * 2 )
L2_err = np.sqrt(np.mean((u - u_ex)**2 ))
return h, L2_err, mesh.nvertices
refinements = range (1 , 7 )
results = [solve_poisson_refined(n) for n in refinements]
h_vals = [r[0 ] for r in results]
err_vals = [r[1 ] for r in results]
nodes = [r[2 ] for r in results]
rates = [np.log(err_vals[i]/err_vals[i-1 ]) / np.log(h_vals[i]/h_vals[i-1 ])
for i in range (1 , len (err_vals))]
print ("Convergence rates:" , [f"{r:.2 f} " for r in rates])
fig, ax = plt.subplots(figsize=(7 , 5 ))
ax.loglog(h_vals, err_vals, 'bo-' , linewidth=2 , markersize=8 , label='FEM P1 error' )
h_ref = np.array([h_vals[0 ], h_vals[-1 ]])
ax.loglog(h_ref, 0.5 * h_ref**2 , 'r--' , label='O(h²) reference' )
ax.set_xlabel("Mesh size h" ); ax.set_ylabel("L² error" )
ax.set_title("FEM Convergence Study — Poisson Equation (P1 Elements)" )
ax.legend(); ax.grid(True , which='both' , alpha=0.3 )
plt.tight_layout()
plt.savefig("fem_convergence.png" , dpi=150 )
plt.show()
Example 2: Thermal Stress Analysis import numpy as np
E = 70e9
nu = 0.33
alpha_T = 23e-6
delta_T = 100.0
eps_thermal = alpha_T * delta_T
print (f"Free thermal strain: {eps_thermal:.4 f} " )
sigma_thermal = -E * alpha_T * delta_T / (1 - nu)
print (f"Thermal stress (fully constrained): {sigma_thermal/1e6 :.1 f} MPa" )
sigma_y = 270e6
print (f"Yield strength: {sigma_y/1e6 :.0 f} MPa" )
print (f"Safety factor: {sigma_y / abs (sigma_thermal):.2 f} " )
Last updated: 2026-03-17 | Maintainer: @xjtulyc
Issues: GitHub Issues
Mehr aus diesem Repository
Verwandte Berufe SOC
Basierend auf der SOC-Berufsklassifikation
