.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/eigen_solve.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_eigen_solve.py: .. _ref_pymath_eigen_solve: Use PyAnsys Math to solve eigenproblems --------------------------------------- This example uses a verification manual input file, but you can use your own sparse or dense matrices. .. GENERATED FROM PYTHON SOURCE LINES 35-38 Perform required imports and start PyAnsys ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Perform required imports. .. GENERATED FROM PYTHON SOURCE LINES 38-50 .. code-block:: Python import time from ansys.mapdl.core.examples import vmfiles import matplotlib.pylab as plt import numpy as np import ansys.math.core.math as pymath # Start PyAnsys Math as a service. mm = pymath.AnsMath() .. rst-class:: sphx-glr-script-out .. code-block:: none /home/runner/work/pyansys-math/pyansys-math/.venv/lib/python3.10/site-packages/ansys/tools/common/cyberchannel.py:187: UserWarning: Starting gRPC client without TLS on 127.0.0.1:21000. This is INSECURE. Consider using a secure connection. warn(f"Starting gRPC client without TLS on {target}. This is INSECURE. Consider using a secure connection.") .. GENERATED FROM PYTHON SOURCE LINES 51-55 Get matrices ~~~~~~~~~~~~ Run the input file from Verification Manual 153 and then get the stiff (``k``) and mass (``m``) matrices from the FULL file. .. GENERATED FROM PYTHON SOURCE LINES 55-62 .. code-block:: Python out = mm._mapdl.input(vmfiles["vm153"]) fullfile = mm._mapdl.jobname + ".full" k = mm.stiff(fname=fullfile) m = mm.mass(fname=fullfile) .. GENERATED FROM PYTHON SOURCE LINES 63-66 Display size of the matrices ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Display size of the ``k`` and ``m`` matrices. .. GENERATED FROM PYTHON SOURCE LINES 66-70 .. code-block:: Python print(m.shape) print(k.shape) .. rst-class:: sphx-glr-script-out .. code-block:: none (63, 63) (63, 63) .. GENERATED FROM PYTHON SOURCE LINES 71-76 Allocate an array to store eigenshapes ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Allocate an array to store the eigenshapes, where ``nev``` is the number of eigenvalues requested, .. GENERATED FROM PYTHON SOURCE LINES 76-80 .. code-block:: Python nev = 10 a = mm.mat(k.nrow, nev) a .. rst-class:: sphx-glr-script-out .. code-block:: none AnsMath dense matrix (63, 10 .. GENERATED FROM PYTHON SOURCE LINES 81-88 Perform modal analysis ~~~~~~~~~~~~~~~~~~~~~~ Perform the modal analysis. The algorithm is automatically chosen with respect to the properties of the matrices (such as scalar, storage, or symmetry). .. GENERATED FROM PYTHON SOURCE LINES 88-95 .. code-block:: Python print("Calling PyAnsys Math to solve the eigenproblem...") t1 = time.time() ev = mm.eigs(nev, k, m, phi=a) print(f"Elapsed time to solve this problem: {time.time() - t1}") .. rst-class:: sphx-glr-script-out .. code-block:: none Calling PyAnsys Math to solve the eigenproblem... Elapsed time to solve this problem: 0.03886890411376953 .. GENERATED FROM PYTHON SOURCE LINES 96-99 Print eigenfrequencies ~~~~~~~~~~~~~~~~~~~~~~ Print the vector of eigenfrequencies. .. GENERATED FROM PYTHON SOURCE LINES 99-102 .. code-block:: Python print(ev) .. rst-class:: sphx-glr-script-out .. code-block:: none APGGOY : Size : 10 2.113e+02 4.865e+02 7.701e+02 1.065e+03 1.375e+03 < 5 1.517e+03 1.517e+03 1.709e+03 1.747e+03 1.747e+03 < 10 .. GENERATED FROM PYTHON SOURCE LINES 103-109 Verify the accuracy of eigenresults ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Check the residual error for the first eigenresult: :math:`R_1=||(K-\lambda_1.M).\phi_1||_2` First, compute :math:`\lambda_1 = \omega_1^2 = (2.\pi.f_1)^2` .. GENERATED FROM PYTHON SOURCE LINES 109-117 .. code-block:: Python # Eigenfrequency (Hz) i = 0 f = ev[0] omega = 2 * np.pi * f lam = omega * omega .. GENERATED FROM PYTHON SOURCE LINES 118-120 Then get the first eigenshape :math:`\phi_1` and compute both :math:`K.\phi_1` and :math:`M.\phi_1`. .. GENERATED FROM PYTHON SOURCE LINES 120-131 .. code-block:: Python # shape phi = a[0] # APDL command: *MULT,K,,Phi,,KPhi kphi = k.dot(phi) # APDL command: *MULT,M,,Phi,,MPhi mphi = m.dot(phi) .. GENERATED FROM PYTHON SOURCE LINES 132-134 Next, compute the :math:`||K.\phi_1||_2` quantity and normalize the residual value. .. GENERATED FROM PYTHON SOURCE LINES 134-143 .. code-block:: Python # APDL command: *MULT,K,,Phi,,KPhi kphi = k.dot(phi) # APDL command: *NRM,KPhi,NRM2,KPhiNrm kphinrm = kphi.norm() .. GENERATED FROM PYTHON SOURCE LINES 144-147 Add these two vectors, using the :math:`\lambda_1` scalar factor, and compute the normalized residual value: :math:`\frac{R_1}{||K.\phi_1||_2}` .. GENERATED FROM PYTHON SOURCE LINES 147-156 .. code-block:: Python # APDL command: *AXPY,-lambda,,MPhi,1,,KPhi mphi *= lam kphi -= mphi # Compute residual res = kphi.norm() / kphinrm print(res) .. rst-class:: sphx-glr-script-out .. code-block:: none 1.9136679485631705e-10 .. GENERATED FROM PYTHON SOURCE LINES 157-158 Compute this residual for all eigenmodes .. GENERATED FROM PYTHON SOURCE LINES 158-198 .. code-block:: Python def get_res(i): """Compute the residual for a given eigenmode.""" # Eigenfrequency (Hz) f = ev[i] # omega = 2.pi.Frequency omega = 2 * np.pi * f # lambda = omega^2 lam = omega * omega # i-th eigenshape phi = a[i] # K.Phi kphi = k.dot(phi) # M.Phi mphi = m.dot(phi) # Normalize scalar value kphinrm = kphi.norm() # (K-\lambda.M).Phi mphi *= lam kphi -= mphi # Return the residual return kphi.norm() / kphinrm pymath_acc = np.zeros(nev) for i in range(nev): f = ev[i] pymath_acc[i] = get_res(i) print(f"[{i}] : Freq = {f}\t - Residual = {pymath_acc[i]}") .. rst-class:: sphx-glr-script-out .. code-block:: none [0] : Freq = 211.2534447881588 - Residual = 1.9136679485631705e-10 [1] : Freq = 486.50419022721263 - Residual = 1.3403406950622093e-10 [2] : Freq = 770.104536570519 - Residual = 6.739462952718192e-11 [3] : Freq = 1064.5432295081987 - Residual = 3.631667679775435e-11 [4] : Freq = 1374.535393913082 - Residual = 2.426773924787909e-11 [5] : Freq = 1516.8545419017319 - Residual = 3.024157834515717e-11 [6] : Freq = 1516.8545419017403 - Residual = 1.2726936409880203e-10 [7] : Freq = 1708.9655656955329 - Residual = 1.4727869229348076e-09 [8] : Freq = 1747.1764513661403 - Residual = 2.311798816791578e-09 [9] : Freq = 1747.176451366141 - Residual = 2.5094281080563296e-09 .. GENERATED FROM PYTHON SOURCE LINES 199-202 Plot accuracy of eigenresults ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Plot tahe accuracy of the eigenresults. .. GENERATED FROM PYTHON SOURCE LINES 202-214 .. code-block:: Python fig = plt.figure(figsize=(12, 10)) ax = plt.axes() x = np.linspace(1, nev, nev) plt.title("PyAnsys Math Residual Error (%)") plt.yscale("log") plt.ylim([10e-13, 10e-7]) plt.xlabel("Frequency #") plt.ylabel("Errors (%)") ax.bar(x, pymath_acc, label="PyAnsys Math Results", color="orange") plt.show() .. image-sg:: /examples/images/sphx_glr_eigen_solve_001.png :alt: PyAnsys Math Residual Error (%) :srcset: /examples/images/sphx_glr_eigen_solve_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 215-218 Stop PyAnsys Math ~~~~~~~~~~~~~~~~~ Stop PyAnsys Math. .. GENERATED FROM PYTHON SOURCE LINES 218-220 .. code-block:: Python mm._mapdl.exit() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.505 seconds) .. _sphx_glr_download_examples_eigen_solve.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: eigen_solve.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: eigen_solve.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: eigen_solve.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_