I Really thanks for your helps.According to your suggestions and provide solutions, we successfully solve the rigid problem of calculation of rigid body temperature.
1. add the following two lines
>compute thermo all temp/asphere
>thermo_modify temp thermo
2.download the updated sources for compute_erotate_rigid.cpp and
compute_ke_rigid.cpp from the LAMMPS-ICMS repository here:
http://git.lammps.org/git/?p=lammps-icms.git;a=tree;f=src/RIGID;hb=HEAD
We use fix it rigid/NPT/small and rigid/npt to calculate but there is an error message (Temperature compute degrees of freedom < 0 (…/compute_temp.cpp:101)****).As following,
group RIGID molecule 1:3000
compute RIGID RIGID temp
velocity all create 300.0 333
neighbor 5.0 bin
timestep 0.5
thermo 1000
compute thermo all temp/asphere
thermo_modify temp thermo
fix 1 RIGID rigid/npt/small molecule temp 300.0 300.0 10.0 iso 1 1 1 #not ok
If there is a wrong rigid body in calculation when the rigid/NPT?
Best Regards
Mark Li
Dear Axel and Trung,
I Really thanks for your helps.According to your suggestions and provide
solutions, we successfully solve the rigid problem of calculation of rigid
body temperature.
1. add the following two lines
compute thermo all temp/asphere
thermo_modify temp thermo
2.download the updated sources for compute_erotate_rigid.cpp and
compute_ke_rigid.cpp from the LAMMPS-ICMS repository here:
http://git.lammps.org/git/?p=lammps-icms.git;a=tree;f=src/RIGID;hb=HEAD
We use fix it rigid/NPT/small and rigid/npt to calculate but there is an
error message (Temperature compute degrees of freedom < 0
(../compute_temp.cpp:101)).As following,
group RIGID molecule 1:3000
compute RIGID RIGID temp
velocity all create 300.0 333
neighbor 5.0 bin
timestep 0.5
thermo 1000
compute thermo all temp/asphere
thermo_modify temp thermo
fix 1 RIGID rigid/npt/small molecule temp 300.0 300.0 10.0 iso 1
1 1 #not ok
If there is a wrong rigid body in calculation when the rigid/NPT?
no. the error message comes from a "compute temp" command.
the whole purpose of step 1. was to switch away from compute temp to
compute temp/asphere, as compute temp only considers point particles,
why you have aspherical extended particles with 3 additional degrees
of freedom per particles that compute temp does not consider, but fix
rigid will (correctly) remove, if they are taken away due to forming a
rigid body.
you have the command right there.
axel.
Dear Axel
Sorry,I didn’t describe the problem more accurately and cause you misunderstood.
We use fix it rigid/NPT/small and rigid/npt to calculate but there is an error message
(Temperature compute degrees of freedom < 0 (…/compute_temp.cpp:101)****).
If there is a wrong rigid body in calculation when the rigid/NPT?
As follow
units real
boundary p p p
atom_style hybrid ellipsoid molecular
pair_style hybrid/overlay gayberne ***
pair_coeff ***
group RIGID molecule 1:3000
compute RIGID RIGID temp
velocity all create 300.0 333
neighbor 5.0 bin
timestep 0.5
thermo 1000
>compute thermo all temp/asphere
>thermo_modify temp thermo
>fix 1 RIGID rigid/nve/small molecule langevin 300.0 300.0 5.0 2211 #NVE is ok
>fix 1 RIGID rigid/npt/small molecule temp 300.0 300.0 10.0 iso 1 1 1 #NPT not ok
B.S.
Mark
Dear Axel
* Sorry,I didn't describe the problem more accurately and cause you
misunderstood.*
what you have below is a *different* input from what you quoted previously.
don't blame me to misunderstand, when you are giving me *wrong*
information.
*We use fix it rigid/NPT/small and rigid/npt to calculate but there is an
error message *
*(Temperature compute degrees of freedom < 0** (../compute_temp.cpp:101)*
*).*
*If there is a wrong rigid body in calculation when the rigid/NPT? *
no. but you *must not* use compute temp for extended bodies.
*As follow*
units real
boundary p p p
atom_style hybrid *ellipsoid *molecular
pair_style hybrid/overlay gayberne ***
pair_coeff ***
>group RIGID molecule 1:3000
>compute RIGID RIGID temp
>velocity all create 300.0 333
>neighbor 5.0 bin
>timestep 0.5
>thermo 1000
*>compute thermo all temp/asphere *
*>thermo_modify temp thermo*
*>fix 1 RIGID rigid/nve/small molecule langevin 300.0 300.0 5.0
2211 #NVE is ok *
*>fix 1 RIGID rigid/npt/small molecule temp 300.0 300.0 10.0 iso
1 1 1 #NPT not ok*
the langevin option in rigid fixes and the temperature management in the
rigid/npt and rigid/nvt fixes work differently.
in the first case, the temperature for the rotational and translational
degrees of freedom are computed internally from the internal information
about the rigid objects and then adjusted. in the second case, an external
temperature compute is used. by default, same as for the thermo output,
this compute assumes point particles and thus you'll get the error
messages, when you apply the rigid fix to objects with additional degrees
of freedom.
so the obviously missing command here is:
fix_modify 1 temp thermo
axel.
Dear Axel,
Thanks,I according to your suggestions for rigid/NPT and rigid/NPT/small calculation is available.
As follow, but results of temp ,ke, and ke_rot are NaN Value
neighbor 2.0 bin
timestep 0.5
thermo 1000
compute thermo all temp/asphere
thermo_modify temp thermo
fix 1 RIGID rigid/npt/small molecule temp 300.0 300.0 10.0 iso 1 1 1 #NPT ok
fix_modify 1 temp thermo ------->> Add in #
compute rigid_ke RIGID ke/rigid 1
compute rigid_rot RIGID erotate/rigid 1
variable total_ke equal “c_rigid_ke+c_rigid_rot” #translational Ke +rotational Ke
variable dof equal 63000 # # each dimer has 6 dof
variable boltz equal 1.987e-3 # Boltzmann constant in units real: 1.987e-3 kcal/mol/K
variable rigid_temp equal "2.0v_total_ke/v_dof/v_boltz"
outputs:
thermo_style custom step c_thermo dt f_1 v_rigid_temp c_rigid_ke c_rigid_rot
0 284.90993 0.01 82869.875 243.63538 2720.3939 1636.5375
1000 -nan 0.43862013 -nan nan -nan -nan
2000 -nan 0.43862013 -nan nan -nan -nan
3000 -nan 0.43862013 -nan nan -nan -nan
4000 -nan 0.43862013 -nan nan -nan -nan
5000 -nan 0.43862013 -nan nan -nan -nan
Something wrong?
Regards
Mark