forked from lijiext/lammps
Merge branch 'qtb_example' of github.com:stanmoore1/lammps into pair-potential-file-unit-convert
This commit is contained in:
Axel Kohlmeyer
commit
d926274911
|
|
@ -35,7 +35,7 @@ Syntax
|
|||
v_name = variable with name that calculates value of aparam
|
||||
|
||||
* zero or more keyword/value pairs may be appended
|
||||
* keyword = *scale* or *reset*
|
||||
* keyword = *scale* or *reset* or *mass*
|
||||
|
||||
.. parsed-literal::
|
||||
|
||||
|
|
@ -45,6 +45,9 @@ Syntax
|
|||
*reset* value = *no* or *yes*
|
||||
*no* = values will remain altered at the end of a run
|
||||
*yes* = reset altered values to their original values at the end of a run
|
||||
*mass* value = *no* or *yes*
|
||||
*no* = mass is not altered by changes in diameter
|
||||
*yes* = mass is altered by changes in diameter
|
||||
|
||||
Examples
|
||||
""""""""
|
||||
|
|
@ -53,7 +56,7 @@ Examples
|
|||
|
||||
fix 1 all adapt 1 pair soft a 1 1 v_prefactor
|
||||
fix 1 all adapt 1 pair soft a 2* 3 v_prefactor
|
||||
fix 1 all adapt 1 pair lj/cut epsilon * * v_scale1 coul/cut scale 3 3 v_scale2 scale yes reset yes
|
||||
fix 1 all adapt 1 pair lj/cut epsilon * * v_scale1 pair coul/cut scale 3 3 v_scale2 scale yes reset yes
|
||||
fix 1 all adapt 10 atom diameter v_size
|
||||
|
||||
variable ramp_up equal "ramp(0.01,0.5)"
|
||||
|
|
@ -86,12 +89,13 @@ the end of a simulation. Even if *reset* is specified as *yes*\ , a
|
|||
restart file written during a simulation will contain the modified
|
||||
settings.
|
||||
|
||||
If the *scale* keyword is set to *no*\ , then the value of the altered
|
||||
parameter will be whatever the variable generates. If the *scale*
|
||||
keyword is set to *yes*\ , then the value of the altered parameter will
|
||||
be the initial value of that parameter multiplied by whatever the
|
||||
variable generates. I.e. the variable is now a "scale factor" applied
|
||||
in (presumably) a time-varying fashion to the parameter.
|
||||
If the *scale* keyword is set to *no*\ , which is the default, then
|
||||
the value of the altered parameter will be whatever the variable
|
||||
generates. If the *scale* keyword is set to *yes*\ , then the value
|
||||
of the altered parameter will be the initial value of that parameter
|
||||
multiplied by whatever the variable generates. I.e. the variable is
|
||||
now a "scale factor" applied in (presumably) a time-varying fashion to
|
||||
the parameter.
|
||||
|
||||
Note that whether scale is *no* or *yes*\ , internally, the parameters
|
||||
themselves are actually altered by this fix. Make sure you use the
|
||||
|
|
@ -107,16 +111,17 @@ style supports it. Note that the :doc:`pair_style <pair_style>` and
|
|||
to specify these parameters initially; the fix adapt command simply
|
||||
overrides the parameters.
|
||||
|
||||
The *pstyle* argument is the name of the pair style. If :doc:`pair_style hybrid or hybrid/overlay <pair_hybrid>` is used, *pstyle* should be
|
||||
a sub-style name. If there are multiple sub-styles using the same
|
||||
pair style, then *pstyle* should be specified as "style:N" where N is
|
||||
which instance of the pair style you wish to adapt, e.g. the first,
|
||||
second, etc. For example, *pstyle* could be specified as "soft" or
|
||||
"lubricate" or "lj/cut:1" or "lj/cut:2". The *pparam* argument is the
|
||||
name of the parameter to change. This is the current list of pair
|
||||
styles and parameters that can be varied by this fix. See the doc
|
||||
pages for individual pair styles and their energy formulas for the
|
||||
meaning of these parameters:
|
||||
The *pstyle* argument is the name of the pair style. If
|
||||
:doc:`pair_style hybrid or hybrid/overlay <pair_hybrid>` is used,
|
||||
*pstyle* should be a sub-style name. If there are multiple
|
||||
sub-styles using the same pair style, then *pstyle* should be specified
|
||||
as "style:N" where N is which instance of the pair style you wish to
|
||||
adapt, e.g. the first, second, etc. For example, *pstyle* could be
|
||||
specified as "soft" or "lubricate" or "lj/cut:1" or "lj/cut:2". The
|
||||
*pparam* argument is the name of the parameter to change. This is the
|
||||
current list of pair styles and parameters that can be varied by this
|
||||
fix. See the doc pages for individual pair styles and their energy
|
||||
formulas for the meaning of these parameters:
|
||||
|
||||
+---------------------------------------------------------------------+--------------------------------------------------+-------------+
|
||||
| :doc:`born <pair_born>` | a,b,c | type pairs |
|
||||
|
|
@ -234,31 +239,32 @@ the coefficients for the symmetric J,I interaction to the same values.
|
|||
|
||||
A wild-card asterisk can be used in place of or in conjunction with
|
||||
the I,J arguments to set the coefficients for multiple pairs of atom
|
||||
types. This takes the form "\*" or "\*n" or "n\*" or "m\*n". If N = the
|
||||
number of atom types, then an asterisk with no numeric values means
|
||||
all types from 1 to N. A leading asterisk means all types from 1 to n
|
||||
(inclusive). A trailing asterisk means all types from n to N
|
||||
types. This takes the form "\*" or "\*n" or "n\*" or "m\*n". If N =
|
||||
the number of atom types, then an asterisk with no numeric values
|
||||
means all types from 1 to N. A leading asterisk means all types from
|
||||
1 to n (inclusive). A trailing asterisk means all types from n to N
|
||||
(inclusive). A middle asterisk means all types from m to n
|
||||
(inclusive). Note that only type pairs with I <= J are considered; if
|
||||
asterisks imply type pairs where J < I, they are ignored.
|
||||
|
||||
IMPROTANT NOTE: If :doc:`pair_style hybrid or hybrid/overlay <pair_hybrid>` is being used, then the *pstyle* will
|
||||
be a sub-style name. You must specify I,J arguments that correspond
|
||||
to type pair values defined (via the :doc:`pair_coeff <pair_coeff>`
|
||||
command) for that sub-style.
|
||||
IMPROTANT NOTE: If :doc:`pair_style hybrid or hybrid/overlay
|
||||
<pair_hybrid>` is being used, then the *pstyle* will be a sub-style
|
||||
name. You must specify I,J arguments that correspond to type pair
|
||||
values defined (via the :doc:`pair_coeff <pair_coeff>` command) for
|
||||
that sub-style.
|
||||
|
||||
The *v_name* argument for keyword *pair* is the name of an
|
||||
:doc:`equal-style variable <variable>` which will be evaluated each time
|
||||
this fix is invoked to set the parameter to a new value. It should be
|
||||
specified as v_name, where name is the variable name. Equal-style
|
||||
variables can specify formulas with various mathematical functions,
|
||||
and include :doc:`thermo_style <thermo_style>` command keywords for the
|
||||
simulation box parameters and timestep and elapsed time. Thus it is
|
||||
easy to specify parameters that change as a function of time or span
|
||||
consecutive runs in a continuous fashion. For the latter, see the
|
||||
*start* and *stop* keywords of the :doc:`run <run>` command and the
|
||||
*elaplong* keyword of :doc:`thermo_style custom <thermo_style>` for
|
||||
details.
|
||||
:doc:`equal-style variable <variable>` which will be evaluated each
|
||||
time this fix is invoked to set the parameter to a new value. It
|
||||
should be specified as v_name, where name is the variable name.
|
||||
Equal-style variables can specify formulas with various mathematical
|
||||
functions, and include :doc:`thermo_style <thermo_style>` command
|
||||
keywords for the simulation box parameters and timestep and elapsed
|
||||
time. Thus it is easy to specify parameters that change as a function
|
||||
of time or span consecutive runs in a continuous fashion. For the
|
||||
latter, see the *start* and *stop* keywords of the :doc:`run <run>`
|
||||
command and the *elaplong* keyword of :doc:`thermo_style custom
|
||||
<thermo_style>` for details.
|
||||
|
||||
For example, these commands would change the prefactor coefficient of
|
||||
the :doc:`pair_style soft <pair_soft>` potential from 10.0 to 30.0 in a
|
||||
|
|
@ -319,23 +325,28 @@ The *atom* keyword enables various atom properties to be changed. The
|
|||
current list of atom parameters that can be varied by this fix:
|
||||
|
||||
* charge = charge on particle
|
||||
* diameter, or, diameter/disc = diameter of particle
|
||||
* diameter or diameter/disc = diameter of particle
|
||||
|
||||
The *v_name* argument of the *atom* keyword is the name of an
|
||||
:doc:`equal-style variable <variable>` which will be evaluated each time
|
||||
this fix is invoked to set, or scale the parameter to a new value.
|
||||
It should be specified as v_name, where name is the variable name. See the
|
||||
discussion above describing the formulas associated with equal-style
|
||||
variables. The new value is assigned to the corresponding attribute
|
||||
for all atoms in the fix group.
|
||||
:doc:`equal-style variable <variable>` which will be evaluated each
|
||||
time this fix is invoked to set, or scale the parameter to a new
|
||||
value. It should be specified as v_name, where name is the variable
|
||||
name. See the discussion above describing the formulas associated
|
||||
with equal-style variables. The new value is assigned to the
|
||||
corresponding attribute for all atoms in the fix group.
|
||||
|
||||
If the atom parameter is *diameter* and per-atom density and per-atom
|
||||
mass are defined for particles (e.g. :doc:`atom_style granular <atom_style>`), then the mass of each particle is also
|
||||
changed when the diameter changes. The mass is set from the particle volume
|
||||
for 3d systems (density is assumed to stay constant). For 2d, the default is
|
||||
for LAMMPS to model particles with a radius attribute as spheres.
|
||||
However, if the atom parameter is *diameter/disc*, then the mass is
|
||||
set from the particle area (the density is assumed to be in mass/distance^2 units).
|
||||
mass are defined for particles (e.g. :doc:`atom_style granular
|
||||
<atom_style>`), then the mass of each particle is, by default, also
|
||||
changed when the diameter changes. The mass is set from the particle
|
||||
volume for 3d systems (density is assumed to stay constant). For 2d,
|
||||
the default is for LAMMPS to model particles with a radius attribute
|
||||
as spheres. However, if the atom parameter is *diameter/disc*, then the
|
||||
mass is set from the particle area (the density is assumed to be in
|
||||
mass/distance^2 units). The mass of the particle may also be kept constant
|
||||
if the *mass* keyword is set to *no*. This can be useful to account for
|
||||
diameter changes that do not involve mass changes, e.g., thermal expansion.
|
||||
|
||||
|
||||
For example, these commands would shrink the diameter of all granular
|
||||
particles in the "center" group from 1.0 to 0.1 in a linear fashion
|
||||
|
|
@ -348,13 +359,63 @@ over the course of a 1000-step simulation:
|
|||
|
||||
----------
|
||||
|
||||
This fix can be used in long simulations which are restarted one or
|
||||
more times to continuously adapt simulation parameters, but it must be
|
||||
done carefully. There are two issues to consider. The first is how
|
||||
to adapt the parameters in a continuous manner from one simulation to
|
||||
the next. The second is how, if desired, to reset the parameters to
|
||||
their original values at the end of the last restarted run.
|
||||
|
||||
Note that all the parameters changed by this fix are written into a
|
||||
restart file in their current changed state. A new restarted
|
||||
simulation does not know their original time=0 values, unless the
|
||||
input script explicitly resets the parameters (after the restart file
|
||||
is read), to their original values.
|
||||
|
||||
Also note, that the time-dependent variable(s) used in the restart
|
||||
script should typically be written as a function of time elapsed since
|
||||
the original simulation began.
|
||||
|
||||
With this in mind, if the *scale* keyword is set to *no* (the default)
|
||||
in a restarted simulation, original parameters are not needed. The
|
||||
adapted parameters should seamlessly continue their variation relative
|
||||
to the preceding simulation.
|
||||
|
||||
If the *scale* keyword is set to *yes*, then the input script should
|
||||
typically reset the parameters being adapted to their original values,
|
||||
so that the scaling formula specified by the variable will operate
|
||||
correctly. An exception is if the *atom* keyword is being used with
|
||||
*scale yes*. In this case, information is added to the restart file
|
||||
so that per-atom properties in the new run will automatically be
|
||||
scaled relative to their original values. This will only work if the
|
||||
fix adapt command specified in the restart script has the same ID as
|
||||
the one used in the original script.
|
||||
|
||||
In a restarted run, if the *reset* keyword is set to *yes*, and the
|
||||
run ends in this script (as opposed to just writing more restart
|
||||
files, parameters will be restored to the values they were at the
|
||||
beginning of the run command in the restart script. Which as
|
||||
explained above, may or may not be the original values of the
|
||||
parameters. Again, an exception is if the *atom* keyword is being
|
||||
used with *reset yes* (in all the runs). In that case, the original
|
||||
per-atom parameters are stored in the restart file, and will be
|
||||
restored when the restarted run finally completes.
|
||||
|
||||
----------
|
||||
|
||||
**Restart, fix_modify, output, run start/stop, minimize info:**
|
||||
|
||||
No information about this fix is written to :doc:`binary restart files <restart>`. None of the :doc:`fix_modify <fix_modify>` options
|
||||
are relevant to this fix. No global or per-atom quantities are stored
|
||||
by this fix for access by various :doc:`output commands <Howto_output>`.
|
||||
No parameter of this fix can be used with the *start/stop* keywords of
|
||||
the :doc:`run <run>` command. This fix is not invoked during :doc:`energy minimization <minimize>`.
|
||||
If the *atom* keyword is used and the *scale* or *reset* keyword is
|
||||
set to *yes*, then this fix writes information to a restart file so
|
||||
that in a restarted run scaling can continue in a seamless manner
|
||||
and/or the per-atom values can be restored, as explained above.
|
||||
|
||||
None of the :doc:`fix_modify <fix_modify>` options are relevant to
|
||||
this fix. No global or per-atom quantities are stored by this fix for
|
||||
access by various :doc:`output commands <Howto_output>`. No parameter
|
||||
of this fix can be used with the *start/stop* keywords of the
|
||||
:doc:`run <run>` command. This fix is not invoked during :doc:`energy
|
||||
minimization <minimize>`.
|
||||
|
||||
For :doc:`rRESPA time integration <run_style>`, this fix changes
|
||||
parameters on the outermost rRESPA level.
|
||||
|
|
@ -371,4 +432,4 @@ Related commands
|
|||
Default
|
||||
"""""""
|
||||
|
||||
The option defaults are scale = no, reset = no.
|
||||
The option defaults are scale = no, reset = no, mass = yes.
|
||||
|
|
|
|||
|
|
@ -53,7 +53,7 @@ include alpha_quartz_potential.mod
|
|||
|
||||
## This part equilibrates your crystal to a pressure of ${pressure}(unit pressure) and a temperature of ${temperature}(unit temperatureture) with quantum nuclear effects
|
||||
variable p_damp equal ${delta_t}*1000 #Recommended pressure damping parameter in fix nph
|
||||
fix scapegoat_qtb all nph iso ${pressure} ${pressure} ${p_damp} #NPH does the time integration
|
||||
fix scapegoat_qtb all nph iso ${pressure} ${pressure} ${p_damp} ptemp ${temperature} #NPH does the time integration
|
||||
fix quartz_qtb all qtb temp ${temperature} damp ${damp_qtb} seed 35082 f_max 120.00 N_f 100 #Change f_max (THz) if your Debye frequency is higher
|
||||
thermo_style custom step temp press etotal vol lx ly lz pxx pyy pzz pxy pyz pxz
|
||||
thermo 100
|
||||
|
|
|
|||
|
|
@ -70,7 +70,7 @@ neigh_modify check yes every 1 delay 0 page 100000 one 2000
|
|||
|
||||
## This part equilibrates your crystal to a pressure of ${pressure}(unit pressure) and a temperature of ${temperature}(unit temperatureture) with quantum nuclear effects
|
||||
variable p_damp equal ${delta_t}*1000 #Recommended pressure damping parameter in fix nph
|
||||
fix scapegoat_qtb all nph iso ${pressure} ${pressure} ${p_damp} #NPH does the time integration
|
||||
fix scapegoat_qtb all nph iso ${pressure} ${pressure} ${p_damp} ptemp ${temperature} #NPH does the time integration
|
||||
fix quartz_qtb all qtb temp ${temperature} damp ${damp_qtb} seed 35082 f_max 120.00 N_f 100 #Change f_max (THz) if your Debye frequency is higher
|
||||
thermo_style custom step temp press etotal vol lx ly lz pxx pyy pzz pxy pyz pxz
|
||||
thermo 100
|
||||
|
|
|
|||
|
|
@ -108,6 +108,7 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
|
|||
} else error->all(FLERR,"Illegal fix adapt command");
|
||||
nadapt++;
|
||||
iarg += 6;
|
||||
|
||||
} else if (strcmp(arg[iarg],"bond") == 0 ){
|
||||
if (iarg+5 > narg) error->all(FLERR, "Illegal fix adapt command");
|
||||
adapt[nadapt].which = BOND;
|
||||
|
|
@ -127,6 +128,7 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
|
|||
} else error->all(FLERR,"Illegal fix adapt command");
|
||||
nadapt++;
|
||||
iarg += 5;
|
||||
|
||||
} else if (strcmp(arg[iarg],"kspace") == 0) {
|
||||
if (iarg+2 > narg) error->all(FLERR,"Illegal fix adapt command");
|
||||
adapt[nadapt].which = KSPACE;
|
||||
|
|
@ -137,10 +139,12 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
|
|||
} else error->all(FLERR,"Illegal fix adapt command");
|
||||
nadapt++;
|
||||
iarg += 2;
|
||||
|
||||
} else if (strcmp(arg[iarg],"atom") == 0) {
|
||||
if (iarg+3 > narg) error->all(FLERR,"Illegal fix adapt command");
|
||||
adapt[nadapt].which = ATOM;
|
||||
if (strcmp(arg[iarg+1],"diameter") == 0 || strcmp(arg[iarg+1],"diameter/disc") == 0) {
|
||||
if (strcmp(arg[iarg+1],"diameter") == 0 ||
|
||||
strcmp(arg[iarg+1],"diameter/disc") == 0) {
|
||||
adapt[nadapt].aparam = DIAMETER;
|
||||
diamflag = 1;
|
||||
discflag = 0;
|
||||
|
|
@ -163,6 +167,7 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
|
|||
|
||||
resetflag = 0;
|
||||
scaleflag = 0;
|
||||
massflag = 1;
|
||||
|
||||
while (iarg < narg) {
|
||||
if (strcmp(arg[iarg],"reset") == 0) {
|
||||
|
|
@ -177,9 +182,22 @@ nadapt(0), id_fix_diam(NULL), id_fix_chg(NULL), adapt(NULL)
|
|||
else if (strcmp(arg[iarg+1],"yes") == 0) scaleflag = 1;
|
||||
else error->all(FLERR,"Illegal fix adapt command");
|
||||
iarg += 2;
|
||||
} else if (strcmp(arg[iarg],"mass") == 0) {
|
||||
if (iarg+2 > narg)error->all(FLERR,"Illegal fix adapt command");
|
||||
if (strcmp(arg[iarg+1],"no") == 0) massflag = 0;
|
||||
else if (strcmp(arg[iarg+1],"yes") == 0) massflag = 1;
|
||||
else error->all(FLERR,"Illegal fix adapt command");
|
||||
iarg += 2;
|
||||
} else error->all(FLERR,"Illegal fix adapt command");
|
||||
}
|
||||
|
||||
// if scaleflag set with diameter or charge adaptation,
|
||||
// then previous step scale factors are written to restart file
|
||||
// initialize them here in case one is used and other is never defined
|
||||
|
||||
if (scaleflag && (diamflag || chgflag)) restart_global = 1;
|
||||
previous_diam_scale = previous_chg_scale = 1.0;
|
||||
|
||||
// allocate pair style arrays
|
||||
|
||||
int n = atom->ntypes;
|
||||
|
|
@ -238,6 +256,7 @@ int FixAdapt::setmask()
|
|||
|
||||
void FixAdapt::post_constructor()
|
||||
{
|
||||
if (!resetflag) return;
|
||||
if (!diamflag && !chgflag) return;
|
||||
|
||||
// new id = fix-ID + FIX_STORE_ATTRIBUTE
|
||||
|
|
@ -432,15 +451,19 @@ void FixAdapt::init()
|
|||
error->all(FLERR,"Fix adapt requires atom attribute diameter");
|
||||
if (!atom->rmass_flag)
|
||||
error->all(FLERR,"Fix adapt requires atom attribute mass");
|
||||
if (discflag && domain->dimension!=2)
|
||||
if (discflag && domain->dimension != 2)
|
||||
error->all(FLERR,"Fix adapt requires 2d simulation");
|
||||
if (!restart_reset) previous_diam_scale = 1.0;
|
||||
}
|
||||
if (ad->aparam == CHARGE) {
|
||||
if (!atom->q_flag)
|
||||
error->all(FLERR,"Fix adapt requires atom attribute charge");
|
||||
if (!restart_reset) previous_chg_scale = 1.0;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (restart_reset) restart_reset = 0;
|
||||
|
||||
// make copy of original pair/bond array values
|
||||
|
||||
|
|
@ -570,48 +593,54 @@ void FixAdapt::change_settings()
|
|||
|
||||
} else if (ad->which == ATOM) {
|
||||
|
||||
// reset radius from diameter
|
||||
// also scale rmass to new value
|
||||
|
||||
// reset radius to new value, for both owned and ghost atoms
|
||||
// also reset rmass to new value assuming density remains constant
|
||||
// for scaleflag, previous_diam_scale is the scale factor on previous step
|
||||
|
||||
if (ad->aparam == DIAMETER) {
|
||||
double density;
|
||||
|
||||
// Get initial diameter if `scale` keyword is used
|
||||
|
||||
double *vec = fix_diam->vstore;
|
||||
double density,scale;
|
||||
double *radius = atom->radius;
|
||||
double *rmass = atom->rmass;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
int nall = nlocal + atom->nghost;
|
||||
|
||||
for (i = 0; i < nall; i++)
|
||||
if (scaleflag) scale = value / previous_diam_scale;
|
||||
|
||||
for (i = 0; i < nall; i++) {
|
||||
if (mask[i] & groupbit) {
|
||||
if (discflag) density = rmass[i] / (MY_PI * radius[i]*radius[i]);
|
||||
else density = rmass[i] / (4.0*MY_PI/3.0 *
|
||||
radius[i]*radius[i]*radius[i]);
|
||||
if (scaleflag) radius[i] = value * vec[i];
|
||||
if (massflag) {
|
||||
if (!scaleflag) scale = 0.5*value / radius[i];
|
||||
if (discflag) rmass[i] *= scale*scale;
|
||||
else rmass[i] *= scale*scale*scale;
|
||||
}
|
||||
if (scaleflag) radius[i] *= scale;
|
||||
else radius[i] = 0.5*value;
|
||||
if (discflag) rmass[i] = MY_PI * radius[i]*radius[i] * density;
|
||||
else rmass[i] = 4.0*MY_PI/3.0 *
|
||||
radius[i]*radius[i]*radius[i] * density;
|
||||
}
|
||||
}
|
||||
|
||||
if (scaleflag) previous_diam_scale = value;
|
||||
|
||||
// reset charge to new value, for both owned and ghost atoms
|
||||
// for scaleflag, previous_chg_scale is the scale factor on previous step
|
||||
|
||||
} else if (ad->aparam == CHARGE) {
|
||||
|
||||
// Get initial charge if `scale` keyword is used
|
||||
|
||||
double *vec = fix_chg->vstore;
|
||||
double scale;
|
||||
double *q = atom->q;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
int nall = nlocal + atom->nghost;
|
||||
|
||||
for (i = 0; i < nall; i++)
|
||||
if (scaleflag) scale = value / previous_chg_scale;
|
||||
|
||||
for (i = 0; i < nall; i++) {
|
||||
if (mask[i] & groupbit) {
|
||||
if (scaleflag) q[i] = value * vec[i];
|
||||
if (scaleflag) q[i] *= scale;
|
||||
else q[i] = value;
|
||||
}
|
||||
}
|
||||
|
||||
if (scaleflag) previous_chg_scale = value;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
@ -672,23 +701,24 @@ void FixAdapt::restore_settings()
|
|||
|
||||
} else if (ad->which == ATOM) {
|
||||
if (diamflag) {
|
||||
double density;
|
||||
double density,scale;
|
||||
|
||||
double *vec = fix_diam->vstore;
|
||||
double *radius = atom->radius;
|
||||
double *rmass = atom->rmass;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
|
||||
if (scaleflag) scale = previous_diam_scale;
|
||||
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (mask[i] & groupbit) {
|
||||
if(discflag) density = rmass[i] / (MY_PI * radius[i]*radius[i]);
|
||||
else density = rmass[i] / (4.0*MY_PI/3.0 *
|
||||
radius[i]*radius[i]*radius[i]);
|
||||
if (massflag) {
|
||||
if (!scaleflag) scale = vec[i] / radius[i];
|
||||
if (discflag) rmass[i] *= scale*scale;
|
||||
else rmass[i] *= scale*scale*scale;
|
||||
}
|
||||
radius[i] = vec[i];
|
||||
if(discflag) rmass[i] = MY_PI * radius[i]*radius[i] * density;
|
||||
else rmass[i] = 4.0*MY_PI/3.0 *
|
||||
radius[i]*radius[i]*radius[i] * density;
|
||||
}
|
||||
}
|
||||
if (chgflag) {
|
||||
|
|
@ -717,3 +747,28 @@ void FixAdapt::set_arrays(int i)
|
|||
if (fix_diam) fix_diam->vstore[i] = atom->radius[i];
|
||||
if (fix_chg) fix_chg->vstore[i] = atom->q[i];
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
write scale factors for diameter and charge to restart file
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixAdapt::write_restart(FILE *fp)
|
||||
{
|
||||
int size = 2*sizeof(double);
|
||||
|
||||
fwrite(&size,sizeof(int),1,fp);
|
||||
fwrite(&previous_diam_scale,sizeof(double),1,fp);
|
||||
fwrite(&previous_chg_scale,sizeof(double),1,fp);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
use scale factors from restart file to restart the Fix
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixAdapt::restart(char *buf)
|
||||
{
|
||||
double *dbuf = (double *) buf;
|
||||
|
||||
previous_diam_scale = dbuf[0];
|
||||
previous_chg_scale = dbuf[1];
|
||||
}
|
||||
|
|
|
|||
|
|
@ -40,13 +40,16 @@ class FixAdapt : public Fix {
|
|||
void setup_pre_force_respa(int,int);
|
||||
void pre_force_respa(int,int,int);
|
||||
void set_arrays(int);
|
||||
void write_restart(FILE *);
|
||||
void restart(char *);
|
||||
|
||||
private:
|
||||
int nadapt,resetflag,scaleflag;
|
||||
int nadapt,resetflag,scaleflag,massflag;
|
||||
int anypair, anybond;
|
||||
int nlevels_respa;
|
||||
char *id_fix_diam,*id_fix_chg;
|
||||
class FixStore *fix_diam,*fix_chg;
|
||||
double previous_diam_scale,previous_chg_scale;
|
||||
int discflag;
|
||||
|
||||
struct Adapt {
|
||||
|
|
|
|||
|
|
@ -909,6 +909,7 @@ void Modify::add_fix(int narg, char **arg, int trysuffix)
|
|||
strcmp(style_restart_global[i],fix[ifix]->style) == 0) {
|
||||
fix[ifix]->restart(state_restart_global[i]);
|
||||
used_restart_global[i] = 1;
|
||||
fix[ifix]->restart_reset = 1;
|
||||
if (comm->me == 0) {
|
||||
if (screen)
|
||||
fprintf(screen,"Resetting global fix info from restart file:\n");
|
||||
|
|
|
|||
Loading…
Reference in New Issue