# Eigen solver syntax in dolfin-x

**URL:** <https://fenicsproject.discourse.group/t/eigen-solver-syntax-in-dolfin-x/7369>\
**Category:** Uncategorized\
**Created:** [December 23, 2021, 4:22pm UTC](https://fenicsproject.discourse.group/t/eigen-solver-syntax-in-dolfin-x/7369 "2021-12-23T16:22:06Z")\
**Posts on this page:** 1\
**Showing post:** 5

<div class="post-metadata">

**Author:** ![Shubham\_Saurabh](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/shubham_saurabh/32/3460_2.png) [@Shubham\_Saurabh](https://fenicsproject.discourse.group/u/Shubham_Saurabh)\
**Post date:** [December 27, 2021, 11:55am UTC](https://fenicsproject.discourse.group/t/eigen-solver-syntax-in-dolfin-x/7369/5 "2021-12-27T11:55:49Z")

</div>

Hello sir,  
I am attaching the `dolfin` code for the same problem as you have asked for it. I have to replicate the same using `dolfin-x`.

```auto
from dolfin import *
import numpy as np
import json
import matplotlib.pyplot as plt
import scipy as sp
%matplotlib inline
mesh = BoxMesh(Point(0.,0.,0.),Point(5,0.6,0.3), 25, 3, 2)
E, nu, rho = (2e11), (0.3), (7850)
mu,lmbda = E/2./(1+nu), E*nu/(1+nu)/(1-2*nu)
def sigma(v):
    dim = v.geometric_dimension()
    return 2.0*mu*eps(v) + lmbda*tr(eps(v))*Identity(dim)
def eps(v):
    return sym(grad(v))
V = VectorFunctionSpace(mesh, 'CG', degree=2)
u = TrialFunction(V)
v = TestFunction(V)
left = CompiledSubDomain("near(x[0], 0) && on_boundary")
bc = DirichletBC(V, Constant((0.,0.,0.)), left)
k_form = inner(sigma(u),eps(v))*dx
m_form = rho*dot(u,v)*dx
K = PETScMatrix()
assemble(k_form, tensor=K)
M = PETScMatrix()
assemble(m_form, tensor=M)
bc.zero(M), bc.apply(K)
eigensolver = SLEPcEigenSolver(K, M)
eigensolver.parameters['problem_type'] = 'gen_hermitian'
eigensolver.parameters['spectral_transform'] = 'shift-and-invert'
eigensolver.parameters['spectral_shift'] = 0.

N_eig = 40 # number of eigenvalues
print("Computing {} first eigenvalues...".format(N_eig))
eigensolver.solve(N_eig)

L, B, H = 5., 0.3, 0.6 #Dimensions of the beam
# Exact solution computation
from scipy.optimize import root
from math import cos, cosh
falpha = lambda x: cos(x)*cosh(x)+1
alpha = lambda n: root(falpha, (2*n+1)*pi/2.)['x'][0]
file_results = XDMFFile("modal_analysis.xdmf")
file_results.parameters["flush_output"] = True
file_results.parameters["functions_share_mesh"] = True

# Extraction
eig_v = []
eig_l = []
for i in range(N_eig):
    # Extract eigenpair
    r, c, rx, cx = eigensolver.get_eigenpair(i)

    # 3D eigenfrequency
    freq_3D = sqrt(r)/2/pi

    # Beam eigenfrequency
    if i % 2 == 0: # exact solution should correspond to weak axis bending
        I_bend = H*B**3/12.
    else: #exact solution should correspond to strong axis bending
        I_bend = B*H**3/12.
    freq_beam = alpha(i/2)**2*sqrt(E*I_bend/(rho*B*H*L**4))/2/pi

    print("Solid FE: {0:8.5f} [Hz] Beam theory: {1:8.5f} [Hz]".format(freq_3D, freq_beam))
    
    
    # Initialize function and assign eigenvector
    eigenmode = Function(V,name="Eigenvector "+str(i))
    eigenmode.vector()[:] = rx
    file_results.write(eigenmode, 0)
    eig_l.append(freq_3D)
    eig_v.append(rx)

```

The results are:

```auto
Computing 40 first eigenvalues...
Solid FE: 9.83051 [Hz] Beam theory: 9.78457 [Hz]
Solid FE: 19.45146 [Hz] Beam theory: 19.56914 [Hz]
Solid FE: 60.61639 [Hz] Beam theory: 61.31884 [Hz]
Solid FE: 114.64296 [Hz] Beam theory: 122.63769 [Hz]
Solid FE: 118.82391 [Hz] Beam theory: 171.69454 [Hz]
Solid FE: 165.68699 [Hz] Beam theory: 343.38908 [Hz]

```

---

_[View the full topic](https://fenicsproject.discourse.group/t/eigen-solver-syntax-in-dolfin-x/7369)._