For the magnetic field, the boundary conditions I chose are: \mathbf{B}\times\mathbf{n}|_{\partial \Omega}=0,\quad \nabla\cdot \mathbf{B}|_{\partial \Omega}=0
Regarding the boundary conditions \mathbf{B}\times\mathbf{n}|_{\partial \Omega}=0, I know the correct approach is:
However, I applied a divergence-free condition to the boundary \nabla\cdot \mathbf{B}|_{\partial \Omega}=0. How is this specifically implemented in the FreeFem++ code?
The issue is indeed to have an appropriate variational formulation, and to choose a good finite element space.
Here since you already set {\mathrm{div} B}=0 in \Omega, you don’t really need to impose {\mathrm{div} B}=0 in \partial\Omega because the latter follows from the former.
If you want nevertheless to enforce such condition you could add a term -\nabla\mathrm{div} B in the B equation. Then you have to take a space of the type P1 for B (you need \nabla B in L^2).
On the contrary if you don’t put such term you can take a space with only \nabla\times B\in L^2, such as RT03d or Edge03d, this should be more efficient.
I notice that you don’t have to really impose {\mathrm{div} B}=0 in \Omega because if it holds at t=0, the B equation implies that it will hold for all t>0 (but maybe not exactly at the discrete level).