NASA-CR-205115 JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 101, NO. A2, PAGES 2547-2560, FEBRUARY 1, 1996 Numerical simulations of mass loading in the solar wind interaction with Venus K. Murawski 1 and R. S. Steinolfson Department of Space Science, Southwest Research Institute, San Antonio, Texas Abstract. Numerical simulations are performed in the framework of nonlinear two - dimensional magnetohydrodynamics to investigate the influence of mass loading on the solar wind interaction with Venus. The principal physical features of the interaction of the solar wind with the atmosphere of Venus are presented. The formation of the bow shock, the magnetic barrier, and the magnetotail are some typical features of the interaction. The deceleration of the solar wind due to the mass loading near Venus is an additional feature. The effect of the mass loading is to push the shock farther outward from the planet. The influence of different values of the magnetic field strength on plasma evolution is considered. Introduction ghang et al. [1991] have used PVO data to show that the convected field gasdynamic model of Spreiter et al. [1966] predicts the correct bow shock location Data from the Pioneer Venus Orbiter (PVO), as well if the magnetic barrier is taken as the obstacle around as data from Venera 9 and 10 and other sources, have which the plasma is deflected. They also found that the contributed to a growing understanding of the inter- magnetic barrier is strongest and thinnest (about 200 action of the solar wind with an unmagnetized planet km thick) at the subsolar point and becomes weaker such as Venus. The general picture derived from the and thicker at the terminator. The magnetic barrier is available data is that, for typical solar wind conditions, responsible for transferring most of the solar wind dy- an electrically conducting ionosphere deflects the super- namic pressure to the ionosphere through the enhanced sonic and superalfv_nic solar wind around the planet. A magnetic pressure. bow shock forms upstream of the ionosphere and serves Observations at Venus have revealed that a magnetic to slow and also assists in deflecting the solar wind. tail forms as a consequence of the solar wind interac- The size of the bow shock depends on the sonic and tion with the Yenusian ionosphere. The portions of IMF alfv_nic Math numbers and on the shape and the size lines passing near the surface of Venus are slowed due of the effective obstacle, whereas asymmetries in the to the interaction with newly ionized atmospheric neu- shock shape are determined by the direction of the in- trals. The two ends of the magnetic ropes continue to be terplanetary magnetic field (IMF) [Khuvana and Kivel- pulled downstream by the solar wind [e.g., Spreiter and son, 1994]. ghang et al. [1990] have proven that at Stahara, 1980]. Slavin et al. [1989] have estimated that Venus the effective obstacle undergoes a systematic size the magnetotail extends to between 50 to 150/L, from change of 0.13 R_ (where R_ is the radius of Venus) dur- the center of the planet for a typical solar wind flow ing a solar cycle. The boundary separating the shocked with the AlfvSn speed V_ = 60 km/s and the Alfv_n and magnetized solar wind plasma from the thermal Mach number MA -- 7.2. ionosphere plasma is referred to as the ionopause, and A study by Nagy et al. [1981] showed that Venus the region between the bow shock and the ionosphere is has a dayside neutral exosphere dominated by oxygen referred to as the magnetosheath. The inner portion of at altitudes above ,_ 400 km calculated from the plane- the magnetosheath contains a region of enhanced mag- tary surface. The exospheric oxygen together with other netic pressure referred to as the magnetic barrier [e.g., particles present can be ionized by solar radiation and Luhmann, 1986]. There are physical processes involved charge exchange. Newly ionized particles are electro- in these relatively large-scale interaction regions as well magnetically coupled to the shocked solar plasma of the as in the magnetic tail that the present work has appfi- magnetosheath. As a consequence, near the ionopause cation to. the plasma becomes mass loaded due to interactions with ions, particularly O+ [e.g., Luhmann et al., 1991]. Model computations by Gombosi et al., [1980] using in- put from PVO ionosheath and neutral atmosphere pro- tAlso at Faculty of Mechanics, Politechnic of Lublin, Lublin, Poland. files indicate that approximately 10 %by number of the solar wind protons may undergo charge exchange with Copyright 1996 by the American Geophysical Union. neutrals at ionopause altitudes. As a result of the mass Paper number 95JA002433. loading (ML), the plasma is slowed and compressed. 0148-0227/96/95JA-02433505.00 When the plasma then flows around the planet into the 2547 2548 MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING magnetosheath, momentum conservation requires an in- ent physical processes than occur at magnetized Earth. crease in flow speed into the magnetotail. In addition, alarge quantity of relevant observations are It has been noted by Luhmann [1986] that the above now available for Venus. When interpreted in associa- generM behavior tends to prevail providing the thermal tion with the simulated results, this data set provides a pressure in the upper ionosphere is sufficiently larger test of the physical processes included in the model. than the solar wind pressure and the pressure balance Single-fluid MHD simulations of the global interac- occurs where the ionospheric plasma is collisionless. As tion of the solar wind with Venus without considera- the solar wind pressure increases, relative to that of tion of the ionosphere have been performed by several the ionospheric plasma, the magnetic barrier eventu- investigators. Wu [1992] limited his study to the day- ally breaks down and the height of the dayside iono- side. Tanaka [1993] included the nightside as well but sphere decreases. When the ionopause is too close only considered magnetic field orientations parallel and to the planet, the transterminator transport of iono- perpendicular to the sol_ wind flow and was in a pa- spheric plasma ceases, thereby leading to a nightside rameter regime not directly applicable to average solar - phenomenon known as the disappearing of the iono- wind conditions at Venus. Cable and Steinolfson [1995] sphere [Cravens el al., 1982]. carried out simulations for typical observed solar wind Cloutier et al. [1987] have noted that ML processes conditions at Venus and for an IMF orientation near at Venus should be asymmetric since the solar wind the Parker spiral. McGarv and Pontius [1994] studied electric field depends on the relative directions of the the effects of ML, but they considered only the dayside. flow velocity and the IMF. Phillips et al. [1987] have Moreover, they discussed the case of a cylinder with claimed that the flow and field configuration of the mag- O/Oz : 0 and B along the cylinder axis (z - axis) and netosheath plasma, together with the large gyroradius the perpendicular flow, namely, V l B. of the pickup ions, cause ML to occur preferentially on Our purpose is to extend the model of McGary and one side of the magnetosheath. More efficient pickup Pontius to study the nightside and to discuss various of newly created ions should occur over the hemisphere strengths of the magnetic field. We also considered a that produces an electric field directed outward, away more realistic spherical model. from the ionosphere. Luhrnann et al. [1991] and Phillips The paper is organized as follows. The physical model et al. [1986] have presented observations indicating an used in the present study is described in the next see- increase in the magnetic field intensity and current over tion. Section 3 describes the numerical model which is the hemisphere where the electric field is directed out- used in the present studies. We present and discuss our ward. The enhanced solar UV during solar maximum detailed results in section 4. The paper concludes with increases the ionization rate of neutrals, which moves a short summary. the bow shock out from the planet both in the equa- torial plane and at the terminator [Zhan 9 et al., 1990]. Physical Model Similarly, Alezander and Russell [1985] have shown that the position of the bow shock at the terminator de- We solve the following form of the compressible MHD pends on solar activity. This dependence may be indi- equations as an initial value problem: rect evidence of ML since, when the solar EUV is high, more mass may be added to the shocked solar wind, Op thereby forcing the bow shock to recede from the planet. 0_- + V. (pV) = q0exp [-(h - ho)/Ho], (1) Linker et al. [1989] have revealed that ML can increase or decrease the plasma temperature, depending on the value of the sonic Math number Ms. When 7 = 5/3, 0(pV) 1 Ms > v/_ is required for heating to occur. Numerical -- 0t + V. [(pV)V] = -Vp+ _(V × B) × B, (2) simulations performed show that the effects of ML on plasma temperature and velocity are most pronounced 0B - v ×(v ×B), (3) in a wake region. However, this theory was developed Ot for the case of Io's atmosphere, which is described by MHD equations with the source terms in the mass con- tinuity, Euler, and energy equations. In the case of Io, p B2 ((pV 2 the magnetic field lines are in the north-south direction + + + + lp)V and are perpendicular to the horizontal flow. Recent computer simulations of the three-dimensional 1 global interaction between the solar wind and unmag- +G[B x (V x B)]}= O, (4) netized bodies have proven to be extremely useful for improving our understanding of the associated large- v. B = 0. (5) scale processes [e.g., Wu, 1992; Tanaka, 1993; Cable and Steinolfson, 1995; Gombosi et al. 1994; McGary Here p is the plasma density. V is the velocity. B is and Pontius, 1994]. Venus is of particular interest since the magnetic field, p is the plasma pressure, and q0 is the planet does not have a significant intrinsic mag- the production rate at areference altitude ho. Ho is the netic field, and the interaction of the solar wind with scale height of the hot neutral oxygen population, and the Venusian ionosphere involves fundamentally differ- the ratio of specific heats is 7 : 5/3. MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING 2549 In the present first-order approximation, the newly ner that it can be applied to a single r, 0 plane when created photoions are considered as cold (T = 0) and the phenomenon being studied is cylindrically symmet- motionless (V = 0) particles. A more rigorous approach ric about the poles at 0 = 0 and/9 = 180 °, as is the case would be to add corresponding source terms to the Eu- when the magnetic field is parallel to the solar wind ler (2) and energy (4) equations [e.g., Gombosi et al., 1994]. However, it was shown by Biermann et al. [1967] ....... 0.1, no load that the major effect is that associated with the mass -- 0.1, load continuity equation and that contributions to the mo- ..... 0.5, no load 0.5, load S-I mentum and energy density are small. Therefore in 3.8 --- 10, no load these preliminary calculations source terms in the mo- 10, load mentum and energy equations are neglected. The ML C) II term is given by the right-hand side of the mass continu- ity equation (1). Oxygen photoionization is considered 4-0-,2',8 (1) to be the only mass source with an altitude distribution t'- described by the exponential function. This production -o /'/; rate is based on PVO observations for solar maximum, I'l: 1.8 where the hot neutral oxygen atmosphere had a scale height H0 = 400 km in the subsolar region and a refer- I ] ence production rate of q0 = 3x 105 cm-3s -1 at altitude ";,,;....-.,.......,-'T_7,".,....,........,...... ha = 400 km [Belotserkovskii et al., 1987]. As numerical 0.8 .4 1.3 1.2 1.1 1.0 simulations, which were made for various mass addition f rates, revealed that the boundary layer develops when oxygen ion production rate is qo = 3 × 10_ cm-3s -t 2.8 [McGary, 1993; McGary and Pontius, 1994], we carried 1_ ....... 0.1, no load I/I1 -- 0.1, load on our simulations for this value of qo- I1 I_ ..... 0.5, no load I/ t| _ 0.5, load .r..2.3 If A'_ --- 10, no load Numerical Model //il -- 10,oad These numerical calculations were performed with HEMIS3D, a three - dimensional ideal MHD code. The code utilizes the modified Lax-Wendroff scheme devel- oped by Rubin and Burstein [1967], which yields stable solutions by adding numerical diffusion to the scheme. The dissipative terms (kept as small as possible) are included to help damp out short-wavelength ripples generated by numerical dispersion and numerically in- duced reflections from the simulation boundaries, while 0._ t I I I II II1.I4 J I II I I I I I_I*gI ¢ I I I II I I_.4I J II I I I I II_.gI J I I I leaving the longer-wavelength phenomena minimally af- r fected. The equations are solved in all three spatial di- 2.7 mensions of a spherical (r, 0, ¢) coordinate system in --'_',, -...... 0.1, no load which the O= 0 axis is directed toward the Sun or into -. \', --0.1, load -"{, X" -.... 0.5, no load the solar wind flow. The generally small deviations of ", \',, _ 0.5 load v _,, --- 10, no load the solar wind flow from the radial direction (with re- "0"2.2 ', "Q, _ 10. load spect to the Sun) are neglected so the solar wind flow 03 \ impacting the Venusian ionosphere can be assumed par- II allel to the Sun-Venus line (0 = 0). The spherical co- I_ 1.7 ordinate system is centered on the planet so r = R_ 4.qJ.) represents the planetary surface. t- Although a grid defined in spherical coordinates in- X2_ troduces some computational complexity, it has a num- 1.2 ber of advantages for studying flow past planetary ob- stacles such as Venus. One obvious advantage is that the boundary at the planetary surface occurs at a fixed 0,7 IIIII' ' I'II'II' 'IIIIItIIItI' ''I' It' IIZ'_I'I"I_ I' 0.9 1.4 1.g 2.4 2.g value of the radial coordinate, which avoids the compli- r- cations of interpolating in a Cartesian coordinate sys- tem. Another distinct advantage is that as the center Figure 1. Plots of the normalized mass density profiles along 0 of Venus is located at the center of the spherical coor- = const lines for 13= 0.1, [3 = 0.5, 13= 10 and for the NL and dinate system, high resolution near Venus is obtained ML simulations. (a) the Sun-Venus (0 = 0) line, (b) the without the expense of high resolution everywhere in terminator (0 = 90*) line, and (c) the 0 = 180" line. Mass the simulation. The code is constructed in such a man- loading shifts the bow shock farther away from Venus. 2550 MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING flow. In this case, the simulation code uses a 73 x 74 grid very close to those presented in this paper. For the cor- (r, 0 grid dimensions), which corresponds to A8 = 2.5, responding discussion, see also Cable and Steinolfson with a radial extent of the simulations of 3Rv. The [1995]. largest Ar -- 0.083Rv. We have checked that the re- The time step limit for methods using explicit tempo- sults for larger extend of the numerical boundaries were ral differencing for advective problems is the well-known 80.0 6.0 ,o°° iL 60.0 10, load /1//_i / ,.0 ,o.o // tii _'_ L ......... 0.1, no load _"_', C).-20.0 ]- ..... 0.5, no load '_, [ _ 0.5, load "_%.. l---- 10. no load "_ > 2.0 I_-- 0.1, load _ _l,,\ [ _------- 10, load 0.0 ,"_5.................... ,.3 1.1 ' .... 0.9 0.01.4 1.3 1.2 1.1 1.0 r- r- 7.0 -- 0.1, load 15.0 ..... 0.5, no load I j_ ....... 00..15,, nlooadload -- 10, no load .--. b) --10, ,oad ,oo 6.0 _._ . "so[- i ........-0._, .o ,oad _ 5.0 ,4 "',. b) '_'5.0 I- ! -- 0.1, load / """ J" I ..... 0.5, no load ! _ 0.5, load ]- / --- 10, no load "q o '% ,., ,.g 2., 2.o 1.5 2.0 2.5 3.0 F F- 5.0 . .-...'7.'. 4o ....... 0.1, no load .I ..-o --'-. -- 0.1, load f.%._...... ..... 0.5, no load 4.0 A ...'" "'', _ 0.5, load ¢....- ,,,--,,30 ",, --- 1100,, nlooadload C3 (,¢'./'" C) oO w- 3.0 yf....... ,--°0 II •' j_.-- II ID //.." /S c)_l"" __ 2o _,__2_.0 ////:' _ ......... 0.I, no load _(- > I/ /.,f -- 6:i', io_i--- o-,0 1.0 [/ _ ..... 0.5, no load ,,,, jr -- 0.5, load : _ --- 10, no load f" _ 10, load 0.q_ ........................................ 8.'9ii,,,,,,ii,,,,,,,,,I,*,,,,,1¢iI,,ii,,,,,,I,,,, I._ z, zo . 1.5 2.0 2.5 3.0 F Figure 2. As in Figure 1 but for the normalized velocity Figure 3. As in Figure 1 but for the normalized pressure profiles, profiles. MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING 2551 Plate 1. Normalized mass densities (p/p_) for the MHD non loaded (top) and mass loaded (bottom) simulations for the ease of a parallel (nllv) magnetic fleld and the plasma fl = 0.5. Here, poo is the density of the solar wind. The effect of mass loading (ML) is to increase the mass density in the overall region, except the tail region in which the mass density is reduced. poo denotes the solar wind mass density. CFL condition, At __AI/[V], where the cell size is AI in more detail by Cable and Steinolfson [1995]. and V is the fastest physical speed in the system. The We used a simple model for the boundary conditions resulting computer code [Cable and Steinolfson, 1995] that is physically motivated and performs adequately is second-order accurate in space and time. well. Venus is modeled as a hard, highly electrically Defining a spherical boundary for the problem is conducting sphere. The inner boundary is taken at the easily accomplished. However, solving the equations planetary surface (1 P_) where the radial components of of magnetohydrodynamics in spherical coordinates re- both the magnetic field and velocity axe set to zero. All quires some modification of the usual Lax-Wendroff other physical variables are computed by extrapolation scheme. In particular, the singular points at 0 = of values along a radial line. At the dayside of the outer 0 and 0 = 180° are handled through application of radial boundary, all quantities are set to the solar wind L'Hospital's rule to extrapolate values from grid points values. The values at the nightside of the outer radial next to the singular points. Since the value of each vari- boundary are determined by extrapolation along the able at the singular points may vary depending on the local flow direction. Moreover, we have never encoun- angle _bat which the singular points are approached, a tered transient conditions where the flow was coming flux-weighted average of the values at all angles is used in from the outer radial boundary. A typical compu- for the final updated value. This procedure is described tation begins with the introduction of the desired solar 2552 MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING Plate 2. Normalized plasma flows (V/aoo) for the MtID non loaded (top) and mass loaded (bottom) simulations for the case of a parallel (BllV) magnetic field and the plasma/3 = 0.5. Here, ado is the sound speed of the solar wind. The effect of ML is to decrease the plasma flow. wind values in the dayside within the numerical box. is initialized in a similar way, except that B is skewed The initial IMF is taken from the potential solution for with respect to the solar wind flow by the desired angle. flow over a sphere so that V. B = 0 in the initial state. In the present paper, this angle is 0°. This is clearly The continued satisfaction of V •B = 0 is assured by not a frequently occurring case [Phillips et al., 1986], solving a Poisson equation and updating the magnetic but it is a reasonable first-order approximation which field. The numerical solution then continues until the will serve as a reference for comparison with results for interaction achieves an approximate steady state. cone angles of 45° and 90°. These processes will be The density and pressure are initially constant through- realistically modeled in future calculations. out the simulation region. The initial velocity field is the flow of an incompressible fluid over a sphere, di- rected along the 0 = 0 line: Results and Discussion v,(,, o, t= o)= -roll - cosO, (6) We present all numerical results for the oxygen ion production rate q0 = 3 x 105 cm-as -1 [McGary and Pontius, 1994] and discuss several values of the plasma V,O', O,qLt = 0) = V0[1 + _(-_)s] sin 0, (7) /3 = 87rpo/B_, namely/3 = 0.1,/3 = 0.5,/3 = 10, and a typical value of/3 = 1.5. The case of/3 = 0.1 (/3 = 10) corresponds to the strong (weak) magnetic field case. v4,(,,o,C,t = o)= o, (8) The physical solar wind values used are: electron den- where V0 is the solar wind speed. The magnetic field sity _tt = 20 cm-s, temperature T = 105 °K, sound MUI AWgAI(II DgTt INOLFgOgIIqM! ULAq IOONPgMAggLOADII C 255 Plate 3. Normalized plasma pressure (p/(pooaL)) for the MHD non loaded (top) and mass loaded (bottom) simulations for the case of a parallel (BItV) magnetic field and the plasma /_= 0.5. The effect of ML is to decrease the plasma pressure. speed 55 km/s, flow velocity 370 kin/s, sonic Maeh num- [e.g., Linker et al., 1989; McGary and Pontius, 1994]. ber 6.73, and thermal pressure 6.06. 10-1° dyn/cm 2. However, a more detailed close-up reveals that in the Other physical parameters are chosen the same as in taft and terminator regions the mass density is reduced the work by Cable and Steinolfson [1995]. by the ML effects (Figure 1). Similar reduction is ob- As a first step to apply the methodology to the prob- served in the case of/3 =- 0.5 along the Sun-Venus line, lem of ML in the atmosphere of Venus, we consider an at r __ 1.07 R_ (Figure la). The density reduction axisymmetric ease with the magnetic field lines parallel does not occur at the terminator for /_ = 0.1, while to the solar wind flow. Consequently, B is parallel to this effect takes place for/_ = 0.5 and /_ --- 10 (Fig- the Sun-Venus line and the entire solution is rotation- ure lb). The density depletion is stronger for a higher ally symmetric around this line, a/0$ = 0. We used value of/_. Just at the planetary surface, at r < 1.3 the numerical technique described in the previous sec- R_, the mass density increase is observed (Figure lb). tion to obtain a solution of the set of (1)-(5). The flow In the magnetotail, the plasma density reduction occurs is from left to right and all physical quantities are either at the planetary surface for all values of the plasma 13. symmetric or antisymmetric about the lower boundary. Again, the rate of density depletion is proportional to As shown in Plate 1, the mass density is increased in ft. The case offl -- 0.1 differs from the other cases just the case of mass loaded (ML) solutions. This finding as at the planetary surface p is increased, but farther is consistent with a simple analysis of mass continu- on there is a region of mass depletion which is followed ity equation (1) and with previously published results by a region of mass enhancement (Figure lc). To our 2554 MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING Plate 4. Normalized plasma temperature (T/T_) for the MHD non loaded (top) and mass loaded (bottom) simulations for the case of a parallel (BttV) magnetic field and the plasma /3= 0.5. The effect of ML is to decrease the plasma temperature. knowledge, this is the first report of mass depletion due the reduction of a velocity magnitude because conserva- to the ML. The density depletion is a consequence of tion of momentum and energy requires that the plasma the nonlinear term V. (pV), which is smaller due to flow is decelerated as a consequence of the mass density the plasma flow reduction by the ML effects. Indeed, reduction. Plate 2presents results both for the ML and the largest density (velocity) increase (decrease) occurs no loading (NL) simulations. By examining the two so- in the magentotail at r = 3 Re for/3 = 0.1. The nor- lutions, one can see a fundamental difference between realized density jumps from PNL -- 1.25 to PML _--1.5 them. The effect of ML is to decrease the flow. A de- (Figure lc), while the magnitude of plasma flow is re- creasing rate is highest in the tail region where slow duced from VIvL _- 4.9 to VML -- 4 (Figure 3c). Con- flows in the case of NL are much more reduced by ML. sequently, pNLVNL = 6.125 > PMLVML = 6. Therefore The decreasing rate is larger for a smaller value of/3 at the effect associated with this nonlinear term overcomes the terminator and in the tail region. The plasma starts the density increasing effect exerted by the exponential decelerating because of ML well ahead of the shock. In term q0exp [-(h- ho)/Ho]. The density reduction in the subsolar region, this rate is largest for/3 = 0.5. A the tail may also be due to the changed shape of the large flow deceleration occurs also at the terminator, bow shock. A different bow shock shape may result in dose to the planetary surface where the plasma flow less mass being carried into the tail, regardless of other drops from V_vL _- 5.7 to VML _--4.2 (Figure 2b). In- factors. spection of Plate 2 reveals the global structure of the It follows from Euler equation (2) that we can expect shock and the stagnation region behind it. The flow MURAWSKI AND STEINOLFSON: SIMULATIONS OF MASS LOADING 2555 decelerates at the shock as approaching solar wind ab- consistent with PVO observations [Mihalov et al., 1982] sorbs an increasing number of ions. At the planetary and the gasdynamic results by McGary [1993]. surface, _' - P_, there is no perpendicular flow as a con- From (4) we can derive an equation which governs sequence of the boundary conditions imposed. Velocity the evolution of the NL plasma pressure p shears in the radial direction are present at the bow 8p shock region both for the ML and NL flows, which is _- + V. (pV) = (I- qc)pV. V. (9) 20.0 ....... 0.1, no load ----__'S,.-2-=C -- 0.1, load 0.8 ..... 0.5, no load 0.5, load 15.0 --- 10, no load 10, load _" 0 II o _ 10.0 I1) C" .._ 0.4 m °°' t-- 5.0 0.2 ....... 0.5, no load _, 0.5, load _,_ ---- 10, no load _1 0.0 .4 1,.3......... 1,.2......... 1,.1......... 1,.0..... 0.01 5 1.4 1.3 1.2 1.1 1.0 F F ....... 0.1, no load .........0.I, no load -- 0.1, load -- 0.I, load 8.8 ..... 0.5, no load 2.2 I I ..... 0.5, no load 0.5, load I 0.5, load --- 10, no load --- 10, no load 0 _ 10, load 0 10, load (3) (J) 1.7 fl II4.8 0 0 (]) [-- C" 1.2 I-- 2.8 CD u) 0.7 i b) .......................................... 0';_3.9 1.4 1.g 2.4 2.9 . 1.5 2.0 2.5 3.0 F r- i C/_ ....... 0.1, no load 0.4 ......... 0.1, no toad ..J -- 0.1, load -- 0.I, load _ ..... 0.5, no load ..... 0.5, no load / 12.8 .............. _ 0.5, load -- 0.5, load "'-... --- 10, no load ,_ --- 10, no load / , 0 "'-- _ 10, load 0 .--- 1o,,oad / / --..... -,....... 0o°._ ..f",," o,/ t-.- " "" II 1./° II s.8 o _ _ -¢-" 0.2 soJ'°*J° C) f- --iil: -4--' 03 0.1 _ 4.a _ °'a.o...;..;....i..'..;....;....i....i.....i....;..... o..o/'._o_,,_,_mm__ .'T,'_.,'.;._..'..... 1.5 2.0 2.5 3.0 F F Figure 4. As in Figure 1 but for the normalized temperature Figure 5. As in Figure 1 but for the normalized magnitude of the profiles, magnetic field profiles. 2556 MURAWSKAINDSTEINOLFSOSNIM: ULATIONOSFMASSLOADING As the plasma flow velocity V is decreased by the ML 2.20 effects we can expect a pressure reduction. The simula- > tion results shown in Plate 3 are qualitatively consistent 052.10 with our expectations. However, Plate 3 shows that the f- plasma pressure is increased near the terminator region, _'2.00 in the region 0 < 90°. This increase is caused by a c- nonlinear interaction which is described by the V. (pV) •0 1.90 and pV.V terms. Our results differ from the conclusion .,J 03 made by McGary and Pontius [1994], who claimed that 0 _.8o CL the dynamics of the plasma flow are not significantly :\ .._+... ..........++..o.+_°..._o affected by ML. A reason for the difference between the _:_1.70 _°..- present results and those of McGary and Pontius has 4J E_ to do with the different models applied; MeGary and c- 1.60 *" °no load Pontius performed gasdynamic (with B = 0) simula- load E tions for the cylinder flow, whereas our model concerns 15°I.6llllt ,,IlI_IIL.LIoI'IIIIIIII'III 7.0 11.0 MHD and a plasma sphere. beta The above scenario may differ for different Values of /9. It occurs in the case of a sufficiently small value of/3 that the magnetic field is so strong that the fast shock Figure 6. Terminator positions for different values of the plasma becomes an intermediate one [Steinolfson and Hund- 13and for the non loaded and mass loaded simulations. The bow hausen, 1990]. For example, the ML increases the gas shock at the terminator is shifter farther away farther away from the surface of Venus. pressure along the whole Sun-Venus line in the case of /3 : 0.1 and/3 : 10 (Figure 3a) but reduces it in the case of/3 = 0.5,/3 : 1, and/3 : 1.5 at the planetary case of/9 = 0.5 for the ML curve. Now, however, ML surface. In the bow shock region, an enhancement of cools the plasma at the planetary surface but warms it p due to the ML is observed for all values of/3. The in the bow shock region. The temperature increase at strongest enhancement seems to occur for/3 : 0.5. At the bow shock was observed recently by Gombosi et al. the terminator (Figure 3b), the gas pressure is increased [1994] in their simulations of the interaction of an ex- by ML. The rate of decrease is highest at the planetary panding cometary atmosphere with the solar wind. At surface for/3 : 0.1, while at larger distances the influ- the terminator, the plasma is warmed by ML (Figure ence of ML on p is negligible. The pressure reduction at 4b). The warming is more significant at the planetary the planetary surface is caused by the right-hand side surface and is higher for a higher value of/9. The cooling term of (9) which becomes small at the surface because occurs at the magnetotail (Figure 4e). The cooling rate of the flow reduction. The bow shock is closer to the is higher at the surface of Venus and for a smaller value planetary surface for a higher value of/3, in accordance of 19. This is a consequence of the pressure decrease with our expectations. At the magnetotail the pressure at the planetary surface (Figure 3c). The mass den- is reduced (Figure 3c). The reduction rate is highest at sity is increased for/3 = 0.1 and decreased for/3 = 0.5 the planetary surface for the smallest value of/3 and is and /3 -- 10, but the decreasing rate is smaller than a monotonic function of/3. the pressure decreasing rate (Figure lc). Consequently, The next point of comparison between the ML and T ._ p/p is much decreased at the surface of Venus. The NL results is the temperature distribution around Venus. cooling rate is much more significant than the warming As for the ideal gas the plasma temperature T is pro- effect. The bow shock position is closest to the plane- portional to pip so both p and p play a major role in tary surface for the smallest presented values of/3. The determining its behavior. We already learned that for bow shock distance from Venus grows with/3. /3 : 0.5 p (p) is increased (decreased) by the ML effects Plate 5 shows the magnitude of the magnetic field. everywhere except in the terminator and tail (nose) re- In the subsolar region there is no jump in the magni- gions. Therefore, we may expect the general behavior of tude of the magnetic field across the bow shock. At the T would be to decrease its magnitude by the ML effects. flanks, however, the magnetic field increases through Indeed, Plate 4 indicates that the ML significantly cools the bow shock. The magnetic tail is a region of reduced the flow relative to the NL case. The cooling occurs as magnetic field. It can be seen that the ML exerts some zero-temperature oxygen ions are introduced into the influence on the magnetic field strength• Small differ- flow and thermalized [McGary and Pontius, 1994]. The ences occur at the sunward side of the shock where the behavior of T is interesting to discuss in more detail as magnetic field is decreased by the ML effects (Figure a function of the radial direction r at different locations 5a). This effect is small for/5 = 0.1 and grows with of 0. In the case of 0 : 0, for/3 = 0.1 and/3 : 10 /3. The decrease in the magnetic field strength nearest the temperature is decreased by ML both at the plane- the planet along the stagnation streamline, is the re- tary surface and at the sunward side of the bow shock sult of the boundary conditions (B,(r : Rv, O,t) : O) (Figure 4a), while in the shock the effect of ML is neg- imposed together with the fact that the parallel shock ligible. The curves which correspond to the ML possess does not change a value of the magnetic field across maxima in the magnetosheath, while no maxima exist the shock. As the planet is conducting there can be for the NL case. A similar maximum occurs also in the no radial magnetic field component at the surface of