Convection of a fluid with strongly temperature and pressure dependent viscosity

Plate tectonics on the Earth is a surface manifestation of convection within the Earth’s mantle, a subject which is as yet improperly understood, and it has motivated the study of various forms of buoyancy-driven thermal convection. The early success of the high Rayleigh number constant viscosity theory was later tempered by the absence of plate motion when the viscosity is more realistically strongly temperature dependent, and the process of subduction represents a continuing principal conundrum in the application of convection theory to the Earth. A similar problem appears to arise if the equally strong pressure dependence of viscosity is considered, since the classical isothermal core convection theory would then imply a strongly variable viscosity in the convective core, which is inconsistent with results from post-glacial rebound studies. In this paper we address the problem of determining the asymptotic structure of high Rayleigh number convection when the viscosity is strongly temperature and pressure dependent, i.e. thermobaroviscous. By a method akin to lid-stripping, we are able to extend numerical computations to extremely high viscosity contrasts, and we show that the convective cells take the form of narrow, vertically-oriented fingers. We are then able to determine the asymptotic structure of the solution, and it agrees well with the numerical results. Beneath a stagnant lid, there is a vigorous convection in the upper part of the cell, and a more sluggish, higher viscosity flow in the lower part of the cell. We then offer some comments on the possible meaning and interpretation of these results for planetary mantle convection.     

Layered convection with an interface at a depth of 1000 km: stability and generation of slab-like downwellings

We investigate the stability of hypothetical layered convection in the mantle and the mechanisms how the downwelling structures originating in the lower layer are generated. The stability is studied by means of numerical simulations of the double-diffusive convection in a 2D spherical model with radially dependent viscosity. We demonstrate that the stability of the layering strongly depends not only on the density contrast between the layers but also on the heating mode and the viscosity profile. In the case of the classical Boussinesq model with an internally heated lower layer, the density contrast of about 4% between the compositionally different materials is needed for the layered flow to be maintained. The inclusion of the adiabatic heating/cooling in the model reduces the temperature contrast between the two layers and, thus, enhances the stability of the layering. In this case, a density contrast of 2-3% is sufficient to preserve the layered convection on a time scale of billions of years. The generation of the downwelling structures in the lower layer occurs via mechanical or thermal coupling scenarios. If the viscosity dependent on depth and average temperature at each depth is considered, the low viscosity zone develops at a boundary between the two convecting layers which suppresses mechanical coupling. Then the downwelling structures originating in the lower layer develop beneath upper layer subductions, thus resembling continuous slab-like structures observed by seismic tomography.     

A matrix-dependent transfer multigrid method for strongly variable viscosity infinite Prandtl number thermal convection

Abstract

We apply a two-dimensional Cartesian finite element treatment to investigate infinite Prandtl number thermal convection with temperature, strain rate and yield stress dependent rheology using parameters in the range estimated for the mantles of the terrestrial planets. To handle the strong viscosity variations that arise from such nonlinear rheology in solving the momentum equation, we exploit a multigrid method based on matrix-dependent intergrid transfer and the Galerkin coarse grid approximation. We observe that the matrix-dependent transfer algorithm provides an exceptionally robust and efficient means for solving convection problems with extreme viscosity gradients. Our algorithm displays a convergence rate per multigrid cycle about five times better than what other published methods (e.g., CITCOM of Moresi and Solomatov, 1995) offer for cases with similar extreme viscosity variation. The algorithm is explained in detail in this paper.

When this method is applied to problems with temperature and strain rate dependent rheologies, we obtain strongly time dependent solutions characterized by episodic avalanching of cold material from the upper boundary layer to the bottom of the convecting domain for a significantly broad range of parameter values. In particular, we observe this behavior for the relatively simple case of temperature dependent Newtonian rheology with a plastic yield stress. The intensity and temporal character of the episodic behavior depends sensitively on the yield stress value. The regions most strongly affected by the yield stress are thickened portions of the cold upper boundary layer which can suddenly become unstable and form downgoing diapirs. These computational results suggest that the finite yield properties of silicate rocks must play a vitally important role in planetary mantle dynamics. Although our example calculations were selected mainly to illustrate the power of our multigrid method, they suggest that many possible exotic behaviors in planetary mantles have yet to be discovered.     

Evaluation of stresses in two geodynamically different areas: Stable foreland and extensional backarc

Areas which are geodynamically different have different behaviors both in their thermal regime and seismic activity. A stable area has a geotherm which can be considered as standard, extensional and compressional areas have, respectively, high and low temperature gradients. The Italian region includes different geodynamical areas and all such situations are present. We consider the Apulian platform as an example of a stable area and the Tuscany-Latium as an example of an extensional area. For both of them the present geotherms are calculated, taking into account, for the Tuscany-Latium, its thermal history. Assuming that each region is subject to a constant strain rate, the stresses are calculated as functions of depth and time. The rheological behavior is assumed to be linear viscoelastic, with viscosity dependent on temperature and elastic parameters dependent on lithology. The geothermal profile and the rheological structure of the lithosphere remarkably affect the processes of stress accumulation which control the distribution of seismic activity. The abrupt decrease of the temperature gradient at the Moho produces considerably higher stress values with respect to the case of uniform gradient, thus favoring subcrustal seismicity. In the case of a standard temperature gradient, subcrustal seismicity is predicted and a gap in seismicity, indicating a soft intracrustal layer, exists if there is a discontinuity in rheology. By contrast, in the case of a high-temperature gradient, subcrustal seismicity is not to be expected, even in the presence of a discontinuity in rheology, since subcrustal temperatures are already too high to permit a sufficient stress accumluation.     

Heat transport by variable viscosity convection and implications for the Earth's thermal evolution

The case is presented that the efficiency of variable viscosity convection in the Earth's mantle to remove heat may depend only very weakly on the internal viscosity or temperature. An extensive numerical study of the heat transport by 2-D steady state convection with free boundaries and temperature dependent viscosity was carried out. The range of Rayleigh numbers (Ra) is 10⁴?10⁷ and the viscosity contrast goes up to 250000. Although an absolute or relative maximum of the Nusselt number (Nu) is obtained at long wavelength in a certain parameter range, at sufficiently high Rayleigh number optimal heat transport is achieved by an aspect ratio close to or below one. The results for convection in a square box are presented in several ways. With the viscosity ratio fixed and the Rayleigh number defined with the viscosity at the mean of top and bottom temperature the increase of Nu with Ra is characterized by a logarithmic gradient β = ?ln(Nu)/? ln(Ra) in the range of 0.23–0.36, similar to constant viscosity convection. More appropriate for a cooling planetary body is a parameterization where the Rayleigh number is defined with the viscosity at the actual average temperature and the surface viscosity is fixed rather than the viscosity ratio. Now the logarithmic gradient β falls below 0.10 when the viscosity ratio exceeds 250, and the velocity of the surface layer becomes almost independent of Ra. In an end-member model for the Earth's thermal evolution it is assumed that the Nusselt number becomes virtually constant at high Rayleigh number. In the context of whole mantle convection this would imply that the present thermal state is still affected by the initial temperature, that only 25–50% of the present-day heat loss is balanced by radiogenic heat production, and the plate velocities were about the same during most of the Earth's history.     

Effect of coastal boundary resolution and mixing upon internal wave generation and propagation in coastal regions

A non-linear two-dimensional vertically stratified cross-sectional model of a constant depth basin without rotation is used to investigate the influence of vertical and horizontal diffusion upon the wind-driven circulation in the basin and the associated temperature field. The influence of horizontal grid resolution, in particular the application of an irregular grid with high resolution in the coastal boundary layer is examined. The calculations show that the initial response to a wind impulse is downwelling at the downwind end of the basin with upwelling and convective mixing at the opposite end. Results from a two-layer analytical model show that the initial response is the excitation of an infinite number of internal seiche modes in order to represent the initial response which is confined to a narrow near coastal region. As time progresses, at the downwind end of the basin a density front propagates away from the boundary, with the intensity of its horizontal gradient and associated vertical velocity determined by both horizontal and vertical viscosity values. Calculations demonstrate the importance of high horizontal grid resolution in resolving this density gradient together with an accurate density advection scheme. The application of an irregular grid in the horizontal with high grid resolution in the nearshore region enables the initial response to be accurately reproduced although physically unrealistic short waves appear as the frontal region propagates onto the coarser grid. Parameterization of horizontal viscosity using a Smagorinsky-type formulation acts as a selective grid size-dependent filter, and removes the short-wave problem although enhanced smoothing can occur if the scaling coefficient in the formulation is too large. Calculations clearly show the advantages of using an irregular grid but also the importance of using a grid size-dependent filter to avoid numerical problems.     

Testing density-dependent groundwater models: two-dimensional steady state unstable convection in infinite,finite and inclined porous layers

This study proposes the use of several problems of unstable steady state convection with variable fluid density in a porous layer of infinite horizontal extent as two-dimensional (2-D) test cases for density-dependent groundwater flow and solute transport simulators. Unlike existing density-dependent model benchmarks, these problems have well-defined stability criteria that are determined analytically. These analytical stability indicators can be compared with numerical model results to test the ability of a code to accurately simulate buoyancy driven flow and diffusion. The basic analytical solution is for a horizontally infinite fluid-filled porous layer in which fluid density decreases with depth. The proposed test problems include unstable convection in an infinite horizontal box, in a finite horizontal box, and in an infinite inclined box. A dimensionless Rayleigh number incorporating properties of the fluid and the porous media determines the stability of the layer in each case. Testing the ability of numerical codes to match both the critical Rayleigh number at which convection occurs and the wavelength of convection cells is an addition to the benchmark problems currently in use. The proposed test problems are modelled in 2-D using the SUTRA [SUTRA––A model for saturated–unsaturated variable-density ground-water flow with solute or energy transport. US Geological Survey Water-Resources Investigations Report, 02-4231, 2002. 250 p] density-dependent groundwater flow and solute transport code. For the case of an infinite horizontal box, SUTRA results show a distinct change from stable to unstable behaviour around the theoretical critical Rayleigh number of 4π² and the simulated wavelength of unstable convection agrees with that predicted by the analytical solution. The effects of finite layer aspect ratio and inclination on stability indicators are also tested and numerical results are in excellent agreement with theoretical stability criteria and with numerical results previously reported in traditional fluid mechanics literature.     

Dielectrophoretic,thermal instability in a spherical shell of fluid

Abstract

This paper is concerned with the dielectrophoretic instability of a spherical shell of fluid. A dielectric fluid, contained in a spherical shell, with rigid boundaries is subjected to a simultaneous radial temperature gradient and radial a.c. electric field. Through the dependence of the dielectric constant on temperature, the fluid experiences a body force somewhat analogous to that of gravity acting on a fluid with density variations. Linear perturbation theory and the assumption of exchange of stabilities lead to an eighth order differential equation in radial dependence of the perturbation temperature. The solution to this equation, satisfying appropriate boundary conditions, yields a critical value of the electrical Rayleigh number and corresponding critical wave number at which convective motion begins. The dependence of each critical number is presented as a function of the gap size and temperature gradient. In the limit of zero shell thickness both the critical Rayleigh number and critical wave number agree with results for the case in the infinite plane problem.     

A simple model of convection in the terrestrial planets

Abstract

In this paper we study analytically the simplest fluid mechanical model which can mimic the convective behavior which is thought to occur in the solid mantles of the terrestrial planets. The convecting materials are polycrystalline rocks, whose creep behavior depends very strongly on temperature and probably also on pressure. As a simple model of this situation, we consider the flow of a Newtonian viscous fluid, whose viscosity depends strongly on temperature (only), and in fact has an infinite viscosity below a certain temperature, and a constant viscosity above this temperature. This model would also be directly relevant to the convection of a melt beneath its own solid phase (e.g. water below ice, though in that case there are other physical complications).

As a consequence of this assumption, there is a vigorous convection zone overlain by a stagnant lid, as also observed in analogous laboratory experiments (Nataf and Richter, 1982). The analysis is then very similar to that of Roberts (1979), but the extension to variable viscosity introduces important differences, most notably that the boundary between the lid and the convecting zone is unknown, and not horizontal. The resulting buoyancy induced stresses near this boundary are much larger than the stresses produced by buoyancy in the side-wall plumes, and mean that the dynamics of this region, and hence also the heat flux, are independent of the rest of the cell. We give a first order approximation for the Nusselt number-Rayleigh number relationship.     

The structure of 2D mantle convection and stress fields: Effects of viscosity distribution

Spatial fields of temperature, velocity, overlithostatic pressure, and horizontal stresses in the Earth’s mantle are studied in two-dimensional (2D) numerical Cartesian models of mantle convection with variable viscosity. The calculations are carried out for three different patterns of the viscosity distribution in the mantle: (a) an isoviscous model, (b) a four-layer viscosity model, and (c) a temperature- and pressure-dependent viscosity model. The pattern of flows, the stresses, and the surface heat flow are strongly controlled by the viscosity distribution. This is connected with the formation of a cold highly viscous layer on the surface, which is analogous to the oceanic lithosphere and impedes the heat transfer. For the Rayleigh number Ra = 10⁷, the Nusselt number, which characterizes the heat transfer, is Nu = 34, 28, and 15 in models with constant, four-layered, and p, T-dependent viscosity, respectively. In all three models, the values of overlithostatic pressure and horizontal stresses σ _xx in a vast central region of the mantle, which occupies the bulk of the entire volume of the computation domain, are approximately similar, varying within ±5 MPa (±50 bar). This follows from the fact that the dimensionless mantle viscosity averaged over volume is almost similar in all these models. In the case of temperature- and pressure-dependent viscosity, the overlithostatic pressure and stress σ _xx fields exhibit much stronger concentration towards the horizontal boundaries of the computation domain compared to the isoviscous model. This effect occurs because the upwellings and downwellings in a highly viscous region experience strong variations in both amplitude and direction of flow velocity near the horizontal boundaries. In the models considered with the parameters used, the stresses in the upper and lower mantle are approximately identical, that is, there is no denser concentration of stresses in the upper or lower mantle. In contrast to the overlithostatic pressure field, the fields of horizontal stresses σ _xx in all models do not exhibit deep roots of highly viscous downwelling flows.     

The effect of stratification and compressibility on anelastic convection in a rotating plane layer

We study the effect of stratification and compressibility on the threshold of convection and the heat transfer by developed convection in the nonlinear regime in the presence of strong background rotation. We consider fluids both with constant thermal conductivity and constant thermal diffusivity. The fluid is confined between two horizontal planes with both boundaries being impermeable and stress-free. An asymptotic analysis is performed in the limits of weak compressibility of the medium and rapid rotation (τ^?1/12???|θ|???1, where τ is the Taylor number and θ is the dimensionless temperature jump across the fluid layer). We find that the properties of compressible convection differ significantly in the two cases considered. Analytically, the correction to the characteristic Rayleigh number resulting from small compressibility of the medium is positive in the case of constant thermal conductivity of the fluid and negative for constant thermal diffusivity. These results are compared with numerical solutions for arbitrary stratification. Furthermore, by generalizing the nonlinear theory of Julien and Knobloch [Fully nonlinear three-dimensional convection in a rapidly rotating layer. Phys. Fluids 1999, 11, 1469–1483] to include the effects of compressibility, we study the Nusselt number in both cases. In the weakly nonlinear regime we report an increase of efficiency of the heat transfer with the compressibility for fluids with constant thermal diffusivity, whereas if the conductivity is constant, the heat transfer by a compressible medium is more efficient than in the Boussinesq case only if the specific heat ratio γ is larger than two.     

Models of solar differential rotation

Abstract

Models of a differentially rotating compressible convection zone are calculated, considering the inertial forces in the poloidal components of the equations of motion. Two driving mechanisms have been considered: latitude dependent heat transport and anisotropic viscosity. In the former case a meridional circulation is induced initially which in turn generates differential rotation, whereas in the latter case differential rotation is directly driven by the anisotropic viscosity, and the meridional circulation is a secondary effect.

In the case of anisotropic viscosity the choice of boundary conditions has a big influence on the results: depending on whether or not the conditions of vanishing pressure perturbation are imposed at the bottom of the convection zone, one obtains differential rotation with a fast (≥ 10 ms^?1) or a slow (～ 1 ms^?1) circulation. In the latter case the rotation law is mainly a function of radius and the rotation rate increases inwards if the viscosity is larger in radial direction than in the horizontal directions.

The models with latitude dependent heat transport exhibit a strong dependence on the Prandtl number. For values of the Prandtl number less than 0.2 the pole-equator temperature difference and the surface velocity of the meridional circulation are compatible with observations. For sufficiently small values of the Prandtl number the convection zone becomes globally unstable like a layer of fluid for which the critical Rayleigh number is exceeded.     

How flat is the lower-mantle temperature gradient?

The temperature gradient in the lower mantle is fundamental in prescribing many transport properties, such as the viscosity, thermal conductivity and electrical conductivity. The adiabatic temperature gradient is commonly employed for estimating these transport properties in the lower mantle. We have carried out a series of high-resolution 3-D anelastic compressible convections in a spherical shell with the PREM seismic model as the background density and bulk modulus and the thermal expansivity decreasing with depth. Our purpose was to assess how close under realistic conditions the horizontally averaged thermal gradient would lie to the adiabatic gradient derived from the convection model. These models all have an endothermic phase change at 660 km depth with a Clapeyron slope of around −3 MPa K⁻¹, uniform internal heating and a viscosity increase of 30 across the phase transition. The global Rayleigh number for basal heating is around 2×10⁶, while an internal heating Rayleigh number as high as 10⁸ has been employed. The pattern of convection is generally partially layered with a jump of the geotherm across the phase change of at most 300 K. In all thermally equilibrated situations the geothermal gradients in the lower mantle are small, around 0.1 K km⁻¹, and are subadiabatic. Such a low gradient would produce a high peak in the lower-mantle viscosity, if the temperature is substituted into a recently proposed rheological law in the lower mantle. Although the endothermic phase transition may only cause partial layering in the present-day mantle, its presence can exert a profound influence on the state of adiabaticity over the entire mantle.     

Nonlinear free thermal convection in a spherical shell:A variable viscosity model

Introduction The velocity field of surface plate motion can be split into a poloidal and a toroidal parts.At the Earth′s surface,the toroidal component is manifested by the existence of transform faults,and the poloidal component by the presence of convergence and divergence,i.e.spreading and subduc-tion zones.They have coupled each other and completely depicted the characteristics of plate tec-tonic motions.The mechanism of poloidal field has been studied fairly clearly which is related to …     

Analytical solutions to the transient,unsaturated transport of water and contaminants through horizontal porous media

We present a range of analytical solutions to the combined transient water and solute transport for horizontal flow. We adopt the concept of a scale and time dependent dispersivity used for contaminant transport in aquifers and apply it to transient, unsaturated horizontal flow to develop similarity solutions for both constant solute concentration and solute flux boundary conditions. Through the use of a specific form of the water profile as used by Brutsaert [Water Resour Res 1968:4;785], the solute profiles can be reduced to a simple quadrature. We also derive a solution for the instantaneous injection of water and solute into a horizontal media for an arbitrary dispersivity. It is found that the solute concentration remains constant in both space and time as the water redistributes, suggesting that the solute does not disperse relative to the water.     

Physical and numerical experiments on layered convection in a density-Stratified fluid

Abstract

The formation and growth of horizontal layered convection cells in a density stratified solution of salt water subject to an impulsively applied lateral temperature gradient is investigated with physical and numerical experiments. Results indicate that lyers are induced by two mechanisms. One is the successive formation of layers due to the presence of the top and bottom boundaries. The other is the spontaneous occurrence of layers when a suitably defined Rayleigh number exceeds a critical value. It is found that well established layers are homogeneous in temperature and salinity and are separated by sharp gradients in density. Lateral heat transfer is of a periodic nature. Numerical experiments were carried out for finite and infinite geometry cases. For the finite geometry case, convection cells are generated successively inward from the horizontal boundaries. For the infinite geometry case, periodic conditions in the vertical direction are assumed. With continuous input of small perturbations, simultaneous occurrence of the convection cells is obtained at supercritical Rayleigh numbers. Criteria for determining the onset of spontaneous cells numerically are explored.     

Evolution of double diffusion convection in a felsic magma chamber

The onset of double diffusion convection (DDC) is modeled in a two-dimensional case in respect to magma chambers. The viscosity model for the melt takes into account the effects of temperature and concentration of the dissolved component (H₂O). The upper boundary of the convecting magma chamber is assumed to be anhydrous and at constant temperature, whereas the lower boundary is treated as being hydrous permeable with a temperature greater than that within the upper boundary. The case of positive compositional and thermal buoyancy of melt is studied assuming a H₂O diffusion coefficient small in comparison with thermal diffusivity. The DDC has been modeled using a system of equations solved by the finite difference method on a square grid. The convective pattern evolution has been studied for fixed boundary conditions as well as for cooling and degassing. Due to the higher viscosity in the upper zone, the upper boundary layer is thicker than the lower one. The variation of water concentration in this zone of the convective cell can be significant. In nature, the high gradient of water concentration can be responsible for the observed variations of water content in minerals crystallized from a granite melt (e.g., biotite). Because of a high Lewis number (= 100), temperature variations in the magma chamber decay much faster than the water concentration. In this case the intensive convection can continue at a constant temperature due to the non-zero water content in the chamber. In principle, the effect can be applied to the formation of magmatic bodies. If the cooling and degassing system reaches a uniform temperature distribution prior to the crystallization temperature, water content throughout the body may still remain variable.     

The effect of a shallow low viscosity zone on the apparent compensation of mid-plate swells

A simple model for mid-plate swells is that of convection in a fluid which has a low viscosity layer lying between a rigid bed and a constant viscosity region. Finite element calculations have been used to determine the effects of the viscosity contrast, the layer thickness and the Rayleigh number on the flow and on the perceived compensation mechanism for the resulting topographic swell. As the viscosity decreases in the low viscosity zone, the effective local Rayleigh number for the top boundary layer of the convecting cell increases. Also, because the lower viscosity facilitates greater velocities in the low viscosity zone, the low viscosity layer produces proportionally greater horizontal flow near the conducting lid, causing the base of the conducting lid to appear like a free boundary. The change in the local Rayleigh number and in the effective boundary condition both cause the top boundary layer to thin. Through a Green's function analysis, we have found that the low viscosity zone damps the response of the surface topography to the temperature anomalies at depth, whereas it causes the gravity and geoid response functions to change sign at depth counteracting the positive contributions from the shallower temperature variations. By increasing the viscosity contrast, the conbined effects of the thinning of the boundary layer and the behaviour of the response functions allow the apparent depth of compensation to become arbitrarily small. Therefore, shallow depths of compensation cannot be used to argue against dynamic support of mid-plate swells. Furthermore, we compared the distribution of the effective compensating densities, which is used to obtain the geoid, to that of Pratt compensation, which is often used to calculate the depth of compensation from geoid and topography data for mid-plate swells. For all of our calculations including those with no low viscosity layer, the effective gravitational mass distribution is more complex than assumed in simple Pratt models, so that the Pratt models are not an appropriate gauge of the compensation mechanism.     

Scale‐dependent Poiseuille flow alternatively explains enhanced dispersion in geothermal environments

Fluid flow exerts a critical impact on the convection of thermal energy in geological media, whereas heat transport in turn affects fluid properties, including fluid dynamic viscosity and density. The interplay of flow and heat transport also affects solute transport. To unravel these complex coupled flow, heat, and solute transport processes, here, we present a theory for the idealized scale‐dependent Poiseuille flow model considering a constant temperature gradient (?T) along a single fracture, where fluid dynamic viscosity connects with temperature via an exponential function. The idealized scale‐dependent model is validated based on the solutions from direct numerical simulations. We find that the hydraulic conductivity (K) of the Poiseuille flow either increases or decreases with scales depending on ?T > 0°C/m or ?T < 0°C/m, respectively. Indeed, the degree of changes in K depends on the magnitude of ?T and fracture length. The scale‐dependent model provides an alternative explanation for the well‐known scale‐dependent transport problem, for example, the dispersion coefficient increases with travel distance when ?T > 0°C/m according to the Taylor dispersion theory, because K (or equivalently flux through fractures) scales with fracture length. The proposed theory unravels intertwined interactions between flow and transport processes, which might shed light on understanding many practical geophysical problems, for example, geothermal energy exploration.     

Numerical simulations of mantle convection: Time-dependent,three-dimensional,compressible, spherical shell

Abstract

We describe nonlinear time-dependent numerical simulations of whole mantle convection for a Newtonian, infinite Prandtl number, anelastic fluid in a three-dimensional spherical shell for conditions that approximate the Earth's mantle. Each dependent variable is expanded in a series of 4,096 spherical harmonics to resolve its horizontal structure and in 61 Chebyshev polynomials to resolve its radial structure. A semiimplicit time-integration scheme is used with a spectral transform method. In grid space there are 61 unequally-spaced Chebyshev radial levels, 96 Legendre colatitudinal levels, and 192 Fourier longitudinal levels. For this preliminary study we consider four scenarios, all having the same radially-dependent reference state and no internal heating. They differ by their radially-dependent linear viscous and thermal diffusivities and by the specified temperatures on their isothermal, impermeable, stress-free boundaries. We have found that the structure of convection changes dramatically as the Rayleigh number increases from 10⁵ to 10⁶ to 10⁷. The differences also depend on how the Rayleigh number is increased. That is, increasing the superadiabatic temperature drop, δT, across the mantle produces a greater effect than decreasing the diffusivities. The simulation with a Rayleigh number of 10⁷ is approximately 10,000 times critical, close to estimates of that for the Earth's mantle. However, although the velocity structure for this highest Rayleigh number scenario may be adequately resolved, its thermodynamic structure requires greater horizontal resolution. The velocity and thermodynamic structures of the scenarios at Rayleigh numbers of 10⁵ and 10⁶ appear to be adequately resolved. The 10⁵ Rayleigh number solution has a small number of broad regions of warm upflow embedded in a network of narrow cold downflow regions; whereas, the higher Rayleigh number solutions (with large δT) have a large number of small hot upflow plumes embedded in a broad weak background of downflow. In addition, as would be expected, these higher Rayleigh number solutions have thinner thermal boundary layers and larger convective velocities, temperatures perturbations, and heat fluxes. These differences emphasize the importance of developing even more realistic models at realistic Rayleigh numbers if one wishes to investigate by numerical simulation the type of convection that occurs in the Earth's mantle.