CFD is used to investigate the performance of a self-evaporating MHD system by utilizing a circumferential ejector. Moreover, liquid sodium is used as the working fluid. The system's performance in terms of velocity augmentation and MHD power density are calculated and presented in terms of vapor/liquid ratio, heat input power, and nozzle ejector configuration. Figure 3 shows the schematics of the proposed system.

Figure 3.

The schematic diagram of ejector

These cases are selected to cover a variety of possible circumferential ejector configurations and to identify the optimal or near-optimal one. Five different cases and four geometrical parameters are taken for the off-design operation. The parameters A, B, C, and D are the radius of the secondary inlet, the difference between the outer radius and inner radius, an axial distance over 200 mm that is added to the primary region, and the distance between the two throats (primary throat and main throat), respectively. In case 1, no distance is added to the primary region (C =0), which means no change in the cross section of the primary region happened. While in case 2, the parameter C is 31 mm, where the smallest diameter of the primary throat is at this value. However, the dimensions of the considered ejector configurations are listed in Table 1.

Assuming a 2D steady flow mixture of the liquid metal and its vapor, the governing equations can be presented as [25]:

The continuity equation for the mixture is:

$\nabla .({\rho }_{\text{m}}{\stackrel{\rightharpoonup }{\nu }}_{\text{m}})=0$ (1)

where ${\stackrel{\rightharpoonup }{\nu }}_{\text{m}}=\frac{{\Sigma }_{k=1}^n{\alpha }_k{\rho }_k{\stackrel{\rightharpoonup }{\nu }}_k}{{\rho }_m}$ is the mass-averaged velocity and ${\rho }_{\text{m}}={\Sigma }_{k=1}^n{\alpha }_k{\rho }_k$ is the mixture density.

α_k is the volume fraction of phase k where n = 2 for this case. Phase 1 is the liquid metal, and phase 2 is its vapor:

${\displaystyle \sum _{k=1}^n{\alpha }_k=1}$ (2)

The momentum equation for the mixture is:

$\nabla .({\rho }_{\text{m}}{\stackrel{\rightharpoonup }{\nu }}_{\text{m}}{\stackrel{\rightharpoonup }{\nu }}_{\text{m}})=-\nabla p+[{\mu }_{\text{m}}(\nabla {\stackrel{\rightharpoonup }{\nu }}_{\text{m}}+\nabla {\stackrel{\rightharpoonup }{\nu }}_{\text{m}}{}^T)]+{\rho }_{\text{m}}\stackrel{\rightharpoonup }{g}+\stackrel{\rightharpoonup }{j}+\stackrel{\rightharpoonup }{B}+({\displaystyle \sum _{k=1}^n{\alpha }_k{\rho }_k{\stackrel{\rightharpoonup }{\nu }}_{\text{dr},k}{\stackrel{\rightharpoonup }{\nu }}_{dr,k})}$ (3)

The term on the left-hand side of eq. (3) is the rate of momentum transfer, and the terms on the right-hand side are the pressure force, the viscous stresses, the gravity body force, the Lorentz force, and the drift force between interphases, respectively. j is the current density, B is the magnetic field intensity, μ_m is the viscosity of the mixture, which is equal to ${\mu }_{\text{m}}={\Sigma }_{k=1}^n{\alpha }_k{\mu }_k$and ${\stackrel{\rightharpoonup }{\nu }}_{\text{dr},k}$is the drift velocity for secondary phase k, which is equal to ${\stackrel{\rightharpoonup }{\nu }}_{\text{dr},k}={\stackrel{\rightharpoonup }{\nu }}_k-{\stackrel{\rightharpoonup }{\nu }}_{\text{m}}$.

The energy equation for the mixture takes the following form:

$\nabla .{\displaystyle \sum _{k=1}^n({\alpha }_k}{\stackrel{\rightharpoonup }{\nu }}_k({\rho }_k{E}_k+p))=\nabla .(k_{\text{eff}}\nabla T)+\nabla .({\stackrel{\rightharpoonup }{\nu }}_{\text{m}}.\stackrel{\rightharpoonup }{\tau })+\stackrel{\rightharpoonup }{j}.\stackrel{\rightharpoonup }{E}$ (4)

The term on the left-hand side of equation (4) is the enthalpy carried by the mixture, and the terms on the right-hand side are the axial conduction, the viscous dissipation, and the energy caused by the Lorentz force, respectively. In this equation, k_eff is the effective conductivity ${\Sigma }_{k=1}^n{\alpha }_k(K_k+K_t)$, K_t is the turbulent thermal conductivity, and E is the electric field.

In addition, the k-ε turbulence model is used in this work. The equations of the turbulence model can be referred to in references [24], [25].

The continuity equation for current is:

$\nabla .\stackrel{\rightharpoonup }{J}=0$ (3)

Ohm's law is:

$\stackrel{\rightharpoonup }{J}=\sigma (\stackrel{\rightharpoonup }{E}+\stackrel{\rightharpoonup }{\nu }\times \stackrel{\rightharpoonup }{B})$ (4)

where σ is the electric conductivity.

The electric field is:

$\stackrel{\rightharpoonup }{E}=-\nabla \phi$ (5)

where φ is the electrical potential which can be calculated by using equation (6):

$\nabla .(\sigma \nabla \phi )=\nabla .(\sigma \stackrel{\rightharpoonup }{\nu }\times \stackrel{\rightharpoonup }{B_{\text{o}}})$ (6)

The geometrical modeling and meshing of the five test cases (Table 1) are conducted using ANSYS-FLUENT. To ensure numerical stability, hexahedral mesh grids were structured using the MultiZone Qua/Tri method with face meshing. Mesh generation nearing the primary and secondary nozzle throat is shown in Figure 4.

Figure 4.

Mesh structure near the nozzle throat

Grid sensitivity analysis is conducted by using the incrementally refined mesh sizes of 32426, 48357, 67626, 104594, and 184043. Figure 5 shows the velocity distribution of liquid sodium along the central x-axis. The difference between the course and the finer mesh size is lower than 1.1%. Therefore, mesh sizes 48357 and 67626 elements were applied throughout the present study.

Figure 5.

Velocity distribution of liquid sodium along the central x-axis at different element size

The governing equations are solved numerically by using the finite volume method by ANSYS-FLUENT. The Mixture and MHD models [25], [28] are utilized for simulating the multi-phase flow under a magnetic field. The Evaporation-Condensation model is embedded to solve the interaction between the phases [25]. SIMPLE scheme is used for the pressure-velocity coupling. The least-square cell-based method is used for computing the cell-centered spatial gradient terms and the second-order discretization scheme for the convective terms in the pressure and density equations. The volume fraction is discretized using QUICK scheme. Table 2 shows the physical properties of liquid Sodium and its vapor at T_ref. =1155 K [8]. Density, thermal conductivity, electrical conductivity, latent heat of vaporization, magnetic permeability, and viscosity are all assumed constant. The k-ε two-equation turbulent model with standard wall functions is employed. Symmetry boundary conditions are applied along the centerline of the ejector.

To validate the new model, Lu et al. [8] geometry is modeled with the same dimensions that they used. Grid sensitivity analysis is conducted by using the incrementally refined five mesh sizes, as shown in Figure 5. The governing equations are solved numerically using the methods presented in the previous section. The outlet boundary condition was pressure-outlet with atmospheric conditions, whereas the inlet boundary condition was velocity-inlet. Despite the strong coupling between the equations introduced in mathematical model section, convergence is achieved with precision 10^-6 of the residual monitor. The solution procedure is described in the following:

Lu et al. [8] case (without the circumferential annular ejector) is resolved using the above methodology. As shown in Figure 6, Figure 7 and Figure 8, our findings coincide very closely with Lu et al. [8]. The discrepancy ranged 0.027-0.75% for the outlet velocity, 0.29 -3.96% for the liquid outlet volume fraction, and 0.52 - 1.13% for the outlet temperature under the same conditions. We used Lu et al. [8] test cases as reference cases. They simulated a conventional Self-Evaporating LMMHD system, which did not contain the ejector device. Contrary to their system, we added a circumferential ejector to investigate the effect of adding this device on the system's performance where the one entrance they used divide to two entrances, primary and secondary.

Figure 6 shows the effect of wall temperature on the outlet velocity. Regardless of whether an ejector is used or not, it is evident that increasing the wall temperature increases the outlet velocity. At 1337 K, the outlet velocity is 7.44 m/s for case 1 and 9.09 m/s for case 2. The higher the wall temperature, the higher the outlet velocity; so, using the ejector increases the outlet velocity in most cases except case 1, where no change in the cross section of the primary region happened. Case 2 shows the most effective case, where the outlet velocity is increased by 18.36% due to the ejector compared with the case without the ejector at 1773K.

Figure 6.

The effect of wall temperature on the outlet velocity with and without ejector

Figure 7 shows that increasing the wall temperature increases the outlet temperature in all cases. At 1337 K, the outlet temperature is 1263.6 K for case 5 and 1210.2 K for case 1. The outlet temperature gradually increases as the wall temperature increases to 1773K. Using the annular ejector with different parameters affects the outlet temperature; contrary to cases 1,2 and 3, the outlet temperature in case 5 increases by 2.46% while it decreases by 6.51% in case 1.

Figure 7.

The effect of wall temperature on the outlet temperature with and without ejector

The wall temperature increases leading to sodium LM evaporation. Figure 8 shows that increasing the wall temperature reduces the outlet volume fraction of liquid in all cases. The outlet volume fraction of liquid at 1337 K wall temperature is 0.376 for case 1 and 0.29 for case 5. The outlet volume fraction of liquid declines steadily as the wall temperature. Figure 8 also shows that the annular ejector affects the outlet volume fraction of liquid which increases by 21.24% in case 1 and decreases by 16% in case 5 compared with the reference case at 1773K wall temperature when the ejector was not used. Temperature increases, reaching 0.25 for case 1 and 0.17 for case 5 at 1773K wall.

Figure 8.

The effect of wall temperature on the outlet volume fraction of liquid with and without ejector

As the previous section demonstrates, increasing the wall temperature can result in a higher outlet velocity, which is beneficial to output power but also results in a lower volume fraction of liquid. Reduced liquid volume fraction (more vapor phase) accelerates the mixture but reduces its electrical conductivity.

The output power was computed for all cases in this section using a constant magnetic field of 0.1 T. All five proposed geometries for the annular ejector and geometry for the central ejector is compared with the reference case (without ejector). The findings show that the output power has increased in all cases when using the ejector. Figure 9 and Table 3. show that the increase of output power due to increased velocity is higher than the loss in power due to the decrease in the volume fraction of the liquid.

As shown in Figure 9, case 2 is the best geometry that gives the maximum output power augmentation. For case 2, output power increases by 14.0% over the reference case (without ejector) and 8.7% over the central ejector case.

Figure 9.

The effect of the annular ejector on the output power

In the next section, the analysis will focus on the performance of the best case, case number 2, where the velocity, temperature, pressure, and volume fraction distributions are presented.

Figure 10, Figure 11, and Figure 12 present the contours of velocity in an annular ejector, the velocity distribution along the x-axis, and the velocity distribution along the y-axis. As shown in these figures, the sodium LM is admitted at a velocity of 3 m/s for primary and secondary nozzles. In the primary nozzle, part of the sodium LM starts converting to vapors which causes acceleration in the primary flow to reach around 90 m/s at the throat of the ejector. Then the velocity decreases due to the momentum exchange between the primary and the secondary flows. The velocity of the primary flow continues to drop to reach the same value as the secondary flow through the diffuser. Finally, bubbly flow leaves the MHD channel with a velocity of around 15 m/s.

Figure 10.

Distribution of velocity [m/s] along X at different Y position

Figure 11.

The contour of velocity [m/s] for case 2

Figure 12.

Distribution of velocity [m/s] along Y at different X position: at X= 0.1, 0.2 0.264 and 0.3(a); at X= 0.4, 0.5, 0.6 and 0.7 (b)

Figure 13, Figure 14 and Figure 15 present the contours of temperature in the annular ejector, the temperature distribution along the x-axis, and the temperature distribution along the y-axis. As shown in these figures, the sodium LM enters with a temperature of 1155 K, less than the sodium's boiling point by 1 K. The heat transfer occurs from the heated wall to the primary at a flux of 7.4×10⁵ W, which causes evaporation of the liquid metal in the annular ejector. At the throat, two stream mix together. The temperature of the primary flow continues dropping to reach the same value as the secondary flow through the diffuser. Finally, the mixture leaves the MHD channel at a temperature of around 1400 K.

Figure 13.

Distribution of temperature [K] along X at different Y position

Figure 14.

The contour of temperature [K] for case 2

Figure 15.

Distribution of temperature [K] along Y at different X position: at X= 0.1, 0.2 0.264 and 0.3 (a); at X= 0.4, 0.5, 0.6 and 0.7 (b)

Figure 16, Figure 17 and Figure 18 present the volume fraction contours of liquid in an annular ejector, the volume fraction of liquid distribution along the X-axis, and the volume fraction of liquid distribution along the Y-axis. These figures show that; the volume fraction of liquid is dominant at the beginning; then part of the primary stream converts to sodium vapor due to heating, while the secondary one remains liquid until two streams mix after the throat. Bubbly flow occurs in the diffuser and leaves the MHD channel with a volume fraction of liquid around 22%.

Figure 16.

Distribution of liquid volume fraction along X at different Y position

Figure 17.

The contour of liquid volume fraction for case 2

Figure 18.

Distribution of liquid volume fraction along Y at different X position: at X= 0.1, 0.2 0.264 and 0.3 (a); at X= 0.4, 0.5, 0.6 and 0.7 (b)

Figure 19 and Figure 20 display the contours of pressure in the annular ejector and the pressure distribution on the X-axis. As demonstrated, the pressure of the primary stream as it passes through the primary region drops dramatically and then falls to a lower value near the main throat, where vacuum and momentum transfer between the two streams occur in the mixing region after the throat. Driven by the pressure gradient, the secondary stream accelerates by the primary stream. Finally, the mixture leaves the MHD channel with atmospheric pressure imposed as a boundary condition.

Figure 19.

Distribution of pressure [Pa] along X at different Y position

Figure 20.

The contour of pressure [Pa] for case 2

The entrainment ratio is the ratio of the primary stream's mass flow rate to the secondary stream's mass flow rate. It measures the amount of the secondary fluid driven by the primary one. The effect of the entrainment ratio on power density was investigated using a magnetic field of 0.1 T, a wall temperature of 1773K, a density of sodium liquid metal of 927 kg/m³, an area of primary section is 7.5×10⁻⁵ m^2, and the area of the secondary section is 6.7×10⁻⁴ m².

Figure 21 and Table 4 show that raising the entrainment ratio reduces the power density. This power reduction is characterized by a decline in the mass flow rate of the secondary stream, i.e., the flow through the generator becomes more diluted, thus less power generation. At the entrainment ratio of 0.04172, the power density is 95.613 kW/m³ for the annular ejector and 68.558 KW/m³ for the central ejector. It decreases gradually as the entrainment ratio increases; at a 0.1113 ratio, the generated power becomes 47.853 kW/m³ for the annular ejector and 61.629 KW/m³ for the central ejector.

Figure 21.

The effect of entrainment ratio on power density

With the increased input heating power and, therefore, the increased wall temperature, more vapor phase is generated to drive the mixture. Speeding up the mixture counterweighted the adverse effect of reducing the volume fraction of the electrically conductive liquid and resulted in a higher net power generation. Figure 22 shows the optimum value of the heat added that delivers the maximum performance in case 2.

Figure 22.

The effect of wall temperature on the outlet parameters, case 2