Submitted as: model description paper 27 May 2021
Submitted as: model description paper  27 May 2021
UBER v1.0: A universal kinetic equation solver for radiation belts
 ^{1}William B. Hanson Center for Space Sciences, Department of Physics, University of Texas at Dallas, Richardson, Texas, USA
 ^{2}Department of Physics and Astronomy, Rice University, Houston, Texas, USA
 ^{3}Department of Earth, Planetary and Space Sciences, University of California at Los Angeles, Los Angeles, California, USA
 ^{1}William B. Hanson Center for Space Sciences, Department of Physics, University of Texas at Dallas, Richardson, Texas, USA
 ^{2}Department of Physics and Astronomy, Rice University, Houston, Texas, USA
 ^{3}Department of Earth, Planetary and Space Sciences, University of California at Los Angeles, Los Angeles, California, USA
Abstract. Recent proceedings in the radiation belt studies have proposed new requirements for numerical methods to solve the kinetic equations involved. In this article, we present a numerical solver that can solve the general form of radiation belt FokkerPlanck equation and Boltzmann equation in arbitrarily provided coordinate systems, and with userspecified boundary geometry, boundary conditions, and equation terms. The solver is based upon the mathematical theory of stochastic differential equations, whose computational accuracy and efficiency are greatly enhanced by specially designed adaptive algorithms and variance reduction technique. The versatility and robustness of the solver is exhibited in three example problems. The solver applies to a wide spectrum of radiation belt modeling problems, including the ones featuring nondiffusive particle transport such as that arises from nonlinear waveparticle interactions.
Liheng Zheng et al.
Status: closed

RC1: 'Comment on gmd2021122', Anonymous Referee #1, 18 Jun 2021
This is a very nice paper describing the UBER 1 code, its potential applications and in particular, the advantages of the stochastic approach to solving convectiondiffusion equations including mixed diffusion terms. I have mostly minor comments that should be easy for authors to address. I would also recommend adding an example that would show how the code behaves in case of strong gradients and not smooth initial distributions.
 I would clarify around line 20 that Js are the canonical momentums; otherwise alpha would go from 1 to 6.
 I thought that we neglect the phases of J simply because we assume isotropy in all phis . I could not precisely follow the argument of electric potential energy. Please elaborate.
 While the paper discusses in detail the advantages of this approach, it omits the disadvantages. I would suggest clearly saying how long the presented runs take to run on a regular PC and on a supercomputer. I would also specify how much wall clock computing would take to calculate one day charging along a given satellite orbit. Also, specify how much slower this approach than more traditional approaches for the multidimensional diffusion equations.
 The presented tests all show examples that are initiated with smooth initial conditions. Please provide examples similar to Aseev et al., 2016 with strong gradients in initial conditions.
Aseev N. A., Y. Y. Shprits, A. Y. Drozdov, A. C. Kellerman, (2016), Numerical applications of the advectivediffusive codes for the inner magnetosphere, Space Weather, 14, 9931010, doi:10.1002/2016SW001484

RC2: 'Comment on gmd2021122', Anonymous Referee #2, 21 Jun 2021
This paper presents a very nice universal solver for transport equations with examples in various situations, utilizing the SDE method. It is always happy to see that significant improvements have been made to the SDE method in terms of speed and handling of boundary conditions, since its first adoption to the radiation belt problem in 2008. The SDE method clearly has its advantage in terms of handling crossdiffusion terms and efficiency the case when only a few solutions in (t,x) are needed. In principle, I think this paper is publishable. I only have a few minor comments and suggestions for the authors.
1, The main problem I have with the paper is that it failed to discuss the main concern of the the SDE method: its speed, especially compared with finite difference methods etc. Clearly, the SDE code with newer techniques is much much more efficient than the traditional one. But how does it compare with simple finite difference methods in the three cases discussed? I think the authors could discuss the comparison in different situations; e.g., when one cares about only a few snapshots, and when needs to know the whole history of evolution. Having a fair discussion about the disadvantage of the method does not make the paper less significant, but instead show future directions where improvements can be made.
2, Lines 4550: These listed numerical codes use either finite difference type methods or SDE/layer methods. None of the methods are limited by the choice of coordinates. Those authors chose a particular choice of coordinates probably because they did not intend to build a general library or because they tended to demonstrate a new method. I think it would nice for the authors to take this into consideration when discussing previous models.
3, Lines 5355: The complicated geometry is a problem for some of the models, mainly because radiation belt people seem to have a preference for finite difference methods. There exist general powerful finite volume/elements methods that can handle complicated boundary geometry. So here “… would be challending for numerical methods” should really be “… would be challenging for finite difference methods.”
4, In Table 1, the UBER library input, can the current code handle timedependent D? I know the method can, but not sure if it has been implemented by the code.
5, Lines 158: No, the grids in the layer method need NOT to be uniform. It was simply chosen for simplicity for estimating error and demonstration purposes.
6, Lines 155159: Yes, in SDE methods, one can design this kind of method to obtain the global solutions of f and its history. However, the key to implementing this in SDE is actually about finding an appropriate interpolation method. “The only operation … would be interpolation” sounds like finding such an interpolation method is easy, while in fact it is NOT. For the described choice of irregular domain or nonuniform domain, the interpolation method still needs to be of high order to reduce systematic error from interpolation, and to preserve positivity of the solution. For example, if a simple 4^{th} order polynomial interpolation method is used, one might introduce oscillation of f, and hence negative value s of f, from interpolation, and hence violates one of the key advantage of the SDE method.That is why Tao et al., 2016 (doi:10.1002/2015JA022064) introduced those higherorder positivitypreserving methods to be used with layer method.
7: Lines 310315: There actually exists fluxlimited LaxWendroff methods that can avoid introducing unphysical negative solutions. And in the case of forming steep gradients, it is wellknown that one should add fluxlimiters to LaxWendroff type methods.
Overall, I think this is a nice paper, and it would be an important contribution to the community of radiation belt modeling. I simply hope that the authors can give a more balanced description of the UBER library, especially when comparing with other methods.

RC3: 'Comment on gmd2021122', Anonymous Referee #3, 24 Jun 2021
The paper presents a new numerical solver, UBER, which can solve the general form of radiation belt FokkerPlanck equation and Boltzmann equation in arbitrarily provided coordinate systems, and with userspecified boundary geometry, boundary conditions, and equation terms. The mathematical theory and numerical techniques of UBER are clearly described, and three example problems are presented which well demonstrates the capability and reliability of the solver. Overall, this is a nicely written paper with convincing results showing that UBER could make important contribution to radiation belt modeling. The reviewer only has a few minor comments before recommending the paper for publication in GMD.
 Could the authors clarify or discuss further if the UBER solver can handle boundary conditions at varying locations or not? If we look at the radial diffusion model as a simple example, in the traditional solvers based on finite different method, it is often driven by datadriven outer boundary conditions at a fixed L*. However, the L* location of the satellite data providing the outer boundary is actually varying in time, which could lead to data gaps in the outer boundary condition and uncertainties in the model results. Discussions on if and how this type of boundary conditions can be implemented in UBER will further increase the significance of the work.
 The reviewer is also curious how the coordinate conversion between adiabatic invariants and energy and pitch angle can be implemented in UBER. As discussed in Subbotin and Shprits [JGR, 2012], several 3D diffusion radiation belt models “utilized two grids to solve the FokkerPlanck equation; one grid, which keeps the first and second adiabatic invariants constant, was used for the computation of the radial diffusion, and the other grid, orthogonal in energy and pitch angle at each fixed radial distance, was used for the computation of energy diffusion, pitch angle diffusion, and mixed energy and pitch angle diffusion.” At each time step, the results were converted and interpolated between the two grids, which could lead to uncertainties in the model results depending on how the conversion and interpolation are performed. Is this coordinate conversion still needed in UBER? If not, why?
 It will greatly enhance the significance of the paper if discussions are included in comparing the efficiency, stability, and accuracy between the UBER code and traditional finite difference codes. It is demonstrated the UBER is more efficient than previous SDE codes, but how does it compare with the traditional solvers based on the finite difference method? Is it still much less efficient if global distributions of radiation belt electrons in L, energy, and pitch angle are targeted? How about if we only need to solve for the distribution locally at certain L, energy, and pitch angle? The authors could perhaps use example problem 3 to compare the efficiency, stability, and the accuracy between the two different types of solvers and then expand the discussion to higherdimension models such as 3D diffusion models.
 AC1: 'Author comments on gmd2021122', Liheng Zheng, 19 Aug 2021
Status: closed

RC1: 'Comment on gmd2021122', Anonymous Referee #1, 18 Jun 2021
This is a very nice paper describing the UBER 1 code, its potential applications and in particular, the advantages of the stochastic approach to solving convectiondiffusion equations including mixed diffusion terms. I have mostly minor comments that should be easy for authors to address. I would also recommend adding an example that would show how the code behaves in case of strong gradients and not smooth initial distributions.
 I would clarify around line 20 that Js are the canonical momentums; otherwise alpha would go from 1 to 6.
 I thought that we neglect the phases of J simply because we assume isotropy in all phis . I could not precisely follow the argument of electric potential energy. Please elaborate.
 While the paper discusses in detail the advantages of this approach, it omits the disadvantages. I would suggest clearly saying how long the presented runs take to run on a regular PC and on a supercomputer. I would also specify how much wall clock computing would take to calculate one day charging along a given satellite orbit. Also, specify how much slower this approach than more traditional approaches for the multidimensional diffusion equations.
 The presented tests all show examples that are initiated with smooth initial conditions. Please provide examples similar to Aseev et al., 2016 with strong gradients in initial conditions.
Aseev N. A., Y. Y. Shprits, A. Y. Drozdov, A. C. Kellerman, (2016), Numerical applications of the advectivediffusive codes for the inner magnetosphere, Space Weather, 14, 9931010, doi:10.1002/2016SW001484

RC2: 'Comment on gmd2021122', Anonymous Referee #2, 21 Jun 2021
This paper presents a very nice universal solver for transport equations with examples in various situations, utilizing the SDE method. It is always happy to see that significant improvements have been made to the SDE method in terms of speed and handling of boundary conditions, since its first adoption to the radiation belt problem in 2008. The SDE method clearly has its advantage in terms of handling crossdiffusion terms and efficiency the case when only a few solutions in (t,x) are needed. In principle, I think this paper is publishable. I only have a few minor comments and suggestions for the authors.
1, The main problem I have with the paper is that it failed to discuss the main concern of the the SDE method: its speed, especially compared with finite difference methods etc. Clearly, the SDE code with newer techniques is much much more efficient than the traditional one. But how does it compare with simple finite difference methods in the three cases discussed? I think the authors could discuss the comparison in different situations; e.g., when one cares about only a few snapshots, and when needs to know the whole history of evolution. Having a fair discussion about the disadvantage of the method does not make the paper less significant, but instead show future directions where improvements can be made.
2, Lines 4550: These listed numerical codes use either finite difference type methods or SDE/layer methods. None of the methods are limited by the choice of coordinates. Those authors chose a particular choice of coordinates probably because they did not intend to build a general library or because they tended to demonstrate a new method. I think it would nice for the authors to take this into consideration when discussing previous models.
3, Lines 5355: The complicated geometry is a problem for some of the models, mainly because radiation belt people seem to have a preference for finite difference methods. There exist general powerful finite volume/elements methods that can handle complicated boundary geometry. So here “… would be challending for numerical methods” should really be “… would be challenging for finite difference methods.”
4, In Table 1, the UBER library input, can the current code handle timedependent D? I know the method can, but not sure if it has been implemented by the code.
5, Lines 158: No, the grids in the layer method need NOT to be uniform. It was simply chosen for simplicity for estimating error and demonstration purposes.
6, Lines 155159: Yes, in SDE methods, one can design this kind of method to obtain the global solutions of f and its history. However, the key to implementing this in SDE is actually about finding an appropriate interpolation method. “The only operation … would be interpolation” sounds like finding such an interpolation method is easy, while in fact it is NOT. For the described choice of irregular domain or nonuniform domain, the interpolation method still needs to be of high order to reduce systematic error from interpolation, and to preserve positivity of the solution. For example, if a simple 4^{th} order polynomial interpolation method is used, one might introduce oscillation of f, and hence negative value s of f, from interpolation, and hence violates one of the key advantage of the SDE method.That is why Tao et al., 2016 (doi:10.1002/2015JA022064) introduced those higherorder positivitypreserving methods to be used with layer method.
7: Lines 310315: There actually exists fluxlimited LaxWendroff methods that can avoid introducing unphysical negative solutions. And in the case of forming steep gradients, it is wellknown that one should add fluxlimiters to LaxWendroff type methods.
Overall, I think this is a nice paper, and it would be an important contribution to the community of radiation belt modeling. I simply hope that the authors can give a more balanced description of the UBER library, especially when comparing with other methods.

RC3: 'Comment on gmd2021122', Anonymous Referee #3, 24 Jun 2021
The paper presents a new numerical solver, UBER, which can solve the general form of radiation belt FokkerPlanck equation and Boltzmann equation in arbitrarily provided coordinate systems, and with userspecified boundary geometry, boundary conditions, and equation terms. The mathematical theory and numerical techniques of UBER are clearly described, and three example problems are presented which well demonstrates the capability and reliability of the solver. Overall, this is a nicely written paper with convincing results showing that UBER could make important contribution to radiation belt modeling. The reviewer only has a few minor comments before recommending the paper for publication in GMD.
 Could the authors clarify or discuss further if the UBER solver can handle boundary conditions at varying locations or not? If we look at the radial diffusion model as a simple example, in the traditional solvers based on finite different method, it is often driven by datadriven outer boundary conditions at a fixed L*. However, the L* location of the satellite data providing the outer boundary is actually varying in time, which could lead to data gaps in the outer boundary condition and uncertainties in the model results. Discussions on if and how this type of boundary conditions can be implemented in UBER will further increase the significance of the work.
 The reviewer is also curious how the coordinate conversion between adiabatic invariants and energy and pitch angle can be implemented in UBER. As discussed in Subbotin and Shprits [JGR, 2012], several 3D diffusion radiation belt models “utilized two grids to solve the FokkerPlanck equation; one grid, which keeps the first and second adiabatic invariants constant, was used for the computation of the radial diffusion, and the other grid, orthogonal in energy and pitch angle at each fixed radial distance, was used for the computation of energy diffusion, pitch angle diffusion, and mixed energy and pitch angle diffusion.” At each time step, the results were converted and interpolated between the two grids, which could lead to uncertainties in the model results depending on how the conversion and interpolation are performed. Is this coordinate conversion still needed in UBER? If not, why?
 It will greatly enhance the significance of the paper if discussions are included in comparing the efficiency, stability, and accuracy between the UBER code and traditional finite difference codes. It is demonstrated the UBER is more efficient than previous SDE codes, but how does it compare with the traditional solvers based on the finite difference method? Is it still much less efficient if global distributions of radiation belt electrons in L, energy, and pitch angle are targeted? How about if we only need to solve for the distribution locally at certain L, energy, and pitch angle? The authors could perhaps use example problem 3 to compare the efficiency, stability, and the accuracy between the two different types of solvers and then expand the discussion to higherdimension models such as 3D diffusion models.
 AC1: 'Author comments on gmd2021122', Liheng Zheng, 19 Aug 2021
Liheng Zheng et al.
Data sets
UBER v1.0: A universal kinetic equation solver for radiation belts Liheng Zheng, Lunjin Chen, Anthony A. Chan, Peng Wang, Zhiyang Xia, Xu Liu https://doi.org/10.5281/zenodo.4686050
Model code and software
zhenglh/UBER: First release of UBER Liheng Zheng https://doi.org/10.5281/zenodo.4671646
Liheng Zheng et al.
Viewed
HTML  XML  Total  BibTeX  EndNote  

221  85  14  320  3  2 
 HTML: 221
 PDF: 85
 XML: 14
 Total: 320
 BibTeX: 3
 EndNote: 2
Viewed (geographical distribution)
Country  #  Views  % 

Total:  0 
HTML:  0 
PDF:  0 
XML:  0 
 1