diff --git a/doc/src/fix_electrode_conp.rst b/doc/src/fix_electrode_conp.rst index 02aa73ad54..246b212b0e 100644 --- a/doc/src/fix_electrode_conp.rst +++ b/doc/src/fix_electrode_conp.rst @@ -31,26 +31,21 @@ Syntax * ID, group-ID are documented in fix command * mode = electrode/conp or electrode/conq or electrode/thermo -* potential = electrode potential (number, or equal-style variable) -* charge = electrode charge -* eta = reciprocal width of electrode charge smearing +* potential = electrode potential (can be an equal-style variable, see below) +* charge = electrode charge (can be an equal-style variable, see below) +* eta = reciprocal width of electrode charge Gaussian distribution * T_v = temperature of thermo-potentiostat * tau_v = time constant of thermo-potentiostat -* rng_v = integer used to initialize random number generator +* rng_v = integer used to initialize thermo-potentiostat random number generator +* Optional keywords and values: .. parsed-literal:: - *algo mat_inv/mat_cg value/cg value* + *algo mat_inv/mat_cg tol/cg tol* specify the algorithm used to compute the electrode charges - *symm(etry) on/off* - turn on/off charge neutrality constraint for the electrodes *couple group-ID value* group-ID = group of atoms treated as additional electrode value = electric potential or charge on this electrode - *etypes* - construct optimized neighbor lists (types of the electrode must be exclusive to them) - *ffield on/off* - turn on/off finite-field implementation *write_mat filename* write elastance matrix to file *write_inv filename* @@ -59,6 +54,12 @@ Syntax read elastance matrix from file *read_inv filename* read inverted matrix from file + *symm(etry) on/off* + turn on/off charge neutrality constraint for the electrodes + *ffield on/off* + turn on/off finite-field implementation + *etypes on/off* + turn on/off type-based optimized neighbor lists (electrode and electrolyte types may not overlap) Examples """""""" @@ -66,32 +67,81 @@ Examples .. code-block:: LAMMPS fix fxconp bot electrode/conp -1.0 1.805 couple top 1.0 couple ref 0.0 write_inv inv.csv symm on - fix fxconp electrodes electrode/conq 0.0 1.805 + fix fxconp electrodes electrode/conq 0.0 1.805 algo cg 1e-5 fix fxconp bot electrode/thermo -1.0 1.805 temp 298 100 couple top 1.0 Description """"""""""" -fix electrode/conp mode implements a constant potential method (CPM) -(:ref:`Siepmann `, :ref:`Reed `). Charges of groups specified -via group-ID and optionally with the `couple` keyword are adapted to meet their respective -potential at every time step. An arbitrary number of electrodes can be set but -the respective groups may not overlap. Electrode charges have a Gaussian charge -distribution with reciprocal width eta. +The *electrode* fixes implement the constant potential method (CPM) +(:ref:`Siepmann `, :ref:`Reed `), and modern variants, to +accurately model electrified, conductive electrodes. This is primarily useful +for studying electrode-electrolyte interfaces, especially at high potential +differences or ionicities, with non-planar electrodes such as nanostructures or +nanopores, and to study dynamic phenomena such as charging or discharging time +scales or conductivity or ionic diffusivities. -Three algorithms are available to minimize the energy: -via matrix inversion (*algo mat_inv*) :ref:`(Wang) `, or with the conjugate gradient -method either using a matrix for the electrode-electrode interaction (*algo -mat_cg value*) or computing the interaction directly (*algo cg value*). This -method requires a threshold value to set the accuracy. +Each *electrode* fix allows users to set additional electrostatic relationships +between the specified groups which model useful electrostatic configurations: -fix electrode/conq enforces a total charge specified in the input on each electrode. The energy is -minimized w.r.t. the charge distribution within the electrode. +* *electrode/conp* sets potentials or potential differences between electrodes -fix electrode/thermo implements a thermo-potentiostat :ref:`(Deissenbeck) -`. Temperature and time constant of the thermo-potentiostat need -to be specified using the temp keyword. Currently, only two electrodes are possible with -this style. This fix is not compatible with the conjugate gradient algorithms. + * (resulting in changing electrode total charges) + +* *electrode/conq* sets the total charge on each electrode + + * (resulting in changing electrode potentials) + +* *electrode/thermo* sets a thermopotentiostat :ref:`(Deissenbeck)` between two electrodes + + * (resulting in changing charges and potentials with appropriate average potential difference and thermal variance) + +The first group-ID provided to each fix specifies the first electrode group, and +more group(s) are added using the *couple* keyword for each additional group. +While *electrode/thermo* only accepts two groups, *electrode/conp* and +*electrode/conq* accept any number of groups, up to LAMMPS's internal +restrictions (see Restrictions below). Electrode groups must not overlap, i.e. +the fix will issue an error if any particle is detected to belong to at least +two electrode groups. + +CPM involves updating charges on groups of electrode particles, per time step, +so that the system's total energy is minimized with respect to those charges. +From basic electrostatics, this is equivalent to making each group conductive, +or imposing an equal electrostatic potential on every particle in the same group +(hence the name CPM). The charges are usually modelled as a Gaussian +distribution to make the charge-charge interaction matrix invertible +(:ref:`Gingrich `). The keyword *eta* specifies the distribution's +width in units of inverse length. + +Three algorithms are available to minimize the energy, varying in how matrices +are pre-calculated before a run to provide computational speedup. These +algorithms can be selected using the keyword *algo*: + +* *algo mat_inv* pre-calculates the capacitance matrix + and obtains the charge configuration in one matrix-vector calculation per time step + +* *algo mat_cg* pre-calculates the elastance matrix (inverse of capacitance matrix) + and obtains the charge configuration using a conjugate gradient solver + in multiple matrix-vector calculations per time step + +* *algo cg* does not perform any pre-calculation and obtains the charge configuration + using a conjugate gradient solver and multiple calculations of the electric potential per time step. + +For both *cg* methods, the command must specify the conjugate gradient +tolerance. *fix electrode/thermo* currently only supports the *mat_inv* +algorithm. + +For *fix electrode/conp*, the keyword *symm* can be set *on* (or *off*) to turn +on (or turn off) the constraint requiring total electrode charge to be zero. +With *symm off*, the potentials specified are absolute potentials, but the +charge configurations satisfying them may add up to an overall non-zero, varying +charge for the electrodes (and thus the simulation box). With *symm on*, the +total charge over all electrode groups is constrained to zero. + +The *symm* keyword is not allowed for *electrode/conq*, for which overall +neutrality is explicitly obeyed or violated by the user input, and for +*electrode/thermo*, for which overall neutrality is always automatically +imposed. For all three fixes, any potential (or charge for *conq*) can be specified as an equal-style variable prefixed with "v\_". For example, the following code will @@ -105,15 +155,16 @@ course of the simulation: Note that these fixes only parse their supplied variable name when starting a run, and so these fixes will accept equal-style variables defined *after* the -fix definition, including variables dependent on the fix's own output. For an -advanced example of this see the in.conq2 input file in the directory +fix definition, including variables dependent on the fix's own output. This is +useful, for example, in the fix's internal finite-field commands (see below). +For an advanced example of this see the in.conq2 input file in the directory examples/PACKAGES/electrode/graph-il. -This fix necessitates the use of a long range solver that calculates and provides the matrix -of electrode-electrode interactions and a vector of electrode-electrolyte -interactions. The Kspace styles *ewald/electrode*, *pppm/electrode* and -*pppm/electrode/intel* are created specifically for this task -:ref:`(Ahrens-Iwers) `. +This fix necessitates the use of a long range solver that calculates and +provides the matrix of electrode-electrode interactions and a vector of +electrode-electrolyte interactions. The Kspace styles *ewald/electrode*, +*pppm/electrode* and *pppm/electrode/intel* are created specifically for this +task :ref:`(Ahrens-Iwers) `. For systems with non-periodic boundaries in one or two directions dipole corrections are available with the :doc:`kspace_modify `. For @@ -135,6 +186,49 @@ moderate mesh size but requires more memory. kspace_modify amat onestep/twostep +For all versions of the fix, the keyword-value *ffield on* enables the +finite-field mode (:ref:`Dufils `, :ref:`Tee `), which uses an +electric field across a periodic cell instead of non-periodic boundary +conditions to impose a potential difference between the two electrodes bounding +the cell. The fix (with name *fix-ID*) detects which of the two electrodes is +"on top" (has the larger maximum *z*-coordinate among all particles). Assuming +the first electrode group is on top, it then issues the following commands +internally: + +.. code-block:: LAMMPS + + variable fix-ID_ffield_zfield equal (f_fix-ID[2]-f_fix-ID[1])/lz + efield fix-ID_efield all efield 0.0 0.0 v_fix-ID_ffield_zfield + +which implements the required electric field as the potential difference divided +by cell length. The internal commands use variable so that the electric field +will correctly vary with changing potentials in the correct way (for example +with equal-style potential difference or with *fix electrode/conq*). This +keyword requires two electrodes and will issue an error with any other number of +electrodes. + +For all versions of the fix, the keyword-value *etypes on* enables type-based +optimized neighbor lists. With this feature enabled, LAMMPS provides the fix +with an occasional neighborlist restricted to electrode-electrode interactions +for calculating the electrode matrix, and a perpetual neighborlist restricted to +electrode-electrolyte interactions for calculating the electrode potentials, +using particle types to list only desired interactions, and typically resulting +in 5--10\% less computational time. Without this feature the fix will simply +use the active pair style's neighborlist. This feature cannot be enabled if any +electrode particle has the same type as any electrolyte particle (which would be +unusual in a typical simulation) and the fix will issue an error in that case. + +.. + (if we merge the overlap_etypes branch) + This feature will provide minimal benefit if any electrode particle has the same type as any + electrolyte particle, since it will be impossible for LAMMPS to list only electrode-electrolyte + neighbor pairs and discard other neighbor pairs from the provided perpetual neighborlist. + + +Restart, fix_modify, output, run start/stop, minimize info +""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" + +This fix currently does not write any information to restart files. The *fix_modify tf* option enables the Thomas-Fermi metallicity model (:ref:`Scalfi `) and allows parameters to be set for each atom type. @@ -144,7 +238,8 @@ The *fix_modify tf* option enables the Thomas-Fermi metallicity model fix_modify ID tf type length voronoi -If this option is used parameters must be set for all atom types of the electrode. +If this option is used parameters must be set for all atom types of the +electrode. The *fix_modify timer* option turns on (off) additional timer outputs in the log file, for code developers to track optimization. @@ -156,8 +251,8 @@ file, for code developers to track optimization. ---------- These fixes compute a global (extensive) scalar, a global (intensive) vector, -and a global array, which can be accessed by various -:doc:`output commands `. +and a global array, which can be accessed by various :doc:`output commands +`. The global scalar outputs the energy added to the system by this fix, which is the negative of the total charge on each electrode multiplied by that @@ -167,48 +262,46 @@ The global vector outputs the potential on each electrode (and thus has *N* entries if the fix manages *N* electrode groups), in :doc:`units ` of electric field multiplied by distance (thus volts for *real* and *metal* units). The electrode groups' ordering follows the order in which they were input in the -fix command using *couple*. The global vector output is useful for -*fix electrode/conq* and *fix electrode/thermo*, -where potential is dynamically updated based on electrolyte configuration -instead of being directly set. +fix command using *couple*. The global vector output is useful for *fix +electrode/conq* and *fix electrode/thermo*, where potential is dynamically +updated based on electrolyte configuration instead of being directly set. The global array has *N* rows and *2N+1* columns, where the fix manages *N* electrode groups managed by the fix. For the *I*-th row of the array, the elements are: * array[I][1] = total charge that group *I* would have had *if it were at 0 V - applied potential* -* array[I][2 to *N* + 1] = the *N* entries of the *I*-th row of the electrode - capacitance matrix (definition follows) -* array[I][*N* + 2 to *2N* + 1] = the *N* entries of the *I*-th row of the electrode - elastance matrix (the inverse of the electrode capacitance matrix) + applied potential* * array[I][2 to *N* + 1] = the *N* entries of the *I*-th + row of the electrode capacitance matrix (definition follows) * array[I][*N* + + 2 to *2N* + 1] = the *N* entries of the *I*-th row of the electrode elastance + matrix (the inverse of the electrode capacitance matrix) The :math:`N \times N` electrode capacitance matrix, denoted :math:`\mathbf{C}` -in the following equation, summarizes how the total charge induced on each electrode -(:math:`\mathbf{Q}` as an *N*-vector) is related to the potential on each -electrode, :math:`\mathbf{V}`, and the charge-at-0V :math:`\mathbf{Q}_{0V}` +in the following equation, summarizes how the total charge induced on each +electrode (:math:`\mathbf{Q}` as an *N*-vector) is related to the potential on +each electrode, :math:`\mathbf{V}`, and the charge-at-0V :math:`\mathbf{Q}_{0V}` (which is influenced by the local electrolyte structure): .. math:: \mathbf{Q} = \mathbf{Q}_{0V} + \mathbf{C} \cdot \mathbf{V} -The charge-at-0V, electrode capacitance and elastance matrices are internally used to -calculate the potentials required to induce the specified total electrode charges -in *fix electrode/conq* and *fix electrode/thermo*. With the *symm on* option, -the electrode capacitance matrix would be singular, and thus its last row is -replaced with *N* copies of its top-left entry (:math:`\mathbf{C}_{11}`) for -invertibility. +The charge-at-0V, electrode capacitance and elastance matrices are internally +used to calculate the potentials required to induce the specified total +electrode charges in *fix electrode/conq* and *fix electrode/thermo*. With the +*symm on* option, the electrode capacitance matrix would be singular, and thus +its last row is replaced with *N* copies of its top-left entry +(:math:`\mathbf{C}_{11}`) for invertibility. The global array output is mainly useful for quickly determining the 'vacuum capacitance' of the system (capacitance with only electrodes, no electrolyte), and can also be used for advanced simulations setting the potential as some -function of the charge-at-0V (such as in the in.conq2 example mentioned above). +function of the charge-at-0V (such as the `in.conq2` example mentioned above). -Please cite :ref:`(Ahrens-Iwers2022) ` in any publication that uses -this implementation. -Please cite also the publication on the combination of the CPM with pppm if you -use *pppm/electrode* :ref:`(Ahrens-Iwers) `. +Please cite :ref:`(Ahrens-Iwers2022) ` in any publication that +uses this implementation. Please cite also the publication on the combination +of the CPM with PPPM if you use *pppm/electrode* :ref:`(Ahrens-Iwers) +`. ---------- @@ -218,6 +311,14 @@ Restrictions For algorithms that use a matrix for the electrode-electrode interactions, positions of electrode particles have to be immobilized at all times. +With *ffield off* (i.e. the default), the box geometry is expected to be +*z*-non-periodic (i.e. *boundary p p f*), and this fix will issue an error if +the box is *z*-periodic. With *ffield on*, the box geometry is expected to be +*z*-periodic, and this fix will issue an error if the box is *z*-non-periodic. + +TODO: will fix check if *kspace_modify slab* is enabled or does it silently give +wrong results? + The parallelization for the fix works best if electrode atoms are evenly distributed across processors. For a system with two electrodes at the bottom and top of the cell this can be achieved with *processors * * 2*, or with the @@ -230,6 +331,22 @@ line which avoids an error if the script is run on an odd number of processors (such as on just one processor for testing). +The fix creates an additional group named *[fix-ID]_group* which is the union of +all electrode groups supplied to LAMMPS. This additional group counts towards +LAMMPS's limitation on the total number of groups (currently 32), which may not +allow scripts that use that many groups to run with this fix. + +The matrix-based algorithms (*algo mat_inv* and *algo mat_cg*) currently store +an interaction matrix (either elastance or capacitance) of *N* by *N* doubles +for each MPI process. This memory requirement may be prohibitive for large +electrode groups. The fix will issue a warning if it expects to use more than +0.5 GiB of memory. + +Default +""""""" +The default keyword-option settings are *algo mat_inv*, *symm off*, *etypes off* +and *ffield off*. + ---------- .. include:: accel_styles.rst @@ -244,14 +361,14 @@ as on just one processor for testing). **(Reed)** Reed *et al.*, J. Chem. Phys. 126, 084704 (2007). -.. _Wang5: - -**(Wang)** Wang *et al.*, J. Chem. Phys. 141, 184102 (2014). - .. _Deissenbeck: **(Deissenbeck)** Deissenbeck *et al.*, Phys. Rev. Letters 126, 136803 (2021). +.. _Gingrich: + +**(Gingrich)** Gingrich, `MSc thesis` ` (2010). + .. _Ahrens-Iwers: **(Ahrens-Iwers)** Ahrens-Iwers and Meissner, J. Chem. Phys. 155, 104104 (2021). @@ -260,6 +377,14 @@ as on just one processor for testing). **(Hu)** Hu, J. Chem. Theory Comput. 10, 5254 (2014). +.. _Dufils: + +**(Dufils)** Dufils *et al.*, Phys. Rev. Letters 123, 195501 (2019). + +.. _Tee: + +**(Tee)** Tee and Searles, J. Chem. Phys. 156, 184101 (2022). + .. _Scalfi: **(Scalfi)** Scalfi *et al.*, J. Chem. Phys., 153, 174704 (2020). @@ -267,4 +392,3 @@ as on just one processor for testing). .. _Ahrens-Iwers2: **(Ahrens-Iwers2022)** Ahrens-Iwers *et al.*, J. Chem. Phys. 157, 084801 (2022). - diff --git a/examples/PACKAGES/electrode/au-aq/in.ffield b/examples/PACKAGES/electrode/au-aq/in.ffield index 35b053ea3a..de594f8c9a 100644 --- a/examples/PACKAGES/electrode/au-aq/in.ffield +++ b/examples/PACKAGES/electrode/au-aq/in.ffield @@ -4,7 +4,7 @@ boundary p p p # ffield uses periodic z-boundary and no slab include settings.mod # styles, groups, computes and fixes -fix conp bot electrode/conp -1.0 1.805132 couple top 1.0 symm on ffield yes etypes +fix conp bot electrode/conp -1.0 1.805132 couple top 1.0 symm on ffield yes etypes on thermo 50 thermo_style custom step temp c_ctemp epair etotal c_qtop c_qbot c_qztop c_qzbot diff --git a/examples/PACKAGES/electrode/au-aq/in.tf b/examples/PACKAGES/electrode/au-aq/in.tf index 48a226456b..23beb357d8 100644 --- a/examples/PACKAGES/electrode/au-aq/in.tf +++ b/examples/PACKAGES/electrode/au-aq/in.tf @@ -6,7 +6,7 @@ boundary p p p # ffield uses periodic z-boundary and no slab include settings.mod # styles, groups, computes and fixes -fix conp bot electrode/conp -1.0 1.805132 couple top 1.0 symm on ffield yes etypes +fix conp bot electrode/conp -1.0 1.805132 couple top 1.0 symm on ffield yes etypes on fix_modify conp tf 6 1.0 18.1715745 fix_modify conp tf 7 1.0 18.1715745 diff --git a/examples/PACKAGES/electrode/graph-il/in.conq b/examples/PACKAGES/electrode/graph-il/in.conq index ea4685fa36..9706ba3b9b 100644 --- a/examples/PACKAGES/electrode/graph-il/in.conq +++ b/examples/PACKAGES/electrode/graph-il/in.conq @@ -5,7 +5,7 @@ boundary p p f # slab calculation include settings.mod # styles, groups, computes and fixes kspace_modify slab 3.0 -fix conq bot electrode/conq -1.0 1.979 couple top 1.0 etypes # conq doesn't take symm option +fix conq bot electrode/conq -1.0 1.979 couple top 1.0 etypes on # conq doesn't take symm option thermo 50 thermo_style custom step temp c_ctemp epair etotal c_qbot c_qtop f_conq[1] f_conq[2] diff --git a/examples/PACKAGES/electrode/graph-il/in.conq2 b/examples/PACKAGES/electrode/graph-il/in.conq2 index ad706bc50b..d860bf6c1a 100644 --- a/examples/PACKAGES/electrode/graph-il/in.conq2 +++ b/examples/PACKAGES/electrode/graph-il/in.conq2 @@ -5,7 +5,7 @@ boundary p p f # slab calculation include settings.mod # styles, groups, computes and fixes kspace_modify slab 3.0 -fix conp bot electrode/conp v_vbot 1.979 couple top v_vtop etypes +fix conp bot electrode/conp v_vbot 1.979 couple top v_vtop etypes on variable qex_bot equal -1.0-f_conp[1][1] # difference between desired and 0V charge variable qex_top equal 1.0-f_conp[2][1] # difference between desired and 0V charge diff --git a/examples/PACKAGES/electrode/graph-il/in.etypes b/examples/PACKAGES/electrode/graph-il/in.etypes index b9a0c65979..990c8bb415 100644 --- a/examples/PACKAGES/electrode/graph-il/in.etypes +++ b/examples/PACKAGES/electrode/graph-il/in.etypes @@ -5,7 +5,7 @@ boundary p p f # slab calculation include settings.mod # styles, groups, computes and fixes kspace_modify slab 3.0 -fix conp bot electrode/conp -1.0 1.979 couple top 1.0 symm on etypes +fix conp bot electrode/conp -1.0 1.979 couple top 1.0 symm on etypes on thermo 50 thermo_style custom step temp c_ctemp epair etotal c_qbot c_qtop diff --git a/examples/PACKAGES/electrode/graph-il/in.ffield b/examples/PACKAGES/electrode/graph-il/in.ffield index 66614899ec..0771e91e62 100644 --- a/examples/PACKAGES/electrode/graph-il/in.ffield +++ b/examples/PACKAGES/electrode/graph-il/in.ffield @@ -4,7 +4,7 @@ boundary p p p # ffield uses periodic z-boundary and no slab include settings.mod # styles, groups, computes and fixes -fix conp bot electrode/conp -1.0 1.979 couple top 1.0 symm on etypes ffield yes +fix conp bot electrode/conp -1.0 1.979 couple top 1.0 symm on etypes on ffield yes thermo 50 thermo_style custom step temp c_ctemp epair etotal c_qbot c_qtop diff --git a/examples/PACKAGES/electrode/graph-il/in.ramp b/examples/PACKAGES/electrode/graph-il/in.ramp index 12461bf455..e4b897937d 100644 --- a/examples/PACKAGES/electrode/graph-il/in.ramp +++ b/examples/PACKAGES/electrode/graph-il/in.ramp @@ -6,7 +6,7 @@ include settings.mod # styles, groups, computes and fixes kspace_modify slab 3.0 variable v equal ramp(2,4) -fix conp bot electrode/conp 0.0 1.979 couple top v_v symm on etypes +fix conp bot electrode/conp 0.0 1.979 couple top v_v symm on etypes on thermo 50 thermo_style custom step temp c_ctemp epair etotal c_qbot c_qtop v_v diff --git a/examples/PACKAGES/electrode/graph-il/in.thermo b/examples/PACKAGES/electrode/graph-il/in.thermo index 6ae29e2de1..253dafc138 100644 --- a/examples/PACKAGES/electrode/graph-il/in.thermo +++ b/examples/PACKAGES/electrode/graph-il/in.thermo @@ -6,9 +6,9 @@ include settings.mod # styles, groups, computes and fixes kspace_modify slab 3.0 unfix nvt # remove NVT thermostat included from "settings.mod" -fix conpthermo bot electrode/thermo -1.0 1.979 couple top 1.0 etypes temp 500 100 7 # temp tau rng +fix conpthermo bot electrode/thermo -1.0 1.979 couple top 1.0 etypes on temp 500 100 7 # temp tau rng # to compare to regular constant potential, switch previous line to this: -#fix conp bot electrode/conp -1.0 1.979 couple top 1.0 etypes symm on +#fix conp bot electrode/conp -1.0 1.979 couple top 1.0 etypes on symm on fix nve electrolyte nve # note ionic liquid does not reach 500K immediately diff --git a/src/ELECTRODE/fix_electrode_conp.cpp b/src/ELECTRODE/fix_electrode_conp.cpp index 5157cf0286..0136217122 100644 --- a/src/ELECTRODE/fix_electrode_conp.cpp +++ b/src/ELECTRODE/fix_electrode_conp.cpp @@ -136,16 +136,6 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) : group_psi_var_styles.push_back(VarStyle::CONST); group_psi.push_back(utils::numeric(FLERR, arg[iarg], false, lmp)); } - } else if ((strncmp(arg[iarg], "symm", 4) == 0)) { - if (iarg + 2 > narg) error->all(FLERR, "Need yes/no command after symm keyword"); - char *symm_arg = arg[++iarg]; - if ((strcmp(symm_arg, "yes") == 0) || (strcmp(symm_arg, "on") == 0)) { - symm = true; - } else if ((strcmp(symm_arg, "no") == 0) || (strcmp(symm_arg, "off") == 0)) { - symm = false; - } else { - error->all(FLERR, "Invalid argument after symm keyword"); - } } else if ((strcmp(arg[iarg], "algo") == 0)) { if (!default_algo) error->one(FLERR, fmt::format("Algorithm can be set once, only")); default_algo = false; @@ -195,8 +185,6 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) : } else { error->all(FLERR, "Illegal fix electrode/conp command with read"); } - } else if ((strcmp(arg[iarg], "etypes") == 0)) { - etypes_neighlists = true; } else if ((strcmp(arg[iarg], "temp") == 0)) { if (iarg + 4 > narg) error->all(FLERR, "Need three arguments after temp command"); if (strcmp(this->style, "electrode/thermo") != 0) @@ -204,16 +192,13 @@ FixElectrodeConp::FixElectrodeConp(LAMMPS *lmp, int narg, char **arg) : thermo_temp = force->boltz / force->qe2f * utils::numeric(FLERR, arg[++iarg], false, lmp); thermo_time = utils::numeric(FLERR, arg[++iarg], false, lmp); thermo_init = utils::inumeric(FLERR, arg[++iarg], false, lmp); + // toggle parameters + } else if ((strcmp(arg[iarg], "etypes") == 0)) { + etypes_neighlists = utils::logical(FLERR, arg[++iarg], false, lmp); + } else if ((strncmp(arg[iarg], "symm", 4) == 0)) { + symm = utils::logical(FLERR, arg[++iarg], false, lmp); } else if ((strcmp(arg[iarg], "ffield") == 0)) { - if (iarg + 2 > narg) error->all(FLERR, "Need yes/no command after ffield keyword"); - char *ffield_arg = arg[++iarg]; - if ((strcmp(ffield_arg, "yes") == 0) || (strcmp(ffield_arg, "on") == 0)) { - ffield = true; - } else if ((strcmp(ffield_arg, "no") == 0) || (strcmp(ffield_arg, "off") == 0)) { - ffield = false; - } else { - error->all(FLERR, "Invalid argument after ffield keyword"); - } + ffield = utils::logical(FLERR, arg[++iarg], false, lmp); } else { error->all(FLERR, "Illegal fix electrode/conp command"); }