Questions about fix aveforce in LAMMPS

Dear all,
I am a beginner and now I have met some problem in my work.
My topic is about diffusion between Cu and Al. A simple model is shown in the figure below.

Now I want to study the effect of pressure on diffusion. I try to use “fix aveforce” command to apply pressure on the move layer. But the movement of move layer is not obvious. It only moves downward for a short displacement, and then does not continue. Part of the code is as follows:
variable p equal 1000 #the set pressure
variable s equal 40.440.4 #the surface area of move layer
variable N equal 1323 # the number of atoms in the move layer
variable ev equal 624150913 #newton to ev/angstrom
variable maxforce equal ($p
s*1.0e-14*{ev})/$N #required force
variable samplerun equal 30000

#The external force changes with time and finally reaches the required force
variable extforce equal {maxforce}*step/{samplerun}
fix externalforce move aveforce 0.0 0.0 -${extforce}

fix boundary_f fix setforce 0 0 0 # fix the boundary
fix 1 move nvt temp 300 300 0.1# I also tried “nve” , the result is not good either

the whole system relaxes in nve ensemble. May I ask where is the error and how should I correct it.
Best wishes
Happy New Year!

It is very difficult to assess all details of your simulation with seeing only part of your input.
That input is also garbled due to the typesetting of the forum software, you need to surround your input with triple backquotes (```) to stop the forum software from doing that.

The probable misconception here is how stress/pressure increases when a material is compressed.
Your input excerpt is showing a linear increase in force on the “move” layer during the simulation, and - to me - that seems to have the expected effect. Initially, it will be easy to compress, but then quickly will require a lot of force to have a larger amount of compression. A simple way to assess how large that would be, is to do a trial run where you do a linear compression with fix move (instead of fix aveforce + nvt) and then monitor the resulting pressure in the compression direction.

Mind you that the effects you are looking for are rather indirect and thus you may be required to increase the force very slowly (or in steps instead of continually) and then determine the desired properties from averaging. While it is certainly possible to set up simulations like macroscopic experiments, the evaluation and relation of microscopic data to macroscopic observables is not always as straightforward. You may need to spend some more time with a text book on statistical thermodynamics for some additional insight.

Dear Dr. Kholmeyer,
thank you very much for your reply. I’m sorry I didn’t explain my work clearly.
I want to apply a constant force on the “move layer” and find out the effect of pressure on diffusion. Let me simplify my code and show it.
“fix extforce equal -7.7e-4” # the constant force
“fix externalforce move aveforce 0.0 0.0 ${extforce}”
“fix boundary_f fix setforce 0 0 0” # fix the boundary
“fix 1 thermoup nvt temp 300 300 0.1”
“fix 2 cu nvt temp 300 300 0.1”
“fix 3 al nvt temp 300 300 0.1”
“fix 4 thermolow nvt temp 300 300 0.1”
“fix 5 move nvt temp 300 300 0.1”
The above is my code for applying a constant external force. Is my “fix command” for “move layer” correct? After applying the constant force, the structure of “move layer” can not be maintained as “fixed layer”. And the atoms of "thermoup layer”(the green part in the picture) will diffuse into the “move layer”. If I fix “move layer”, it won’t move. How can I solve this problem and get the correct result? Thank you!
Best wishes,
Happy New Year

Sorry, but this doesn’t really help much.
Please provide a complete input, explain what you expect to see happening and what happens instead.

Have you studied diffusion before? How do you plan to quantify it? What is your base value for a system with no external pressure?

Also, you appear to be mixing technical problems and scientific questions. Let us sort out the technical issues first and one at a time.

Please also note that I won’t have the time to dive into the science of your project. That is where you need to talk to your adviser.

Dear Dr. Kholmeyer
I’m sorry it took so long to reply. Below is my input.
units metal
boundary p p p
dimension 3
atom_style atomic
timestep 0.001
neighbor 2.0 bin
neigh_modify delay 1 check yes

region fix block 0 40.4 0 40.4 0 10 units box
region thermo1 block 0 40.4 0 40.4 10 20 units box
region al block 0 40.4 0 40.4 20 60.4 units box
region cu block 0 40.4 0 40.4 60.4 100.8 units box
region thermo2 block 0 40.4 0 40.4 100.8 110.8 units box
region move block 0 40.4 0 40.4 110.8 120.8 units box
region box block 0 40.4 0 40.4 0 120.8 units box
create_box 6 box
lattice fcc 4.05
create_atoms 1 region fix
create_atoms 2 region thermo1
create_atoms 3 region al
lattice fcc 3.61
create_atoms 4 region cu
create_atoms 5 region thermo2
create_atoms 6 region move
mass 1 27
mass 2 27
mass 3 27
mass 4 64
mass 5 64
mass 6 64

group fix region fix
group thermo1 region thermo1
group al region al
group cu region cu
group thermo2 region thermo2
group move region move

pair_style eam/alloy
pair_coeff * * AlCu.eam.alloy Al Al Al Cu Cu Cu
run_style verlet
minimize 1.0e-6 1.0e-7 1000 1000
thermo 1000
dump process1 all xyz 1000
compute thermo1temp thermo1 temp
compute cutemp cu temp
compute altemp al temp
compute thermo2temp thermo2 temp
compute s1 al stress/atom NULL
compute s11 al reduce sum c_s1[3]
variable press1 equal -(c_s11)/vol
compute s2 cu stress/atom NULL
compute s22 cu reduce sum c_s2[3]
variable press2 equal -(c_s22)/vol
compute s3 move stress/atom NULL
compute s33 move reduce sum c_s3[3]
variable press3 equal -(c_s33)/vol

velocity all create 300 4928459 dist gaussian
fix 1 thermo2 nvt temp 300 300 0.1
fix 2 cu nvt temp 300 300 0.1
fix 3 al nvt temp 300 300 0.1
fix 4 thermo1 nvt temp 300 300 0.1
thermo_style custom step temp v_press1 v_press2 v_press3
run 10000
reset_timestep 0

compute msd1 al msd
variable twopoint1 equal c_msd1[4]/6/(stepdt+1.0e-6)
fix msd1 al vector 1000 c_msd1[4]
variable fitslope1 equal slope(f_msd1)/6/(1000
compute msd2 cu msd
variable twopoint2 equal c_msd2[4]/6/(stepdt+1.0e-6)
fix msd2 cu vector 1000 c_msd2[4]
variable fitslope2 equal slope(f_msd2)/6/(1000

variable p equal 100 #Mpa
variable s equal 40.4*40.4e-20
variable N equal count(move)
variable ev equal 624150913 #ev/Angstrom to newton
variable extforce equal 10*(v_p*v_s*v_ev)/v_N
fix externalforce move addforce 0.0 0.0 $(-v_extforce) region move
fix 5 move nvt temp 300 300 0.1
thermo_style custom step temp press etotal c_thermo2temp c_thermo1temp c_cutemp c_altemp c_msd1[4] c_msd2[4] v_twopoint1 v_twopoint2 v_fitslope1 v_fitslope2 v_press1 v_press2 v_press3
run 100000

The result I want is that after applying an external force the move layer moves down and achieve a compression effect, promoting the diffusion of Cu and Al.
Best wishes,

Sorry, but this is not helping much either. In part because it is not escaped with triple backquotes “```” and thus subject to typesetting with the markup language of the editor here, in part because it is very convoluted and contains lots of unrelated commands (some of which are redundant, too, and thus needlessly verbose since the same could be implemented with fewer commands), but mostly because you are not addressing the rather fundamental concerns I raised in my first response. It looks to me that you are misunderstanding the relation between force and resulting compression for the system at hand.

Thank you for your reply. I have re-edited my input. What is the better command to achieve the compaction process? Thank you!
Best wishes

I don’t have any advice beyond what I have already provided.