subject: compact temp/com

Hi,
I try to simulate a collision between a nano_particle and a fixed substrate. When i applied the “compute temp/com” to calculate the temperature of the" upper particle", this temperature was higher than the actual temperature. According to the command, it subtract the
center of mass velocity. Therefore, i think that the " upper particle" temperature should be around the initial temperature before having a collision between them.
How should i calculate the " upper particle" temperature which subtract its “center of mass velocity”

Partial input script is following:

group upper id >= 3970 #upper particle
region 1 block INF INF INF INF INF -5 units box
region 2 block INF INF INF INF -5 0 units box
region 3 block INF INF INF INF 0 15 units box

group bulk_down region 1
group bulk_middle region 2
group bulk_upper region 3

compute T_middle bulk_middle temp
compute T_upper upper temp/com

velocity upper create 300.0 123456 temp T_upper

velocity bulk_middle create 300 123456 temp T_middle

fix 1 all nve
fix 2 bulk_down setforce 0.0 0.0 0.0
fix 3 bulk_middle temp/rescale 10 293.0 293.0 0.001 1.0
fix_modify 3 temp T_middle

velocity bulk_upper set 0.0 0.0 0.0 units box
velocity upper set 0.0 0.0 -0.02 units box

thermo_style custom step temp press pe ke etotal
thermo_modify lost warn temp T_upper

Hi,
       I try to simulate a collision between a nano_particle and a fixed
substrate. When i applied the "compute temp/com" to calculate the
temperature of the" upper particle", this temperature was higher than the
actual temperature.

higher than the "actual" temperature of what? the entire system?

According to the command, it subtract the
center of mass velocity. Therefore, i think that the " upper particle"
temperature should be around the initial temperature before having a
collision between them.

again, imprecise formulation makes it impossible to give any advice.

      How should i calculate the " upper particle" temperature which
subtract its "center of mass velocity"

Partial input script is following:

partial inputs are almost always useless, since typically the most
important information is in the part that is *not* quoted.

axel.

Hi Dr. Axel Kohlmeyer,
      Thanks your replies.
For useful help, follow is my complete input script, and simulation results. My queation is Why the temperature of bulk_upper rapidly increased?
Like the last letter said, the "compute temp" command cannot work. If you have any good idea, pls tell me.

Hi Dr. Axel Kohlmeyer,
      Thanks your replies.
For useful help, follow is my complete input script, and simulation results. My queation is Why the temperature of bulk_upper rapidly increased?
Like the last letter said, the "compute temp" command cannot work. If you have any good idea, pls tell me.

this is a bit of a complex input and tries to do too many things in
one go, which makes it difficult to debug.
you should approach your simulation in stages. first your should do a
proper equilibration simulation with your nanoparticle at rest. for
this you should remove any input from your input that is not
essential. best start with an empty file and only copy over what
absolutely necessary. if you cannot get the system equilibrated, you
need to figure out why. that should also solve the "temperature
mystery", or make it easier to find a solution. please find a few more
comments in between the lines:

#########
units real
dimension 3
newton on
boundary p p p
atom_style hybrid charge atomic

this is redundant and pointless. atom_style charge is a superset of
atomic and thus hybrid is not needed.

communicate multi vel yes

this is not needed for your system as well.

neighbor 0.3 nsq

such a small neighbor skin only works well in solids and with reduced
(=lj) units.

neigh_modify delay 0 every 1 check yes page 100000 one 10000

neither setting page nor setting one should be needed here.

timestep 0.5

# --------------------------

read_data data.TIO_bulk
  orthogonal box = (-40 -40 -400) to (40 40 600)
  1 by 1 by 1 MPI processor grid
  4335 atoms

these are a rather small number of atoms for the dimensions of the box.

#---------------Matsui-Akaogi potential----------------------
## Ti type 1, O type 2

pair_style buck/long/coul/long cut long 10.0

pair_coeff 1 1 7.177e+5 0.154 1.210e+2
pair_coeff 1 2 3.911e+5 0.194 2.904e+2
pair_coeff 2 2 2.717e+5 0.234 6.969e+2

pair_modify table 0

why turn off tabulation? you don't use a very accurate kspace solver,
so this will only make you slower, but not more accurate.

kspace_style ewald 1.0e-4

delete_atoms overlap 0.3 all all

why this? you provide initial input from a data file. if you need to
delete atoms, why add them in the first place?

Ewald initialization ...
  G vector (1/distance) = 0.244615
  estimated absolute RMS force accuracy = 0.0434719
  estimated relative force accuracy = 0.000130914
  KSpace vectors: actual max1d max3d = 112466 148 13099036
Deleted 0 atoms, new total = 4335

#-------------------------------define group-------------------
group upper id >= 3970 #upper particle
366 atoms in group upper

region 1 block INF INF INF INF INF -5 units box
region 2 block INF INF INF INF -5 0 units box
region 3 block INF INF INF INF 0 19 units box

group bulk_down region 1
467 atoms in group bulk_down
group bulk_middle region 2
759 atoms in group bulk_middle
group bulk_upper region 3
2743 atoms in group bulk_upper

###compute
# computer energy
compute new all temp
compute 1 all pe/atom
compute 2 upper ke
compute 4 all stress/atom pair
compute T_middle bulk_middle temp
compute bulk_upper bulk_upper temp/com
compute T_upper upper temp/com

variable velocity_x equal vcm(upper,x)
variable velocity_y equal vcm(upper,y)
variable velocity_z equal vcm(upper,z)

variable position_x equal xcm(upper,x)
variable position_y equal xcm(upper,y)
variable position_z equal xcm(upper,z)

#initial velocity and temp

velocity upper create 300.0 123456 temp T_upper
velocity bulk_upper create 300.0 123456 temp bulk_upper

velocity bulk_middle create 300 123456 temp T_middle

why so complicated? why not initialize all atoms to 300k and then set
the x/y/z velocity for the "down" group to 0.0/0.0/0.0

###fix

fix 1 all nve
fix 2 bulk_down setforce 0.0 0.0 0.0
#fix 3 upper setforce 0.0 0.0 0.0
fix 3 bulk_middle temp/rescale 10 293.0 293.0 0.001 1.0

ouch. temp/rescale is a very bad option to control temperature,
particularly with such a small margin.

fix_modify 3 temp T_middle

####
#velocity bulk_upper set 0.0 0.0 0.0 units box
velocity upper set 0.0 0.0 -0.002 units box

before moving your nano particle, you need to make certain, that you
can correctly run a static system. it doesn't look like it. rather
that all your problems were hidden by using temp/rescale.

axel.

Dear Dr. Axel,
       Thank you for your suggestions!
        Accordingly, i modified my input script. Firstly simualte at rest at the condition of nvt, then move the particle at the condition of nve.
The "bulk_upper" temperature can be normally calculated. However, the system temperature at nve was lower than that at nvt.
Someone suggests that using temp/rescale at nve controls the system temperature, but i want to know the variation of the "upper particle" temperature.
        The reason for the variation of the temperature from nvt to nve has always puzzled me. And having no good idea for solving this problem.
      Thank you again.

Yao
在 2013-08-31 23:10:17,"Axel Kohlmeyer" <[email protected]> 写道:

Dear Dr. Axel,
       Thank you for your suggestions!
        Accordingly, i modified my input script. Firstly simualte at rest at the condition of nvt, then move the particle at the condition of nve.
The "bulk_upper" temperature can be normally calculated. However, the system temperature at nve was lower than that at nvt.
Someone suggests that using temp/rescale at nve controls the system temperature, but i want to know the variation of the "upper particle" temperature.

ah, the highly reputable source of "somebody"...

please let me remind you that it is close to impossible to give
recommendations without having sufficiently detailed and specific
information.

        The reason for the variation of the temperature from nvt to nve has always puzzled me. And having no good idea for solving this problem.

what variation for what kind of system? if you do the simulation
correctly and have a properly equilibrated system, there should be no
difference in temperature from simply switching from fix nvt to fix
nve. if you disagree, please provide a (simple!) example supporting
your opinion.

axel.

ok , i post to the list.
Thank you.

TiO_compact_5.txt (459 KB)

ok , i post to the list.
Thank you.

i see nothing unexpected in the output
(btw: for posting to the list, you could have used a significantly
larger interval than printing out thermo info every 10 steps).

in the first part you equilibrate your system to 650K

then you re-set the velocities of the left/right block to the desired
values. this will - of course - result in losing all current velocity
information except for the new center of mass velocity. consequently
the temperature for each of these blocks goes to OK when the center of
mass velocity is removed and for the total system, the temperature is
the kinetic energy you set with the velocity command.

since you did equilibrate your system, it will relax and convert about
half of the stored potential energy back into kinetic energy and you
reach the temperature that you are getting.

as i said, this all looks like LAMMPS is behaving correct. if this is
not what you want, then the error is in how you write your input. most
likely, you want to *add* the velocity not overwrite it. also, please
note, that your thermostat relaxation time is extremely short, as is
your total equilibration time. you system is quite small, so it may
not be as drastic, but for larger systems, you should correct that.
also, it is often more effective, to start equilibration with using
fix langevin+nve for a bit to enforce good equipartitioning of the
kinetic energy, before switching to fix nvt to equilibrate into a
proper NVT system.

axel.

Hi Dr. Axel,
     Thank you very much. After discussing with you, my question is much clearer. I want to *add* the velocity of two blocks when they have a equilibration
temperature not overwrite it.
      Resetting the "velocity" command may break the existed equilibration. It means that the "velocity" command must not be used again.
Then, How to add velocity to the blocks without breaking. I search for commands that both "fix move" and "velocity" can set the velocity.
Beside the "velocity" command, only the "fix move" command can be used. May be "fix move" can add the velocity without breaking.

Yao
在 2013-09-02 12:08:34,"Axel Kohlmeyer" <[email protected]> 写道:

Hi Dr. Axel,
     Thank you very much. After discussing with you, my question is much clearer. I want to *add* the velocity of two blocks when they have a equilibration
temperature not overwrite it.
      Resetting the "velocity" command may break the existed equilibration. It means that the "velocity" command must not be used again.
Then, How to add velocity to the blocks without breaking. I search for commands that both "fix move" and "velocity" can set the velocity.
Beside the "velocity" command, only the "fix move" command can be used. May be "fix move" can add the velocity without breaking.

please read the docs of the velocity command more carefully. you can
just use an atom style variable like this

variable addl atom vx+0.002
variable addr atom vz-0.002

velocity left set v_addl NULL NULL units box
velocity right set v_addr NULL NULL units box

and thus add a COM velocity to each group while preserving the
pre-existing equilibration.

axel.

typo aliert: this is what i meant to write.

variable addl atom vx+0.002
variable addr atom vx-0.002

velocity left set v_addl NULL NULL units box
velocity right set v_addr NULL NULL units box

axel.

Thank you Dr. Axel. I will try them and post the results.

Hi, Dr. Axel,
        Results for two particles show that the temperature of two particles was normal,and they are contact with each other.
However, their velocity always fluctuated and could not stop. i think whether i can fix one of them?
i use the "fix_setforce" command to fix one part of a particle, like in cutting simulation. Then i compute the temperature
of other part of the fixed particle and the unfixed one.However, their tempature become adnormal, one is very high and other is very low.
Should i fix the particle like in cutting simulation, the down is fixed by "fix_setforce 0 0 0", and the middle layer is fixed to be constant_temperature .
and the upper layer is not thermostat?

Yao

Hi, Dr. Axel,
        Results for two particles show that the temperature of two particles was normal,and they are contact with each other.
However, their velocity always fluctuated and could not stop. i think whether i can fix one of them?

it is not clear to me, what you mean by that. remember, i am not a
psychic, so i don't know all the things that you are thinking about.
also, i am not doing your research, so i cannot know many of the
details that are obvious to you.

i use the "fix_setforce" command to fix one part of a particle, like in cutting simulation. Then i compute the temperature

fix setforce by itself does not immobilize a particle. you also need
to set its velocity to zero.

of other part of the fixed particle and the unfixed one.However, their tempature become adnormal, one is very high and other is very low.
Should i fix the particle like in cutting simulation, the down is fixed by "fix_setforce 0 0 0", and the middle layer is fixed to be constant_temperature .
and the upper layer is not thermostat?

it is impossible to give any advice based on such vague formulations.
it also is not clear to me any more what this is exactly about.
so my suggestion is the following. you post a *new* message with a
*new* subject line, where you - again - *clearly* state the problem
that you worry about and provide *sufficient* information to
*reproduce* it, while leaving out everything that is not necessary.
this obviously requires some more effort from you, but that is the
whole point. if you don't make it easy to give specific help based on
specific and easy to understand information, you will end up without
help or get vague answers like "if it doesn't work, then you must be
doing something wrong", that mirror the vagueness of your description
and information.

axel.