Fix_ttm for laser ablative

LAMMPS LAMMPS Development
,
1

Dear All,

I have modified the fix_ttm_mod.cpp for laser ablative simulation. The following codes were modified.

image
image639×325 16.1 KB

After simulation, I found the atoms at the top center region (blue) move along negative x direction, instead of explosion (right figure).

image
image828×306 55.3 KB

Does anyone know the reasons? Thanks. The simulation files (input script and ttm_mod parameters) can be seen below

in.ttm_al (1.8 KB)

Al.ttm_mod (903 Bytes)

2

Note: This response was drafted with the assistance of an AI tool (Claude) after this thread sat without any replies for over five months. While every effort has been made to reason carefully from the LAMMPS source code and the attached files, AI-assisted answers carry a real risk of error — please treat the points below as a starting point for debugging rather than authoritative guidance, and do verify against the LAMMPS documentation and the original Pisarev & Starikov papers.

Looking at the attached in.ttm_al, Al.ttm_mod, and the modified fix_ttm_mod.cpp, there are a few issues that together explain the unexpected −x motion.

1. The Gaussian modification includes ix, which reverses the x-direction force

The intended change is a Gaussian beam profile in the transverse (y-z) plane. However, the modification also includes ix in the Gaussian:

*exp((-1.0)*(pow(ix-nxgrid/2,2.0) + pow(iy-nygrid/2,2.0) + pow(iz-nzgrid/2,2.0))/dia)

The Gaussian is centred at (nxgrid/2, nygrid/2, nzgrid/2) = (25, 20, 20). With dia=12 (σ ≈ 2.4 grid cells) and surface_l=2, the surface is 23 cells away from the Gaussian peak, giving a factor of exp(−23²/12) = exp(−44) ≈ 0 at ix=2. Almost all laser energy ends up deposited at ix=25 (the right boundary of the coupling zone, rsurface), not at the surface. This reverses the electron temperature gradient in x, so every atom in the coupling zone feels an electron pressure force in the −x direction — which is exactly what you observe.

The fix is to remove ix from the Gaussian. The depth profile is already handled by the Beer-Lambert term; the Gaussian belongs only in y and z:

double dia = 12.0;
if (ix >= t_surface_l)
  mult_factor = (intensity/(dx*skin_layer_d))
                * exp((-1.0)*(pow(iy-nygrid/2.0,2.0) + pow(iz-nzgrid/2.0,2.0))/dia)
                * exp((-1.0)*(ix_d - surface_d)/skin_layer_d);

With this corrected form the Beer-Lambert term drives a +x compressive wave into the material (standard TTM ablation shock), while the y-z Gaussian creates an outward radial electron pressure gradient — pushing atoms away from the beam axis, which is the lateral “explosion” effect you are after.

2. skin_layer is in grid-cell units, not Angstrom

fix ttm/mod reads skin_layer as an integer and uses it directly as a number of grid cells in the Beer-Lambert exponent:

exp((-1.0)*(ix_d - surface_d)/skin_layer_d)   // skin_layer in grid cells

With nxgrid=50 and a coupling zone of 23 cells (lsurface=2 to rsurface=25), setting skin_layer=50 means the intensity varies only from 1.0 to exp(−23/50)≈0.63 across the entire slab — essentially uniform volumetric heating rather than surface-concentrated absorption. For realistic laser heating of Al the optical skin depth is roughly 7 nm; divide by your grid-cell size in Angstrom to get the correct integer value (typically 2–5 cells for a ~100–200 Angstrom slab).

3. Working around the periodic-boundary constraint

fix ttm/mod requires fully periodic boundaries. True ablation (atoms escaping the box) is not directly possible, but it can be approximated by:

  1. Adding a large vacuum region above the free surface (right side of the slab, beyond rsurface).
  2. Using fix evaporate with nevery 1 and N larger than the maximum number of atoms that could reach the vacuum in a single step, to remove expelled atoms every timestep before they wrap around:
fix evap vapor evaporate 10000 1 vacuum_region 12345

The vacuum region must be large enough that no atom — even the fastest — can traverse it within one nevery interval and re-enter the slab from the far side.

4. Consider movsur 1 for a proper ablation run

With movsur=0 the electron–ion coupling region is fixed at cells lsurface–rsurface throughout the run. Once the surface starts receding (atoms removed by fix evaporate), vacated cells near rsurface remain coupled and apply Langevin forces to empty space. Setting movsur=1 in Al.ttm_mod lets the right surface boundary track the actual material surface, which is more physical once ablation begins.