Hi LAMMPS developers and users,
I am confused by thermostatting a non-equilibrium system. As a test system, I applied shear on a pure solvent system, and checked the validity of ‘compute temp/deform’ after ‘fix deform’. The following are the relevant codes.
''velocity all create 0.5 235820 dist gaussian
fix 1 all nve
fix shear all deform 1 xy erate 0.05
compute stream all temp/deform
compute myTemp all temp’’
I also thermostated system with dpd/tstat at kbT = 0.5, with a friction coefficient = 5.0
As an outcome, I found ‘c_stream’ is even higher than ‘c_myTemp’, which I suppose it should be the other way around when streaming velocity is subtracted. The following are some relevant output.
Step myTemp stream
6000 0.50292292 0.75206524
6500 0.50477264 0.75377328
7000 0.50288556 0.75090999
7500 0.49943975 0.74711878
Can someone kindly points out if I have made some mistakes here? Thanks a lot in advance!
Sincerely,
Ling
As a follow up investigation, I turned off dpd/tstat as I would just want to compare the temperature calculation invoked by ‘compute temp’ and ‘compute temp/deform’, and thus I assume I would be able to match these two together if I add the streaming velocity back.
However, the temperature calculated after streaming velocity has been taken out*(compute temp/deform)* is actually higher than a sum of thermal velocity and streaming velocity (compute temp), and the difference is almost a constant = 0.25.
I wonder if anyone else had this ever happened before? or simply I made a silly mistake? Please let me know anyways, I am very confused.
Thanks very much!
Sincerely,
Ling
Compute temp/deform is doing something simple. It simply subtracts
the deformation velocity (which may or may not be the velocity
the atoms are moving at) and uses the remaineder as the thermal velocity.
If your atoms are not moving (mostly) with the box deformation
then you will get a funny result. You should verify that your atoms
are flowing with the box deformation. And you should look
at the fix deform remap option and the discussion of what it does
for liquids and solids. With remap on, it could appear (in viz) that
your atoms are moving with the box, but that is not the case. Their
velocity could still be 0.0.
Steve
Hi Steve,
Thanks for your reply. I had the ‘remap v’ on for a system I just tested, and I think the results just confirmed what you have stated by ‘the velocity are still zero’. I am posting it so as to share with everyone.
To recap my problem, initially I didn’t understand at first why the temperature (Td) after subtracting streaming velocity (v_s) from thermal velocity (v_i) could be possibly higher than the temperature (Tk) without subtracted, which I presume the former be always smaller.
To figure out, I tested a simple LJ system, without initial velocity assigned ( see following for the input script). And the results from time 0 can explain the situation. The initial thermal energy is zero, thus compute temp gave a ‘zero’ output, while compute temp/deform computes temperature by alpha* (0 - v_s )^2 . That’s why Td at time zero is a positive number.
Step Tk(by compute temp) Td (by compute temp/deform)
0 0 0.7659
500 0.2656 0.3379
1000 0.3641 0.4708
… … …
By definition, Tk = alpha * (v_i)^2; Td = alpha * (v_i - v_s)^2. So the conclusion is that there is no absolute comparison between Tk and Td, but just based on the correlation with v_i and v_s.
Thanks for the help. It seems to a trivial problem, but to me it’s important to figure out why. Please point out anything wrong as it will be of great help.
Sincerely,
Ling
Please see following for the input script.
units lj
atom_style atomic
neighbor 0.6 bin
neigh_modify delay 0 every 10 check no
read_data coord.n4000_long
fix 1 all nve
fix shear all deform 1 xy erate 0.05 remap v
pair_style lj/cut 3.0
pair_coeff * * 1.0 1.0 3.0
compute myTemp all temp
compute stream all temp/deform
thermo 500
thermo_style custom step temp c_stream c_myTemp
timestep 0.001
run 10000
Another way to say it is that you should only use compute temp/deform
(or trust its values) when your fluid is in fact flowing in a manner consistent
with the box deformation rate. You can check that by monitoring its
velocity profile, e.g. via fix ave/spatial.
Steve
Good to know! Thanks for the help, Steve.
Regards,
Ling