Formation Mechanism of Laser-driven Magnetized “Pillars of Creation”

The Pillars of Creation, one of the most recognized objects in the sky, are believed to be associated with the formation of young stars. However, so far, the formation and maintenance mechanism of the pillars are still not fully understood due to the complexity of the nonlinear radiation magnetohydrodynamics (RMHD). Here, assuming laboratory laser-driven conditions, we studied the self-consistent dynamics of pillar structures in magnetic fields by means of two-dimensional and three-dimensional (3D) RMHD simulations, and the results support our proposed experimental scheme. We find that only when the magnetic pressure and ablation pressure are comparable, the magnetic field can significantly alter the plasma hydrodynamics. For medium-magnetized cases (β initial ≈ 3.5), the initial magnetic fields undergo compression and amplification. This amplification results in the magnetic pressure inside the pillar becoming large enough to support the sides of the pillar against radial collapse due to pressure from the surrounding hot plasma. This effect is particularly pronounced for the parallel component (B y ), which is consistent with observational results. In contrast, a strong perpendicular (B x , B z ) magnetic field (β initial < 1) almost retains its initial distribution and significantly suppresses the expansion of blown-off gas plasma, leading to the inability to form pillar-like structures. The 3D simulations suggest that the bending at the head of “Column I” in the Pillars of Creation may be due to nonparallel magnetic fields. After similarity scaling transformation, our results can be applied to explain the formation and maintenance mechanism of the pillars, and can also provide useful information for future experimental designs.


INTRODUCTION
The interstellar medium is photo-ionized by ultraviolet (UV) radiation from nearby massive stars, as documented by Oort & Spitzer (1955) and Pound et al. (2003).Most of the UV radiation in this process comes from two stars of spectral type O5 V and O5.5 V in the nearby young cluster NGC 6611, as noted by Pound (1998).This photoionization is a widespread phenomenon in the universe, resulting in the formation of over-pressurized bubbles of ionized gas known as H II regions (Churchwell 2002).Over a period of several million years, this process may give rise to the formation of thousands of stars due to the UV radiation ionization (White et al. 1999;Kalari et al. 2018;Robertson et al. 2010).The Pillars of Creation (POC,as shown in Fig. 1(b), also called 'elephant trunks'), located in Eagle Nebula, are typical representative objects of H II region and the reason why they are named so is that gas and dust are creating new stars, and they are also eroded by the UV radiation of the nearby stars.Recently, the James Webb Space Telescope captured a new high-resolution image of the POC, and the detail for the formation and development of young stars can be seen at the edges of these pillars (Adkins 2022;Dennis 2022).
Although a lot of observation (Hester et al. 1996;Indebetouw et al. 2007;Pattle et al. 2018) and numerical results (Lim & Mellema 2003;Mackey & Lim 2010, 2011;Hennebelle & Inutsuka 2019) have been acquired, the formation and maintenance mechanism of these pillars are still not fully understood due to the nonlinear RMHD process.The prevailing view is that the primary physical process for the pillar formation is related to the "rocket effect" (Spitzer Jr 1954), which is similar to the compression process of indirect-drive Inertial Confinement Fusion (ICF) fuel (Atzeni & Meyer-ter Vehn 2004).The UV radiation energy is absorbed near the cloud surface, which can lead to continued ablating of the surface and generating ionized blow-off flows away from the cloud surface.Because the ablation pressure generated by the "rocket effect" is several orders of magnitude higher than the initial pressure, strong shocks are driven into the cloud, compressing and heating it.There exist about three models for the formation of these pillars, namely the instability model, cometary model, and shielding model (Remington et al. 2006).The instability model suggested by Spitzer Jr. 1954 (Spitzer Jr 1954) mainly pertains to the Rayleigh-Taylor instability (RTI) that continues to grow at the ablation front due to preexisting density fluctuations or fluctuations in UV radiation.The pillars are analogous to the "spikes" of a heavy fluid penetrating through a light fluid, as noted by previous work (Frieman 1954).This analogy is also applicable to the ablative RTI during the compression process of ICF (Ye et al. 2010;Wang et al. 2014Wang et al. , 2017)).The cometary model is that the pillars are formed by stand-alone dense clumps.A portion of the clumped material is pushed away from the radiation source by ablation pressure, forming a comet-like structure (Lefloch & Lazareff 1994).The shielding model assumes that there are preexisting dense clumps within the low-density surrounding cloud, and the ablation pressure cannot push these dense clumps.These clumps roughly stay at their initial positions, while the rest of the cloud is moved a long distance by the ablation pressure, forming pillar structures (Williams et al. 2001).We can find that the latter two models are very similar, and the main difference is whether there is a background cloud.In recent years, the shielding model has received more attention, and a lot of numerical simulations are carried out based on this model (Lim & Mellema 2003;Pound et al. 2007;Mackey & Lim 2011).
Another confusing point about Eagle Nebula is the "missing stiffness" (Ryutov et al. 2005;Remington et al. 2006).According to the observation results (Levenson et al. 2000;Ryutov & Remington 2002), the material inside the clouds is very cold (10 − 30 K), and the number density of the unshocked material is about 10 3 − 10 4 cm −3 .However, the evaporated blow-off gas from the surface is very hot (10 4 K), and the ablation pressure is about two orders of magnitude higher than the gas pressure inside the cloud, which causes the clouds or pillars to collapse to a higher density (10 6 cm −3 ), contrary to observations.Some models, such as the quasi-homogeneous magnetic field model or the magnetostatic turbulence model (Ryutov et al. 2005), have been proposed to solve this paradox, and magnetic fields are believed to play a key role in these models.The near-infrared extinction observations (Sugitani et al. 2007) of the POC show magnetic fields are aligned along the pillars and the first high-resolution, submillimeter-wavelength polarimetric observations (Pattle et al. 2018) also show that the magnetic field runs along the length of the pillars and the strength is on the order 100 µGs.These results suggest that the magnetic field is important in the formation of pillar-like structures and it also may be propping up the pillars of creation (Pattle et al. 2018).Due to the limitations of astronomical observational facilities and the complex background, the finer magnetic field structure still remains indistinct.
Over the past two decades, laboratory astrophysics has become an important tool in astrophysics research (Remington et al. 1999(Remington et al. , 2006;;Lebedev et al. 2019), and it can help us to investigate the individual specific factor on the formation of pillars under the conditions of scaling laws (Ryutov et al. 1999(Ryutov et al. , 2001)).A novel experiment scheme was proposed to create a long-duration X-ray source to model the UV radiation produced by O stars (Casner et al. 2015).Then they use this X-ray source to simulate the shielding model in OMEGA and NIF facility (Kane et al. 2015;Pound 2017).The experiment results show a strong shock is launched into the CH foam, and the high-density clumps still roughly stay at the initial positions, forming a short pillar-like structure.
However, the effect of magnetic fields on the formation of a pillar-like structure is still a subject of debate, and it also has not been studied in laboratory experiments.In this paper, we use 2D and 3D RMHD simulations to study the detailed hydrodynamics process for the formation of the POC within the magnetic field.And these simulation results also provide information for future experimental designs.To mimic the UV radiation in Eagle Nebula, a long-duration X-ray source can be achieved through an array of small laser-driven cavities (Casner et al. 2015).A low-density foam with preexisting spherical dense clumps is illuminated by the X-rays, and three materials (pure plastic (CH), doping aluminum (Al), and doping gold (Au)) foam targets are used in our RMHD simulations to study the radiation cooling effect.The simulation results show that a strong shock can be launched in low-Z material (CH and Al), but the pillar-like structure collapses radially inward under the pressure of the surrounding hot plasma, unable to form a long collimated pillar structure.However, the existence of a magnetic field, especially the field parallel to the pillar structure (B y component in this work), can effectively support the sides of pillars against collapsing, which is consistent with the observation results (Sugitani et al. 2007;Pattle et al. 2018).The perpendicular components (B x , B z in Fig. 1(a)) of the magnetic field inside the foam are compressed and amplified to provide the balance with the ablation pressure.According to the 3D RMHD simulations, an asymmetrical pillar may be formed when the direction of the initial magnetic field is not perfectly parallel to the direction of shock propagation, which is similar to the 'Column I' in the POC.These findings are confirmed through RMHD simulations and the similarities between laboratory and astrophysical systems (Ryutov et al. 2001) suggest that our results can be applied to explore the formation mechanism of pillars of creation.

EXPERIMENTAL SCHEME AND SIMULATION SETUP
Figure 1.Schematic view of the experimental setup.A long-duration X-ray source is achieved through an array of small laser-driven cavities.These cavities are irradiated by a series of laser beams (as shown by "1" to "4" laser in panel (a)) in time (one after another), which can mimic the UV radiation from the nearby massive stars (as shown in panel (b)).The quasi-static magnetic field can be generated through the pulse-power system (for weak magnetic fields) or laser-driven capacitor-coil target (for strong magnetic fields), and a uniform magnetic field is generated within the Helmholtz coil.We could change the direction of the coil to study the hydrodynamic evolution under different magnetic field components.The X-ray illuminates a low-density foam with preexisting spherical dense clumps, and three materials (CH foam, doping Al foam, doping Au foam) are used in our radiation simulation, which can be used to mimic the formation mechanism of the POC.The corresponding astrophysical system is shown in panel (b).Image credits: NASA, ESA, and the Hubble Heritage Team (STScI/AURA).The initial density perturbation distribution is plotted in panel (c), and the colormap is corresponding to the logarithmic mass density in g/cm 3 .Panel (d) plots the initial radial distribution of foam density at x=0 (the position at the white dot line in panel (c)) As shown in Fig. 1(a), We propose a new experimental scheme to study the formation of pillar structure within magnetic fields.The quasi-static magnetic field can be generated using either a pulse-power system (for weak magnetic fields) or a laser-driven capacitor-coil target (for strong magnetic fields).A uniform magnetic field is produced within the Helmholtz coil.And the hydrodynamic evolution under different magnetic field components can be achieved by changing the direction of the coil.When simulating the formation of pillars in the Eagle Nebula, it is more appropriate for high-power laser facilities to illuminate the target with X-rays from hohlraums rather than directly irradiating the target with lasers.Direct laser irradiation would introduce complications that do not exist in astrophysical systems, such as laser-plasma interaction instability and non-local electron heat transport (Kruer 2000;Marocchino et al. 2013).A long-duration X-ray source is achieved through an array of small laser-driven cavities.These cavities are irradiated by a series of laser beams in time ("1" to "4, one after another), and the radiation temperature (T r ) of X-rays is about 60 − 100 eV.Such a platform can be used to mimic the formation of the pillars in the Eagle Nebula and the corresponding astrophysical system is shown in Fig. 1(b), and the similarity scaling laws between the astrophysical systems and laboratory systems are given in Section 4.
Based on such experimental setup, we conducted numerical radiation simulations.The 2D and 3D RMHD simulations are carried out by using the FLASH code (Fryxell et al. 2000), which includes many high-energy-density physics modeling capabilities, such as laser energy deposition, multi-temperature (T e ̸ = T i ̸ = T rad ), electron thermal conduction, and radiation transport, etc.In 2D RMHD simulations, the simulation domain size is set as (x, y) = 4000 µm×8000 µm with adaptive mesh refinement applied and the highest resolution 8µm achieved.And the simulation geometry is Cartesian geometry.For 3D RMHD simulations, the simulation domain size is set as (x, y, z) = 4000 µm × 8000 µm × 4000 µm with the highest resolution 40 µm achieved.The X-ray radiation with T r = 60 eV is set at the bottom boundary to simulate the long-duration X-ray source by an array of small laser-driven cavities, as shown in Fig. 2(g).The Courant-Friedrichs-Lewy (CFL) is set as 0.4 and a third-order interpolation is used to reach a balance between accuracy and stability in simulation.The equation of state (EoS) and opacity of different materials (CH, Al, Au) are calculated from the code BADGER (Heltemes & Moses 2012) and IONMIX (MacFarlane 1989) respectively.For foam targets doped with Al or Au elements, our targets will comprise 70% CH and 30% high-Z metal element doping.Such the composition of the target is expected to significantly decrease the electrical conductivity (Alhamidi et al. 2022), and the magnetic fields generated by the coils can diffuse smoothly throughout the entire foam target.
The initial density of foam target is set as ρ foam = 30 mg/cm 3 .For the preexisting spherical dense clumps, considering the distribution of target manufacturing in the experiment and astrophysical perturbations (Mackey & Lim 2011), we use a triple-density sphere as a perturbation in the simulation.The initial density perturbation distribution is plotted in Fig. 1(c).The distribution consists of three distinct regions.The first region, located within a radius of 100 µm, is a highly dense core with a density of ρ core = 1.05 g/cm 3 .The second region, with a radius ranging from 100 µm to 300 µm, is the middle part of the target where the density is set to 0.5 g/cm 3 .The third region, with a radius ranging from 300 µm to 800 µm, is the outer part of the target where the density is set to 0.1 g/cm 3 .And the initial radial distribution for the density of the target at x = 0 is plotted in Fig. 1(d) Due to the configuration requirement of the code "FLASH", a low-density (1 × 10 −6 g/cm 3 ) Helium (He) background plasma is set in other places in simulations.All the initial temperatures are set to be uniform as room temperature 290 K.The boundary conditions for fluid and radiation transport are all set to open.

Radiation cooling for different materials
In this section, we investigate the radiation-cooling effect on plasma dynamics by examining three materials (CH, Al, Au).According to previous work Neufeld et al. (1995), the radiation cooling time, which is the ratio of thermal energy density to radiation power per unit volume, ranges from 10 2 − 10 3 years for the POC.As for the hydrodynamic motion of POC, it has a characteristic time scale of 10 4 − 10 5 years, indicating that the shocked material cools down quickly due to radiation cooling.In laboratory experiments and simulations, we also need to consider the radiation-cooling effect for the formation mechanism of pillar structure.As mentioned in previous work (Tikhonchuk et al. 2008), the radiative cooling time is inversely proportional to the nuclear charge Z when the temperature and number density are approximately equal.To investigate this effect, we selected three materials: pure CH plastic foam, Al-doped foam, and Au-doped foam.These materials represent low-Z, middle-Z, and high-Z materials, respectively.
As shown in Fig. 2(a), a high-density clump with a radius of 800 µm is located at y = 2 mm.X-rays inject from the bottom boundary of the simulation box and X-ray radiation is absorbed in the low-density foam surface by the photo-ionization process (Fig. 2(a) to 2(e)).The ionized blow-off gas plasma begins to expand away from the foam surface, but the expansion velocities are different for the three materials.The CH gas plasma expands fastest, and the velocity is about −4.7 × 10 7 cm/s.The expansion velocities of Al and Au gas plasma are −3.0 × 10 7 cm/s and −6 × 10 6 cm/s respectively.Then the pressure from the "rocket effect" compresses the foam, and according to the rocket model, the compression velocity is proportional to the mass ablation rate and blow-off gas velocity (Atzeni & Meyer-ter Vehn 2004).We can find that the compression velocity is fastest in CH plasma, and a strong shock is launched into the foam for CH and Al material, compressing and heating the foam.For the Au material, there is no shock formation within the foam due to strong radiation and low-speed blow-off gas.
Although there is no shock formation within the Au material, the preexisting dense clump within the foam is also ablated by X-rays, as indicated by the black arrow in Fig. 2(e).According to the distribution of the radiation temperature, we find there is a significant difference in the deposition mode of X-ray energy between Au and CH (or Al) material.For CH or Al material, the X-ray energy is deposited in a very thin layer near the surface, and the temperature of undisturbed foam is much lower than the ablated gas plasma.While for the Au material, due to the strong radiation characteristics of Au materials, the energy of X-ray has been significantly deposited to the inner side, and the high-density clump is also preheated.However, as the density increases, X-rays cannot continue to deposit energy inward.As shown in Fig. 2(i), the X-ray does not deposit energy into the interior of the high-density clump.This can be explained by the photon mean free path of radiation.The typical absorption coefficient (Johnston & Dawson 1973) in plasma can be expressed as where ln Λ is Coulomb logarithm, k B T ≈ 80 eV are the Au plasma temperature (in 2(l)), Z ≈ 25, N e are the average ionization and electron number density, ν and ν p are radiation frequency and plasma frequency.The estimated radiation photon (energy = 60 eV) mean free path l f ree ≈ 1/K ∝ 1/N 2 e .As shown in 2(f), the electron number densities for the foam and clump of Au are 10 20 cm −3 and 4 × 10 21 cm −3 respectively, and the corresponding radiation photon mean free paths are 7.7 cm and 50 µm.According to the above results, for strong radiation materials, a relatively higher-density foam is required to allow the X-ray energy to be deposited in a thin layer, then a shock can be launched into the foam, but it may need more time for the shock traversing the target.
As mentioned above, radiation cooling is significant to the evolution of the POC.For the laboratory system, the radiative cooling time is the ratio of the plasma thermal energy and radiated power.For simplicity, here we can ignore the line emission and account for the bremsstrahlung emission only (Tikhonchuk et al. 2008).The bremsstrahlung emission power per unit energy and volume is scaled as J br (ν, T ) ∝ T − 1 2 n e n i Z 2 e −hν/kT (ergs −1 cm −3 Hz −1 ) (Halverson 1972), where ν, T represent the photon frequency and temperature; and n e , n i , Z are the number density of electron, ion, and the average ionic charge.We can integrate emission over the whole spectrum to get the power per volume scales as P br ∝ T 1/2 n e n i Z 2 ≈ 1.7 × 10 −32 Zn 2 e T 1/2 W/cm 3 (Richardson 2019).The energy density of plasma can be defined as (Nicolai et al. 2006) 2 n e T ≈ 2.4 × 10 −19 n e T J/cm 3 .Then the radiation cooling time can be calculated as τ rad = E p /P br = 1.4 × 10 22 T 1/2 /(Zn e ) ns.Here, we mainly focus on the CH and Al material, where a strong shock is launched in these materials, and the related parameters are T CH = 100 ev, ZCH = 2.5, n e,CH ≈ 6 × 10 21 cm −3 , T Al = 70 ev, ZAl = 10, n e,Al ≈ 6 × 10 21 cm −3 (see Fig. 2(f) and 2(l)).Then the typical radiation cooling times can be obtained as τ rad,CH ≈ 9.3 ns, τ rad,Al ≈ 1.9 ns respectively.The hydrodynamic time is about τ hd = 120 ns for the shock wave traversing all foam, and the radiation cooling parameter is χ = τ rad /τ hd .For the POC, the radiation cooling parameter is in the range ∼ 10 −3 − 10 −2 , and the radiation cooling parameter χ Al ≈ 1.5 × 10 −2 of Al material is closer to the astrophysical system, so for the experiment, CH foam doped with Al is a better choice.In the following section, we mainly focus on the Al material.The formation of the pillar-like structure is mainly due to the interaction of the shock wave with the pre-existing highdensity clumps, which are similar to the shock-bubble interaction in fluid mechanism (Ranjan et al. 2011;Niederhaus et al. 2008).The Atwood number is defined as A = (ρ clump −ρ foam )/(ρ clump +ρ foam ), where ρ clump , ρ foam represent the density of unshocked clumps and foam, and here we can obtain A > 0. The specific heats γ clump , γ foam are the same due to the same material, then we can get the sound speed in unshocked clumps to be smaller than foam (c s,clump < c s,f oam ).Thus the incident shock wave is refracted while crossing the curved upstream clump surface, owing to the change in sound speeds, and here the refraction is convergent so that the transmitted shock wave is concave curvature.As shown in Fig. 3(a), an incident shock moves inward and the clump is compressed by this shock (compression phase).And according to the relationship:

The formation of pillar-like structure
where M i is the Mach number of the incident shock wave, we can obtain that R tt 0 = 1.The initial density ratio R 0 = ρ clump /ρ foam > R tt 0 = 1, thus the reflected wave is a shock wave (Niederhaus et al. 2008).Fig. 3(c) and 3(d) are corresponding to this compression phase.The shock structure is also similar to Fig. 3(a), and a weak reflected shock wave moves back-forward to the shocked low-density foam.
When the incident shock pass through all clump, due to the high Atwood number (A max = (1.05− 0.03)/(1.05+ 0.03) ≈ 0.94), portions of the shock front sweeping around the clump periphery are diffracted (Niederhaus et al. 2008).They turn towards the axis so that the discontinuous surface is almost perpendicular to the interface.These diffracted shock waves may then converge with the transmitted shock wave at the downstream pole, leading to shock focusing, as shown in Fig. 3(b).After this shock passes through the clump, a rarefaction wave travels back to the shocked materials, and then the whole clump is accelerated (acceleration phase).Figure 3(e) is corresponding to the acceleration phase.The shock passes through all clumps and is focused towards the axis, and a short pillar-like structure is formed after the high-density clump core, which is similar to the previous experiment result ( Kane et al. 2015).The velocity in the y-axis near the head (where the position of the high-density clump core) of this pillar is about 5 × 10 5 cm/s, while at the rare surface of the pillar, the velocity is as high as 6 × 10 6 cm/s.So this heavy clump core can be viewed as nearly stationary.In Fig. 3(g), we find the pillar-like structure collapses radially inward due to the surrounding hot plasma, which leads the pillar heads to be disconnected from the tails.Figure 3(h) plots the mass densities at x = 0 for t = 0, 30, 40, 60, 100 ns respectively, and both compression and acceleration phase are shown.In the compression phase, the maximum density at the shock position is about 3.5 − 4 times higher than the initial density, which means a strong shock is launched into the clump.
The aspect ratio of length to width for the laboratory pillar is much smaller than the POC in Eagle Nebula.If there is no mechanism to support the sides of pillars against collapsing radially due to the surrounding hot plasma, the pillar may be destroyed in the end.The formation of the POC is related to many factors, such as multiple pre-existing disturbances (Lim & Mellema 2003;Pound et al. 2007), magnetic fields (Mackey & Lim 2011;Pattle et al. 2018), etc.In the following sections, we mainly focus on the effect of different components of magnetic fields.

The effect of magnetic field
The previous numerical simulation results (Mackey & Lim 2011) suggest that magnetic fields are dynamically important in the formation of pillar structures.Sugitani et al. 2007(Sugitani et al. 2007) used near-infrared polarimetry to measure the magnetic field in the Eagle Nebula, and they found a difference in the direction of magnetic fields between the pillars and the surrounding photo-ionized clouds.K. Pattle et al. 2018(Pattle et al. 2018) used the 850 µm polarized light to give a high-resolution map of the magnetic field in the dense gas of the photo-ionized pillars.They show that the magnetic fields run along the length of the pillars, perpendicular to and decoupled from the field in surrounding clouds, and this configuration can also support the pillars.These observation results show that there exist large-scale ordered magnetic fields near the pillar, but the morphology of the magnetic field in Eagle Nebula is complicated, and the finer magnetic field structure is still unclear.
In this section, we study the effect of different magnetic field components on forming a pillar-like structure.The initial thermal pressure of ionized blow-off gas plasma is about P T ≈ 3×10 11 dyn/cm 2 , then we define the initial plasma beta as β initial = P T /P B = P T /(B 2 /8π).The simulation in this study includes cases of magnetic field components perpendicular (B x , B z ) and parallel (B y ) to the pillars, as shown in Table 1.The cases can be categorized into three types: weakly magnetized (cases 1 to 3), medium magnetized (cases 4 to 6 and cases 10 to 11), and strongly magnetized (cases 7 to 9).Cases 10 to 11 correspond to 3D simulations.The magnetic field is frozen in the plasma flux due to the large magnetic Reynolds number (Crutcher 2012).c) correspond to the logarithmic mass density in g/cm 3 for case 1, case 2, and case 3 at 5 ns respectively, and the color bar is given at the bottom of these panels.The white dot line shows the front position of the blow-off gas plasma for the case without magnetic fields (as shown in 2(b)).Panels (d) to (f) are the corresponding magnetic pressure for these cases at 5 ns in logarithmic values, and the black arrows represent the magnetic field lines.The distributions of mass density at 120 ns for these cases are shown in panel (g) to (i), and the corresponding magnetic pressures are given in (j) to (l) in logarithmic values.

Weakly magnetized cases
The weakly magnetized cases of different magnetic components are shown in figure 4. The white dot line in 4(a) to 4(c) represents the front position of the blow-off gas plasma for the case without magnetic fields.We find that the magnetic field components (case 1, case 3) perpendicular to the pillar will suppress the movement of the blow-off low-density plasma.At the same time, the magnetic field is also significantly compressed and amplified in these cases.The peak strength of the magnetic field is 60 × 10 4 Gs, which is 3 times higher than the initial magnetic field.However, the magnetic pressure (1 × 10 10 dyn/cm 2 ) is still much smaller than the thermal pressure value, so the hydrodynamic effect is still dominant.For case 2 (4(b)), the front position of the blow-off gas plasma is similar to the case without magnetic fields, and there is no significant change in the distribution of the magnetic field in the early stage (4(e)).
At 120 ns, the morphologies of pillar-like structure (for case 1 to case 3) are almost the same as the case without magnetic fields (3(g)), but we can see that the field orientation is significantly changed by the dynamics of the photonionization process.The magnetic Alfvén velocity within the unshocked material can be calculated as v A = B/ √ 4πρ = 20 × 10 4 Gs/ 4π × 0.03 g/cm 3 ≈ 5.3 × 10 5 cm/s, which is much smaller than the shock speed.Then, the formation of a pillar-like structure is much faster than the time of the magnetic field relaxing to a lower energy configuration.For case 1 (4(j)), the magnetic fields roughly keep the initial orientation near the shock front, while within the pillar-like structure, we find that the magnetic field is compressed and amplified to 100×10 4 Gs, and the direction of the magnetic field becomes parallel to the pillar, which is similar to the high-resolution observations (Sugitani et al. 2007;Pattle et al. 2018).For case 2, at the position of the shock front in 120 ns (4(k)), a perpendicular component appears.This is because when X-rays interact with a spherical high-density clump, a plasma is generated that expands approximately spherically, leading to a significant change in the direction of the magnetic field in case 2. Within the pillar-like structure, the magnetic fields still keep parallel to the pillar.For the B z component case, magnetic field lines cannot be drawn directly in 2D.According to the principle of symmetry, the evolution should be similar to the B x case, and as shown in Fig. 4(l), the distribution of magnetic pressure is also similar to Fig. 4(j).

Medium magnetized cases
The medium magnetized cases of different magnetic components are shown in figure 5.The magnetic components perpendicular to the pillar (case 4,6) can suppress the expansion of low-density plasma more effectively than the weakly magnetized cases.As shown in 5(a), the B x component effectively suppresses the expansion of low-density plasma, and the magnetic field is compressed and amplified to 350T (on the right-side of 5(a)).The amplified magnetic pressure is already close to the thermal pressure, which plays an important role in the subsequent evolution of plasma dynamics.For the initial components parallel to the pillar (case 5), low-density plasma motion along the magnetic field lines is barely affected and there is also little change in the magnetic field distribution in the early stage (5(b)).
Here, we mainly focus on the formation of a pillar-like structure in the late stage.The distributions of logarithmic mass density at 120 ns for different cases are plotted in 5(c) to 5(e), and the corresponding magnetic field lines are given in 5(c) and 5(d) through black arrow lines.For case 4 (5(c)), a very narrow pillar-like structure forms.The magnetic field inside the cloud is compressed and amplified, so the magnetic pressure can provide the balance with the ablation pressure, and then the shock wave propagation speed drops than without the magnetic field case.We also find a noticeable perturbation in the undisturbed low-density foam behind the shock front, and the perturbation thickness is about 2 mm.This perturbation may be due to Alfvén waves, which are perpendicular to the magnetic field.The mean shock velocity is about 4 − 5 × 10 6 cm/s, and the magnetic Alfvén velocity (after magnetic field compression and amplification) is about 9 × 10 6 cm/s, so the disturbance of Alfvén waves can transfer faster.The initial magnetic field lines are perpendicular to the pillar and at 120 ns, most of the field lines still remain perpendicular.This indicates that with the increase of the magnetic field, especially when the magnetic pressure and thermal pressure are comparable, the dynamic behavior of the plasma is difficult to significantly change the topology of the magnetic field.
For case 5, a long and stable pillar is formed, and the morphology of the pillars is most similar to the POC.As shown in Fig. 5(d), the field lines are perpendicular to the boundary of the pillar, and within the pillar, the magnetic field still remains in the initial direction, but the strength has been amplified to 300 T. The amplified magnetic pressure prevents the shock waves on both sides of the high-density clump from converging toward the axis (as shown in the right side of 5(h)), and the magnetic field strength is large enough to magnetically support the sides of pillars against collapsing radially, so a long and stable pillar can be formed.For the B z component case, a triangle-like structure is formed at the head of the pillar, and a vortex appears after the triangle head.This vortex instability may be due to the large Larmor radius instability within the vertical magnetic field (Ripin et al. 1987).This instability is mainly driven by gravitation or acceleration drifts such that v g ≫ v di , where v g = g/ω ci , v di = T i /ZeL n B are the gravity drift and diamagnetic drift respectively, and g, ω ci , T i , Z, e, L n , B are effective gravity, ion gyro-frequency, ion temperature, ion  b) correspond to the logarithmic mass density in g/cm 3 for case 4, case 5 at 5 ns, and the right sides correspond to the distribution of magnetic fields.Panels (c) to (d) plot the logarithmic mass density for cases 4 to 6 at 120 ns.The mass densities of low-density foam at x = −1.5 mm for different cases are shown in panel (f).The left sides of panels (g) and (h) plot the magnetic pressure for case 4, and case 5 at 120 ns, and the right sides are the corresponding thermal pressure.The left sides of panels (i) to (l) are the radiation temperature in eV for case 4 to case 6 and case without magnetic fields, and the right sides of these panels are the corresponding electron temperature.
charge state, the elementary charge, plasma-density-gradient scale length, and strength of the magnetic field.Here, the condition of v g ≈ 4700 ≫ v di ≈ 5 is met, so the large Larmor radius instability can grow up at the boundary of pillar, and it is similar to the experiment results (Tang et al. 2020).The mass densities of low-density foam at x = −1.5 mm for different cases are shown in Fig. 5(f).We find that the mass density distribution of case 4 and case 6 are similar.The peak number density of case 5 is similar to the case without fields, and it is about 3 − 4 times higher than the peak density in case 4 or 6.
The distributions of electron temperature and radiation temperature for cases 4 to 6 at 120 ns are shown in Figs.5(i) to 5(k), and 5(l) plots the case without magnetic fields.We find that in the cases with magnetic fields, the head, middle, and tail of the pillar maintain a low-temperature structure, and the radiation temperature inside the pillar is also very low, so the tail will not be separated from the head.In Fig. 5(l), for the case without magnetic fields, only the head of the pillar keeps a low temperature, while the tail has been effectively heated, whose temperature is about 30 − 40 eV.
The gas emission in recombination radiation (e.g.H α or X-ray) is one of the most important observational features for the POC (Linsky et al. 2007;Sofue 2020), and the most luminescent part is located on the pillar's head, followed by a long shadow tail (as shown in Fig. 1(b)).In the experiments, the bremsstrahlung emission, as we mentioned above, is the main radiation mechanism for the high-temperature low-Z plasma (Nicolai et al. 2006), and the radiation maps for medium magnetized cases are shown in figure 6. 6(a) corresponds to the case without magnetic fields, and the distinct radiation features appear on both the head and tail of the pillar, which is quite different from the observation results.For cases with magnetic fields, their radiation signature is very similar to the observations, where the most intense radiation is located on the pillar's head, followed by a long shadow tail.The case in which the B y component dominates (case 5, where the magnetic field lines are parallel to the pillar) is most similar to the observed result.
Strong magnetized cases If we continue to increase the magnetic field strength, when the initial magnetic pressure is larger than the thermal pressure, the plasma dynamics are significantly affected by the magnetic field.As shown in Fig. 7(a) and (c), the movement of blow-off gas flows is significantly suppressed when the magnetic field component is perpendicular to the pillar, and almost no obvious forward propagating shock is formed.The initial magnetic field is compressed and amplified to about 350 T. Compared with the weak magnetized cases and medium magnetized cases, the magnetic field is amplified to a smaller extent, only 1.1 times the initial field strength.And the magnetic pressure inside the pillar is also much higher, but this does not cause rapid lateral diffusion and destruction of the pillar structure.For strong magnetized cases, the topology and distribution of the magnetic field are not significantly affected because the thermal pressure of the blow-off gas is lower than the initial magnetic pressure.And total pressure of residual magnetic pressure outside the pillar and the thermal pressure can maintain the stability of the pillar structure.As shown in Fig. 7(e), pillar-like structures cannot be formed even at later stages (120 ns), and only an arc-like structure is formed.For case 9 (7(b)), since the plasma is not hindered along the direction of the magnetic field lines, even if the initial plasma β is less than 1, the motion of the plasma is not significantly different from the case without magnetic fields.Similar to the medium case, the magnetic pressure can support the sides of pillars against collapsing radially under pressure from the surrounding hot ablated plasma (7(d)), but the head of the pillar forms a more pointed structure.
In brief, for the weakly magnetized cases, the plasma thermal pressure clearly dominates and the morphology of pillar-like structure is similar to the case without magnetic fields.The magnetic field orientation is obviously changed by the plasma dynamics, and we find the magnetic field lines all become parallel to the pillar within the pillar structure in a later stage for different cases.In medium magnetized cases, the initial magnetic fields are compressed and amplified from 150 T to 300 T, resulting in an increase in magnetic pressure by 4 − 5 times.This magnetic pressure is close to the thermal pressure.And it is sufficient to support the sides of pillars against radial collapse under pressure from the surrounding hot plasma, thereby maintaining the structure of the pillar and increasing its longevity.According to the radiation maps of medium magnetized cases, when the magnetic field lines are parallel to the pillar, the radiation feature is most similar to the observed results, and the field orientation within the pillar is also consistent with the observation results (Sugitani et al. 2007;Pattle et al. 2018).For the strong magnetized cases, the formation of pillars is significantly suppressed when the magnetic field component (B x , B z ) is parallel to the shock front.In order to avoid dimensional effects due to 2D simulation, we also performed full 3D simulations.As shown in Fig. 8(a) and (b), the setup of case 10 (3D simulation) is the same as case 5 (2D simulation).The distribution of mass density (8(a)) shows a nearly axisymmetric and long pillar structure is formed, and the 2D simulation result (5(d)) is also similar to this result, which confirms the reliability of the 2D simulation.The magnetic field strength before the pillar is reduced (as shown by arrow lines in 8(a)), and the magnetic field within the pillar is compressed and amplified to 280 T. The magnetic field lines are perpendicular to the pillar boundary outside the pillar.The distribution of radiation shows the most luminescent is located on the head of this pillar structure, and it is also a pillar-like structure.

3D simulation results
However, it is difficult for the magnetic field to be perfectly parallel or perpendicular to the pillar structure in Eagle Nebula.So we refer to Ryutov's model (Ryutov et al. 2005), making a natural assumption that the parallel (B y ) and perpendicular (B x ) components are initially comparable (B x = B y = 106 T), and the total initial strength of the magnetic field remains at 150 T. As shown in Fig. 8(c), the perpendicular component of the magnetic field inside the foam is compressed by the shock wave, and magnetic fields left in the shocked blow-off gas is mainly the parallel component, which is well consistent with Ryutov's model (Ryutov et al. 2005).But the mass density distribution map shows that the pillar structure is asymmetrical, not a collimated long pillar structure, which gradually bends from the head to the tail, forming an arc-shaped structure.And the morphology of this asymmetrical pillar is similar to the head of 'Column I' in the POC (as shown by the green circle in 1(b)).The radiation distribution also has a similar asymmetric structure.
As shown in Fig. 8(d) and 8(e), we find that the distribution of magnetic pressure is asymmetric during the formation of the pillar, and the magnetic pressure is larger on the right side of the pillar due to the shock compressed.Then this would lead to an asymmetry of the Lorentz force, as shown in the white arrows in 8(d).The right side of the highdensity clump is subjected to greater downward resistance, so compared to the fluid on the left side of the pillar, the movement of the fluid on the right side will be more hindered, then the pillar eventually develops into an asymmetric structure.And at 120 ns, We can clearly see that the asymmetrical structure of the magnetic pressure still exists, and this asymmetrical force leads to the generation of an asymmetrical pillar structure in case 11.

DISCUSSION AND CONCLUSION
The model of a quasi-homogeneous magnetic field is proposed based on the assumption that there is a large-scale magnetic field in the cloud, and the large-scale magnetic field is also considered to play an important role in the process of star formation (Blandford & Payne 1982).In this work, we consider the effect of magnetic fields on the formation of laser-driven scaled pillars of creation.A low-density foam with preexisting triple-density spherical dense clumps (maximum 1 g/cm 3 ) are used in our simulations, and such target configuration makes manufacturing in experiments more convenient.We also use our radiation simulation to reproduce the POC experiments in the OMEGA EP facility, and the simulation results are very similar to their experimental results, which show the radiation hydrodynamic simulation is reliable in predicting experimental results.Regarding the generation of the magnetic field in the experiment, we plan to use the pulse-power system (Hu et al. 2022) for weakly magnetized cases.For medium and strong magnetized cases, we plan to generate magnetic fields using the laser-driven capacitor coil target method.And in recent years, a large number of magnetic field generation experiments (Zhang et al. 2018) have been conducted on the SG-II laser facility.
The astrophysical POCs occur on spatial and temporal scales that are typically 10-20 orders of magnitude greater than those of laboratory experiments that are designed to simulate them.As a result, ensuring similarity between the astrophysical phenomenon and its laboratory counterpart becomes a critical issue.According to previous works (Ryutov et al. 1999(Ryutov et al. , 2001;;Ryutov & Remington 2002), if the key dimensionless numbers are kept much larger than unity in two systems, such as magnetic Reynolds number, Reynolds number, and Péclet number, both systems will behave as ideal compressible hydrodynamic fluids.And the control equations will remain invariant in the two systems when the following transformation conditions are satisfied, r lab = ar ast , ρ lab = bρ ast , P lab = cP ast , t lab = a b/ct ast , v lab = c/bv ast , B lab = √ cB ast , where a, b, and c are arbitrary positive numbers, and they also called free transformation parameters.The subscripts "ast" and "lab" represent the astrophysical pillars and those of the laboratory.
For the astrophysical POC, the estimated magnetic diffusion time is about 10 13 years (Kane et al. 2015), so the magnetic Reynolds number is much larger than unity (note that the age of pillars is about 10 4 − 10 5 years).Viscosity and thermal conductivity are negligible, because of the very large spatial scales involved, so that the Reynolds number and Péclet number are also much larger than unity.In the experiments, these key dimensionless parameters can be calculated as follows.The Reynolds number is R e = vL/ν, where v, L, ν are the velocity, system size, and viscosity respectively, and the viscosity can be calculated as ν e n i Z(Z + 1)e 4 ln Λ] is the thermal diffusivity (Brandenburg & Subramanian 2005).As shown in table 2, the Magnetic Reynolds number, Reynolds number, and Peclet number all satisfy the condition of being far greater than the unity, so the scaling laws are valid, and as mentioned above, the radiation cooling parameter is also closer to the astrophysical system.Then, the same set of equations can accurately describe both the astrophysical and laboratory scenarios.Table 2 gives the detailed comparison between the scaled-up laboratory system and specifically the POC in Eagle Nebula, where the transformation parameters a = 3.3 × 10 −20 , b = 3 × 10 18 , c = 2.0 × 10 18 are taken.According to astrophysical observation results (Pattle et al. 2018), the typical magnetic field strength within the pillar is about 100 − 300 µGs, and the directions of the magnetic field are parallel to the cylindrical structure.The range of magnetic field strength corresponding to the laboratory system is approximately 100-300 T after conversion.In this work, this is very similar to the case 5 situation, and our laboratory radiation simulation also indirectly indicates that the pillar structure is mostly collimated and stable when there is a parallel magnetic field component, which corresponds to astronomical observations.The scaling laws confirm that our laboratory simulation results can be applied to explore the formation mechanism of POC.
In summary, we propose a new experimental scheme to study the formation and support mechanism of POC.Three materials (CH, Al, Au) are used in our simulations, and the Al-doped CH foam should be used as the first choice for experimental targets according to our simulation results, whose radiative ratio is closer to the POC and a strong shock can be launched by the ablation process.And the 2D and 3D RMHD simulations all show that the magnetic fields, especially the B y component, play an important role in the formation of POC.For the medium magnetized cases, the magnetic pressure can effectively support the sides of pillars against collapsing, and a long collimated pillar structure is formed.The radiation maps are also similar to the observed results.The bending at the head of 'Column I' in POC may be due to the non-parallel magnetic field lines according to our 3D RMHD simulations.The similarities between laboratory and astrophysical systems suggest that our results can be applied to explore the formation mechanism of pillars of creation.In addition, the mechanism of magnetic topology changes within the pillars (Ryutov et al. 2005) and the multiple pre-existing disturbances (Lim & Mellema 2003;Pound et al. 2007) are also the key factors in the formation of POC, and we will plan to carry out further related research in the future.

Figure 2 .
Figure 2. Radiation cooling for different materials.Panel (a) to (e), these color maps correspond to the logarithmic mass density in g/cm 3 for CH, Al, and Au at different times.The electron number densities at x = 0 (the position at white dot line in (a) to (c)) and t = 5 ns for different materials are shown in panel (f).Panel (g) to (i) are the corresponding radiation temperature distribution, and panels (j), (k) are the electron temperature distribution.Panel (l) plots the electron temperature at x = 0 (the position at white dot line in (j) and (k)) and t = 10 ns of different materials.

Figure 3 .
Figure 3.The formation of pillar-like structure without magnetic fields.Panel (a) and (b) are the schematic views of the interaction between shock and clump, (a) during initial shock wave transit, and (b) after initial shock wave transit.Panel (c) to (f), these color maps correspond to the logarithmic mass density in g/cm 3 for Al at 40 ns, 60 ns, 100 ns, and 120 ns respectively, and an apparent rarefaction wave travels back to the front surface after the shock has traversed the high-density clumps in panel (f).The distribution of velocity in the y-axis is shown in (g).The mass densities at x = 0 for different times are shown in panel (h).

Figure 4 .
Figure 4. Weakly magnetized cases.Panels (a) to (c) correspond to the logarithmic mass density in g/cm 3 for case 1, case 2, and case 3 at 5 ns respectively, and the color bar is given at the bottom of these panels.The white dot line shows the front position of the blow-off gas plasma for the case without magnetic fields (as shown in 2(b)).Panels (d) to (f) are the corresponding magnetic pressure for these cases at 5 ns in logarithmic values, and the black arrows represent the magnetic field lines.The distributions of mass density at 120 ns for these cases are shown in panel (g) to (i), and the corresponding magnetic pressures are given in (j) to (l) in logarithmic values.

Figure 5 .
Figure5.Medium magnetized cases.The left sides of panels (a) and (b) correspond to the logarithmic mass density in g/cm 3 for case 4, case 5 at 5 ns, and the right sides correspond to the distribution of magnetic fields.Panels (c) to (d) plot the logarithmic mass density for cases 4 to 6 at 120 ns.The mass densities of low-density foam at x = −1.5 mm for different cases are shown in panel (f).The left sides of panels (g) and (h) plot the magnetic pressure for case 4, and case 5 at 120 ns, and the right sides are the corresponding thermal pressure.The left sides of panels (i) to (l) are the radiation temperature in eV for case 4 to case 6 and case without magnetic fields, and the right sides of these panels are the corresponding electron temperature.

Figure 6 .
Figure 6.Radiation for medium magnetized cases.The distributions of radiation power for different cases at 120 ns are shown in panels (a) to (d).

Figure 7 .
Figure 7. Strong magnetized cases.Panels (a) to (c) correspond to the logarithmic mass density in g/cm 3 for case 7 to case 9 at 5 ns.The left sides of panels (d) and (e) are the distribution of mass density at 120 ns, and the right side is the corresponding magnetic fields.

Figure 8 .
Figure 8. 3D simulation results.Panels (a) and (b) correspond to the distribution of mass density in g/cm 3 and radiation power for case 10 at 120 ns.Panel (c) is the distribution of mass density for case 11, and the color bar for mass density is given under panel (a).The arrows in panels (a) and (c) represent the magnetic field vector and the colors on the arrows represent the strength (the color bar is given under panel (c)).Panel (d) and (e) plot the distribution of magnetic pressure for case 11 at 60 ns and 120 ns respectively in the XY plane, and the white arrows represent the Lorentz force in the unit of dyn/cm 3 .

2 e
2018), where k B , T, n i , m i , Z, ln Λ are Boltzmann constant, temperature, ion number density, ion mass, average ionization and Coulomb logarithm; the magnetic Reynolds number of R m = vL/η, where η = c 2 m 1/Ze 2 lnΛ/[4(k B T ) 3/2 ] is the magnetic diffusivity (Brandenburg & Subramanian 2005); and the Peclet number P e = vL/κ, where κ = (k B T ) 5/2 /[m 1/2 This work is supported by the National Key R&D Program of China, Grants No. 2022YFA1603200 and No. 2022YFA1603204; Science Challenge Project, No. TZ2018005; National Natural Science Foundation of China, Grant No. 11825502 and 11921006; the Strategic Priority Research Program of Chinese Academy of Sciences Grant No. XDA25050900.L. F. W. is supported by the National Natural Science Foundation of China (Grant No. 11975053).The simulations are carried out on the supercomputers in China.