This project looks at the application of general-purpose CFD to two plasma propulsion engines. Verification against experimental data for a 200 kW engine was undertaken, achieving similar thrust levels and jet exhaust divergence angles. A scaled-up 30 MW engine, capable of interplanetary travel, showed similar operating behaviour to the test engine. The increased seventy-fold thrust was achieved by increased mass efflux and operating conditions. The results show that general-purpose CFD can be successfully applied to plasma propulsion in a qualitative manner.
August 2020
1.Introduction
The concept of spaceflight has improved significantly since the first rocket was sent into space. However, propulsion systems used in rockets have changed little. While different fuels have been investigated and rockets made more efficient, the nature of modern-day chemical propulsion systems are not significantly different from their predecessors.
Plasma propulsion systems are viewed as more viable and efficient options for interplanetary space missions. A plasma engine that reached many headlines was the Variable Specific Impulse Magnetoplasma Rocket (VASIMR) due to its possible high efficiency, thrust and specific impulse resulting in a great reduction in the time of flight.
The aim of this project is to use general purpose Computational Fluid Dynamics (CFD) to build a test case according to the experimentally measured data from the 200 kW VASIMR engine. After optimization, this test case will be scaled for a 30 MW VASIMR engine. Applications of general purpose CFD to plasma engines have not been observed in published literature up until August 2020.
2. Theoretical Background
2.1. Variable Specific Impulse Magnetoplasma Rocket (VASIMR)
The main function that makes the VASIMR (a plasma propulsion system) different is that it has higher exhaust velocities and specific impulse compared to other chemical propulsion systems. AdAstra Company states that a 200 MW VASIMR engine is capable of reaching Mars in 39 days [6]. VASIMIR can vary its thrust and specific impulse, optimizing fuel usage and allowing greater payloads. It is capable of processing more power than other plasma propulsion systems, making it more suitable for interplanetary space travel. 
Figure 1 – VASIMR Schematic
A gas such as hydrogen, xenon, or argon is injected into a tube surrounded by magnets. The helicon coupler launches helical waves at the gas, resulting in an electron being knocked off the gas atom, converting it to plasma. Magnetic fields contain the plasma, shielding the walls of the engine. The Ion Cyclotron Heating (ICH) coupler uses ICH waves to exert a force on the ions, causing acceleration and higher temperatures up to 10 million Kelvin. A magnetic nozzle converts the orbital motion of these ions into linear momentum resulting in high exhaust velocities.
A 200 kW VASIMR engine was tested experimentally. The thrust efficiency reached 54% and specific impulse reached 3500 seconds.[7] Models operating at higher power inputs have not been experimentally tested. However, AdAstra states that a 30 MW VASIMR engine at 60% efficiency human mission launched from low earth orbit, would reach Mars in 66 days.[8]
2.2 Computational Fluid Dynamics (CFD)
The ever-growing pace in computer power has resulted in the development of more detailed prediction models for flowfield behaviour such as CFD. CFD methods offer a wide range of possible applications due to the different physics models that can be applied in complex geometries to a project. For example, there is flexibility of application to steady/transient, laminar/turbulent and reacting/non-reacting fluid behavior.
The above approaches are based on the Navier-Stokes partial differential equations representing conservation of mass, momentum, energy and chemical species for a reacting flow. Additionally, heat and mass diffusion are represented. The general form which is numerically solved in CFD programs is as follows:
The dependent variable () can represent several quantities such as mass fraction of a chemical species for conservation of a species; enthalpy for conservation of energy; velocity component for conservation of momentum in the component direction. Appropriate diffusion constants( ) and source terms () need to be defined for the conservation equations. The rate of change of the dependent variable is denoted by and the velocity field by .
The equation above, provides the theoretical basis to general-purpose CFD programs. Such programs use a general set of instructions for solving the above system of equations in their numerical forms alongside appropriate boundary and initial conditions. The numerical forms of these conservation equations are solved in a finite volume CFD approach where the geometry of the model is represented using finite grids and meshes. This project uses a general-purpose CFD software which implements such a methodology.
3. Modelling Methodology
A 200 kW VASIMR engine was modelled and verified according to available experimental data for this engine. A scaled-up model was then built and tested for the 30 MW VASIMR engine with geometry being scaled up seven times as the power input was scaled up by around 150 times. According to the equation for thrust:
Where:
= change in momentum over time
= density
= area of nozzle (circle)
= velocity squared
The velocity and density parameters need to be comparable to the 200 kW engine so that the mass flow, which was the main factor for the scaling of thrust, is proportional to area. The radius was scaled up seven times as square rooting the power input scaling gave a value of around 12. The additional scaling (x1.7) was achieved by increasing operating conditions in the 30 MW model.
This project uses Ansys general-purpose CFD software, version 2019-R3, including: SpaceClaim for CAD modelling; Fluent Mesh for generating the computational grid; Fluent as the CFD solver and post-processor.
3.1 VASIMR geometry
The geometry of the 200 kW VASIMR engine was constructed using SpaceClaim by approximations using key diameters and lengths stated in experimental data and the diagrams showing the engine.
Figure 2- A diagram representing showing the different components in the VASIMR engine
Figure 3 – side on view of geometry with dimensions labelled. Black font represents 200 kW geometry and blue font represents the geometry of the 30 MW engine
Figure 4 – Slanted view of the geometry showing the different components
As presented in Figures 3 and 4, given that fluid behaviour was of interest, the geometry represents only the inner components of the VASIMR engine. Additional components were introduced using sub-models in the solver such as a fan representing the magnets that provide rotational motion containing the plasma at the core of the engine, away from wall surfaces. Figure 3 shows that the 30 MW engine geometry was simply scaled up from the 200 KW engine geometry by seven times.
3.2 Computational Mesh
Both geometries were imported into Fluent Mesh, where a volume mesh was generated.
Figure 5 – An overall view of the mesh 
|
30MW engine |
200kW engine |
|
|
Local mesh sizing on walls of fan core (cm) |
2.5 |
0.25 |
|
Surface mesh minimum (cm) |
10 |
1 |
|
Surface mesh maximum (cm) |
30 |
2.5 |
|
Volume mesh layers |
3 |
3 |
|
Number of cells (2sf) |
160000 |
450000 |
|
Mesh type |
Polyhedral mesh & prism cell layers |
Polyhedral mesh & prism cell layers |
Figure 8 – A table of the values used to generate the meshes for both engines
Figures 6 and 7show the fan was refined to achieve a good mesh resolution by implementing local sizing. Figure eight illustrates that the number of cells for the 200 KW engine was lower than for the 30 MW engine. However, the same number of cells was tested for the 30 MW engine as for the 200 kW engine and the results were similar to using a lower number of cells. Thus, a lower number of cells were used allowing a greater efficiency in the optimizations computations required.
3.3. Solver
The plasma was represented as a fluid, so several models were used to ensure the gas had similar behavior to plasma in the VASIMR engine. A steady state solution was used. An inlet boundary represented plasma generation. A fan model illustrated plasma confinement and introduced rotational velocity. Additionally, a user function for heating was applied in the fan zone representing RF heating. At the nozzle, the gas expanded to represent the exhaust.
3.3.1 Energy Model
This was activated to allow for the calculation of the energy conservation equation in the flow field solution.
3.3.2 Standard k-ε Viscous Model
A two-equation turbulence model represented turbulence generation and dissipation, depending on transport equations for turbulent kinetic energy (k) and dissipation rate (ε). For k, the transport equation is obtained from it’s exact equation, however for ε, the transport equation can be found through physical reasoning. Turbulent viscosity (), an input to the turbulent diffusive terms in the conservation equations, is computed using:
Where = model constant = 0.09.
3.3.3 Heater 3D fan zone
A 3D fan zone was used to introduce rotational velocity. 3D fan zones are fluid cell zones simulating the effect of an axial fan by applying a distributed momentum source in a blade-swept volume. Radial momentum was not activated as the flow needed to be confined rather than swept out. The following are expressions for the momentum sources in the axial and tangential directions:
Where:
= pressure rise across fan for given axial flow rate Q
= thickness of fan swept volume in axial direction
= fan power
r = local radial distance
= fan operating angular velocity
= radius fan hub
= radius of a point on fan blade based on inflection point ratio
= radius fan blade tip
= constant based on inputted parameters
The above parameters were entered into the solver. A relatively small nominal value of 1000 Pa was chosen for the pressure jump to ensure that the increase in pressure was not heavily impacting the flow expansion. The inflection point ensured the core of the rotation was confined close to the engine axis representing plasma confinement away from engine wall surfaces.
3.3.4 Heater Energy Source Term
A user function for heating was introduced in the fan zone representing ICH waves which then brought forward the concept of heating and higher temperatures. This directly impacted the enthalpy conservation in the CFD solver allowing for gas expansion in the nozzle resulting in an increase in the flow velocities.
3.3.5 Material Properties
Air was chosen to represent the fluid medium. An incompressible form of ideal gas was used for representation of density calculations as a function of primarily temperature.
3.3.6 Boundary Conditions
Gas was introduced at the inlet with a velocity inlet boundary representing the formation of plasma. Wall Boundary Conditions were activated via the Wall Function method in the Fluent, where the solution at the near wall cell layer illustrated the flow turbulent boundary layer. At the domain exit, pressure outflow and nozzle walls were chosen to be far away, so flow expansion would not be affected.
4. Results
4.1 200 kW VASIMR engine
This section presents simulation results for the 200 kW engine. Figure 8 shows the scaled residuals (representing the imbalance in conservation for all equations) which had dropped several orders of magnitude, showing good convergence is achieved for the Steady State calculations.
Figure 8 – Scaled Residuals Plot
Figure 9 illustrates gas entering the inlet at relatively low velocities. In the ICH coupler, the gas is then rotated and heated in the fan core, causing flow velocity magnitude to increase via gas expansion. The gas is subsequently ejected from the nozzle. The rotational component of the velocity field introduced by the fan is illustrated in the Path-line plots of Figure 10.
Figure 9 – Velocity Magnitude Computed (Maximum Range) Contour
Figure 10 – Pathlines of Velocity Magnitude
Figure 11 shows the temperature rise around the core of the fan where heating was applied correlating to the high velocities shown in Figure 9 and the density drop resulting in gas expansion due to heating is illustrated in Figure 12. 
Figure 11 – Temperature Contour
Figure 12 – Density Contour
Plots of figures 13-16 are area-averages of parameters taken at a section of the nozzle therefore, through the core heater section together with the surrounding flow. The jet momentum represented the thrust settling at around 10N (Figure 13).. Average z-velocity reached around 12 m/s (Figure 14), whilst tangential velocity was at 3.5 m/s (Figure 15).
Figure 13 – Integral of Jet Momentum Plot
Figure 14 – Area-weighted Average of z-velocity Plot
Figure 15- Area-Weighted Average of tangential-velocity Plot
The average temperature had been found at approximately 400 K, but on the other hand, the temperature had reached up to a maximum of 530.7 K in the fan core and this is highlighted in Figure 10.
Figure 16 – Area-Weighted Average of temperature Plot
4.2 30 MW VASIMR engine
The following section presents the results obtained for the 30 MW engine simulation. As shown in figure 17, the scaled residuals dropped several orders of magnitude, suggesting good convergence for all steady state calculations.
Figure 17 – Scaled Residuals Plot
At the inlet, the gas entered with the same low velocity magnitude as the 200 kW engine. In the ICH coupler, velocity increased significantly, and flow expanded from the nozzle. Rotational velocity was introduced via the fan model(Figure 19).
Figure 18 – Velocity Magnitude Contour
Figure 19 – Pathlines of Velocity Magnitude
In the fan’s core, the temperature had risen (Figure 20) and was accompanied by a drastic decrease in density (Figure 21).
Figure 20 – Temperature contour
Figure 21 – Density Contour
Figure 22 illustrates jet momentum reaching around 700 N while z-velocity settling at approximately 21.5 m/s shown by Figure 23. The tangential velocity, seen in Figure 24, reached around 6 m/s however, it remained there for a lower number of iterations.
Figure 22 – Integral of Jet momentum Plot
Figure 23 – Area-weighted average of z-velocity Plot
Figure 24 – Area-Weighted Average of tangential-velocity Plot
Temperature settled at around 800 K, however maximum temperatures reached up to 2000 K in the core of the engine, illustrated by Figure 20.
Figure 25 – Area-weighted Average of temperature plot
5. Discussion
5. 1 Comparison of 200kW VASIMR results against experimental results
A thrust of 10 N was produced experimentally using a 200 kW VASIMR engine[10]. A similar jet momentum was seen in the simulation of a 200 kW engine, showing successful calibration. Additionally a half angle of 24 -30° was observed in the exhaust experimental results[7]. To measure a half angle for the simulation results, a velocity contour limited to 20 m/s was used, shown by Figure 26. An angle of approximately 22° was measured. Considering representation of plasma cannot be exact, a percentage error of less than 10%, is considered acceptable.
Figure 26 – Velocity magnitude contour limited to 20m/s with half angle indicated
5.2 Comparison of 30MW to 200kW
Comparing figures 9 and 18, it’s apparent that the velocity fields are similar between the two engines. Additionally, the ratio of tangential to axial velocities of the jet efflux between the two engines are similar. These observations show that the operation principle is similar between the 30 MW and test engine.
On the contrary, it is observed from the density contours (figures 12 and 21) that the expansion taking place in the core of the 30 MW engine is greater than the test. The greater expansion is achieved by higher operating temperatures of the 30 MW engine (figures 11 and 20). This highlights that the increase in thrust from 10 N to 700 N is achieved by a combination of increased mass flow rate of fluid and engine operating conditions (i.e. increased core temperature).
5.3 Developments
The results showed that the application of general purpose CFD was successful to qualitatively assess the behavior of plasma and to model the different sections in the engine. However, many plasma properties were unable to be modelled.
The representation of gas being converted to plasma in the helicon coupler was not possible as the plasma had to be represented as a fluid. The heating and rotational velocity were successfully introduced in the ICH coupler. Additionally, the plasma couldn’t be confined only to the core of the engine but, in the ICH, heating and rotational velocity were confined to try and represent the plasma behavior as closely as possible.
Density of the flow was represented via the ideal gas equation. However, plasma density doesn’t follow the ideal gas equation but rather the Boltzmann distribution. Plasma density can be modelled using the following:
Where the number of electrons () is directly proportional to the exponential () of potential energy () over temperature (). Inputting this equation as a user defined function, it was attempted to model density using the above equation. However, the solver was unable to handle this representation of density in a stable manner.
To represent plasma properties as a fluid using plasma physics, a Magnetohydrodynamics (MHD) model could be used in Fluent. However, this was beyond the scope of the current project as it required an extra license and it would not be possible to find all the input parameters required from the published papers on VASIMR.
The geometry for the model was not completely accurate as some values had to be scaled up from available figures resulting in errors occurring in the model. Mesh sensitivity was investigated, and the results produced were similar. Thus, a slightly lower mesh size was used to increase efficiency whilst tuning the models for optimum performance.
6.0 Conclusions
General-purpose CFD has been applied qualitatively to a plasma propulsion system. Verification of the 200 kW model against the experimental test results showed encouraging similarities in jet efflux angle and effective engine thrust. The test case was successfully scaled up to a 30 MW engine with 700 N thrust, representative of operations to reach Mars with relatively short journey times.
The scaled-up model showed similar operation behaviour to the test engine. The seventy-fold increase in thrust was achieved by a combination of increased mass flow and engine operating conditions.
Adapting general purpose CFD for the property of density for plasma is challenging as the solver was unable to perform in a stable manner. Specialist models, such as MHD, are required for accurate modelling of plasma.
7. Appendix: User routine for plasma density
8. Bibliography
Research
[1] Choueiri, Edgar Y. ‘New Dawn for Electric Rockets’, Scientific American, February 2009,
https://alfven.princeton.edu/publications/choueiri-sciam-2009. [Accessed 21 Aug.2020]
[2] Crew, Bec. ‘This Plasma Engine Could Get Humans to Mars on 100 Million Times Less Fuel’. Science Alert. Last modified 28 October, 2015.
https://www.sciencealert.com/this-plasma-engine-could-get-humans-to-mars-on-100-million-times-less-fuel. [Accessed 21 Aug.2020]
[3] Cassady, Leonard D. Longmier, Benjamin W. Olsen, Chris S. Ballenger, Maxwell G. McCaskill, Greg E. Ilin, Andrew V. Carter, Mark D. et al. \’VASIMR Performance Results’. 2010. American Institute of Aeronautics and Astronautics: 2010-6772
http://www.adastrarocket.com/AIAA-2010-6772-196_small.pdf. [Accessed 21 Aug.2020]
[4]AdAstra Rocket. ‘NEP-VASIMR Human Mission to Mars’. Youtube, Aug 14, 2013.
https://www.youtube.com/watch?time_continue=2&v=IdcOWVm_Ov4&feature=emb_title. [Accessed 21 Aug.2020]
[5] Patankar, Suhas V. 1980. Numerical Heat Transfer and Fluid Flow. 1st ed. Washington; Hemisphere Publishing Cooperation.
[6] Marquardt, Timothy A. ‘An Overview of the VASIMR Thruster Concept’. 2014. Auburn University: AERO 6530
http://www.objective-europa.com/wp-content/uploads/2013/08/VASIMR_Overview_OE_5.pdf. [Accessed 21 Aug.2020]
Additional Resources
1.Ansys. ‘Fluent Magnetohydrodynamics (MHD) Module Manual’. 2013
https://www.afs.enea.it/project/neptunius/docs/fluent/html/mhd/main_pre.htm. [Accessed 30 Aug. 2020]
2.Ansys. ‘Fluent Theory Guide’. 2019.
https://www.afs.enea.it/project/neptunius/docs/fluent/html/th/main_pre.htm. [Accessed 30 Aug. 2020]
3.Hutchinson, Ian. ‘Fluid Description of Plasma’. MIT OpenCourseWare. 2003.
https://ocw.mit.edu/courses/nuclear-engineering/22-611j-introduction-to-plasma-physics-i-fall-2006/readings/chap4.pdf. [Accessed 30 Aug.2020]
4.Peraire, J. Widnall, S. ‘Lecture L14 – Variable Mass Systems: The Rocket Equation’ MIT OpenCourseWare. 2008.
https://ocw.mit.edu/courses/aeronautics-and-astronautics/16-07-dynamics-fall-2009/lecture-notes/MIT16_07F09_Lec14.pdf. [Accessed 30 Aug. 2020]
About The Author
Mona is a 17 year old student studying Physics, Maths, Further Maths and Chemistry for her A levels. She is especially interested in Physics and Engineering and so is hoping to read Physics at university. In her free time, she loves to play the piano and violin, row and read about new scientific discoveries