I was trying to model liquid-liquid interface between two immiscible liquids. In order to model it
I decided on using the repulsive part of the LJ potential by changing the cutoff radius to ~2^(1/6) for interactions
between the two liquids. For all other interactions with wall and between the same liquid molecules I have
used a cutoff radius of 2.2.
I am not sure what I am doing wrong but the interface seems to breakdown on visualizing it and the two types of liquid
end up mixing.
Would appreciate any help/suggestions. Have attached the code below.
Thank you.
BR/Joseph
# 3-d LJ flow simulation
units lj
dimension 3
boundary p s p
atom_style atomic
neighbor 0.3 bin
neigh_modify delay 1
# create geometry
region simulation_box block 0.0 24.136 0.0 54.82 0.0 7.22 units box
create_box 4 simulation_box
region fluid_1 block INF 11.0 2.39 52.5 INF INF units box
region fluid_2 block 12.0 INF 2.39 52.5 INF INF units box
region btm_wall block INF INF INF 0.95 INF INF units box
region top_wall block INF INF 53.86 INF INF INF units box
lattice fcc 0.81 orient x 1 1 -2 orient y 1 1 1 orient z 1 -1 0
create_atoms 1 region fluid_1 units box
lattice fcc 0.81 orient x 1 1 -2 orient y 1 1 1 orient z 1 -1 0
create_atoms 2 region fluid_2 units box
lattice fcc 3.24 orient x -1 1 0 orient y 1 1 1 orient z 1 1 -2
create_atoms 3 region btm_wall units box
create_atoms 4 region top_wall units box
mass 1 1.0
mass 2 1.0
mass 3 1.0
mass 4 1.0
# LJ potentials
pair_style lj/cut 2.2
pair_coeff 1 1 1.0 1.0
pair_coeff 2 2 1.0 1.0
pair_coeff 3 3 1.0 1.0
pair_coeff 4 4 1.0 1.0
pair_coeff 1 3*4 0.6 0.75
pair_coeff 2 3*4 0.6 0.75
pair_coeff 1 2 0.6 0.75 1.16
pair_modify shift yes
# define groups
group lower region btm_wall
group upper region top_wall
group boundary union lower upper
group flow subtract all boundary
set group flow type 1
set group lower type 2
set group upper type 3
# initial velocities
compute mobile flow temp
velocity flow create 1.10 482748 temp mobile units box
thermo_modify temp mobile
# thermostating
compute mobile_thermo all temp/partial 0 0 1
fix 3 all langevin 1.10 1.10 1.0 699483
fix 1 flow nve
fix 2 lower nve
fix_modify 3 temp mobile_thermo
thermo_modify temp mobile_thermo
# Couette flow
velocity lower set 0.0 0.0 0.0 units box
velocity upper set 0.0 0.0 0.0 units box
#fix 8 upper move wiggle 10.0 0.0 0.0 40.0 units box #<<<<< CHANGE AMPLITUDE AND TIMEPERIOD
fix 4 boundary setforce 0.0 0.0 0.0
# Run
timestep 0.002
thermo 1000
dump 10 all atom 100 dump.flow
dump 12 all xyz 100 dump_xyz.flow
run 6000