Hi LAMMPS developers and users,
In my work, there are two groups of atoms. Atoms of group 1 constitute the solid boundaries and atoms of group 2 are the SPH fluid particles. I want to compute the force exerted on each solid particle by all the SPH fluid particles (of course within the cutoff distance). I have tried
compute 2 all pair/local force
, but was given this error:
Pair style does not support compute pair/local
I checked my pair styles. They are as follows:
pair_style hybrid/overlay sph/rhosum 1 sph/taitwater soft ${cut} zero ${cut}
pair_coeff 1 1 zero
pair_coeff 1 2 soft 3.e-6
pair_coeff * 2 sph/taitwater ${rho} ${c} ${alpha} ${h}
pair_coeff 2 2 sph/rhosum ${h}
And then I found the style that did not support pair/local is sph/taitwater. But what I want to compute is merely the per-atom force from the fluid on the solid, i.e., the force resulting from the pair_style soft in my work. I don’t need the force between the fluid particles. That is to say, when I compute pair/local, I don’t need to incorporate the troublesome sph/taitwater, which will also save much computation time. Is there any way to exclude sph/taitwater when computing pair/local? Or, may there be some better ways to calculate the force from the fluid on each solid particle?
By the way, the error is the same for the command property/local, which does not support sph/taitwater, too.
Another idea: use compute group/group. But the limitation is its output is the total force, not the per-atom force.
Any comments or advice will be appreciated!
My LAMMPS version: Aug 2023
System: Windows 10
The restart file have been uploaded.
relax_kB0.5kR0.01vis0.05.restart (2.0 MB)
Minimum working example:
dimension 3
atom_style hybrid bond sph
units si
newton on
boundary p p p
read_restart relax_kB0.5kR0.01vis0.05.restart
reset_timestep 0
variable R equal 0.025
variable L equal 0.6283185307180
variable dL equal 0.006283185307180
variable h equal 2.*${dL}
variable cut equal 1.*${dL}
variable skin equal 0.3*${h}
variable c equal 0.05
variable rho equal 1000
variable alpha equal 1.
variable nu equal ${alpha}*${h}*${c}/10.
variable kB equal 0.5
variable f0 equal 5.e-4
variable Nsteps equal 50
variable Dt equal 0.002 # timestep
variable Nout equal ${Nsteps} # output interval
variable N equal 10.
variable Dz equal ${L}/${N} # length of section
group wall type 1
group fluid type 2
variable index equal ramp(8,9)
variable z0 equal (v_index)*${Dz}%${L}
variable z1 equal (v_index+1.)*${Dz}%${L}
variable z2 equal (v_index+2.)*${Dz}%${L}
pair_style hybrid/overlay sph/rhosum 1 sph/taitwater soft ${cut} zero ${cut}
pair_coeff 1 1 zero
pair_coeff 1 2 soft 3.e-6
pair_coeff * 2 sph/taitwater ${rho} ${c} ${alpha} ${h}
pair_coeff 2 2 sph/rhosum ${h}
bond_style harmonic
bond_coeff 1 ${kB} ${dL}
neighbor ${skin} bin
variable PP atom ((v_z0<v_z1)|^((z>=v_z0)|^(z<=v_z1)))*gmask(wall)
#section that is squeezed (push) z0<z<z1
variable SS atom ((v_z1<v_z2)|^((z>v_z1)|^(z<=v_z2)))*gmask(wall)
#section that senses (sense) z1<z<z2
variable RR atom ((v_z0>v_z2)|^((z>v_z0)|^(z<v_z2)))*gmask(wall)
#the remaining section (rigid)
group push dynamic wall var PP every 1
group sense dynamic wall var SS every 1
variable rDist atom sqrt(x^2+y^2)
variable cosQ atom x/v_rDist
variable sinQ atom y/v_rDist
variable squeeze_factor atom (v_rDist/${R})^2.
variable FP atom -${f0}*gmask(push)
variable Fx atom v_FP*v_cosQ*v_squeeze_factor
variable Fy atom v_FP*v_sinQ*v_squeeze_factor
###the following only because spring/self cannot be dynamic####
fix x0 wall store/state 0 x # store the initial value of x
fix y0 wall store/state 0 y
fix z0 wall store/state 0 z
fix zu wall store/state 0 zu
variable varX atom (f_x0-x)*gmask(wall) # change of x
variable varY atom (f_y0-y)*gmask(wall) # change of y
variable KR equal 0. # rigid section
variable FRx atom v_RR*${KR}*v_varX
variable FRy atom v_RR*${KR}*v_varY
fix G fluid gravity -0.005 vector 0 1 0
fix 10 all sph
fix 20 wall viscous 0.01
fix 30 wall addforce v_Fx v_Fy 0.
fix 40 wall addforce v_FRx v_FRy 0.
fix 50 wall setforce NULL NULL 0.
# compute 1 all property/local patom1 patom2
compute 2 all pair/local force
dump 1 all local 5 pair_local.dump c_2
thermo 1
timestep ${Dt}
run 10