@testset ExtendedTestSet "Laplace Eigen" begin
  function residual(cell, u_el)
    @unpack X, N, ∇N_X, JxW = cell
    ∇u_q = ∇N_X' * u_el[1, :]
    R_q = (∇N_X * ∇u_q)'
    return JxW * R_q[:]
  end
  
  function stiffness(cell, _)
    @unpack X, N, ∇N_X, JxW = cell
    K_q = ∇N_X * ∇N_X'
    return JxW * K_q
  end
  
  function mass(cell, _)
    @unpack X, N, ∇N_X, JxW = cell
    M_q = N * N'
    return JxW * M_q
  end

  q_degree = 2
  mesh = read_mesh("./laplace_eigen/laplace_eigen.g", Int[])
  bcs = EssentialBC[]
  fspaces, dof, asm = container_setup(mesh, [1], q_degree, 1; is_dynamic=true)
  λ, X = simple_eigen_solver(mesh, fspaces, dof, asm, bcs, residual, stiffness, mass)
  simple_eigen_post_processor("./laplace_eigen/laplace_eigen.g", λ, X, "u")
  # @test exodiff("./laplace_eigen/laplace_eigen.gold", "./laplace_eigen/laplace_eigen.e")
  # exo_1 = ExodusDatabase("./laplace_eigen/laplace_eigen.gold", "r")
  # times_1 = read_times(exo_1)
  # close(exo_1)
  # exo_2 = ExodusDatabase("./laplace_eigen/laplace_eigen.e", "r")
  # times_2 = read_times(exo_2)
  # close(exo_2)
  # display(times_1)
  # display(times_2)
  # @test isapprox(times_1, times_2, atol=1e-6)
  
  rm("./laplace_eigen/laplace_eigen.e")
end