I am learning to use Elmer for magnetostatics in 3D. As an exercise, I setup a simple axisymmetric problem: a cylindrical magnet centered at the origin, which I can compare to FEMM.
One of my goals is to obtain a divergencefree potential. I am using the WhitneyAVSolver with Lagrange gauge. I have set the Lagrange gauge penalization coefficient [variable k in eq. (17.54) of the Elmer Models Manual] to different values, and expected the vector potential to change. I tried k = 1, 2, 10 and, finally 1e8. No change (except for small numerical perturbations) is observed. I have attached a figure comparing values of the zcomponent of the potential (on a grid in the x = 0 plane). The different lines correspond to different y coordinates. Notice that the change in Az from k = 1 to k=1e8 is nearly none (Az is between 1e5 and 1, the difference is concentrated around 1e9).
There are two other related issues. The Az in the plot changes with the coordinate z, an indication that the divergence is not zero. In addition, Az should be zero on this plane (horizontal and vertical components of the vector should be zero at any plane containing the symmetry axis z).
In summary: changing the Lagrange gauge penalization coefficient of the WhitneyAVSolver does not seem to have any meaningful effect on the vector potential.
I have attached a zip file with the sif file and the geometry. What am I doing wrong?
Thanks in advance.
3D Magnetic Vector Potential Lagrange Gauge Penalization
3D Magnetic Vector Potential Lagrange Gauge Penalization
 Attachments

 singlemag_air_sphere.zip
 (4.63 KiB) Downloaded 20 times

 Site Admin
 Posts: 3697
 Joined: 22 Aug 2009, 11:57
 Antispam: Yes
 Location: Espoo, Finland
 Contact:
Re: 3D Magnetic Vector Potential Lagrange Gauge Penalization
Hi
Haven't used it much myself. Does it affects the results if you turn off the Piola transform?
Is there a special reason that you want to gauge your vector potential?
Peter
Haven't used it much myself. Does it affects the results if you turn off the Piola transform?
Is there a special reason that you want to gauge your vector potential?
Peter

 Posts: 641
 Joined: 25 Jan 2019, 01:28
 Antispam: Yes
Re: 3D Magnetic Vector Potential Lagrange Gauge Penalization
According to your geo file the model is built in mm, in the solver the units should be meters. Have you tried
Coordinate Scaling = 0.001
I think this geo file produces the same geometry
Coordinate Scaling = 0.001
I think this geo file produces the same geometry
 Attachments

 mygeo.geo
 (374 Bytes) Downloaded 18 times
Re: 3D Magnetic Vector Potential Lagrange Gauge Penalization
Hi,
In the first place I would compare solutions obtained without a gauge condition (Use Lagrange Gauge = False) and with the parameter value k = 1 to see how the Lagrange multiplier formulation affects the vector potential. In exact arithmetic I would expect that the vector potential solution is the same for all positive constants k, while the scalar variable solution depends on k (if k differs from 1, consider what follows from the change of scalar variable V as V' = k V; the effect seems to be that the penalization term in the equation for the vector potential is just scaled by a factor 1/k in comparison with solving the problem with k = 1). An extreme value of k might cause perturbations due to illconditioning, but I think seeing no changes in the vector potential solution is indeed expected.
 Mika
In the first place I would compare solutions obtained without a gauge condition (Use Lagrange Gauge = False) and with the parameter value k = 1 to see how the Lagrange multiplier formulation affects the vector potential. In exact arithmetic I would expect that the vector potential solution is the same for all positive constants k, while the scalar variable solution depends on k (if k differs from 1, consider what follows from the change of scalar variable V as V' = k V; the effect seems to be that the penalization term in the equation for the vector potential is just scaled by a factor 1/k in comparison with solving the problem with k = 1). An extreme value of k might cause perturbations due to illconditioning, but I think seeing no changes in the vector potential solution is indeed expected.
 Mika
Re: 3D Magnetic Vector Potential Lagrange Gauge Penalization
First of all, thank you for the observations and questions.
(1) Peter,
(a) I have run cases without Piola transformation and have not seen any change.
(b) I used the gauge because of the information in the models manual: "possibility to seek a unique vector potential solution to the magnetostatics equations is included here as a special case via using the scalar variable to impose the divergencefree constraint on A"
(2) Kevinarden, my geo file is in millimeters and I change the node coordinates in mesh.nodes to meters ( I did it in Emacs). Your geo file produces the same geometry, but I have not used it. I like its simplicity (compare to the one I generated from a brep cad file)
(3) Mika, I have run cases (Piola transformation off) for k from 0.1 to 1e8. Following your suggestion, I also run one case without the Lagrange gauge.
I am posting one plot comparing results with Lagrange gauge penalty coefficient k=1 to the no Lagrange gauge case.
There is a slight increase in the difference (compared to k=1, k=1e8 posted before), but they are not numerically meaningful.
I have used the divergence solver to compute the divergence of the vector potential (divA), which is the quantity that should be reduceded by the Lagrange gauge. I am using, as measure, the square root of the sum of the squares of divA at all points in the rectangular region on the YZ plane indicated in the plot (a kind of total energy of divA). As the divergence decreases, I would expect this measure to decrease.
However, as you can see in the plot, there is no change in the total "energy" of divA. The red dot is for the case with no Lagrange gauge. Notice the small range of the vertical axis on this plot.
The analysis suggests that, at least for this particular problem, Lagrange gauge has no effect on the solution: no change in the vector potential and no reduction in its divergence.
(1) Peter,
(a) I have run cases without Piola transformation and have not seen any change.
(b) I used the gauge because of the information in the models manual: "possibility to seek a unique vector potential solution to the magnetostatics equations is included here as a special case via using the scalar variable to impose the divergencefree constraint on A"
(2) Kevinarden, my geo file is in millimeters and I change the node coordinates in mesh.nodes to meters ( I did it in Emacs). Your geo file produces the same geometry, but I have not used it. I like its simplicity (compare to the one I generated from a brep cad file)
(3) Mika, I have run cases (Piola transformation off) for k from 0.1 to 1e8. Following your suggestion, I also run one case without the Lagrange gauge.
I am posting one plot comparing results with Lagrange gauge penalty coefficient k=1 to the no Lagrange gauge case.
There is a slight increase in the difference (compared to k=1, k=1e8 posted before), but they are not numerically meaningful.
I have used the divergence solver to compute the divergence of the vector potential (divA), which is the quantity that should be reduceded by the Lagrange gauge. I am using, as measure, the square root of the sum of the squares of divA at all points in the rectangular region on the YZ plane indicated in the plot (a kind of total energy of divA). As the divergence decreases, I would expect this measure to decrease.
However, as you can see in the plot, there is no change in the total "energy" of divA. The red dot is for the case with no Lagrange gauge. Notice the small range of the vertical axis on this plot.
The analysis suggests that, at least for this particular problem, Lagrange gauge has no effect on the solution: no change in the vector potential and no reduction in its divergence.

 Site Admin
 Posts: 3697
 Joined: 22 Aug 2009, 11:57
 Antispam: Yes
 Location: Espoo, Finland
 Contact:
Re: 3D Magnetic Vector Potential Lagrange Gauge Penalization
Hi
The fact that no gauging can work uses the properties of the Krylov subspace of the iterative solvers, see
https://en.wikipedia.org/wiki/Krylov_subspace
Upon successful completion of the iterative solution for the particular Krylov method is therefore well defined. However, you cannot use direct solvers without an explicit gauge.
This feature could be useful for example when a) source is not divergence free and the augmented lagrangian relaxes the this condition or b) there are nonlinear iterations that may eventually pollute the vector potential making the nullspace solution grow or c) you reallyreally want to compare the vector potential between some reference solution with the same gauging. These are deeply rooted issues and I'm not a mathematician myself....
Peter
The fact that no gauging can work uses the properties of the Krylov subspace of the iterative solvers, see
https://en.wikipedia.org/wiki/Krylov_subspace
Upon successful completion of the iterative solution for the particular Krylov method is therefore well defined. However, you cannot use direct solvers without an explicit gauge.
This feature could be useful for example when a) source is not divergence free and the augmented lagrangian relaxes the this condition or b) there are nonlinear iterations that may eventually pollute the vector potential making the nullspace solution grow or c) you reallyreally want to compare the vector potential between some reference solution with the same gauging. These are deeply rooted issues and I'm not a mathematician myself....
Peter
Re: 3D Magnetic Vector Potential Lagrange Gauge Penalization
I would interpret this that by chance the standard solution doesn't have a significant component which could be represented in terms of a gradient field.
On the other hand, how do you compute div A? As A is approximated in terms of curlconforming basis functions, I think div cannot be applied directly to A, so div A must be considered in weak sense. Might it be easiest to think that the approximation of the Laplacian of the scalar variable gives an indication of the size of div A, cf. the equation div A  k div grad V = 0? Postprocessing A to create a nodal approximation of A and then generating an approximation of div A might not be an ideal way. It seems that the standard divergence solver (DivergenceSolver) of Elmer cannot operate directly on a field approximated in terms of curlconforming basis functions.
 Mika