{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import pyvista\npyvista.OFF_SCREEN = True\npyvista.set_jupyter_backend('static')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Clamped-clamped beam \u2014 first three bending natural frequencies\n\nClassical Bernoulli-beam eigenvalue problem. A prismatic beam\nclamped at both ends has bending natural frequencies (Rao 2017\n\u00a78.5 Table 8.1; Timoshenko 1974 \u00a75.3)\n\n\\begin{align}f_n\n \\;=\\;\n \\frac{(\\beta_n\\, L)^{2}}{2\\pi\\, L^{2}}\n \\sqrt{\\frac{E I}{\\rho A}},\\end{align}\n\nwith the dimensionless eigenvalues $\\beta_n L$ defined as\nthe positive roots of\n\n\\begin{align}1 - \\cos(\\beta L)\\cosh(\\beta L) \\;=\\; 0.\\end{align}\n\nThe first three roots are\n\n\\begin{align}\\beta_1 L = 4.730041, \\quad\n \\beta_2 L = 7.853205, \\quad\n \\beta_3 L = 10.995608.\\end{align}\n\nThese are slightly higher than the cantilever's $\\beta_n\nL \\in \\{1.875, 4.694, 7.855\\}$ (Rao Table 8.1) because clamping\nboth ends raises the bending stiffness \u2014 at any given mode the\nCC beam is a **factor of ~6.3** higher in $\\omega_1^2$\nthan the same-geometry cantilever (\u03b2\u2081\u00b2-ratio\n$(4.730 / 1.875)^{2} \\approx 6.36$).\n\n## Implementation\n\nDrives the existing\n:class:`~femorph_solver.validation.problems.ClampedClampedBeamModes`\nproblem class on a HEX8 enhanced-strain prism mesh; the FE\nspectrum is filtered for transverse-bending modes (UY-dominant)\nand matched to the published $\\beta_n L$ values by\nfrequency-proximity. Higher modes pick up a small Timoshenko\nshear / rotary-inertia correction that grows with mode\nnumber \u2014 visible in the table below as the FE result drifts\n~0.5 % at mode 1 to ~4 % at mode 3.\n\n## References\n\n* Rao, S. S. (2017) *Mechanical Vibrations*, 6th ed., Pearson,\n \u00a78.5 Table 8.1 (clamped-clamped beam eigenvalue roots).\n* Timoshenko, S. P. (1974) *Vibration Problems in Engineering*,\n 4th ed., Wiley, \u00a75.3.\n* Cook, R. D., Malkus, D. S., Plesha, M. E., Witt, R. J.\n (2002) *Concepts and Applications of Finite Element\n Analysis*, 4th ed., Wiley, \u00a711 (modal analysis).\n* Simo, J. C. and Rifai, M. S. (1990) \"A class of mixed\n assumed strain methods \u2026\" *International Journal for\n Numerical Methods in Engineering* 29 (8), 1595\u20131638.\n\n## Vendor cross-references\n\n================================ ===================== ===================================\nSource Reported f_1 [Hz] Problem ID / location\n================================ ===================== ===================================\nClosed form (Euler-Bernoulli) 259.42 Rao 2017 \u00a78.5 Table 8.1\nTimoshenko (1974) \u00a75.3 259.42 clamped-clamped char. eq.\nMAPDL Verification Manual \u2248 259 VM-89 family (clamped variant)\nAbaqus Verification Manual \u2248 259 AVM 1.6.x clamped-beam-NF family\nNAFEMS FV-3 259.42 clamped-clamped beam fundamental + higher modes\n================================ ===================== ===================================\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from __future__ import annotations\n\nimport numpy as np\nimport pyvista as pv\n\nfrom femorph_solver.validation.problems import ClampedClampedBeamModes" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Build the model from the validation problem class\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "problem = ClampedClampedBeamModes()\nm = problem.build_model()\nprint(f\"clamped-clamped beam mesh: {m.grid.n_points} nodes, {m.grid.n_cells} HEX8 cells\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Closed-form references\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "print()\nprint(f\"{'n':>3} {'\u03b2_n L':>10} {'f_pub [Hz]':>12}\")\nprint(f\"{'-' * 3:>3} {'-' * 10:>10} {'-' * 12:>12}\")\nfor n in (1, 2, 3):\n pv_ = [pv0 for pv0 in problem.published_values if pv0.name == f\"f{n}_bending\"][0]\n beta_str = [\"4.7300\", \"7.8532\", \"10.9956\"][n - 1]\n print(f\"{n:>3} {beta_str:>10} {pv_.value:12.3f}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Modal solve + frequency-proximity match\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "res = m.solve_modal(n_modes=problem.default_n_modes)\nfreqs = np.asarray(res.frequency, dtype=np.float64)\nshapes = np.asarray(res.mode_shapes).reshape(-1, 3, freqs.size)\n\nprint()\nprint(f\"{'n':>3} {'k':>3} {'f_FE [Hz]':>12} {'f_pub [Hz]':>12} {'err':>9} {'tol':>5}\")\nprint(f\"{'-' * 3:>3} {'-' * 3:>3} {'-' * 12:>12} {'-' * 12:>12} {'-' * 9:>9} {'-' * 5:>5}\")\n\nclaimed = np.zeros(freqs.size, dtype=bool)\nfor n in (1, 2, 3):\n pv_ = [pv0 for pv0 in problem.published_values if pv0.name == f\"f{n}_bending\"][0]\n f_pub = pv_.value\n candidates = [k for k in range(freqs.size) if not claimed[k]]\n k_best = min(candidates, key=lambda k: abs(freqs[k] - f_pub))\n f_fe = float(freqs[k_best])\n err_pct = (f_fe - f_pub) / f_pub * 100.0\n tol_pct = pv_.tolerance * 100.0\n print(f\"{n:>3} {k_best:>3} {f_fe:12.3f} {f_pub:12.3f} {err_pct:+8.3f}% \u2264{tol_pct:.1f}%\")\n claimed[k_best] = True\n assert abs(err_pct) < tol_pct, f\"f_{n} deviation {err_pct:.2f}% exceeds {tol_pct}%\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Render the first three mode shapes side-by-side\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Re-acquire the matching mode indices for rendering.\nclaimed = np.zeros(freqs.size, dtype=bool)\npanels: list[tuple[int, int, float]] = []\nfor n in (1, 2, 3):\n pv_ = [p for p in problem.published_values if p.name == f\"f{n}_bending\"][0]\n candidates = [k for k in range(freqs.size) if not claimed[k]]\n k_best = min(candidates, key=lambda k: abs(freqs[k] - pv_.value))\n panels.append((n, k_best, float(freqs[k_best])))\n claimed[k_best] = True\n\nplotter = pv.Plotter(shape=(1, 3), off_screen=True, window_size=(1080, 320))\nfor n, k, f in panels:\n plotter.subplot(0, n - 1)\n phi = shapes[:, :, k]\n warp = m.grid.copy()\n warp.points = m.grid.points + (0.05 / np.max(np.abs(phi))) * phi\n warp[\"|phi|\"] = np.linalg.norm(phi, axis=1)\n plotter.add_mesh(m.grid, style=\"wireframe\", color=\"grey\", opacity=0.3)\n plotter.add_mesh(warp, scalars=\"|phi|\", cmap=\"plasma\", show_edges=False, show_scalar_bar=False)\n plotter.add_text(f\"mode {n}: f = {f:.0f} Hz\", position=\"upper_edge\", font_size=10)\nplotter.link_views()\nplotter.view_xy()\nplotter.show()" ] } ], "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.12.3" } }, "nbformat": 4, "nbformat_minor": 0 }