the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Scalarpotential mapping of the steadystate magnetosheath model
Yasuhito Narita
Daniel Schmid
Simon Toepfer
The steadystate magnetosheath model has various applications for studying the plasma and magnetic field profile around the planetary magnetospheres. In particular, the magnetosheath model is analytically obtained by solving the Laplace equation for parabolic boundaries (bow shock and magnetopause). We address the question, how can we utilize the magnetosheath model by transforming into a more general, empirical, nonparabolic magnetosheath geometry? To achieve the goal, we develop the scalarpotential mapping method which provides a semianalytic estimate of steadystate flow velocity and magnetic field in the empirical magnetosheath domain. The method makes use of a coordinate transformation from the empirical magnetosheath domain into the parabolic magnetosheath domain and evaluates a set of variables (shell variable and connector variable) to utilize the solutions of the Laplace equation obtained for the parabolic magnetosheath domain. Our model uses two invariants of transformation: the zenith angle within the magnetosheath with respect to the direction to the Sun and the ratio of the distance to the magnetopause and the thickness of magnetosheath along the magnetopausenormal direction. The use of magnetopausenormal direction makes a marked difference from the earlier model construction using the radial direction as reference. The plasma flow and magnetic field can be determined as a function of the upstream condition (flow velocity or magnetic field) in a wide range of zenith angles. The scalarpotential mapping method is computationally inexpensive, using analytic expressions as much as possible, and is applicable to various planetary magnetosheath domains.
 Article
(4373 KB)  Fulltext XML
 BibTeX
 EndNote
Steadystate plasma flow and magnetic field can be regarded as a realization of potential field in the planetary magnetosheath region when the vorticity and the electric current are ignored. In such a case, the potential is obtained by solving the Laplace equation, which was elegantly and analytically solved by Kobel and Flückiger (1994) for a parabolic shape of magnetosheath (hereafter KF). The KF potential was further extended to the stream function in the magnetosheath by Guicking et al. (2012). The KF solution made a series of breakthroughs in magnetosheath research. One of the most successful applications is the ability to track the plasma parcel along the streamline in the modeled magnetosheath. The tracking method was extensively used to observationally study the mirror mode growth (e.g., Tatrallyay et al., 2002; Génot et al., 2011) and the streamwise turbulence evolution in the magnetosheath (Guicking et al., 2012). Predictive models of plasma flow and magnetic field serve as a useful tool when combined with the numerical simulation or the observational data.
The KF potential is obtained using the assumption that the planetary bow shock and magnetopause have a parabolic shape sharing the same focal point. Empirical models of the bow shock and magnetopause (fitted to the spacecraft data), on the other hand, are not necessarily parabolically or cofocally shaped. For example, the empirical Earth bow shock model by Farris et al. (1991) and Cairns et al. (1995) has a parabolic shape, but the focal point differs from that of the KF solution. The empirical magnetopause model by Shue et al. (1997) applies a powerlaw scaling to the parabolic shape such that the magnetic field lines appear stretched in the tail region. The gap between the KF parabolic magnetosheath and the empirical magnetosheath needs to be filled when applying the KF potential in the empirical magnetosheath.
Naively speaking, one wishes to find a conformal mapping (anglepreserving mapping) from the KF parabolic magnetosheath onto a nonparabolic empirical magnetosheath shape such as the analytic extension of magnetopause shape (Narita et al., 2023). However, no general mathematical algorithm is known so far to obtain the conformal mapping when the spatial domain is not properly bounded. The problem lies in the fact that the magnetosheath is bounded only by two sides, i.e., the standing shock and the magnetopause in the radial direction to the planet, but it is not bounded along the flow in the tail region. The algorithms of numerical conformal mapping are so far proposed for spatially bounded domains (Papamichael and Whiteman, 1973; Chakravarthy and Anderson, 1979; Fornberg, 1980; Karageorghis et al., 1996) or domains with a closed shape of internal boundaries (Wei et al., 2014).
Here we address the question, how can we utilize the KF magnetosheath model by transforming into a more general, empirical, nonparabolic magnetosheath geometry? To achieve the goal, we develop a mapping method which provides a semianalytic estimate of steadystate flow velocity and magnetic field in the empirical magnetosheath domain. Our scalarpotential mapping method is computationally inexpensive by using the analytic expression as much as possible, and it is applicable to various planetary magnetosheath domains.
This work is organized in the following fashion. After reviewing the magnetosheath model constructed by Kobel and Flückiger (1994) (Sect. 2) and discussing different mapping methods (Sect. 3), the detailed procedure of the magnetopausenormal mapping is presented (Sect. 4) and discussed (Sect. 5), with concluding remarks (Sect. 6).
2.1 Parabolic coordinates
In the KF parabolic coordinates, the shell variable v (isocontour lines enveloping the magnetosphere) and the connector u (isocontour lines connecting from the bow shock to the magnetopause) play an important role in computing the flow velocity and magnetic field in the magnetosheath. These variables are explicitly evaluated using Cartesian coordinates and the radial distance from the focal point as
where r_{0} is the distance to the focus at x_{0}:
The focus is along the x axis and is defined as
x_{k}, y_{k}, and z_{k} are the Cartesian representation of the KF magnetosheath model (i.e., with the prefixed bow shock and magnetopause shapes) obtained by projecting the position vector onto the unit vectors e_{x}, e_{y}, and e_{z}:
To complete the variable set for computing the potentials and the stream function, the azimuthal angle ϕ is introduced as
2.2 Velocity potential
In the frame of potential field theory, the flow velocity U is obtained either from the velocity potential (a scalar potential) Φ^{(vel)} or from the stream function (also a scalar potential) Ψ as
The symbol e_{ϕ} is the unit vector in the azimuthal direction around the symmetry axis (Suntoplanet direction). Kobel and Flückiger (1994) and Guicking et al. (2012) obtained the analytic expression of the velocity potential Φ^{(vel)} using the shell variable v and the connector variable u.
where U_{x} is the upstream flow velocity, v_{mp} is the shell variable at the magnetopause, v_{bs} is the shell variable at the bow shock, v is the shell variable, u is the connector variable, and ${\mathrm{\Phi}}_{\mathrm{0}}^{\left(\mathrm{vel}\right)}$ is a free parameter (integration constant) which is set to zero without loss of generality. The boundary shell values v_{mp} and v_{bs} contain the information on the standoff distances (R_{mp} and R_{bs}) in the subsolar region, and they are defined by Kobel and Flückiger (1994) as
2.3 Stream function
Guicking et al. (2012) transformed the KF potential and obtained analytically the stream function Ψ as a function of the shell variable and the connector variable:
Hereafter, one may set ${U}_{x}=\mathrm{1}$ so that the velocity potential Φ^{(vel)} is normalized to the upstream velocity. Isocontour lines of the stream function represent the streamline.
2.4 Magnetic scalar potential
The magnetic field in the magnetosheath is derived from the scalar potential in the same fashion as the flow velocity; that is,
The magnetic potential is a function of the shell variable v and the connector u (Kobel and Flückiger, 1994):
where ${B}_{x}^{\left(\mathrm{up}\right)}$ is the sunward component of the upstream magnetic field (corresponding to the GSEX in nearEarth space), and ${B}_{y}^{\left(\mathrm{up}\right)}$ and ${B}_{z}^{\left(\mathrm{up}\right)}$ are two components of the upstream magnetic field perpendicular to the x direction. ϕ is the azimuthal angle of the position around the symmetry axis (the y direction is given by the angle ϕ=0). The integration constant is chosen as ${\mathrm{\Phi}}_{\mathrm{0}}^{\left(\mathrm{mag}\right)}=\mathrm{0}$. The magnetic potential cannot be further transformed into the form of a stream function since the magnetic field distribution is essentially threedimensional in the magnetosheath.
3.1 Mapping problem
Our task is to find the shell variable v and the connector u in the empirical magnetosheath by finding a suitable mapping of the position vector from the empirical magnetosheath (denoted by r) onto the KF parabolic system (denoted by r^{(k)}). The evaluated v and u are then readily used to obtain the scalar potentials and the stream function. The flow velocity and the magnetic field in the empirical magnetosheath are obtained by computing the gradient of the respective potential.
A practically useful mapping procedure to utilize the KF potential is proposed by Soucek and Escoubet (2012) by using the radial direction as a reference. While the radial mapping can reasonably (i.e., with a relatively high accuracy) transform the KF potential into the empirical magnetosheath domain on the dayside, the mapping quality becomes degraded in the flank region due to the conversion effect associated with the nonorthogonal grid construction. Our approach differs from the radial mapping by using the magnetopausenormal direction as a reference. The azimuthal coordinate ϕ is still orthogonal to the u and v coordinates. We briefly compare between the two mapping methods here.
3.2 Radial direction as reference
Soucek and Escoubet (2012) presented in their pioneering work an algorithm of radial mapping by transforming the KF magnetosheath model into a general, empirical magnetosheath shape by referring to the radial direction from the planet and scaling the radial position in the magnetosheath to the KF model. While the radial mapping can reasonably (i.e., with a relatively high accuracy) transform the KF potential into the empirical magnetosheath domain on the dayside, the mapping quality becomes degraded in the flank region due to the strongly nonorthogonal grids. Figure 1 displays a comparison of the radial grids between the KF magnetosheath model and the empirical magnetosheath model. The grids span the radial direction to the planet and transfinite interpolation between the bow shock and the magnetopause. The radial mapping has a drawback in a stronger grid nonorthogonality effect, which causes an artificial converging flow pattern in the flank region (velocity potential shown in Fig. 2) when the scalar potential is directly transformed. In the Soucek–Escoubet method, the problem of the flow conversion effect was avoided by solving the magnetohydrodynamics (MHD) Rankine–Hugoniot relation and tracking the streamline iteratively between the KF parabolic magnetosheath model and the empirical magnetosheath model.
Although the method introduced by Génot et al. (2011) and later adapted by Soucek and Escoubet (2012) is computationally less expensive than global magnetosheath simulations, the density and velocity fields from the bow shock to a given point in the magnetosheath still need to be computed along the streamline in an incremental way. Moreover, the Rankine–Hugoniot relations need to be solved before calculating iteratively the streamline, the flow velocity vector, to track the plasma density flow velocity along the streamline. Naturally, the uncertainty in this calculation depends on the step size (larger uncertainties for larger step sizes), and the errors accumulate along the streamline. The method introduced in this work overcomes this issue by constructing a magnetopausenormal grid system such that computational efforts are improved (no need to solve the Rankine–Hugoniot relations, and the error does not accumulate in the flank region).
3.3 Magnetopausenormal direction as reference
Our mapping method differs from the radial mapping method in that the magnetopausenormal direction is used as a reference to the magnetopause. Our method guarantees the grid orthogonality around the magnetopause both on the dayside and in the flank region. The magnetopausenormal grids are shown in Fig. 3 for the KF magnetosheath model (with parabolic boundaries) and the empirical magnetosheath model (with nonparabolic boundaries). Even though the exact conformal mapping is not available, the magnetopausenormal mapping method retains the grid orthogonality around the magnetopause. This feature (orthogonality around the magnetopause) plays a crucial role in mapping the scalar potentials. An example of the scalarpotential mapping by referring to the magnetopausenormal direction (our final results) is shown in Sect. 4.8.
4.1 Overview of the procedure
The magnetopausenormal mapping is performed with two transformations. In the first transformation, the position vector is mapped from the empirical magnetosheath r onto the KF magnetosheath model r_{k}. This is based on the assumption that the distance to the magnetopause along the magnetopausenormal direction when scaled to the magnetosheath thickness (defined as the distance from the magnetopause to the bow shock along the magnetopausenormal direction) remains constant. The azimuthal angle ϕ is the same between the empirical magnetosheath and the KF model. The first transformation is divided into computing the distance to the magnetopause (step 1), the thickness of the empirical magnetosheath (step 2), the thickness of the KF magnetosheath (step 3), and the mapping of the position vector onto the KF model (step 4).
In the second transformation, the mapped position vector is used to compute the shell variable v and the connector variable u (step 5) and to obtain the potentials and the stream function in the empirical magnetosheath using Eqs. (10), (13), and (15) (step 6). Here again, the azimuthal angle ϕ is treated as the same.
Figure 4 illustrates the mapping procedure and graphically explains the variables that need to be determined to perform the mapping such as the zenith angle θ_{mp} associated with the minimum distance to the magnetopause the distance from the planet to the bow shock r_{bs}, the distance from the planet to the magnetopause r_{mp}, the relative distance to the magnetopause α_{emp}, and the magnetosheath thickness ${\mathit{\alpha}}_{\mathrm{emp}}^{\left(\mathrm{bs}\right)}$. The position vector r, the bow shock standoff distance R_{bs}, the bow shock shape, the magnetopause standoff distance R_{mp}, and the magnetopause shape are assumed to be known in our mapping.
4.2 Setup
We begin with a position vector in the empirical magnetosheath domain, and we express the position vector as $\mathit{r}=x{\mathit{e}}_{x}+y{\mathit{e}}_{y}+z{\mathit{e}}_{z}$. Hereafter, we present the mapping procedure in the twodimensional plane spanning the x and y directions for simplicity, but the computation in three dimensions is straightforward by representing the y component of position vector in the cylindrical fashion as ρcos ϕ and the z component as ρsin ϕ using the distance ρ to the x axis. The boundaries (bow shock and magnetopause) are specified by the users and do not need to be parabolic. In this paper, we use the following bow shock and magnetopause models.

The empirical bow shock position expressed in GSE (Geocentric Solar Ecliptic) coordinates proposed and discussed by Farris et al. (1991) and Cairns et al. (1995):
$$\begin{array}{}\text{(16)}& x={R}_{\mathrm{bs}}{b}_{\mathrm{emp}}\phantom{\rule{0.125em}{0ex}}{y}^{\mathrm{2}},\end{array}$$where R_{bs} is the bow shock standoff distance, and b_{emp} is the empirical flaring parameter. We note here that the original Farris empirical bow shock model is not a paraboloid model; it is an ellipsoid model (with an eccentricity of 0.81), describing the bow shock on the dayside. It is not a proper representation of the farflank bow shock. Also, the Cairns paraboloid bow shock model does not properly represent the farflank bow shock. The distant bow shock shape approaches that of a hyperboloid.

The empirical magnetopause position by Shue et al. (1997):
$$\begin{array}{}\text{(17)}& {x}^{\mathrm{2}}+{y}^{\mathrm{2}}{\displaystyle \frac{\mathrm{4}{R}_{\mathrm{mp}}^{\mathrm{4}}}{\mathrm{4}{R}_{\mathrm{mp}}^{\mathrm{2}}{y}^{\mathrm{2}}}}=\mathrm{0},\end{array}$$in the Cartesian representation and
$$\begin{array}{}\text{(18)}& {r}_{\mathrm{mp}}={R}_{\mathrm{mp}}\sqrt{{\displaystyle \frac{\mathrm{2}}{\mathrm{1}+\mathrm{cos}\mathit{\theta}}}}\end{array}$$in the polar representation.
We use a specific exponent for the Shue model (with an alpha exponent of 0.5) in an effort to show that the analytic model is “simple”. The solar wind conditions for which this exponent is applicable is not often encountered (e.g., interplanetary magnetic field has the B_{z} component larger than +8 nT, with specific values of solar wind dynamic pressure).
In our setup, the radial distance from the planet to the bow shock is expressed as
By introducing the zenith angle θ and inserting x=r_{bs}cos θ and y=r_{bs}sin θ in Eq. (16), we obtain the equation for the radial distance to the empirical bow shock:
Equation (20) can be algebraically solved, and we take the positive value of the solution as presented in Eq. (19).
The radial distance to the magnetopause is given conveniently by Eq. (18). The Shue model reproduces the magnetopause standoff distance R_{mp} in the subsolar direction (θ=0), and the cylindrical distance asymptotes to 2R_{mp} in the tail. It is worth noting here that one needs to compute the radial distance from the planet to the bow shock or magnetopause as a function of the zenith angle when using different shapes.
4.3 Step 1: measuring the distance to magnetopause
In the first step, the distance from the given position in the magnetosheath to the nearest magnetopause is computed (see Fig. 5). We express the position vector along the magnetopausenormal direction as
where r_{mp} is the magnetopause position nearest to the position vector, and e_{mp} is the unit vector in the magnetopausenormal direction. The unit vector points away from the planet and satisfies the condition
The symbol α_{emp} is the distance to the magnetopause along the magnetopausenormal direction e_{mp} in the empirical magnetosheath.
The nearest magnetopause position is obtained by searching for the zenith angle θ_{mp} for the minimum distance from the sample position to the magnetopause. The distance D is defined as
The search for the minimum distance is implemented in a bruteforce fashion as a function of μ_{mp}=cos θ_{mp} in our study.
Using the minimum distance to the magnetopause r_{mp} and the zenith angle θ_{mp}, we are ready to compute the magnetopausenormal direction and the distance α_{emp}. To obtain the magnetopausenormal direction, we define the magnetopause shape function f_{mp} as
and we compute the normal direction by the gradient of f_{mp} as
The magnetopausenormal direction is obtained by normalizing the gradient vector $({\partial}_{x}{f}_{\mathrm{mp}},{\partial}_{y}{f}_{\mathrm{mp}})$ and representing it with the basis vectors (e_{x} and e_{y}) as
evaluated at the magnetopause (x=r_{mp}cos θ_{mp} and y=r_{mp}sin θ_{mp}). The magnetopausenormal vector e_{mp} has a unit length, and the sign ($\mathrm{sgn}=\pm \mathrm{1}$) is chosen such that the normal vector is pointing outward (Eq. 22). The distance α_{emp} to the magnetopause along the normal direction is obtained from Eq. (21) as
Equation (28) is constructed to be robust against the singular behavior on the dayside (${\mathit{e}}_{\mathrm{mp}}\cdot {\mathit{e}}_{y}=\mathrm{0}$) and in distant tail (${\mathit{e}}_{\mathrm{mp}}\cdot {\mathit{e}}_{x}=\mathrm{0}$).
4.4 Step 2: computing the thickness of empirical magnetosheath
In the second step, the magnetosheath thickness is computed using the position vector and the magnetopause normal direction (Fig. 6). For our mapping purpose, the distance α_{emp} is normalized to the magnetosheath thickness ${\mathit{\alpha}}_{\mathrm{emp}}^{\left(\mathrm{bs}\right)}$ such that the relative distance ${\mathit{\alpha}}_{\mathrm{emp}}/{\mathit{\alpha}}_{\mathrm{emp}}^{\left(\mathrm{bs}\right)}$ serves as an invariant of the mapping from the empirical magnetosheath onto the KF magnetosheath. To achieve this, we combine Eq. (16) with Eq. (21) and analytically determine the thickness from the bow shock to the magnetopause in the empirical magnetosheath. That is, the thickness ${\mathit{\alpha}}_{\mathrm{emp}}^{\left(\mathrm{bs}\right)}$ is obtained by rewriting the bow shock quadratic equation (Eq. 16) for the position vector using the variable ${\mathit{\alpha}}_{\mathrm{emp}}^{\left(\mathrm{bs}\right)}$ (Eq. 21) extended to the bow shock location. The equation is again quadratic, and the solution is algebraically obtained as:
where d_{α} is an auxiliary variable defined as
In the subsolar direction (y_{mp}=0), the thickness is simply given as
Equation (29) becomes singular in the subsolar direction, and Eq. (31) needs to be set separately.
4.5 Step 3: computing the magnetosheath thickness in the KF system
In the third step, the magnetosheath thickness is computed in the KF model (Fig. 7). We repeat the procedures of steps 1 and 2 for the KF system and determine the KF magnetosheath thickness as reference. We treat the zenith angle θ_{mp} and the relative distance ${\mathit{\alpha}}_{\mathrm{emp}}/{\mathit{\alpha}}_{\mathrm{emp}}^{\left(\mathrm{bs}\right)}$ as invariants of the mapping between the empirical magnetosheath and the KF system. The KF bow shock location is given as
where the KF bowshock flaring parameter b_{k} is prefixed as (Kobel and Flückiger, 1994)
The radial distance from the planet to the KF bow shock is
The KF magnetopause is defined in Kobel and Flückiger (1994) as
From Eq. (35), the radial distance from the planet to the KF magnetopause is computed as
To obtain the magnetopausenormal direction in the KF system, we compute the gradient of the magnetopause shape function:
The gradient is analytically given as
The magnetopausenormal direction ${\mathit{e}}_{\mathrm{mp}}^{\left(\mathrm{k}\right)}$ is then obtained by applying Eqs. (38) and (39) to Eq. (27), which reads as
The thickness in the KF system ${\mathit{\alpha}}_{\mathrm{k}}^{\left(\mathrm{bs}\right)}$ is determined by combining the bow shock shape (Eq. 32) with the position vector at the bow shock:
Equation (32) becomes again a quadratic equation with respect to the thickness ${\mathit{\alpha}}_{\mathrm{k}}^{\left(\mathrm{bs}\right)}$, and the solution reads
where the auxiliary variable ${d}_{\mathit{\alpha}}^{\left(\mathrm{k}\right)}$ is defined as
4.6 Step 4: mapping the position vector onto the KF system
In the fourth step, the mapping of the position vector is performed from the empirical magnetosheath onto the KF system (Fig. 8). An assumption is made such that the relative distance to the magnetopause along the magnetopausenormal direction is the same between the two systems. The distance from the magnetosheath position vector to the magnetopause along the magnetopausenormal direction in the KF system α_{k} is then determined by the relative distance in the empirical magnetosheath α_{emp}, the thickness of the empirical magnetosheath ${\mathit{\alpha}}_{\mathrm{emp}}^{\mathrm{bs}}$, and magnetosheath thickness in the KF system ${\mathit{\alpha}}_{\mathrm{k}}^{\left(\mathrm{bs}\right)}$ as
The mapped position vector is then computed as
using the nearest magnetopause position ${\mathit{r}}_{\mathrm{mp}}^{\left(\mathrm{k}\right)}$ (Eq. 36), the magnetosheathtomagnetopause distance α_{k} (Eq. 44), and the magnetopausenormal direction ${\mathit{e}}_{\mathrm{mp}}^{\left(\mathrm{k}\right)}$ (Eq. 40).
4.7 Step 5: evaluating the shell and connector variables
In the fifth step, the shell variable v and the connector variable u are computed from the mapped position vector r^{(k)} using Eqs. (1) and (2), respectively. The variables v and u are the same as the parabolic coordinates used in the KF potential with a focus at ${x}_{\mathrm{0}}={R}_{\mathrm{mp}}/\mathrm{2}$. In our algorithm, the focus is explicitly given in Eqs. (1), (2), and (3). The azimuthal angle around the symmetry axis ϕ is treated in the same way as in the KF model.
Figure 10 compares the isocontours of the shell v and the connector u represented in the KF system (left panel) and the empirical magnetosheath (right panel) for a bow shock standoff distance of 12.8 R_{E} (Génot et al., 2011), a bow shock flaring of $\mathrm{0.0223}\phantom{\rule{0.25em}{0ex}}{R}_{\mathrm{E}}^{\mathrm{1}}$ (Farris et al., 1991; Cairns et al., 1995), and a magnetopause standoff distance 9.8 R_{E} (Génot et al., 2011). The shell variable v is characterized by the lines with the curvature center on the right side in the panel and contains the parabolic bow shock (at v=v_{bs}) and magnetopause (at v=v_{mp}) marked by thick lines. The connector variable u has the curvature center on the left side in the panel, and the isocontour lines are orthogonal to the bow shock and magnetopause. The computation of the u and v variables and their gradient and curl is performed in Cartesian coordinates so that the connection represented by the Christoffel symbol vanishes in the computation. Computation in the Cartesian domain is also beneficial to the practical application because spacecraft trajectories are often represented in Cartesian coordinates.
4.8 Step 6: computing the potentials and stream function
The scalar potentials (velocity potential and magnetic potential) and the stream function are obtained from the shell v and the connector u using Eqs. (10), (13), and (15). The velocity potential (normalized to the upstream flow) is displayed in Fig. 11 left panel, and the stream function is displayed in the right panel. The isocontours of the velocity potential represent the lines of the same flow velocity. The isocontours of the stream function represent the streamlines in the magnetosheath. The flow is deflected around the nose of the magnetopause (the subsolar point), and the streamlines are tangential to the magnetopause. Due to the grid orthogonality around the magnetopause, the streamlines are constructed as tangential to the magnetopause shape, which qualifies the magnetopausenormal mapping method as a useful tool for the magnetosheath model.
The magnetic potential and the derived magnetic field are displayed in Fig. 12. The magnetic potential and the magnetic field (the gradient of the potential multiplied by the minus sign) depend on the upstream field. Figure 12 shows an example with an upstream field angle of 135° to the x axis (i.e., 45° to the upstream flow direction). The magnetic field is computed using the central difference scheme. Near the boundaries (bow shock or magnetopause), the mesh resolution is enhanced so that the mesh points do not cross the boundary when computing with the central difference scheme. The upstream field is deflected on the positive y side (right panel, lower half plane) and is draping the magnetopause on the negative y side (right panel, upper half plane).
5.1 Extension of the mapping approach
It is straightforward to extend our method to different shapes of the Shue magnetopause model. The following form of gradient can be used for a general value of the Shue exponent α_{mp}:
where f_{mp} in Eqs. (46)–(47) is defined as
Mesh pattern applied to different values of the Shue exponent for an opentype magnetopause α_{mp}=0.6 (Fig. 13a) and a closedtype magnetopause α_{mp}=0.4 (Fig. 13b).
Also, our method can be extended to a threedimensional, nonaxisymmetric geometry of magnetosheath (e.g., Dimmock and Nykyri, 2013). To achieve this goal, a suitable set of the variables needs to be found for the mapping: the shell variable, the connector variable, and the azimuthal angle around the subsolar axis (solar wind direction intersecting the planetary magnetic dipole). Namely, one needs to give a nonaxisymmetric bow shock shape and a nonaxisymmetric magnetopause shape; compute the magnetopausenormal direction and construct grids in the magnetosheath; measure the distance to the magnetopause and the bow shock; scale and relate the distance to the KF model; and evaluate the scalar potentials through the u, v, and ϕ variables.
5.2 Other approaches
It is possible to obtain the steadystate magnetosheath potential in different ways.

First, one may numerically solve the Laplace equation for a given set of boundary shapes (bow shock and magnetopause). Various numerical solvers are known for solving the Laplace equation such as the Jacobi method, the Gauss–Seidel method, and the successive overrelaxation (SOR) method. These Laplace solvers are numerically more expensive than the mapping method, but the computation in 3D is feasible with the contemporary computational resources. On the other hand, the magnetosheath is not bounded but extends in the tail direction. The challenge here is thus to construct a properly bounded area for the Laplace equation.

Second, one may expand the magnetosheath magnetic field in different orthogonal functions. The KF solution makes use of the Bessel functions (Kobel and Flückiger, 1994). For flexible magnetopause and bow shock boundary models, a magnetosheath magnetic field model is constructed by making use of Legendre polynomials (Romashets and Vandas, 2019).

Third, one may introduce a suitable conformal mapping by limiting the magnetosheath modeling to a complex plane (twodimensional domain). The harmonic functions such as the KF solution are transformed from the parabolic magnetosheath shape into the nonaxisymmetric magnetosheath shape. The problem here is that finding the conformal mapping is not an easy task, because the magnetosheath is not a spatiallyclosed domain and one has to set the boundary in the magnetosheath to complete the domain bounded by the bow shock, the subsolar axis, and the magnetopause.

Fourth, one may solve the Rankine–Hugoniot relations and track the streamline stepwise by referring to the KF solution, as is done in Soucek and Escoubet (2012). This method computes the magnetosheath flow and magnetic field along the streamline, and the computing for the entire magnetosheath domain is numerically expensive.
Our potential mapping method may be regarded as an updated version of the radial mapping method (Soucek and Escoubet, 2012) by retaining the orthogonality near the magnetopause in the flank to tail region and also by computing the field through the potential mapping. Velocity potential, stream function, and magnetic potential are evaluated in the empirical magnetosheath.
The advantages of our methods are as follows:

The method makes extensive use of the exact solution of the Laplace equation (the Kobel–Flückiger potential and the Guicking stream function). The plasma flow and magnetic field can be determined semianalytically in a wide range of zenith angles in the magnetosheath when the solar wind conditions and the boundary shapes are given.

The method is applicable to a nonparabolic shape of magnetosheath domain, opening the door to develop a tool to assist numerical simulations and spacecraft observations of not only the Earth but also the planetary magnetosheath domain.

The method is computationally inexpensive. In particular, if the shapes of bow shock and magnetopause are analytically given, most of the computational steps in the potential mapping method have an analytic expression.
As stated in Sect. 1, one ideally needs to find a conformal mapping from the KF magnetosheath model onto the empirical magnetosheath. While the conformal mapping is known both for the empirical bow shock and the empirical magnetopause, the conformal mapping of the entire magnetosheath domain still remains a challenge. There are two problems with this approach. First, the closing boundary (the u contours) connecting between the bow shock and the magnetopause is not known, and moreover, the uniqueness of finding such a boundary is not guaranteed. Second, the gradients along u are not the same between the empirical bow shock and the empirical magnetopause such that a naive transfinite interpolation ends up with highly nonorthogonal grids in the magnetosheath.
Our method of computing the plasma flow and magnetic field should be compared against the observations and simulations. For example, THEMIS and ARTEMIS spacecraft (Angelopoulos, 2008) and MMS spacecraft (Burch et al., 2016) are providing a huge amount of data on both sides of the bow shock in the equatorial plane; Cluster spacecraft (Escoubet et al., 2001) are collecting data in polar orbit; ACE spacecraft data (Stone et al., 1998) may be used as an upstream monitor; and Earth Flyby data of planetary missions (such as Cassini, BepiColombo) cover the fardistance tail region. In reality, nonaxisymmetric structure arises in the magnetosheath. There is no restriction regarding the choice of the magnetopause model. The magnetopausenormal direction needs to be computed either analytically using the gradient of the magnetopause function as ∇f_{mp} or numerically for a userdefined magnetopause shape.
The IDL (Interactive Data Language) codes are available upon request to the authors.
YN, DS, and ST developed the idea of the potential mapping method; checked mathematics; and wrote the manuscript. YN prepared the figures. All authors listed have made a substantial, direct, and intellectual contribution to the work and approved it for publication.
The contact author has declared that none of the authors has any competing interests.
Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.
This openaccess publication was funded by Technische Universität Braunschweig.
This paper was edited by Minna Palmroth and reviewed by Steven Petrinec, Octav Marghitu, and two anonymous referees.
Angelopoulos, V. The THEMIS mission, Space Sci. Rev, 141, 5, https://doi.org/10.1007/s1121400893361, 2008. a
Burch, J. L., Moore, T. E., Torbert, R. B., and Giles, B. L., Magnetospheric Multiscale overview and science objectives, Space Sci. Rev., 199, 5–21, https://doi/org/10.1007/s1121401501649, 2016. a
Cairns, I. H., Fairfield, D. H., Anderson, R. R., Carlton, V. E. H., Paularenas, K. I., and Lazarus, A.: Unusual locations of Earth's bow shock on September 24–25, 1987: Mach number effects, J. Geophys. Res., 100, 47–62, https://doi.org/10.1029/94JA01978, 1995. a, b, c
Chakravarthy, S. and Anderson, D.: Numerical conformal mapping, Math. Comp. 33, 953–969, 1979. a
Dimmock, A. P., and Nykyri, K: The statistical mapping of magnetosheath plasma properties based on THEMIS measurements in the magnetosheath interplanetary medium reference frame, J. Geophys. Res.Space, 118, 4963–4976, https://doi.org/10.1002/jgra.50465, 2013. a
Escoubet, C. P., Fehringer, M., and Goldstein, M.: Introduction The Cluster mission, Ann. Geophys., 19, 1197–1200, https://doi.org/10.5194/angeo1911972001, 2001. a
Farris, M. H., Petrinec, S. M., and Russell, C. T.: The thickness of the magnetosheath – Constraints on the polytropic index, Geophys. Res. Lett., 18, 1821–1824, https://doi.org/10.1029/91GL02090, 1991. a, b, c
Fornberg, B.: A numerical method for conformal mappings, SIAM J. Sci. Stat. Comput., 1, 386–400, https://doi.org/10.1137/0901027, 1980. a
Génot, V., Broussillou, L., Budnik, E., Hellinger, P., Trávníček, P. M., Lucek, E., and Dandouras, I.: Timing mirror structures observed by Cluster with a magnetosheath flow model, Ann. Geophys., 29, 1849–1860, https://doi.org/10.5194/angeo2918492011, 2011. a, b, c, d
Guicking, L., Glassmeier, K.H., Auster, H.U., Narita, Y., and Kleindienst, G.: Lowfrequency magnetic field fluctuations in Earth's plasma environment observed by THEMIS, Ann. Geophys., 30, 1271–1283, https://doi.org/10.5194/angeo3012712012, 2012. a, b, c, d
Karageorghis, A., Stylianopoulos, N. S., and Zachariades, H. A.: A numerical conformal mapping method for harmonic mixed boundary value problems, J. Sci. Comp., 11, 167–178, https://doi.org/10.1007/BF02088814, 1996. a
Kobel, E., and Flückiger, E. O.: A model of the steady state magnetic field in the magnetosheath, J. Geophys. Res., 99, 23617–23622, https://doi.org/10.1029/94JA01778, 1994. a, b, c, d, e, f, g, h
Narita, Y., Toepfer, S., and Schmid, D.: Magnetopause as conformal mapping, Ann. Geophys., 41, 87–91, https://doi.org/10.5194/angeo41872023, 2023. a
Papamichael, N., and Whiteman, J. R.: A numerical conformal transformation method for harmonic mixed boundary value problems in polygonal domains, J. Appl. Math., 24, 304–316, https://doi.org/10.1007/BF01595198, 1973. a
Romashets, E. P., and Vandas, M.: Analytic modeling of magnetic field in the magnetosheath and outer magnetosphere, J. Geophys. Res.Space, 124, 2697–2710, https://doi.org/10.1029/2018JA026006, 2019. a
Schmid, D., Narita, Y., Plaschke, F., Volwerk, M., Nakamura, R., and Baumjohann, W.: Magnetosheath plasma flow model around Mercury, Ann. Geophys., 39, 563–570, https://doi.org/10.5194/angeo395632021, 2021.
Shue, J.H., Chao, J. K., Fu, H. C., Russell, C. T., Song, P., Khurana, K. K., and Singer,H. J.: A new functional form to study the solar wind control of the magnetopause size and shape, J. Geophys. Res.Space, 102, 9497–9511, https://doi.org/10.1029/97JA00196, 1997. a, b
Soucek, J. and Escoubet, C. P.: Predictive model of magnetosheath plasma flow and its validation against Cluster and THEMIS data, Ann. Geophys., 30, 973–982, https://doi.org/10.5194/angeo309732012, 2012. a, b, c, d, e
Stone, E. C., Frandsen, A. M., Mewaldt, R. A., Christian, E. R., Margolies, D., Ormes, J. F., and Snow, F. The Advanced Composition Explorer, Space Sci. Rev., 86, 1–22, https://doi.org/10.1023/A:1005082526237, 1998. a
Tat́rallyay, M. and Erdős, G.: The evolution of mirror mode fluctuations in the terrestrial magnetosheath, Planet. Space Sci., 50, 593–599, https://doi.org/10.1016/S00320633(02)000387, 2002. a
Toepfer, S., Narita, Y., and Schmid, D.: Reconstruction of the interplanetary magnetic field from the magnetosheath data: A steadystate approach, Front. Phys., 10, 1050859, https://doi.org/10.3389/fphy.2022.1050859, 2022.
Wei, L. K., Murid, A. H. M., and Hoe, Y. S.: Conformal mapping and periodic cubic spline interpolation, Mathematika, 30, 8–20, 2014. a