Minimize with force constraint

Hi,

I try to place C atoms on the octahedral position of an Fe single crystal. Therefore after creating the C atom, I want to run a minimization only for the C atoms, to get the position perfectly right.

The code I use:

units metal
atom_style atomic
boundary p p p
lattice bcc 2.85531 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
read_data md.data
mass 1 55.845000
mass 2 12.0107
pair_style eam/fs
pair_coeff * * Fe-C_Hepburn_Ackland.eam.fs Fe C
thermo 1
thermo_style custom step temp ke pe etotal press pxx pyy pzz pxy pxz pyz ly lx lz vol
thermo_modify lost warn norm yes
group iron type 1
fix 1 iron setforce 0.0 NULL 0.0
minimize 0.0 1.0e-10 100000 1000000
unfix 1

Unfortunatley, what happens is that the Fe atoms change their position as well as the C atoms.

Before: [x,y,z, distance initial C position]

79.94868 69.955095 62.81682 0 <- C atom
79.94868 68.52744 62.81682 1.427655 <- Fe atom
79.94868 71.38275 62.81682 1.427655 <- Fe atom
78.521025 69.955095 64.244475 2.01900906338976 <- Fe atom
81.376335 69.955095 64.244475 2.01900906338976 <- Fe atom
78.521025 69.955095 61.389165 2.01900906338976 <- Fe atom
81.376335 69.955095 61.389165 2.01900906338976 <- Fe atom

After: [x,y,z, distance initial C position]

79.94868 69.955095 62.81682 2.1316282072803e-14 <- C atom, movement as expected
79.94868 68.1896640213388 62.81682 1.76543097866119 <- Fe atom moves even though the force was set to zero
79.94868 71.7205259786612 62.81682 1.76543097866119 <- Fe atom moves even though the force was set to zero
78.521025 69.955095 64.244475 2.01900906338976 <- Fe atom
81.376335 69.955095 64.244475 2.01900906338976 <- Fe atom
78.521025 69.955095 61.389165 2.01900906338976 <- Fe atom
81.376335 69.955095 61.389165 2.01900906338976 <- Fe atom

Even worse, this effect appeared only for about 30% of the octahedral position in the sample. So I am kind of confused, I thought that the combination of minimize plus force constraint, would result in limiting the degress of freedom for the minimization, meaining in this case only the carbon atoms are allowed to move.

Tested with the stable version from 30 Oct. 2014.

  • Jan

Hi,

I try to place C atoms on the octahedral position of an Fe single crystal. Therefore after creating the C atom, I want to run a minimization only for the C atoms, to get the position perfectly right.

The code I use:
units metal
atom_style atomic
boundary p p p
lattice bcc 2.85531 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
read_data md.data
mass 1 55.845000
mass 2 12.0107
pair_style eam/fs
pair_coeff * * Fe-C_Hepburn_Ackland.eam.fs Fe C
thermo 1
thermo_style custom step temp ke pe etotal press pxx pyy pzz pxy pxz pyz ly lx lz vol
thermo_modify lost warn norm yes
group iron type 1
fix 1 iron setforce 0.0 NULL 0.0
minimize 0.0 1.0e-10 100000 1000000
unfix 1

Unfortunatley, what happens is that the Fe atoms change their position as well as the C atoms.

Before: [x,y,z, distance initial C position]
79.94868 69.955095 62.81682 0 <- C atom
79.94868 68.52744 62.81682 1.427655 <- Fe atom
79.94868 71.38275 62.81682 1.427655 <- Fe atom
78.521025 69.955095 64.244475 2.01900906338976 <- Fe atom
81.376335 69.955095 64.244475 2.01900906338976 <- Fe atom
78.521025 69.955095 61.389165 2.01900906338976 <- Fe atom
81.376335 69.955095 61.389165 2.01900906338976 <- Fe atom

After: [x,y,z, distance initial C position]
79.94868 69.955095 62.81682 2.1316282072803e-14 <- C atom, movement as expected
79.94868 68.1896640213388 62.81682 1.76543097866119 <- Fe atom moves even though the force was set to zero
79.94868 71.7205259786612 62.81682 1.76543097866119 <- Fe atom moves even though the force was set to zero
78.521025 69.955095 64.244475 2.01900906338976 <- Fe atom
81.376335 69.955095 64.244475 2.01900906338976 <- Fe atom
78.521025 69.955095 61.389165 2.01900906338976 <- Fe atom
81.376335 69.955095 61.389165 2.01900906338976 <- Fe atom

Even worse, this effect appeared only for about 30% of the octahedral position in the sample. So I am kind of confused, I thought that the combination of minimize plus force constraint, would result in limiting the degress of freedom for the minimization, meaining in this case only the carbon atoms are allowed to move.

Fe atoms are free to move in y-direction in your input and that is what they do. There is a difference between 0.0 and NULL when using fix setforce.

Axel