# Path: sucemfem/Consts.py
# Z0 = mu0*c0                              


import sys
import numpy as N
import os
import dolfin as dol
import postproc_eigres
from dolfin import dot, curl, inner, dx
from scipy.sparse.linalg.eigen.arpack import speigs
from sucemfem.Consts import eps0, mu0, c0
from material_properties import MaterialProperties
lu = dol.LUSolver(smat)
lu.parameters["reuse_factorization"] = True
lu.parameters["report"] = False
bb = dol.Vector(M.size(0))
xx = dol.Vector(M.size(0))
def sigma_solve(b):
    bb[:] = b
    lu.solve(xx, bb)
    return xx[:]
M_matvec = lambda x: M*x

arpack_eigs,arpack_v = speigs.ARPACK_gen_eigs(M_matvec, sigma_solve, M.size(0), sigma, 51, ncv=91)

# Create eigensolver
# esolver = dol.SLEPcEigenSolver(S,M)
# esolver.parameters["spectrum"] = "smallest real"
# esolver.parameters["solver"] = "arnoldi"
# esolver.parameters["spectral_shift"] = sigma
# esolver.parameters["spectral_transform"] = "shift-and-invert"
# esolver.solve()

# eigs = [esolver.get_eigenvalue(i)[0] for i in range(4)]
# filtered_eigs = []
# for i in range(S.size(0)):
#     ev_r, ev_i = esolver.get_eigenvalue(i)
#     if ev_r > 0.001:
#         filtered_eigs.append(ev_r)

res = N.array(sorted(arpack_eigs)[0:10])
print N.sqrt(res)

# Path: sucemfem/Consts.py
# Z0 = mu0*c0                              


# Path: sandbox/eig/eig_materials.py
import sys
import numpy as N
import os
import dolfin as dol
import postproc_eigres
from dolfin import dot, curl, inner, dx
from scipy.sparse.linalg.eigen.arpack import speigs
from sucemfem.Consts import eps0, mu0, c0
from material_properties import MaterialProperties

lu = dol.LUSolver(smat)
lu.parameters["reuse_factorization"] = True
lu.parameters["report"] = False
bb = dol.Vector(M.size(0))
xx = dol.Vector(M.size(0))
def sigma_solve(b):
    bb[:] = b
    lu.solve(xx, bb)
    return xx[:]
M_matvec = lambda x: M*x

arpack_eigs,arpack_v = speigs.ARPACK_gen_eigs(M_matvec, sigma_solve, M.size(0), sigma, 51, ncv=91)

# Create eigensolver
# esolver = dol.SLEPcEigenSolver(S,M)
# esolver.parameters["spectrum"] = "smallest real"
# esolver.parameters["solver"] = "arnoldi"
# esolver.parameters["spectral_shift"] = sigma
# esolver.parameters["spectral_transform"] = "shift-and-invert"
# esolver.solve()

# eigs = [esolver.get_eigenvalue(i)[0] for i in range(4)]
# filtered_eigs = []
# for i in range(S.size(0)):
#     ev_r, ev_i = esolver.get_eigenvalue(i)
#     if ev_r > 0.001:
#         filtered_eigs.append(ev_r)

res = N.array(sorted(arpack_eigs)[0:10])
print N.sqrt(res)

print c0*N.sqrt(res)/2/N.pi/1e6