[lammps-users] LAMMPS 25Jun07 restart region wierdness?

David_Farrell · June 25, 2007, 4:26pm

Hi all,

I recently upgraded to the 25Jun07 release of LAMMPS and continued to try and get my restart files up and running.

I noticed that for some reason after restarting, 2 odd things occur (I don’t think these are necessarily new issues though)

when re-assigning regions, the numbering is not the same for a region based off an intersect command:

original:

Groups

region allatoms block INF INF INF INF INF INF
region interior block 2 18 2 18 2 INF
region exterior intersect 2 allatoms interior side out

group interior region interior
123099 atoms in group interior
group exterior region exterior
69701 atoms in group exterior

restarted (box size is 20x20x60 lattice spacings or 80x80x~240 angstrom):

Groups

region allatoms block INF INF INF INF INF INF units box
region interior block 8 72 8 72 8 INF units box
region exterior intersect 2 allatoms interior side out units box

group interior region interior
123099 atoms in group interior
group exterior region exterior
72340 atoms in group exterior

However, the dump file still indicates there is 192800 atoms (the same for both before and after restarting)… perhaps something is up with the way I am using the intersect command in the restart file?

fix temp/rescale doesn’t seem to have any effect after restarting. I would imagine its effect should be more pronounced since it is only applied to the exterior region, and after the restart there are more atoms in the exterior (from the above):

original (prior to this an equilibriation period after which a restart file is written and used as the basis for the restarted run)

fix 1 all nve
fix 2 all temp/rescale 1 300.0 300.0 0.5 1.0 region exterior

set PKA velocity to correspond to ~ 10 keV

set atom 192421 vx 109.5108
1 settings made for vx
set atom 192421 vy 301.1547
1 settings made for vy
set atom 192421 vz -2600.8815
1 settings made for vz

Memory usage per processor = 12.7724 Mbytes
Step Temp PotEng TotEng Press Lx Ly Lz
100 630.12875 -1381475.3 -1365771.7 -8151.1659 80 80 241.00353
110 316.36918 -1381457.8 -1373573.5 -13522.539 80 80 241.00353
120 315.64542 -1381439.8 -1373573.5 -13497.301 80 80 241.00353
130 314.84976 -1381420 -1373573.5 -13450.724 80 80 241.00353
140 313.96317 -1381397.9 -1373573.5 -13426.089 80 80 241.00792
150 313.10068 -1381376.4 -1373573.5 -13416.585 80 80 241.00792
160 312.28905 -1381356.2 -1373573.5 -13399.863 80 80 241.00792
170 311.45752 -1381335.4 -1373573.5 -13372.601 80 80 241.00792
180 310.64974 -1381315.3 -1373573.5 -13356.811 80 80 241.00792
190 309.88068 -1381296.1 -1373573.5 -13339.892 80 80 241.01255
200 309.12292 -1381277.3 -1373573.5 -13322.892 80 80 241.01255

Restarted:

fix 1 all nve
fix 2 all temp/rescale 1 300.0 300.0 0.5 1.0 region exterior

set PKA velocity to correspond to ~ 10 keV

set atom 192421 vx 109.5108
1 settings made for vx
set atom 192421 vy 301.1547
1 settings made for vy
set atom 192421 vz -2600.8815
1 settings made for vz

Step Temp PotEng TotEng Press Lx Ly Lz
100 630.12875 -1381475.3 -1365771.7 -8151.1659 80 80 241.00353
110 629.53009 -1381455.6 -1365767 -8108.4845 80 80 241.00353
120 628.46098 -1381429 -1365766.9 -8052.8284 80 80 241.00353
130 627.42029 -1381403 -1365767 -8033.313 80 80 241.00353
140 626.40931 -1381377.8 -1365767 -8002.8135 80 80 241.008
150 625.4768 -1381354.6 -1365767 -7972.4492 80 80 241.008
160 624.53008 -1381331 -1365767 -7941.7513 80 80 241.008
170 623.39271 -1381302.7 -1365767 -7895.743 80 80 241.008
180 621.9946 -1381267.8 -1365767 -7864.2409 80 80 241.01241
190 620.34451 -1381222 -1365762.3 -7838.5596 80 80 241.01241
200 617.74516 -1381157.3 -1365762.3 -7807.4073 80 80 241.01241

Any Ideas? If the whole log file is needed, I can attach those/post them too (they aren’t too big)

Thanks,

Dave

David E. Farrell

Graduate Student

Mechanical Engineering

Northwestern University

email: d-farrell2@…469…

David_Farrell · June 25, 2007, 5:47pm

OK - found out this was actually more than likely due to the differences in geometry and specification of the region, so I added a bit of a buffer region to the restarted specification and it worked the way I wanted… still odd that the interior+exterior doesn’t equal the total in my mind though.

Dave

David E. Farrell

Graduate Student

Mechanical Engineering

Northwestern University

email: d-farrell2@…435…

sjplimp · June 26, 2007, 1:07pm

group interior region interior
123099 atoms in group interior
group exterior region exterior
69701 atoms in group exterior

You likely don't want to do this in your restart script. The group
assignments and names are in the restart file. So whatever atoms
were originally assigned to group "exterior" are still part of it.

If you re-issue a group command then new atoms may be added
to "exterioir". Assuming atoms have moved since the original
assignment, then likely some new atoms are now inside the region
and you are adding them. E.g. one that was "interior" before but
has now crossed into exterior, will be in both groups.

Steve

David_Farrell · June 26, 2007, 2:12pm

group interior region interior
123099 atoms in group interior
group exterior region exterior
69701 atoms in group exterior

You likely don’t want to do this in your restart script. The group
assignments and names are in the restart file. So whatever atoms
were originally assigned to group “exterior” are still part of it.

If you re-issue a group command then new atoms may be added
to “exterioir”. Assuming atoms have moved since the original
assignment, then likely some new atoms are now inside the region
and you are adding them. E.g. one that was “interior” before but
has now crossed into exterior, will be in both groups.

I had tried to not include the definition but LAMMPS repeatably would return (and still does): ERROR: Fix temp/rescale region ID does not exist

when I attempted to re-specify the fix for my thermostat. This lead me to believe that the names were either not stored the same way, or that I was somehow not accessing them correctly, and the documentation on the particulars of restart files seems a bit sparse.

Any ideas on the below, because this one seems to be a show stopper, the above I found a fix for (just re-specify the regions with a bit of fiddling to make sure the regions are the same before and after, ok for my current simulation but probably not in general):

Basically the issue is that the result of the fix after I restart and from the continued simulation are nowhere near the same thing. They should be pretty close to identical, but it appears that for some reason, the restarted simulation’s thermostat isn’t working quite right (that or the one in the continued simulation is incorrect…)

fix temp/rescale doesn’t seem to have any effect after restarting. I
would imagine its effect should be more pronounced since it is only applied
to the exterior region, and after the restart there are more atoms in the
exterior (from the above):

original (prior to this an equilibriation period after which a restart file
is written and used as the basis for the restarted run)

fix 1 all nve
fix 2 all temp/rescale 1 300.0 300.0 0.5 1.0 region exterior

set PKA velocity to correspond to ~ 10 keV

set atom 192421 vx 109.5108
1 settings made for vx
set atom 192421 vy 301.1547
1 settings made for vy
set atom 192421 vz -2600.8815
1 settings made for vz

Memory usage per processor = 12.7724 Mbytes
Step Temp PotEng TotEng Press Lx Ly Lz
100 630.12875 -1381475.3 -1365771.7 -8151.1659 80
80 241.00353
110 316.36918 -1381457.8 -1373573.5 -13522.539 80
80 241.00353
120 315.64542 -1381439.8 -1373573.5 -13497.301 80
80 241.00353
130 314.84976 -1381420 -1373573.5 -13450.724 80
80 241.00353
140 313.96317 -1381397.9 -1373573.5 -13426.089 80
80 241.00792
150 313.10068 -1381376.4 -1373573.5 -13416.585 80
80 241.00792
160 312.28905 -1381356.2 -1373573.5 -13399.863 80
80 241.00792
170 311.45752 -1381335.4 -1373573.5 -13372.601 80
80 241.00792
180 310.64974 -1381315.3 -1373573.5 -13356.811 80
80 241.00792
190 309.88068 -1381296.1 -1373573.5 -13339.892 80
80 241.01255
200 309.12292 -1381277.3 -1373573.5 -13322.892 80
80 241.01255

Restarted:

fix 1 all nve
fix 2 all temp/rescale 1 300.0 300.0 0.5 1.0 region exterior

set PKA velocity to correspond to ~ 10 keV

set atom 192421 vx 109.5108
1 settings made for vx
set atom 192421 vy 301.1547
1 settings made for vy
set atom 192421 vz -2600.8815
1 settings made for vz

Step Temp PotEng TotEng Press Lx Ly Lz
100 630.12875 -1381475.3 -1365771.7 -8151.1659 80
80 241.00353
110 629.53009 -1381455.6 -1365767 -8108.4845 80
80 241.00353
120 628.46098 -1381429 -1365766.9 -8052.8284 80
80 241.00353
130 627.42029 -1381403 -1365767 -8033.313 80
80 241.00353
140 626.40931 -1381377.8 -1365767 -8002.8135 80
80 241.008
150 625.4768 -1381354.6 -1365767 -7972.4492 80
80 241.008
160 624.53008 -1381331 -1365767 -7941.7513 80
80 241.008
170 623.39271 -1381302.7 -1365767 -7895.743 80
80 241.008
180 621.9946 -1381267.8 -1365767 -7864.2409 80
80 241.01241
190 620.34451 -1381222 -1365762.3 -7838.5596 80
80 241.01241
200 617.74516 -1381157.3 -1365762.3 -7807.4073 80
80 241.01241

Any Ideas? If the whole log file is needed, I can attach those/post them too
(they aren’t too big)

Thanks,

Dave

David E. Farrell

Graduate Student

Mechanical Engineering

Northwestern University

email: d-farrell2@…435…

This SF.net email is sponsored by DB2 Express
Download DB2 Express C - the FREE version of DB2 express and take
control of your XML. No limits. Just data. Click to get it now.
http://sourceforge.net/powerbar/db2/

lammps-users mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/lammps-users

David E. Farrell

Graduate Student

Mechanical Engineering

Northwestern University

email: d-farrell2@…435…

sjplimp · June 27, 2007, 1:31pm

Redefine the regions, but not the groups.

Steve

David_Farrell · June 27, 2007, 4:05pm

I did do that, and it seemed to work fine - but it had no effect on the thermostat in the restarted case or the case where I ran the simulation all the way through.

Dave

sjplimp · June 28, 2007, 8:06pm

Looking at your scripts, it appears the thermostat is not working
in the restart case. This could be b/c the atom you have made "hot"
is not inside the region. The region could be different in the 2 cases
b/c it is defined relative to the size of the box at the time of
defiition. In the first run, it is the size at the start of the simulation.
In the 2nd run, it is the size at the time of the restart. With
shrink-wrapped boundaries that could be different.

So you need to figure out why the hot atom isn't included. You can
check if that's right by putting a print statement in fix_temp_rescale.cpp
and seeing it that atom is being thermostatted or included in the T calc.

As an aside, I don't think you want to rescale the T of a high-energy
cascase of atoms or of the very hot atom. That would seem
to defeat the point of creating the cascade.

Steve

David_Farrell · June 28, 2007, 9:04pm

Hmm… actually, I had originally intended the ‘hot’ atom (PKA, whatever) to not be included, along with not thermostating the layers of atoms adjacent to it in the periodic cell (so the ‘upper’ layer) … I will have to take another look to be sure I didn’t mess that up.

But you are correct - the PKA should not be thermostated in any way… if it is, that is a mistake on my part.

But yes, that is the behavior I observed, seemed as though the thermostat wasn’t working… I will give your suggestion a try.

Dave

David_Farrell · June 28, 2007, 9:06pm

Just to clarify - the idea was a ‘thermostat skin’ of 2 layers of atoms on each side of the periodic cell, except for the ‘top’ … not sure if that makes any more sense though… essentially all but 1 side of the box has 2 layers of atoms which are to be rescaled, providing a thermostat for the system.

Dave

David_Farrell · July 3, 2007, 10:22pm

In my most recent attempts to try and sort out what is going on with my simulations, I ran across something I found to be odd:

As Steve recommended, I put some print statements into FixTempRescale::end_of_step so that it reads as follows below.

The odd thing I found was that when I ran the 2-part simulation straight through, it ouput all of the atoms numbers which are in the specified exterior region… however, when I ran the second part (re-specifying the fixes and regions) from a restart file of the first, it didn’t output anything (and even the flags to tell me when it made it to the end_of_step routine didn’t print). So, it sounds to me like it just isn’t using the thermostat stuff at all in the second case.

So now I have a couple questions… can someone clarify how the temp/rescale fix should be used (i.e. can the region in the region keyword be a subgroup of the group-IDs given? the docs are a bit confusing on this point) and any ideas why the restarted simulation may not be making it to the rescaling calls? I am using a fix nve on the whole system, and a fix temp/rescale on a ‘skin’ of atoms on 5 sides of my simulation box… is this the right way to set this up?

If anyone wants me to prepare a more detailed writeup of the problem and gather the files I have been testing with, I can do that… as this seems like it could be a bunch of things all causing problems. I am still debugging this to see if I can get some more information on what exactly is happening.

Thanks,

Dave

void FixTempRescale::end_of_step()
{
double t_current = temperature->compute_scalar();
double delta = update->ntimestep - update->beginstep;
delta /= update->endstep - update->beginstep;
double t_target = t_start + delta * (t_end-t_start);

// Debug
ofstream fFile;
fFile.open(“thermo_atoms.txt”);

if (fabs(t_current-t_target) > t_window) {
t_target = t_current - fraction*(t_current-t_target);
double factor = sqrt(t_target/t_current);

cout << "inside scaling routine " << endl;

double **x = atom->x;
double **v = atom->v;
int *mask = atom->mask;
int nlocal = atom->nlocal;

if (iregion == -1) {
energy += (t_current-t_target) * efactor;
for (int i = 0; i < nlocal; i++) {
if (mask[i] & groupbit) {
v[i][0] *= factor;
v[i][1] *= factor;
v[i][2] *= factor;
// Debug
fFile << i << endl;

}
}

} else {
efactor = (0.5 * force->boltz * temperature->dof);
energy += (t_current-t_target) * efactor;
for (int i = 0; i < nlocal; i++) {
if (mask[i] & groupbit &&
domain->regions[iregion]->match(x[i][0],x[i][1],x[i][2])) {
v[i][0] *= factor;
v[i][1] *= factor;
v[i][2] *= factor;
// Debug
fFile << i << endl;
}
}
}
} else energy = 0.0;

fFile.close();
}

/* ---------------------------------------------------------------------- */

David E. Farrell
Graduate Student
Mechanical Engineering
Northwestern University
email: d-farrell2@…435…

sjplimp · July 3, 2007, 10:47pm

i.e. can the region in the region keyword be a subgroup of the group-IDs given?

Regions and groups have nothing to do with each other. Regions are
geometric volumes. Groups are sets of atoms. You can give the same
ID to
a region and a group (which I think you might have done), but they
still have nothing to do with each other.

If you use fix temp/rescale with the region keyword, it will thermostat
the collection of atoms that meet two criteria: they are in the fix group AND
their coordinates are in the region. The latter test is done every timestep
the fix attempts to do something. It also happens in the temperature that
the rescale fix computes. Nothing about the fix is in the restart
file, except the assignement of atoms to groups. So as I said before,
you
should not redefine groups in your restart script.

If the fix isn't doing anything, it's probably b/c it has no atoms to act on.
It should not be possible that the end_of_step() routine is not being
entered at all during the restarted simulation.

Hope that helps,
Steve

David_Farrell · July 10, 2007, 7:23pm

OK, I have spent the last couple weeks iterating on this problem, and running some tests to get an idea of what is happening.

So, indeed the end_of_step() routine is being entered, but in the restart case, the thermostat region temperature is within the tolerances.

However, I noticed by examining the dump files and the logs (with output of the temperature of the thermostat region) or the straight-run case, for some reason the kinetic energy of the atom I excite to start the cascade drops by about 1/2 (happens in the first timestep after the KE assigned), and the temperature of the thermostat region drops to about 65 K (far lower than it should) … this is in the second phase, noticed about 10 steps into the simulation (but probably happens sooner, I printed the temp from the thermostat region once every 10 timesteps).

The timestep makes it unlikely that my energetic atom would have moved far enough to cause a big enough interaction to reduce the KE that much, but am currently looking into the trajectories to see what those look like.

Its odd that I haven’t been able to sort out what is going on, but now I am wondering if either of my simulations are set up properly (I am no longer naming my regions and groups the same, and everything else seems OK though)… so I may start over from scratch on a smaller system and see what happens.

Dave

sjplimp · July 11, 2007, 1:50pm

My suggestion would be work with a small system and be sure
your different kinds of atoms (substrate, thermostatted, hot, etc)
are assigned the way you expect and are doing what you want.
For small systems it's easy to add a few print statements to the
code and make sure the various fixes, groups, regions, all are
working on the atoms you expect.

Then move to the big system.

Steve