Maxwell Stress Tensor and force calculation

Hello,
I achieved to calculate the Force density related to Maxwell Stress Tensor (MST) with two different ways.
The first way is to use the Div function :

sigma_cf = 1/(mu0 * mur)  * (OuterProduct(B, B) - 0.5 * (B*B) * Id(2) )

fes = HDivDiv(mesh, order=3,discontinuous=True)
sigmaa = GridFunction(fes)
sigmaa.Set(sigma_cf)

Draw(div(sigmaa), mesh)

The other way is to calculate with mass-matrix :

sigma = 1/(mu0 * mur)  * (OuterProduct(B, B) - 0.5 * (B*B) * Id(2) )

fes = VectorH1(mesh, order=3)
u,v = fes.TnT()
force = LinearForm(fes)
force += SymbolicLFI( InnerProduct(sigma, grad(v)))
force.Assemble()

mass = BilinearForm(fes)
mass += SymbolicBFI(u*v, VOL)
mass.Assemble()

gfforce = GridFunction(fes)
gfforce.vec.data = mass.mat.Inverse(inverse="sparsecholesky", freedofs=fes.FreeDofs()) * force.vec

Draw(gfforce, mesh)

My question is why do I have different results, and also, which is the good way to calculate the force density ?

Thank you

1 Like