A second-order self-adjusting steepness based remapping method for arbitrary quadrilateral meshes
AA second-order self-adjusting steepness based remapping method for arbitrary quadrilateral meshes
Zhiwei He (corresponding author: [email protected])
Institute of Applied Physics and Computational Mathematics, Beijing 100094, China
Abstract
In this paper, based on the idea of self-adjusting steepness based schemes[5], a two-dimensional calculation method of steepness parameter is proposed, and thus a two-dimensional self-adjusting steepness based limiter is constructed. With the application of such limiter to the over-intersection based remapping framework, a low dissipation remapping method has been proposed that can be applied to the existing ALE method. Keyword: ALE, remapping, limiter,
The classical ALE method includes three steps [1]: (1) Lagrange step - In this step, the physical governing equations and the computing mesh are updated, thereby deriving the corresponding Lagrangian step numerical solution and the Lagrangian step computing mesh. (2) Rezone step - In this step, the nodes of the Lagrangian step computing mesh are relocated to a better position to form an optimized computing mesh. (3) Remapping step - In this step, the numerical solution obtained by the Lagrangian step would be shifted from the Lagrangian step computing mesh to the rezoned computing mesh. This article focuses on the remapping step. Remapping [2,3,4] does not involve physics, and does not include the time evolution of physical quantities. It only maps the physical quantities on the old grid to the new grid, so it is a relatively independent static process. However, the quality of a remapping algorithm directly affects the actual solution of the entire ALE method. Therefore, the study of remapping algorithms has always been a hot but difficult problem. This article further focuses on how to preserves sharpness (or steepness) of various discontinuous structures during the remapping process, thus to reduce unnecessary umerical dissipation. Recently, the author [5] proposed a class of so-called self-adjusting steepness based schemes. While maintains nominal high order, such schemes can perform anti-diffusion operations on the above discontinuities, thus the resolution of these discontinuities can be improved as much as possible. In this article, the author extends this method to a two-dimensional case. Furthermore, applied to the standard overlay-intersection-based remapping method with linear reconstruction, the author finally obtained a second-order self-adjusting steepness based remapping method suitable for any quadrilateral grid. The final numerical results show that while maintaining the nominal second-order accuracy of the smooth region, the resolution of the discontinuity has been significantly improved. The framework of this paper is as follows. Section 2 introduces the standard overlay-intersection-based remapping method; Section 3, the author will introduce the newly proposed two-dimensional self-adjusting steepness based remapping method in detail; Section 4, the above new method have been tested using a large number of numerical examples on different types of grids; Section 5 is the conclusion.
2. Standard overlay-intersection-based remapping method
Consider two sets of computing meshes: The original (old, Lagrangian) computational mesh whose cell is c , volume is c V , and centroid position is c r ; and the new (rezoned) computational mesh whose cell is c , volume is c V , and centroid position is c r . In this article, it has been assumed that the two sets of meshes have the same convexity, and for mesh cell c , its adjacent cell set is defined as ' C c . Generally, overlay-intersection-based remapping is based on the following formula [3,2] ' ' c c c c (1) This formula provides global remapping between two arbitrary meshes in the same computing zone. If the new mesh is obtained by moving the nodes of the original mesh by a small displacement (such as various mesh smoothing algorithms), we can assume that the intersection operation is just a local operation with cell c and its neighborhood ' C c . Therefore ' C c C c c , and Eq.(1) can be further written as [3,2]: ' ' c C c c c c (2) If two sets of meshes possess the same convexity, the above formula can be written further as follows [3,2] : ' ' ' \ ' c C c c c c c c c (3) Using this expression, for any scalar quantity per unit volume f ( its mean value on the original mesh cell c is c f ) , its remapping formula can be written as: , '' ' c c c c c cc C c f V f V F (4) where , ' ' ' c c c c c c F f dV f dV r r (5) In order to perform the above remapping operation, it is necessary to reconstruct the physical distribution f r based on the cell average value c f on the old mesh. Second order self-adjusting steepness based remapping method
This paper considers piece-wise linear reconstruction: c c f f f r r r (6) However, the above results cannot guarantee the essential non-oscillation property and bound of the final solution. Therefore, the gradient on this cell c needs to be limited. Under the MUSCL framework [26], the construction of multi-dimensional limiters can generally be divided into two types [27,23]: monoslope method and multislope methods. For simple, efficient and better matching of the intersection algorithm in the remapping method, we adopts the monoslope method, i.e., given a cell c , there is only one unique gradient: lim c c f f f r r r (7) where lim f is the final gradient obtained by utilizing various nonlinear limiters. Frequently used limiters of this type are: Barth-Jespersen limiter [15], Venkatakrishnan limiter [16], etc. However, the above-mentioned commonly used limiters will, in turn, severe smear various discontinuities, especially linear discontinuities (such as contact iscontinuities in compressible fluids, diffusive interfaces in multi-material fluids, etc.) And as the calculation time increases, this smearing effect continues. In the previous references [5], the author proposed a new concept called steepness-adjustable limiters. The biggest character of this class of limiters is that they should have a steepness parameter that provides a mechanism to accurately solve both smooth and discontinuous solution by adjusting according to the flow structure. In that paper, the author propose to a simple method to construct such limiters: extend some existing total-variation-diminishing (TVD) limiters into steepness-adjustable limiters. To capture both smooth and discontinuous solutions, the classic harmonic a limiter [26] is modified into a limiter with adjustable steepness: (8) where is the ratio of the upwind increment to the local increment in the standard TVD schemes, and is the steepness parameter. It can be proved theoretically that the steepness-adjustable limiter achieves second-order accuracy when the steepness parameter takes a specific value (remarked as s ). However, when the value of increases (remarked as l ), the slope adjustable limiter introduce the anti-di ff usion. In the reference[5], related works have been applied to the flux-splitting based finite difference method, and the conclusions are verified by the numerical test results. Further details on this topic can be referred from the author’s previous articles [5]. This article plans to extend such method to the two-dimensional case in order to reduce the numerical dissipation during the remapping process. To constructing such type limiter, the following core issue needs to be addressed: how to give a method to calculate the steepness parameter of a mesh cell (i.e.: there is only one steepness parameter corresponding to a mesh cell). Assuming that the steepness parameter of the cell is known, we can directly propose a multidimensional self-adjusting steepness based limiter. The details are as follows: Taking the mesh cell c into consideration, the set of edges is denoted as E c , the total number of edges of this cell is c E . For any edge e , its length is e l , and the adjacent cell sharing this edge is * c . The smoothness indicator corresponding to this dge can be defined as: e c c IS f f (9) The linear weight corresponding to this edge is defined as the ratio of the length of this edge to the sum of the lengths of all the sides of this cell, which is represented as: '' ee ee E c lD l (10) Further, the WENO-JS methodology [25] is used to compute the non-linear weight corresponding to this edge: '' ee ee E c (11) where ee pe DIS (12) where p is the power parameter whose value is usually 2. The final steepness parameter of this cell is defined as: c c s c l (13) where c can be defined as: '' '' ee E cc ee E c D (14) Substituting the above formulas (Eqs.9-12) into Eq.(14), the final equation can be derived as: cc c EE pee E cc Epe se E c s E c e
ISD IS (15) where is a small parameter ,used to avoid the zero value of denominator, usually its value are
20 40
10 10 . Finally, we give the complete process of multi-dimensional second-order self-steepness based remapping method: (1) In the mesh cell c , perform a linear reconstruction to get f , and form c c f f f r r r . (2) Calculate the steepness parameter of this cell: c . (3) At each node of this cell n N c , where N c is the set of nodes of this cell, calculate an allowable maximum n as where n c n c f f f r r is the solution of linear reconstruction without limitation at the node n . min f and max f are the minimum and maximum values of the current cell and its neighboring cell, respectively. (4) Set min c nn N c where n nn c n . (5) Define lim c f f , and the final reconstruction function is lim c c f f f r r r . (6) Execute Eq.(5), and obtain the average value c f on the rezoned mesh cell c .
4. Numerical examples
In the ALE method, rezone/remapping is not mandatory for every computing step. In principle, the only in the occasion that the Lagrangian computational mesh is kinked, the rezone of the mesh and remapping of the physical quantities are needed. However, rezone/remapped per step (i.e., continuous rezone/remapping) can ensure the convexity of the grid and the intersection operation can be performed locally [3]. Therefore, the continuous rezone/remapping has proved to be a more effective fashion, and is adopted in this article. Additionally, all the computations of this paper are in the following mesh { , , , n ni j i j x y } are carried on: , max max max , 1, , , 1, , ; 0 , ; ni j x i i j j n = , n , max max max , 1, , , 1, , ; 0 , ; ni j y i i j j n = , n .1 Tensor Product Mesh In the computing zone , the author uses the following functions [3]: , , 1 x t t t , , , 1 y t t t , sin 42 tt , ; ; t . to generate a series of tensor product meshes: , , , n ni j i j x x t ; , , , n ni j i j y y t where max max / , 0, , n t n n n n , and maxmax maxmax ij i i ii j j jj On the basis of these meshes, the author carried out the following three physical quantity distribution tests [3].
In this example [3], the function f is a smooth function whose expression is: , 1 sin 2 sin 2 f x y x y where max max i j and max n . The author uses the Barth-Jespersen limiter and self-adjusting steepness (SAS) based limiter to calculate this example. Note that, in the SAS limiter, the lower bound of the steepness parameter l is theoretically determined (ensuring second-order accuracy), but its upper bound s is still uncertain. Therefore, the authors performed the calculation with s and s as an example. Figure 1 give the final results after continuous rezone and remapping, where (a) is the initial field, (b) is the result obtained using the Barth-Jespersen limiter, (c) is the result obtained using SAS limiter ( s ), and (d) is the result obtained using SAS limiter ( s ) . Figure 1. Sine test results on the tensor product mesh. (a) Initial field; (b) Barth-Jespersen limiter; (c) self-adjusting steepness based limiter ( s ); (d) self-adjusting steepness based limiter ( s ) . It can be seen from the figure that no matter what the upper bound of the steepness parameter in SAS limiter is, the final result is very close to the result obtained using the Barth-Jespersen limiter. This result shows that the SAS limiter achieve the nominal second-order accuracy calculation for this smooth structure. In this example [3], the function f is a discontinuous function. Its expression is:
0, 0.4 / 0.3, 1, 0.4 / 0.3 y xf x y y x where max max i j and max n . Figure 2 give the final results after continuous rezone and remapping, where (a) is the initial field, (b) is the result obtained using the arth-Jespersen limiter, (c) is the result obtained using SAS limiter ( s ), and (d) is the result obtained using SAS limiter ( s ) . Figure 2. Shock test results on the tensor product mesh. (a) Initial field; (b) Barth-Jespersen limiter; (c) self-adjusting steepness based limiter ( s ); (d) self-adjusting steepness based limiter ( s ) . Comparing the results obtained with the Barth-Jespersen limiter, it can be seen that for discontinuous structures, the SAS limiter has an anti-diffusion mechanism. And as the upper bound of the steepness parameter increases, the anti-diffusion mechanism becomes more obvious, which leads to a significant improvement in the resolution of the discontinuity. In this section, the author will consider the type of grid movement that is more in line with the ALE method. On a uniform mesh, the nodes undergo independent andom disturbances to obtain a random mesh [3]: , n ni j i i x h , n ni j i j y h where n ni j is a random number. At the same time, to mimic the rezone process, the author uses the following smoothing method to obtain the new mesh [3]:
1, , 1, , 1 , , 11, n n n n n ni j i j i j i j i j i jni j x x x x x xx
1, , 1, , 1 , , 11, n n n n n ni j i j i j i j i j i jni j y y y y y yy
Figure 4 highlights the difference between the initial random mesh and the mesh obtained after different computation steps, where max max i j . In Figure 4, (a) is the initial mesh (black solid line) and the mesh obtained after one step of smoothing (red dotted line); (b) is the initial mesh (black solid line) and the mesh obtained after 20 steps of smoothing (red dotted line). Figure 4. The comparison chart of the difference between the initial random mesh and the mesh obtained after various computation steps, where max max i j . (a) The initial mesh (black solid line) and the mesh obtained after one step of smoothing (red dotted line); (b) the initial mesh (black solid line) and the mesh obtained after 20 steps of smoothing (red dotted line). On the basis of these meshes, the author carried out the remapping of the above three physical quantity distribution test cases. .2.1 “Sine” Test The setup of this example is similar to that of 4.1.1. Figure 3 give the final results after continuous rezone and remapping, where (a) is the initial field, (b) is the result obtained using the Barth-Jespersen limiter, (c) is the result obtained using SAS limiter ( s ), and (d) is the result obtained using SAS limiter ( s ) . Figure 3. Sine test results on a random mesh. (a) Initial field; (b) Barth-Jespersen limiter; (c) self-adjusting steepness based limiter ( s ); (d) self-adjusting steepness based limiter ( s ) . It can be seen from the figure that even on a random mesh, the final result obtained by using the SAS limiter is still very close to the result obtained using the Barth-Jespersen limiter. This result once again shows that the SAS limiter better guarantees the nominal second-order accuracy for this smooth structure, independent of the type of mesh. .2.2 Shock Test The setup of this computing example is similar to that of 4.1.3. Figure 4 give the final results after continuous rezone and remapping, where (a) is the initial field, (b) is the result obtained using the Barth-Jespersen limiter, (c) is the result obtained using SAS limiter ( s ), and (d) is the result obtained using SAS limiter ( s ) . Figure 4. Shock test results on a random mesh. (a) Initial field; (b) Barth-Jespersen limiter; (c) self-adjusting steepness based limiter ( s ); (d) self-adjusting steepness based limiter ( s ) . It can be observed from the figure that with the increase in the upper bound of the steepness parameter in the SAS limiter, the resolution of the discontinuity is significantly improved. This further justifies the above series of conclusions. During the remapping process, the physical quantity needs to be reconstructed on he old mesh, and a limiter would be further applied to limit the reconstruction result. However, such methods will bring serious numerical dissipation on some linear discontinuities, such as contact discontinuities in compressible fluids, diffusive interfaces in multi-material fluids, etc. Thus such physical structures are severer smeared. Furthermore, as the calculation time increases, this smearing effect continues. In the previous references [5], the author proposed a new concept called steepness-adjustable limiters. The biggest character of this class of limiters is that they should have a steepness parameter that provides a mechanism to enable the scheme to accurately solve both smooth and discontinuous solution by adjusting the steepness parameter according to the flow structure. In this paper, we extend such limiter to a two-dimensional case in order to reduce the numerical dissipation during the remapping process. After proposing a method to the steepness parameter of a whole mesh cell, we directly give a multidimensional self-adjusting steepness based limiter. Furthermore, applied to the standard overlay-intersection-based remapping method with linear reconstruction, the author finally obtained a second-order self-adjusting steepness based remapping method suitable for any quadrilateral grid. The numerical results further suggest that while maintain the nominal second-order accuracy of the smooth region, the resolution of the discontinuity can be significant improved by simply adjusting upper bound of the steepness parameter.
Reference [1] A. Barlow, P. Maire, W. J.Rider, R. Rieben, M. Shashkov. Arbitrary Lagrangian-Eulerian methods for modeling high-speed compressible multimaterial flows. Journal of Computational Physics 322(2016)603-665. [2] M. Kucharik, M. Shashkov. Conservative multi-material remap for staggered multi-material Arbitrary Lagrangian-Eulerian methods. Journal of Computational Physics 258(2014)268-304. [3] L.G. Margolin, M. Shashkov. Second-order sign-preserving conservative interpolation (remapping) on general grids. Journal of Computational Physics 184(1) 2003) 266–298. [4] J. K. Dukowicz and J. R. Baumgardner. Incremental remapping as a transport/advection algorithm. Journal of Computational Physics, 160(1):318–335, 2000. [5] Zhiwei He, Yucang Ruan, Yaqun Yu, Baolin Tian, Feng Xiao. Self-adjusting steepness based schemes that can preserve discontinuities in compressible flows. submitted. [6] M. Berndt, J. Breil, S. Galera, M. Kucharik, P. Maire, M. Shashkov. Two-step hybrid conservative remapping for multimaterial arbitrary Lagrangian–Eulerian methods. Journal of Computational Physics 230 (2011) 6664–6687. [7] M. Kucharik, M. Shashkov. One-step hybrid remapping algorithm for multi-material arbitrary Lagrangian–Eulerian methods.