Moving particle semi-implicit simulation on the molten Wood’s metal downward relocation process

Moving particle semi-implicit simulation on the molten Wood’s metal downward relocation process 1 first-author https://orcid.org/0000-0002-8278-0375 2 corresp Abdul Waris, awaris@fi.itb.ac.id Abdul Waris, awaris@fi.itb.ac.id awaris@fi.itb.ac.id 2 2 2 Physics Department, Faculty of Mathematics and Natural Science, Institut Teknologi Bandung, Jalan Ganesha 10, Bandung 40132, Indonesia Physics Department, Faculty of Mathematics and Natural Science, Institut Teknologi Bandung, Jalan Ganesha 10, Bandung 40132, Indonesia Nuclear Physics and Biophysics Research Division, Physics Department, Faculty of Mathematics and Natural Science, Institut Teknologi Bandung, Jalan Ganesha 10, Bandung 40132, Indonesia Nuclear Physics and Biophysics Research Division, Physics Department, Faculty of Mathematics and Natural Science, Institut Teknologi Bandung, Jalan Ganesha 10, Bandung 40132, Indonesia In the case of a severe accident involving nuclear reactors, an important aspect that should be considered is the leakage of molten material from the inside of the reactor into the environment. These molten materials damage other reactor components, such as electrical tubes, grid plates, and core catchers. In this study, the moving particle semi-implicit (MPS) method is adopted and improved to analyze the two-dimensional downward relocation process of molten Wood’s metal as a representation of molten material in a nuclear reactor. The molten material impinges the Wood’s metal plate (WMP), which is mounted on a rigid dummy stainless steel in a cylindrical test vessel. The breaching process occurs because of heat transfer between the molten material and WMP. The formed breach areas were in good agreement with the experimental results, and they showed that the molten Wood’s metal spread above the WMP. The solid WMP fraction decreased with time until it reached the termination time of the simulation. The present results show that the MPS method can be applied to simulate and analyze the downward relocation process of molten material in the grid plate of a nuclear reactor. Heat transfer Moving particle semi-implicit Phase change Relocation Severe accident 1. Introduction 1 In recent years, melt leakage behaviors have been increasingly recognized as an important, worldwide nuclear technology concern globally, especially when a severe accident occurs. The molten materials that leak during an accident may relocate to other locations inside the nuclear reactor and damage other components[ 1 1 ]. The damage can be avoided by quenching the molten materials using coolants, such as sodium, to decrease their temperature and fragment them. Consequently, some fragments or freeze materials plug or clog and obstruct the other components [ 2 2 ] . Furutani et al. conducted experiments to determine the erosion rate of a solid plate impinged by a liquid jet using paraffin and Wood’s metal material [ 3 3 ] . These materials were selected because of their low melting points. The conjugated nature of melt breakup and steam generation has been experimentally investigated using the non-boiling conditions of Wood’s metal [ 4 4 ] , which explains why the Kelvin–Helmholtz instability caused jet breakup. Similarly, Wood’s metal has been used to analyze the molten behavior when it is poured into a stainless-steel tube and injected into the coolant through a nozzle [ 5 5 , 6 6 ] . It was found that the frozen layer of metal adhered to the structure, and its shape was thin and wide in the air coolant experiment. However, in the water coolant experiment, the frozen metal was broken and fell down as debris. Additional research into molten behavior has been conducted by the Korea Atomic Energy Research Institute by losing some feed water and coolant without safety injection. These scenarios were referred to as level 1 probabilistic safety assessments [ 7 7 ] . Recently, molten Wood’s metal was used by Sudha et al. (2018) to investigate the downward relocation to a plate with the same material as the molten Wood’s metal [ 8 8 ] . The MPS method has been used to simulate the free surface and multiphase flow of incompressible media [ 9 9 ] , and it was developed by Koshizuka and Oka in 1996. This is a numerical method that does not use meshes, but uses particles as a part of the simulation. The Navier–Stokes equation is mainly used for dynamic particle calculations. The equation is transformed into particle interactions within an effective radius, and it has been used to simulate the penetration rates of molten uranium into steel test sections [ 10 10 ] . It was found that there was a small discrepancy between the result and the experiment because of the assumption of the diffusion coefficient, which was applied in the calculation. Li performed a two-dimensional (2D) MPS calculation to predict the molten behavior of the Fukushima Daiichi Unit 1 [ 11 11 ] . It focused on the lower head tube rupture, resulting in debris being released into the other part of the reactor. Similarly, Chen used it to investigate the melt-freezing process of uranium dioxide in a tube, based on several experiments [ 12 12 ] . In previous studies, the MPS method was used to simulate the process of vertical fluid descent on static fluids such as water and cooking oil using pipes. This has been validated using experimental data [ 13 13 , 14 14 ] . Moreover, the stratification of these fluids was simulated using a 2D MPS code [ 15 15 ] . These previous studies were simulated without heat transfer calculations between particles in the geometry because they were conducted at 300 K. Recently, a molten metal relocation experiment was conducted to investigate the performance in the absence of coolant [ 8 8 ] . The results showed how molten material damages the supporting system when a severe nuclear accident occurs. Although it was not an exact scale model of a nuclear accident, some physical phenomena may be experimentally simulated. However, some behaviors observed when the molten material interacts with WMP were not explained in detail. The physical behavior can be approximated using a numerical simulation. Therefore, this study aims to use the 2D MPS method to analyze the WMP downward relocation process. The focus of the numerical simulation is on the relocation process of molten Wood’s metal and the impact on hole formation of the WMP. The heat transfer process was included to simulate the melting process of the WMP owing to the molten impingement of Wood’s metal. The simulation results were compared and validated by the experimental results. Moreover, some predictions of the behavior of molten Wood’s metal can also be explained numerically. 2 Methods 1 2.1 MPS method 2 The MPS method is a fully mesh-less Lagrangian method for simulating incompressible media, such as solids and liquids [ 9 9 ] . The Lagrangian fluid dynamics are represented by particles that possess their own physical properties. For example, the velocity field is not a function of position, as is the case with Eulerian fluids. MPS has two specific definitions, namely the radius of interaction r e and the weight function w ( r ). The radius r e limits the interaction of the particles, as shown in Fig. 1 Fig. 1 . Moreover, it can be implemented as the physical influence of a particle on others within a certain distance. The weight function is calculated between two particles at a distance r , as shown in Eq. (1), which satisfies the condition that there is no self-interaction in one particle. <math display="block" id="M1"><mrow><mi>w</mi><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mrow><mo>{</mo><mrow><mtable><mtr><mtd><mrow><mn>1</mn><mo>−</mo><mfrac><mi>r</mi><mrow><msub><mi>r</mi><mtext>e</mtext></msub></mrow></mfrac><mn>0</mn><mo><</mo><mi>r</mi><mo>≤</mo><msub><mi>r</mi><mtext>e</mtext></msub></mrow></mtd></mtr><mtr><mtd><mrow><mn>0</mn><msub><mi>r</mi><mtext>e</mtext></msub><mo>≤</mo><mi>r</mi></mrow></mtd></mtr></mtable></mrow></mrow></mrow></math> <math display="block" id="M2"><mrow><msub><mi>n</mi><mi>i</mi></msub><mo>=</mo><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow></mrow></math> The summation of all particle weight functions is defined as the particle number density, shown in Eq. (2). The particle number density was used to evaluate the incompressible fluid density. The integration of the weight function with the particle number density can represent the flow of fluid in the interaction zone, and the number of particles in the interaction zone fluctuates. However, the weight function balances it by the definition that <math id="M3"><mrow><mi>w</mi><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow></mrow></math> becomes low as the distance of the particle approaches r e . Moreover, <math id="M4"><mrow><mi>w</mi><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow></mrow></math> will have an asymptotic value near a zero distance, which prevents fluid particles from passing through other fluids or wall particles [ 16 16 ] . The MPS method has an expression for gradient, divergent, and Laplacian models, shown in Eqs. (3), (4), and (5), respectively. <math display="block" id="M5"><mrow><mover accent="true"><mo>∇</mo><mo>→</mo></mover><msub><mi>ϕ</mi><mi>i</mi></msub><mo>=</mo><mfrac><mi>d</mi><mrow><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mfrac><mrow><msub><mi>ϕ</mi><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>ϕ</mi><mo>^</mo></mover><mi>i</mi></msub></mrow><mrow><msup><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mn>2</mn></msup></mrow></mfrac><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow></mrow></math> <math display="block" id="M6"><mrow><mover accent="true"><mo>∇</mo><mo>→</mo></mover><mo>.</mo><msub><mover accent="true"><mi>φ</mi><mo>→</mo></mover><mi>i</mi></msub><mo>=</mo><mfrac><mi>d</mi><mrow><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mfrac><mrow><msub><mi>φ</mi><mi>j</mi></msub><mo>−</mo><msub><mi>φ</mi><mi>i</mi></msub></mrow><mrow><msup><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mn>2</mn></msup></mrow></mfrac><mrow><mo>(</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>)</mo></mrow><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow></mrow></math> <math display="block" id="M7"><mrow><msup><mo>∇</mo><mn>2</mn></msup><msub><mi>ϕ</mi><mi>i</mi></msub><mo>=</mo><mfrac><mrow><mn>2</mn><mi>d</mi></mrow><mrow><mi>λ</mi><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mrow><mo>(</mo><mrow><msub><mi>ϕ</mi><mi>j</mi></msub><mo>−</mo><msub><mi>ϕ</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow></mrow></math> All three models are used to discretize the differential operator in the MPS governing equation. These three models depend on the weight function. In Eq. (5) of the Laplacian model, the λ value can be approximated by Eq. (6) to be proportional to the analytical solution [ 9 9 ] . <math display="block" id="M8"><mrow><mi>λ</mi><mo>=</mo><mfrac><mrow><msub><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></msub><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow><msup><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mn>2</mn></msup></mrow><mrow><msub><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></msub><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow></mrow></mfrac><mo>≅</mo><mfrac><mrow><msubsup><mstyle mathsize="140%" displaystyle="true"><mo>∫</mo></mstyle><mi>V</mi><mrow/></msubsup><mi>w</mi><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><msup><mi>r</mi><mn>2</mn></msup><mtext>d</mtext><mi>V</mi></mrow><mrow><msubsup><mstyle mathsize="140%" displaystyle="true"><mo>∫</mo></mstyle><mi>V</mi><mrow/></msubsup><mi>w</mi><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mtext>d</mtext><mi>V</mi></mrow></mfrac></mrow></math> The study includes three governing equations: mass, momentum, and energy conservation, shown in Eqs. (7), (8), and (9), respectively. <math display="block" id="M9"><mrow><mfrac><mrow><mi>D</mi><mi>ρ</mi></mrow><mrow><mi>D</mi><mi>t</mi></mrow></mfrac><mo>=</mo><mn>0</mn></mrow></math> <math display="block" id="M10"><mrow><mfrac><mrow><mi>D</mi><mover accent="true"><mi>u</mi><mo>→</mo></mover></mrow><mrow><mi>D</mi><mi>t</mi></mrow></mfrac><mo>=</mo><mo>−</mo><mfrac><mn>1</mn><mi>ρ</mi></mfrac><mo>∇</mo><mi>P</mi><mo>+</mo><mi>v</mi><msup><mo>∇</mo><mn>2</mn></msup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mo>+</mo><mover accent="true"><mi>g</mi><mo>→</mo></mover></mrow></math> <math display="block" id="M11"><mrow><mfrac><mrow><mi>D</mi><mi>h</mi></mrow><mrow><mi>D</mi><mi>t</mi></mrow></mfrac><mo>=</mo><mi>k</mi><msup><mo>∇</mo><mn>2</mn></msup><mi>T</mi></mrow></math> The fluid density is assumed to be constant as the incompressibility of the media is implemented, as shown in Eq.(7) [ 17 17 ] . The Navier–Stokes equation is used to calculate the movement of an incompressible fluid represented by the velocity change. The fluid particle moving mechanism is conducted by an external force, which is represented as the gravitational acceleration. Moreover, the viscosity determines the frictional force and momentum transfer between the particles when a collision occurs. Equation (9) shows that the energy conservation equation is composed of the conduction and heat source terms. Following the Laplacian form of the MPS method, Eq. (9) can be rewritten as <math display="block" id="M12"><mrow><msub><mi>h</mi><mi>i</mi></msub><msup><mrow/><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msup><mo>=</mo><msub><mi>h</mi><mi>i</mi></msub><msup><mrow/><mi>k</mi></msup><mo>+</mo><mi>k</mi><mi>Δ</mi><mi>t</mi><mfrac><mrow><mn>2</mn><mi>d</mi></mrow><mrow><mi>λ</mi><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mrow><mo>(</mo><mrow><msubsup><mi>T</mi><mi>j</mi><mi>k</mi></msubsup><mo>−</mo><mo> </mo><msubsup><mi>T</mi><mi>i</mi><mi>k</mi></msubsup></mrow><mo>)</mo></mrow><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><msup><mrow/><mo>*</mo></msup><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub><msup><mrow/><mo>*</mo></msup></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow><mo>.</mo></mrow></math> Moreover, the enthalpy is a function of temperature, as described in Eq. (11): <math display="block" id="M13"><mrow><mi>h</mi><mrow><mo>(</mo><mi>T</mi><mo>)</mo></mrow><mo>=</mo><mi>ρ</mi><msub><mi>C</mi><mi>p</mi></msub><mrow><mo>(</mo><mrow><mi>T</mi><mo>−</mo><msub><mi>T</mi><mn>0</mn></msub></mrow><mo>)</mo></mrow><mo>+</mo><msub><mi>h</mi><mn>0</mn></msub><mo>,</mo></mrow></math> <math display="block" id="M14"><mrow><msub><mtext>T</mtext><mi>i</mi></msub><msup><mrow/><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msup><mo>=</mo><msub><mtext>T</mtext><mi>i</mi></msub><msup><mrow/><mi>k</mi></msup><mo>+</mo><mfrac><mrow><mtext>k</mtext><mi>Δ</mi><mi>t</mi></mrow><mrow><mtext>ρ</mtext><msub><mi>C</mi><mi>p</mi></msub></mrow></mfrac><mfrac><mrow><mn>2</mn><mi>d</mi></mrow><mrow><mi>λ</mi><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mrow><mo>(</mo><mrow><msubsup><mi>T</mi><mi>j</mi><mi>k</mi></msubsup><mo>−</mo><mo> </mo><msubsup><mi>T</mi><mi>i</mi><mi>k</mi></msubsup></mrow><mo>)</mo></mrow><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><msup><mrow/><mo>*</mo></msup><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub><msup><mrow/><mo>*</mo></msup></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow><mo>.</mo></mrow></math> Each particle temperature is explicitly calculated based on Eq. (12). The calculation is performed after the particle positions are updated, and it does not include the temporal value of each step. The external force and viscosity term are calculated simultaneously with the particle temperature calculation, as shown in Eq. (13). Meanwhile, in Equation (14), the pressure term is calculated implicitly. The temporal position is updated owing to the velocity change, as shown in Eqs. (15) and (16): <math display="block" id="M15"><mrow><msubsup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mi>k</mi><mtext>*</mtext></msubsup><mo>=</mo><msubsup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mi>i</mi><mi>k</mi></msubsup><mo>+</mo><mi>v</mi><mi>Δ</mi><mi>t</mi><mfrac><mrow><mn>2</mn><mi>d</mi></mrow><mrow><mi>λ</mi><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mrow><mo>(</mo><mrow><msubsup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mi>j</mi><mtext>*</mtext></msubsup><mo>−</mo><msubsup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mi>i</mi><mtext>*</mtext></msubsup></mrow><mo>)</mo></mrow><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msubsup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mi>j</mi><mi>k</mi></msubsup><mo>−</mo><msubsup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mi>i</mi><mi>k</mi></msubsup></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow><mo>,</mo></mrow></math> <math display="block" id="M16"><mrow><msup><mo>∇</mo><mn>2</mn></msup><msubsup><mi>P</mi><mi>i</mi><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup><mo>=</mo><mo> </mo><mfrac><mrow><mn>2</mn><mi>d</mi></mrow><mrow><mi>λ</mi><msup><mi>n</mi><mn>0</mn></msup></mrow></mfrac><munder><mstyle mathsize="140%" displaystyle="true"><mo>∑</mo></mstyle><mrow><mi>j</mi><mo>≠</mo><mi>i</mi></mrow></munder><mrow><mo>(</mo><mrow><msubsup><mi>P</mi><mi>j</mi><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup><mo>−</mo><mo> </mo><msubsup><mover accent="true"><mi>P</mi><mo>^</mo></mover><mi>i</mi><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup></mrow><mo>)</mo></mrow><mi>w</mi><mrow><mo>(</mo><mrow><mrow><mo>|</mo><mrow><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>j</mi></msub><msup><mrow/><mo>*</mo></msup><mo>−</mo><msub><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>i</mi></msub><msup><mrow/><mo>*</mo></msup></mrow><mo>|</mo></mrow></mrow><mo>)</mo></mrow><mo>,</mo></mrow></math> <math display="block" id="M17"><mrow><msup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mrow><mo>*</mo><mo>*</mo></mrow></msup><mo>=</mo><msup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mo>*</mo></msup><mo>+</mo><mi>Δ</mi><mi>t</mi><mover accent="true"><mi>g</mi><mo>→</mo></mover><mo>,</mo></mrow></math> <math display="block" id="M18"><mrow><msup><mover accent="true"><mi>r</mi><mo>→</mo></mover><mrow><mo>*</mo><mo>*</mo></mrow></msup><mo>=</mo><msup><mover accent="true"><mi>r</mi><mo>→</mo></mover><mi>k</mi></msup><mo>+</mo><mi>Δ</mi><mi>t</mi><msup><mover accent="true"><mi>u</mi><mo>→</mo></mover><mrow><mo>*</mo><mo>*</mo></mrow></msup><mo>.</mo></mrow></math> Figure 2 Figure 2 shows the flowchart of MPS calculations; it begins with the preparation of the initial conditions. All particles have their initial conditions, such as the position, velocity, pressure, and temperature, and these data are stored in a grid data file. Other data that must be included in the initial condition are the average distance between particles, mass density, compressibility, and kinematic viscosity. The value of gravity is 9.81 m/s 2 for a particle influenced by gravity, such as a fluid. After all preparations, the calculation starts with the initial particle number density. There is a temporary velocity that is influenced by the viscosity and gravity. This temporary velocity triggers the movement of some particles. Eventually, the configuration of the particle changes. In addition, the temperature and pressure are calculated for each particle using Eqs. (12) and (14), respectively. The simulation is terminated after the steady-state condition is reached. 2.2 Verification of heat transfer calculation 2 To examine the accuracy of the MPS method, the heat transfer was verified using a 2D slab model, as shown in Fig.3(a) Fig.3(a) . The MPS calculation was performed assuming that no gravitational acceleration worked on the particles. Both materials are listed in Table 1 Table 1 ; the material specifications are close to those of stainless steel. The analytical solutions of the transient temperature profile were obtained using Eqs. (17) and (18) [ 18 18 ] : 1 1 Specification 1 1 Material 1 1 1 Material 2 1 1 Density (kg/m 3 ) 1 1 8.0 × 10 3 1 1 8.0 × 10 3 1 1 Temperature T n,i (K) 1 1 600 1 1 300 1 1 Thermal conductivity (W/m∙K) 1 1 3 1 1 30 1 1 Specific heat capacity (J/kg∙K) 1 1 500 1 1 500 <math display="block" id="M19"><mrow><msub><mi>T</mi><mn>1</mn></msub><mo>=</mo><msub><mi>T</mi><mrow><mn>1</mn><mo>,</mo><mi>i</mi></mrow></msub><mo>+</mo><mrow><mo>(</mo><mrow><msub><mi>T</mi><mi>s</mi></msub><mo>−</mo><msub><mi>T</mi><mrow><mn>1</mn><mo>,</mo><mi>i</mi></mrow></msub></mrow><mo>)</mo></mrow><mrow><mo>{</mo><mrow><mn>1</mn><mo>−</mo><mtext>erf</mtext><mrow><mo>(</mo><mrow><mfrac><mrow><msub><mi>x</mi><mn>1</mn></msub></mrow><mrow><mn>2</mn><msqrt><mrow><msub><mi>α</mi><mn>1</mn></msub><mi>t</mi></mrow></msqrt></mrow></mfrac></mrow><mo>)</mo></mrow></mrow><mo>}</mo></mrow><mo>,</mo></mrow></math> <math display="block" id="M20"><mrow><msub><mi>T</mi><mn>2</mn></msub><mo>=</mo><msub><mi>T</mi><mrow><mn>2</mn><mo>,</mo><mi>i</mi></mrow></msub><mo>+</mo><mrow><mo>(</mo><mrow><msub><mi>T</mi><mi>s</mi></msub><mo>−</mo><msub><mi>T</mi><mrow><mn>2</mn><mo>,</mo><mi>i</mi></mrow></msub></mrow><mo>)</mo></mrow><mrow><mo>{</mo><mrow><mn>1</mn><mo>−</mo><mtext>erf</mtext><mrow><mo>(</mo><mrow><mfrac><mrow><msub><mi>x</mi><mn>2</mn></msub></mrow><mrow><mn>2</mn><msqrt><mrow><msub><mi>α</mi><mn>2</mn></msub><mi>t</mi></mrow></msqrt></mrow></mfrac></mrow><mo>)</mo></mrow></mrow><mo>}</mo></mrow><mo>.</mo></mrow></math> The thermal diffusion coefficient is shown in Eq. (19). The values T 1,i and T 2,i are the temperatures of the material at the end of the slab, as shown in Fig. 3 Fig. 3 . At the beginning of the simulation, the red and blue particles were 600 and 300 K, respectively. The red and blue particles are materials 1 and 2, respectively. T s is the interface temperature calculated using Eq. (20). The interface temperature depended on the density, heat conductivity, specific heat coefficient, and initial temperature of both slabs. <math display="block" id="M21"><mrow><mi>α</mi><mo>=</mo><mfrac><mi>k</mi><mrow><mi>ρ</mi><msub><mi>c</mi><mtext>p</mtext></msub></mrow></mfrac></mrow></math> <math display="block" id="M22"><mrow><msub><mi>T</mi><mi>s</mi></msub><mo>=</mo><mfrac><mrow><msub><mi>T</mi><mrow><mn>1</mn><mo>,</mo><mi>i</mi></mrow></msub><msqrt><mrow><msub><mi>ρ</mi><mn>1</mn></msub><msub><mi>c</mi><mrow><mi>p</mi><mn>1</mn></mrow></msub><msub><mi>k</mi><mn>1</mn></msub></mrow></msqrt><mo>+</mo><msub><mi>T</mi><mrow><mn>2</mn><mo>,</mo><mi>i</mi></mrow></msub><msqrt><mrow><msub><mi>ρ</mi><mn>2</mn></msub><msub><mi>c</mi><mrow><mi>p</mi><mn>2</mn></mrow></msub><msub><mi>k</mi><mn>2</mn></msub></mrow></msqrt></mrow><mrow><msqrt><mrow><msub><mi>ρ</mi><mn>1</mn></msub><msub><mi>c</mi><mrow><mi>p</mi><mn>1</mn></mrow></msub><msub><mi>k</mi><mn>1</mn></msub></mrow></msqrt><mo>+</mo><msqrt><mrow><msub><mi>ρ</mi><mn>2</mn></msub><msub><mi>c</mi><mrow><mi>p</mi><mn>2</mn></mrow></msub><msub><mi>k</mi><mn>2</mn></msub></mrow></msqrt></mrow></mfrac></mrow></math> Figure 3(b) Figure 3(b) depicts the MPS simulation results of the two-slab heat transfer process at 4 s. The temperature tended to be constant at the edges of each slab. However, at the interface, there was a temperature change due to heat transfer from material 1 to material 2. Figure 4 Figure 4 shows a comparison between the MPS and the analytical solution. The discrepancy between the MPS and analytical calculations was miniscule. The discrepancy was 2 K, which was insignificant compared to the initial condition of the slabs with a temperature difference of 300 K. Therefore, this MPS method could be used for the next heat transfer calculation between two materials with a high temperature difference. 2.3 Analysis of experiment 2 Based on Sudha’s experiment, the initial condition was set in a stainless steel SS304L melt chamber consisting of a cylindrical and conical section [ 8 8 ] . Six ceramic beads and quartz-insulated tape heaters were used to heat the Wood’s metal. Moreover, an auxiliary pressure system could regulate the pressure of the material at 2 atm. The melt chamber system was capable of heating the metal up to 673 K. The chamber was thermally insulated to minimize ambient heat loss during the experiment. Seven K-type thermocouples were placed in various locations to control and monitor the Wood’s metal and vessel test temperatures. The thermocouples were insulated by SS-sheathed MgO. Two thermocouples were placed on the top and bottom surfaces of the first WMP. Then, the molten Wood’s metal was placed inside the melt chamber and preheated until it reached 573 K. It was released through a nozzle 60 mm in diameter from the melt chamber into a cylindrical test vessel. The stainless steel SS304L cylindrical test vessel was 950 mm in diameter and 1100 mm in height. However, the report did not mention the length of the melt-chamber pipe. Therefore, some assumptions were used in the simulation. Two WMPs were mounted on dummy SS inside the cylindrical test vessel. The diameter and thickness of the WMPs were 470 mm and 7 mm, respectively. The distance between the nozzle and the WMPs was 270 mm. The cylindrical test vessel was partially filled with water immediately below the first WMP. This scenario was analogous to a real condition inside a nuclear reactor in which the grid plate had a closed outer boundary. Figure 5 Figure 5 illustrates the experiment. The simulation focused on the red dashed line, which shows the downward relocation of the molten material to the plate. The experiments were conducted under two scenarios. First, the two WMPs were placed inside the cylindrical test vessel without any barrier restricting the outflow of molten material, and this scenario was used for the first to the fourth experiments. This caused the molten material to spill out directly into the water after impinging the WMP. In the second scenario in the fifth to eleventh experiments, there were perspex cylindrical rings to restrict the flow of all molten material on the plate. The sequence of the Wood’s metal impingement process was captured using a 1,000 fps high-speed camera. After conducting the first four experiments, the holes that formed on the WMP on both surfaces were measured. The formed hole areas on the top and bottom surfaces were 170 cm 2 and 20 cm 2 , respectively. The diameter of each surface was approximately 9.44 cm for the top surface and 5.04 cm for the bottom surface. The formed hole diameter data were used to validate the hole formation in the simulation results. 2.4 Simulation of downward relocation process 2 In this study, the 2D MPS simulation was mainly used to analyze the heat transfer between the molten Wood’s metal and the first WMP. Then, we approximated the hole formation behavior when the molten Wood’s metal breached the WMP. Therefore, not all geometries were built in the simulation. It focused only on the first to the fourth experiment without the perspex cylindrical rings mounted above the WMP. It is simulated in an x-y Cartesian coordinate system, as shown in Fig. 6 Fig. 6 . The initial pressure of the molten Wood’s metal inside the melt chamber was set at a zero-bar gauge, as experimentally conducted. Moreover, in the original MPS method, a pressure of 0 Pa was applied at the beginning of the calculation to avoid numerical instability owing to negative pressures [ 15 15 , 16 16 , 18 18 ] . Originally, there were four types of particles in MPS, namely fluid, rigid body, wall, and dummy wall particles. The fluid particles were used to represent the molten Wood’s metal. The rigid body particles represented a more viscous fluid than the previous particle. However, in this simulation, the rigid body particle was not used because it was necessary to present only one type of fluid particle. The dummy wall particle was a type of particle that behaved like a wall without the application of pressure and temperature calculations. It was used to give geometrical boundaries of the wall particle to the void condition. The wall particle was a type of solid particle in which the pressure and temperature calculations were performed, as well as fluid particles. Here, it was necessary to represent the WMP and SS304L as wall-type particles. Therefore, to distinguish between the two, a new stainless steel wall particle, namely a second wall, was declared, and is represented by the magenta color in Fig. 6 Fig. 6 . The particles in red, purple, and blue represent the molten Wood’s metal, WMP, and dummy wall, respectively. The initial temperatures of the molten Wood’s metal, stainless steel, and WMP were set to 573 K, 673 K, and 300 K, respectively. The diameter of all the particles was 1 mm. The total number of particles in the simulation was 15,875. The ratio of the particle number density was set to be 2.1 times larger than the particle diameter, which means that the r e value was 2.1 mm. Based on the experimental results, the phase change due to heat transfer occurred quickly. Therefore, the simulation performed the main critical time when the impingement and phase change of the WMP was 2 s. The time step was 0.001 s, so there were 2000 sequence results in one simulation. They can be represented as the sequences provided by the high-speed camera used in the experiment. However, the result sequences were presented with a time step of 0.1 s to minimize the memory utility. The calculation is assumed to be convergent when the temperature gradient in a limiting radius r e is sufficiently small, and may not conduct the heat transfer between particles, especially for the molten material and the metal plate. Table 2 Table 2 lists the physical properties of the Wood’s metal and stainless steel SS304L used in the simulation. Subsequently, the obtained experimental temperature data of both main thermocouples placed on the top and bottom surfaces of the first WMP were used to validate the simulation results. Both were obtained by digitizing the reported data particularly for the time evolution of the temperature and the prediction time when the hole on the surface started forming. Some information was also provided, including the remaining solid fraction of the WMP and the formed hole due to the impingement of the molten Wood’s metal. Therefore, it is the validation data for the MPS simulation results. a  at melting point 1 1 Parameter 1 1 Wood’s metal 1 1 SS 304 L 1 1 Density (kg/m³) 1 1 Molten 1 1 9.67 × 10 3 1 1 - 1 1 Solid 1 1 9.54 × 10 3 1 1 7.9 × 10 3 1 1 Melting point (K) 1 1 345 1 1 1727 1 1 Specific heat capacity (J/kg∙K) 1 1 168 1 1  500 1 1 Thermal conductivity (W/m∙K) 1 1 19 1 1  21.5 1 1 Latent heat (J/kg) 1 1 3.7 × 10 4 1 1  - 1 1 Viscosity (Pa∙s) a 1 1 1.5 × 10 -3 1 1 8.0 × 10 -3 3 Result and Discussion 1 Figure 7 Figure 7 depicts the simulation result of the downward relocation process at 0.4 s when the molten material impinged the WMP. It can be seen that there were color changes on some WMP surface particles, indicating heat transfer between the molten Wood’s metal and the plate. However, the hole was not observed at that time because the WMP surface particles remained in the original position. Along with the heat transfer between the two materials, the phase change would occur simultaneously with the velocity change, which causes the particle to move. After flowing down to the WMP, all of the molten Wood’s metal particles spread horizontally on the plate. The spread was counted as the leading-edge distance of the molten Wood’s metal particles. The detailed sequences of the molten Wood’s metal impingement on the WPM are shown in Fig. 7 Fig. 7 . Impingement was observed with no temperature change on the plate during 0.2 s of the molten Wood’s metal relocation. The range of velocities of the molten Wood’s metal when it reached the WMP surface was 1.25–2.0 m/s. The breaching process to the bottom surface of the WMP was found at 1.2 s. Figure 8(b) Figure 8(b) showed the breaching process of molten Wood’s metal at 1.3 s when the melted WMP started to flow down owing to the impingement and the melting process. It also indicated that the molten Wood’s metal had spread horizontally on the top surface of the WMP. The molten Wood’s metal temperature at the WMP top surface decreased owing to the interaction with the solid particles. Moreover, the hole formations became more prominent at 1.4 s, as shown in Fig. 8(c) Fig. 8(c) . The melted WMP particles gradually mixed with the molten particles. The particles spread on the top surface tended to move and then flow down through the formed hole, then reach the stainless steel cylindrical test vessel. Therefore, there was a decrease in the local height of the molten Wood’s metal, which spread qualitatively. Moreover, because of the high melting point of the cylindrical test vessel and the low melting point of the Wood’s metal, it could be assumed that there was no heat transfer between the molten material and the vessel. Figure 9 Figure 9 shows a comparison of the temperature profiles obtained from the simulation and experimental results. There was a significant change in the temperature of the WMP top surface at the beginning of the relocation of the molten Wood’s metal. The experimental data showed that the temperature surpassed the WMP melting point at 0.3 s. Meanwhile, in the simulation results, this occurred at 0.4 s, which indicates a phase change and the start of hole formation on the top surface of the WMP. Both experimental and simulation results showed that the maximum recorded temperature was 553 K with a maximum top side temperature at 0.5 s. Furthermore, there was a sharp temperature decrease in the experimental result, while the simulation results showed a gradual decrease in temperature. When the molten Wood’s metal vertically impinged the WMP, it led to the penetration of the bottom surface of the plate. Moreover, the cooling process of the molten Wood’s metal in air accelerated the temperature decrease on the WMP surface. In the experiment, the temperature change on the top surface appeared to have a sharply decreasing slope, while the temperature at the bottom surface increased significantly. The same temperature trend was observed in the simulation results. However, at a specific time step, the temperature tended to have a gentle slope until the end of the simulation. This was caused by the cooling process, which depended on the particle interactions, and was not simultaneously cooled by the real air condition. The bottom surface temperature also exhibited the same pattern, as shown in Figure 8 Figure 8 . The results showed that the WMP melting point was surpassed at 0.6 s and 1.2 s in the experiment and simulation, respectively. Experimentally, the molten Wood’s metal penetrated the plate in an interval time of 0.3 s. Meanwhile, the penetration time of molten Wood’s metal took 0.8 s in the simulation result. The penetration delay could be caused by the independence of the temperature calculation in the pressure module calculation. When the particles underwent a gradual temperature change over a certain time, the viscosity changed, which led to a pressure change around the particles at a certain effective radius. This change is significant when there is a small temperature gradient between the particles and small fluid thermal conductivity, but a large heat capacity is required to increase the temperature of the solid particles. In the simulation, the temperature gradient was high at the beginning, and the penetration occurred after less than 1 s. Therefore, there may be some delay, but the results showed that a similar bottom side maximum temperature was reached at 412 K and 410 K for the experiment and simulation, respectively. The small discrepancy did not significantly impact hole formation during the impingement of molten Wood’s metal on the WMP. The phase change of WMP particles due to the heat transfer process first occurred at 0.4 s, as shown in Figure 10 Figure 10 . The solid and liquid phases changed linearly and oppositely until the end of the simulation. The results indicated that the molten Wood’s metal phase increased, while that of the solid phase of the WMP decreased. The fraction change became stagnant after 1.6 s, along with the flowing down of molten Wood’s metal to the cylinder chamber. Generally, it is assumed that all molten materials reach the reactor catcher, which is located near the main reactor vessel. However, the experimental results showed that the fraction of molten material that reached the core catcher was only 0.720. Similarly, this simulation showed that only 0.750 of the molten material was formed until the end of the simulation time. Therefore, the simulation result has a 5.50% absolute error compared with the experimental results. <math display="block" id="M23"><mrow><mi>T</mi><mi>D</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>=</mo><mn>95</mn><mo>×</mo><mrow><mo>[</mo><mrow><mn>1</mn><mo>−</mo><mi>exp</mi><mrow><mo>(</mo><mrow><mo>−</mo><mn>3</mn><msup><mi>t</mi><mn>5</mn></msup></mrow><mo>)</mo></mrow></mrow><mo>]</mo></mrow><mo>,</mo></mrow></math> <math display="block" id="M24"><mrow><mi>B</mi><mi>D</mi><mrow><mo>(</mo><mi>t</mi><mo>)</mo></mrow><mo>=</mo><mn>50</mn><mo>×</mo><mrow><mo>[</mo><mrow><mn>1</mn><mo>−</mo><mi>exp</mi><mrow><mo>(</mo><mrow><mo>−</mo><mn>0.005</mn><msup><mi>t</mi><mrow><mn>15</mn></mrow></msup></mrow><mo>)</mo></mrow></mrow><mo>]</mo></mrow><mo>.</mo></mrow></math> The hole size or any impingement trace measurement may be performed after the experiment is completed. However, in the simulation case, it could be predicted step-by-step. The erosion and heat transfer process simulation results in previous studies were also validated by approaching the results of the experiment [ 12 12 , 19 19 - 21 21 19-21 ] . However, the validation did not clarify the dynamic behavior of the hole-size increment during the impingement process. Therefore, this simulation result revealed the physical trend of diameter size enlargement during the impingement of molten Wood’s metal on the WMP. Figures 11 Figures 11 and 12 12 depict the diameter of the formed hole due to the melting process in the WMP. The top and bottom holes started to form at 0.4 s and 1.2 s, respectively. The discrete data of hole size enlargement were approached by the fitting pattern, as shown in Eqs. (21) and (22). The coefficient values were 95 and 50 for the top and hole diameter enlargements, respectively. The calculation may result in a maximum of unity if there is no amplification factor from the coefficient. Therefore, the coefficient values were approximated based on the maximum diameter achieved in the experiment. The approximation was similar to the equation reported by Johnson, Mehl, Avrami, and Kolmogorov (JMAK) for the transformation kinetics [ 22 22 , 23 23 ] . In the JMAK equation, k and n characterize the reaction as a function of time, temperature, and other variables that have significant influences on the transformation process. Initially, the JMAK equation was used to explain the solid-state transformation kinetics due to isothermal diffusion. However, recent studies have revealed that the transformation also worked under non-isothermal conditions [ 24 24 – 26 26 24–26 ] . The enlargement of the hole diameter due to the Wood’s metal impingement occurred when there was a high temperature change in a short time, approximately 2 s. The enlargement was also caused by the phase change of the WMP from the solid phase to the liquid phase. Data fitting with Eq. (21) was performed to determine the enlargement behavior on the top surface of the WMP. The calculation result showed a significant increase at 0.6 s, while the simulation result was not the same. However, the pattern of calculation and simulation results showed the same value from 0.7–0.9 s. The saturation result started from 1.1 s along with a significant increase in the result from Eq. (22). Physically, it was explained as an indication that the size of the top surface-formed hole significantly increased at the beginning of the simulation until 0.9 s, as shown in Figure 10 Figure 10 . After that, it changed slightly owing to the horizontal particle movement so that the vertical penetration decreased. Thus, from 1 s to 1.6 s, the top diameter changed from 80 mm to 95 mm. As a result of vertical penetration, the WMP bottom surface would form a hole, although its size was not as large as that on the top surface, as shown in Fig. 12 Fig. 12 . The experiment result obtained previously showed that the formed hole diameters on the top and bottom surfaces were 94.43 mm and 50.48 mm, respectively. Likewise, the simulation results indicated that the diameters of the top and bottom sides were 95 mm and 50 mm, respectively. It can be concluded that the absolute errors for the top and bottom formed hole diameters were 1.2% and 1.8 %, respectively. Figure 13 Figure 13 shows the molten Wood’s metal outermost particle distance when it impinged on the WMP plate. This distance was calculated from particles distributed radially to the x-positive and x-negative. The previous research reported that the molten particles spread quickly above the WMP. The MPS simulation showed that at 1.1 s, all molten particles had spread evenly until the edge of the WMP. After that, some molten particles fell into the cylindrical test vessel. At the first time that molten particles impinged the WMP, their maximum velocity was 1.25 m/s. However, when all particles had spread to the edge of the WMP, the maximum velocity decreased to 0.22 m/s. The decrease was caused by a change of momentum and transfer of energy from the molten particles to WMP particles. 4 Conclusion 1 The improved MPS method has been used to analyze the molten metal behavior during relocation on the plate. The results show that during the spread of molten Wood’s metal, some energy or heat was also transferred to the WMP so there was a phase change from a solid to liquid phase. From 0.4 s to 1.6 s, there was a linear volume fraction change due to the solid–liquid phase change. At the end of the simulation, the total volume fraction of the formed-molten WMP was 0.750, and it had a 5.50% absolute error compared to the experiment. Another consequence of heat transfer between the molten Wood’s metal and WMP was hole formation in the center of the WMP. At the top and bottom surfaces of the WMP, holes with respective diameters of 95 mm and 50 mm were formed. The absolute errors for the top and bottom formed-hole diameters were 1.2% and 1.8%, respectively. The temperature profile was observed in detail in the simulation results, and showed that there was a decrease in the molten Wood’s metal temperature after breaching the plate. All the simulations showed good agreement, and provided an optimal prediction of the behavior of the molten material during the relocation and impingement processes, which was not explained by the experiment. Some empirical formulas have been introduced to explain the hole formation mechanism of the metal plate. For some improvement, a 3D simulation of the downward relocation process of Wood’s metal can be applied and compared with the 2D simulation results. Moreover, an improved MPS can be implemented to analyze the molten material behavior inside a nuclear reactor when a severe accident occurs. References T. Wada , Y. Abe , A. Kaneko et al ., Experimental study on jet breakup behavior with surface solidification. ICONE18-29618 , 635 - 644 , ( 2010 ). doi: 10.1115/ICONE18-29618 http://doi.org/10.1115/ICONE18-29618 Experimental study on jet breakup behavior with surface solidification. ICONE18-29618 A. Jasmin Sudha and K. Velusamy , Evaluation of post accident molten material relocation time to the core catcher in a fast breeder reactor . Prog. Nucl. Energy 75 , 149 - 157 ( 2014 ). doi: 10.1016/j.pnucene.2014.04.018 http://doi.org/10.1016/j.pnucene.2014.04.018 Evaluation of post accident molten material relocation time to the core catcher in a fast breeder reactor A. Furutani , S. Imahori , K. Sato et al ., Erosion behavior of a solid plate by a liquid jet - Effect of molten layer . Nucl. Eng. Des. 132 , 153 - 169 ( 1991 ). doi: 10.1016/0029-5493(91)90263-H http://doi.org/10.1016/0029-5493(91)90263-H Erosion behavior of a solid plate by a liquid jet - Effect of molten layer K. H. Bang , J. M. Kim , D. H. Kim , Experimental study of melt jet breakup in water . J. Nucl. Sci. Technol. 40 , 807 - 813 ( 2003 ). doi: 10.1080/18811248.2003.9715422 http://doi.org/10.1080/18811248.2003.9715422 Experimental study of melt jet breakup in water M. Mizanur Rahman , T. Hino , K. Morita et al ., Experimental investigation of molten metal freezing on to a structure . Exp. Therm. Fluid Sci. 32 , 198 - 213 ( 2007 ). doi: 10.1016/j.expthermflusci.2006.11.009 http://doi.org/10.1016/j.expthermflusci.2006.11.009 Experimental investigation of molten metal freezing on to a structure M. M. Rahman , Y. Ege , K. Morita et al ., Simulation of molten metal freezing behavior on to a structure . Nucl. Eng. Des. 238 , 2706 - 2717 ( 2008 ). doi: 10.1016/j.nucengdes.2008.04.008 http://doi.org/10.1016/j.nucengdes.2008.04.008 Simulation of molten metal freezing behavior on to a structure R. Park , K. Kang , S. Hong et al ., Corium behavior in the lower plenum of the reactor vessel under IVR-ERVC condition: Technical issues . Nucl. Eng. Technol. 44 , 237 - 248 ( 2012 ). doi: 10.5516/NET.03.2012.701 http://doi.org/10.5516/NET.03.2012.701 Corium behavior in the lower plenum of the reactor vessel under IVR-ERVC condition: Technical issues A. J. Sudha , M. Kumaresan , J. Anandan et al ., Experimental simulation of downward molten material relocation by jet ablation of structures and fuel coolant interaction in a fast reactor . Prog. Nucl. Energy 105 , 194 - 201 ( 2018 ). doi: 10.1016/j.pnucene.2018.01.010 http://doi.org/10.1016/j.pnucene.2018.01.010 Experimental simulation of downward molten material relocation by jet ablation of structures and fuel coolant interaction in a fast reactor S. Koshizuka and Y. Oka , Moving-particle semi-implicit method for fragmentation of incompressible fluid . Nucl. Sci. Eng. 123 , 421 - 434 ( 1996 ). doi: 10.13182/NSE96-A24205 http://doi.org/10.13182/NSE96-A24205 Moving-particle semi-implicit method for fragmentation of incompressible fluid A. P. A. Mustari and Y. Oka , Molten uranium eutectic interaction on iron-alloy by MPS method . Nucl. Eng. Des. 278 , 387 - 394 ( 2014 ). doi: 10.1016/j.nucengdes.2014.07.028 http://doi.org/10.1016/j.nucengdes.2014.07.028 Molten uranium eutectic interaction on iron-alloy by MPS method G. Li and Y. Oka , Sensitivity study of melt behavior of Fukushima Daiichi unit 1 type accident with MELCOR code and MPS method . J. Nucl. Sci. Technol. 52 , 109 - 121 ( 2015 ). doi: 10.1080/00223131.2014.933136 http://doi.org/10.1080/00223131.2014.933136 Sensitivity study of melt behavior of Fukushima Daiichi unit 1 type accident with MELCOR code and MPS method R. Chen , Y. Oka , G. Li et al ., Numerical investigation on melt freezing behavior in a tube by MPS method . Nucl. Eng. Des. 273 , 440 - 448 ( 2014 ). doi: 10.1016/j.nucengdes.2014.03.049 http://doi.org/10.1016/j.nucengdes.2014.03.049 Numerical investigation on melt freezing behavior in a tube by MPS method M. Ilham , Y. Yulianto , A. P. A. Mustari , Simulation on relocation of non-compressed fluid flow using moving particle semi-implicit (MPS) method . IOP Conf. Ser. Mater. Sci. Eng. 407 , 12100 ( 2018 ). doi: 10.1088/1757-899x/407/1/012100 http://doi.org/10.1088/1757-899x/407/1/012100 Simulation on relocation of non-compressed fluid flow using moving particle semi-implicit (MPS) method A. N. Hidayati , A. P. A. Mustari , A. Waris et al ., Study on fluid descends vertically on static fluid using moving particle semi-implicit method . J. Phys. Conf. Ser. 1772 , 012007 ( 2021 ). doi: 10.1088/1742-6596/1772/1/012007 http://doi.org/10.1088/1742-6596/1772/1/012007 Study on fluid descends vertically on static fluid using moving particle semi-implicit method Y. Yulianto , A. N. Hidayati , A. P. A. Mustari et al ., Moving particle semi-implicit (MPS) utilization in analyzing the stratification behavior of immiscible liquid . IOP Conf. Ser. Mater. Sci. Eng. 407 , 12189 ( 2018 ). doi: 10.1088/1757-899x/407/1/012189 http://doi.org/10.1088/1757-899x/407/1/012189 Moving particle semi-implicit (MPS) utilization in analyzing the stratification behavior of immiscible liquid S. Koshizuka , K. Shibata , M. Kondo et al ., Chapter 2 - Fundamentals of Fluid Simulation by the MPS Method . ( Academic Press , 25 - 109 , 2018 ). doi: 10.1016/B978-0-12-812779-7.00002-3 http://doi.org/10.1016/B978-0-12-812779-7.00002-3 Chapter 2 - Fundamentals of Fluid Simulation by the MPS Method R. B. Bird , E. W. Stewart , and N. L. Edwin , Transport Phenomena , The 2 nd Ed.( Wiley , 2001 ) Transport Phenomena A. Z. Şahin , Transient heat conduction in semi-infinite solid with spatially decaying exponential heat generation . Int. Commun. Heat Mass Transf. 19 , 349 - 358 ( 1992 ). doi: 10.1016/0735-1933(92)90081-R http://doi.org/10.1016/0735-1933(92)90081-R Transient heat conduction in semi-infinite solid with spatially decaying exponential heat generation M. Kondo and S. Koshizuka , Improvement of stability in moving particle semi-implicit method . Int. J. Numer. Methods Fluids 65 , 638 - 654 ( 2011 ). doi: 10.1002/fld.2207 http://doi.org/10.1002/fld.2207 Improvement of stability in moving particle semi-implicit method T. Kawahara and Y. Oka , Ex-vessel molten core solidification behavior by moving particle semi-implicit method . J. Nucl. Sci. Technol. 49 , 1156 - 1164 ( 2012 ). doi: 10.1080/00223131.2012.740944 http://doi.org/10.1080/00223131.2012.740944 Ex-vessel molten core solidification behavior by moving particle semi-implicit method G. Li , M. Liu , G. Duan et al ., Numerical investigation of erosion and heat transfer characteristics of molten jet impinging onto solid plate with MPS-LES method . Int. J. Heat Mass Transf. 99 , 44 - 52 ( 2016 ). doi: 10.1016/j.ijheatmasstransfer.2016.03.090 http://doi.org/10.1016/j.ijheatmasstransfer.2016.03.090 Numerical investigation of erosion and heat transfer characteristics of molten jet impinging onto solid plate with MPS-LES method F. Liu , F. Sommer , C. Bos et al ., Analysis of solid state phase transformation kinetics: Models and recipes . Int. Mater. Rev. 52 , 193 - 212 ( 2007 ). doi: 10.1179/174328007X160308 http://doi.org/10.1179/174328007X160308 Analysis of solid state phase transformation kinetics: Models and recipes X. Liu , H. Li , M. Zhan , A review on the modeling and simulations of solid-state diffusional phase transformations in metals and alloys . Manuf. Rev. 5 , 10 ( 2018 ). doi: 10.1051/mfreview/2018008 http://doi.org/10.1051/mfreview/2018008 A review on the modeling and simulations of solid-state diffusional phase transformations in metals and alloys R. E. Lyon , An integral method of nonisothermal kinetic analysis . Thermochim. Acta 297 , 117 - 124 ( 1997 ). doi: 10.1016/S0040-6031(97)00158-5 http://doi.org/10.1016/S0040-6031(97)00158-5 An integral method of nonisothermal kinetic analysis K. Fan , F. Liu , W. Yang et al ., Analysis of soft impingement in nonisothermal precipitation . J. Mater. Res. 24 , 3664 - 3673 ( 2009 ). doi: 10.1557/jmr.2009.0434 http://doi.org/10.1557/jmr.2009.0434 Analysis of soft impingement in nonisothermal precipitation S. Cui , J. Wan , X. Zuo et al ., Three-dimensional, non-isothermal phase-field modeling of thermally and stress-induced martensitic transformations in shape memory alloys . Int. J. Solids Struct. 109 , 1 - 11 ( 2017 ). doi: 10.1016/j.ijsolstr.2017.01.001 http://doi.org/10.1016/j.ijsolstr.2017.01.001 Three-dimensional, non-isothermal phase-field modeling of thermally and stress-induced martensitic transformations in shape memory alloys All authors contributed to the study conception and design. Anni Nuril Hidayati performed the numerical calculation, processed the experimental data, performed the analysis, drafted the manuscript and designed the figures. Abdul Waris, Asril Pramutadi Andi Mustari and Dwi Irwanto were involved in planning and supervised the work. Asril Pramutadi Andi Mustari and Nur Asiah Aprianti aided in interpreting the results and analysis. The first draft of the manuscript was written by Anni Nuril Hidayati and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript. The authors wish to thank Prof. S. Koshizuka, Prof. M. Sakai, and Dr. K. Shibata of the University of Tokyo for their helpful comments and for providing the basic MPS code for fluids. 10.1007/s41365-021-00922-x This work was supported by the Pendidikan Magister Menuju Doktor untuk Sarjana Unggul (PMDSU) and a research program from the Ministry of Research, Technology, and Higher Education, Indonesia. 2021-02-06 2021-05-17 2021-05-22 2021-08-23 2021-08-01 2021-09-07 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2021 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-021-00922-x http://dx.doi.org/10.1007/s41365-021-00922-x . 2021 NOMENCLATURE: BD BOTTOM SURFACE DIAMETER CHANGE CP HEAT CAPACITY D NUMBER OF DIMENSIONS FPS FRAMES PER SECOND G GRAVITATIONAL ACCELERATION H ENTHALPY K THERMAL CONDUCTIVITY MPS MOVING PARTICLE SEMI IMPLICIT N PARTICLE NUMBER DENSITY P PRESSURE R RADIUS BETWEEN TWO PARTICLES OR POSITIONS R E EFFECTIVE RADIUS S INTERFACE T TIME T TEMPERATURE TD TOP SURFACE DIAMETER CHANGE U VELOCITY W ( R ) WEIGHT FUNCTION WMP WOOD’S METAL PLATE Χ POSITION AT X-AXIS DIRECTION <math id="M25"><mo>∇</mo></math> GRADIENT <math id="M26"><mrow><mo>∇</mo><mo>.</mo></mrow></math> DIVERGENCE <math id="M27"><mrow><msup><mo>∇</mo><mn>2</mn></msup></mrow></math> LAPLACIAN GREEK LETTERS: Α THERMAL DIFFUSIVITY Δ INCMENT Λ COEFFICIENT IN LAPLACIAN MODEL Ν KINEMATIC VISCOSITY Ρ DENSITY Φ A SCALAR QUANTITY <math id="M28"><mover accent="true"><mi>ϕ</mi><mo>^</mo></mover></math> MINIMUM VALUE OF A SCALAR QUANTITY Φ. A VECTOR SUBSCRIPT: J NEIGHBORING PARTICLE I TARGET PARTICLE 0 INITIAL OR STANDARD VALUE 1 MATERIAL 1 2 MATERIAL 2 SUPERSCRIPT: 0 STANDARD VALUE K TIME STEP * TEMPORAL VALUE Nuclear Science and Techniques 08 32