Hi!

I am working on a program that can calculate the induced eddy current given a periodic changing magnetic field in the frequency domain.

For the domain where eddy durrents can occur I use the following finite element space:

fes = HCurl(mesh, order = 3, nograds = True, complex = True, definedon = “iron”, dirichlet = “gamma_iron”)

I need the dirichlet boundary condition to make sure that the eddy current cannot escape the iron domain.

When using the inverse or the CG solver the system of equations converges and the result is physically correct.

To make sure the result is correct I want to calculate the residual of the problem: r = a*sol - f

r := residual, a := bilinear form, sol := solution, f := linear form

But the residual is not as low as expected although the solver says otherwise.

This only occurs in this finite element space when there is a dirichlet boundary condition.

I tried deleting all non-free dofs from the matrix and vectors using fes.FreeDofs(True) but that did not change the result.

Is there a way to calculate the residual just like the CG solver does?

Here is a short summary of all settings I set:

fes = HCurl(mesh, order = 3, nograds = True, complex = True, definedon = “iron”, dirichlet = “gamma_iron”)

a = BilinearForm(fes, symmetric = True)

a += SymbolicBFI(…)

a.Assemble()

a_inv = a.mat.Inverse(freedofs = fes.FreeDofs())

f = LinearForm(fes)

f += SymbolicLFI(…, definedon = mesh.Materials(“iron”))

f.Assemble()

sol = GridFunction(fes)

sol.vec.data = a.inv*f.vec
res = GridFunction(fes)
res.vec.data = a.mat*sol.vec - f.vec