Changeset e1fe7e for src/Potentials


Ignore:
Timestamp:
Jun 27, 2014, 9:32:55 PM (11 years ago)
Author:
Frederik Heber <heber@…>
Branches:
Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable
Children:
550f2a
Parents:
16227a
git-author:
Frederik Heber <heber@…> (02/27/14 20:15:41)
git-committer:
Frederik Heber <heber@…> (06/27/14 21:32:55)
Message:

FunctionModel now uses list_of_arguments to split sequence of subsets of distances.

  • this fixes ambiguities with the set of distances: Imagine the distances within a water molecule as OH (A) and HH (B). We then may have a sequence of argument_t as AABAAB. And with the current implementation of CompoundPotential::splitUpArgumentsByModels() we would always choose the latter (and more complex) model. Hence, we make two calls to TriplePotential_Angle, instead of calls twice to PairPotential_Harmonic for A, one to PairPotential_Harmonic for B, and once to TriplePotential_Angle for AAB.
  • now, we new list looks like A,A,B,AAB where each tuple of distances can be uniquely associated with a specific potential.
  • changed signatures of EmpiricalPotential::operator(), ::derivative(), ::parameter_derivative(). This involved changing all of the current specific potentials and CompoundPotential.
  • TrainingData must discern between the InputVector_t (just all distances) and the FilteredInputVector_t (tuples of subsets of distances).
  • FunctionApproximation now has list_of_arguments_t as parameter to evaluate() and evaluate_derivative().
  • DOCU: docu change in TrainingData.
Location:
src/Potentials
Files:
27 edited

Legend:

Unmodified
Added
Removed
  • src/Potentials/CompoundPotential.cpp

    r16227a re1fe7e  
    228228{
    229229  arguments_by_model_t partial_args;
    230   // go through each model and have it filter out its arguments
     230  // go through each model and have it filter out its arguments, this already
     231  // returns a list of tuples associated with the specific model
    231232  for(models_t::const_iterator modeliter = models.begin();
    232233      modeliter != models.end(); ++modeliter) {
    233234    FunctionModel::filter_t filterfunction = (*modeliter)->getSpecificFilter();
    234     arguments_t tempargs = filterfunction(arguments);
     235    list_of_arguments_t tempargs = filterfunction(arguments);
    235236    // then split up all the bunches, too.
    236     arguments_t::const_iterator advanceiter = tempargs.begin();
    237     for (arguments_t::const_iterator argiter = tempargs.begin();
    238         argiter != tempargs.end(); argiter = advanceiter) {
    239       advanceiter += (*modeliter)->getSpecificArgumentCount();
     237    for (list_of_arguments_t::const_iterator argiter = tempargs.begin();
     238        argiter != tempargs.end(); ++argiter) {
     239      const arguments_t &args = *argiter;
    240240      partial_args.push_back(
    241241          std::make_pair(
    242242              *modeliter,
    243               arguments_t(argiter, advanceiter)
     243              args
    244244            )
    245245          );
     
    251251
    252252CompoundPotential::arguments_by_model_t CompoundPotential::splitUpArgumentsByModels(
    253     const arguments_t &arguments) const
     253    const list_of_arguments_t &listarguments) const
    254254{
    255255  arguments_by_model_t partial_args;
    256   arguments_t::const_iterator argiter = arguments.begin();
    257256  particletypes_per_model_t::const_iterator typesiter = particletypes_per_model.begin();
    258257  models_t::const_iterator modeliter = models.begin();
    259258
    260   // add constant model (which is always first model) with empty args if present
     259  /// add constant model (which is always first model) with empty args if present
    261260  if (typesiter->empty()) {
    262261    partial_args.push_back(
     
    266265    ++typesiter;
    267266  }
     267
    268268  // then check other models
    269   while (argiter != arguments.end()) {
     269  /// we only have to check whether the current model still matches or whether
     270  /// have to use the next model.
     271  for (list_of_arguments_t::const_iterator argiter = listarguments.begin();
     272      argiter != listarguments.end(); ++argiter) {
     273    const arguments_t &arguments = *argiter;
    270274    if (typesiter+1 != particletypes_per_model.end()) {
    271275      // check whether next argument bunch is for same model or different one
     
    275279
    276280      // we always expect N(N-1)/2 distances for N particle types
    277       arguments_t::const_iterator enditer = argiter+(types.size()*(types.size()-1)/2);
    278       arguments_t::const_iterator nextenditer = argiter+(nexttypes.size()*(nexttypes.size()-1)/2);
    279       arguments_t args(argiter, enditer);
    280       arguments_t nextargs(argiter, nextenditer);
    281       if (types.size() < nexttypes.size()) {
    282         if (areValidArguments(nexttypes, nextargs)) {
     281      // check first from sizes alone
     282      const size_t tuplesize = types.size()*(types.size()-1)/2;
     283      const size_t nexttuplesize = nexttypes.size()*(nexttypes.size()-1)/2;
     284      if ((tuplesize != nexttuplesize)) {
     285        if ((arguments.size() == tuplesize) &&  areValidArguments(types, arguments)) {
     286          // only former still matches, don't increment
     287          partial_args.push_back(
     288              std::make_pair(*modeliter, arguments)
     289              );
     290        } else if ((arguments.size() == nexttuplesize) &&  areValidArguments(nexttypes, arguments)) {
    283291          // latter matches, increment
    284292          ++typesiter;
    285293          partial_args.push_back(
    286               std::make_pair(*(++modeliter), arguments_t(argiter, nextenditer))
     294              std::make_pair(*(++modeliter), arguments)
    287295              );
    288           argiter = nextenditer;
    289         } else if (areValidArguments(types, args)) {
    290           // only former matches, don't increment
    291           partial_args.push_back(
    292               std::make_pair(*modeliter, arguments_t(argiter, enditer))
    293               );
    294           argiter = enditer;
    295         } else
     296        } else {
    296297          ASSERT(0,
    297               "CompoundPotential::splitUpArgumentsByModels() - neither type matches its argument bunch.");
    298       } else {
    299         if (areValidArguments(types, args)) {
    300           // only former matches, don't increment
    301           partial_args.push_back(
    302               std::make_pair(*modeliter, arguments_t(argiter, enditer))
    303               );
    304           argiter = enditer;
    305         } else if (areValidArguments(nexttypes, nextargs)) {
    306           // latter matches, increment
    307           ++typesiter;
    308           partial_args.push_back(
    309               std::make_pair(*(++modeliter), arguments_t(argiter, nextenditer))
    310               );
    311           argiter = nextenditer;
     298              "CompoundPotential::splitUpArgumentsByModels() - neither this model nor next model match (size) with current tuple.");
     299        }
     300      } else { // same size, now we have to check the types individually
     301        size_t encodeValidity = 0;
     302        encodeValidity += 1*areValidArguments(types, arguments);
     303        encodeValidity += 2*areValidArguments(nexttypes, arguments);
     304
     305        switch (encodeValidity) {
     306          case 1:
     307            // only former still matches, don't increment
     308            partial_args.push_back(
     309                std::make_pair(*modeliter, arguments)
     310                );
     311            break;
     312          case 2:
     313            ++typesiter;
     314            partial_args.push_back(
     315                std::make_pair(*(++modeliter), arguments)
     316                );
     317            break;
     318          case 0:
     319          case 3:
     320          default:
     321            ASSERT(0,
     322                "CompoundPotential::splitUpArgumentsByModels() - neither this model nor next model match (type) with current tuple.");
     323            break;
    312324        }
    313325      }
    314326    } else {
    315327      const SerializablePotential::ParticleTypes_t &types = *typesiter;
    316       // we always expect N(N-1)/2 distances for N particle types
    317       arguments_t::const_iterator enditer = argiter+(types.size()*(types.size()-1)/2);
    318       partial_args.push_back(
    319           std::make_pair(*modeliter, arguments_t(argiter, enditer))
    320           );
    321       argiter = enditer;
     328      if (areValidArguments(types, arguments)) {
     329        // only former matches, don't increment
     330        partial_args.push_back(
     331            std::make_pair(*modeliter, arguments)
     332            );
     333      } else {
     334        ASSERT(0,
     335            "CompoundPotential::splitUpArgumentsByModels() - last model does not match with current tuple.");
     336      }
    322337    }
    323338  }
     
    326341}
    327342
    328 CompoundPotential::results_t CompoundPotential::operator()(const arguments_t &arguments) const
     343CompoundPotential::results_t CompoundPotential::operator()(
     344    const list_of_arguments_t &listarguments) const
    329345{
    330346  /// first, we have to split up the given arguments
    331347  arguments_by_model_t partial_args =
    332       splitUpArgumentsByModels(arguments);
     348      splitUpArgumentsByModels(listarguments);
    333349  // print split up argument list for debugging
    334350  if (DoLog(4)) {
     
    348364      iter != partial_args.end(); ++iter) {
    349365    partial_results.push_back(
    350         (*iter->first)(iter->second)
     366        (*iter->first)(
     367            FunctionModel::list_of_arguments_t(1, iter->second))
    351368    );
    352369  }
     
    366383}
    367384
    368 CompoundPotential::results_t CompoundPotential::parameter_derivative(const arguments_t &arguments, const size_t index) const
     385CompoundPotential::results_t CompoundPotential::parameter_derivative(
     386    const list_of_arguments_t &listarguments,
     387    const size_t index) const
    369388{
    370389  // first, we have to split up the given arguments
    371390  arguments_by_model_t partial_args =
    372       splitUpArgumentsByModels(arguments);
     391      splitUpArgumentsByModels(listarguments);
    373392  // then, with each bunch of arguments, we call the specific model
    374393  // get parameter dimensions per model
     
    399418      const size_t indexbase = (iter == dimensions.begin()) ? 0 : *(iter-1);
    400419      CompoundPotential::results_t results =
    401           model->parameter_derivative(argiter->second, index-indexbase);
     420          model->parameter_derivative(
     421              FunctionModel::list_of_arguments_t(1, argiter->second), index-indexbase);
    402422
    403423      // either set results or add
     
    462482  // create initial returnfunction
    463483  FunctionModel::filter_t returnfunction =
    464       boost::bind(&Helpers::returnEmptyArguments);
     484      boost::bind(&Helpers::returnEmptyListArguments);
    465485
    466486  // every following fragments combines its arguments with the initial function
     
    468488      modeliter != models.end(); ++modeliter) {
    469489    returnfunction =
    470           boost::bind(&Extractors::concatenateArguments,
     490          boost::bind(&Extractors::concatenateListOfArguments,
    471491              boost::bind(returnfunction, _1),
    472492              boost::bind((*modeliter)->getSpecificFilter(), _1)
  • src/Potentials/CompoundPotential.hpp

    r16227a re1fe7e  
    9393  /** Evaluates the harmonic potential function for the given arguments.
    9494   *
    95    * @param arguments single distance
     95   * @param listarguments list of tuples of distances
    9696   * @return value of the potential function
    9797   */
    98   results_t operator()(const arguments_t &arguments) const;
     98  results_t operator()(const list_of_arguments_t &listarguments) const;
    9999
    100100  /** Evaluates the derivative of the function with the given \a arguments
    101101   * with respect to a specific parameter indicated by \a index.
    102102   *
    103    * \param arguments set of arguments as input variables to the function
     103   * \param arguments list of sets of arguments as input variables to the function
    104104   * \param index derivative of which parameter
    105105   * \return result vector containing the derivative with respect to the given
    106106   *         input
    107107   */
    108   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     108  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    109109
    110110  /** States whether lower and upper boundaries should be used to constraint
     
    153153  /** Helper function to split up concatenated arguments for operator() calls.
    154154   *
    155    * \param arguments arguments to split up
     155   * \param listarguments list of tuples of distances to assign to a model each
    156156   * \return vector of partial arguments with associated model
    157157   */
    158   arguments_by_model_t splitUpArgumentsByModels(const arguments_t &arguments) const;
     158  arguments_by_model_t splitUpArgumentsByModels(const list_of_arguments_t &listarguments) const;
    159159
    160160  /** Helper function to split up total list of arguments for operator() calls.
  • src/Potentials/EmpiricalPotential.hpp

    r16227a re1fe7e  
    6060   * parameters.
    6161   *
    62    * \param arguments set of arguments as input variables to the function
     62   * \param listarguments list of sets of arguments as input variables to the function
    6363   * \return result of the function
    6464   */
    65   virtual results_t operator()(const arguments_t &arguments) const=0;
     65  virtual results_t operator()(const list_of_arguments_t &listarguments) const=0;
    6666
    6767  /** Evaluates the derivative of the function with the given \a arguments
    6868   * for each component.
    6969   *
    70    * \param arguments set of arguments as input variables to the function
     70   * \param listarguments list of sets of arguments as input variables to the function
    7171   * \return result vector containing the derivative with respect to each
    7272   *         input comonent of the function
    7373   */
    74   virtual derivative_components_t derivative(const arguments_t &arguments) const=0;
     74  virtual derivative_components_t derivative(const list_of_arguments_t &listarguments) const=0;
    7575
    7676  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/ConstantPotential.cpp

    r16227a re1fe7e  
    112112ConstantPotential::results_t
    113113ConstantPotential::operator()(
    114     const arguments_t &arguments
     114    const list_of_arguments_t &listarguments
    115115    ) const
    116116{
    117   ASSERT( arguments.size() == 0,
    118       "ConstantPotential::operator() - requires no argument.");
    119   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    120       arguments, getParticleTypes()),
    121       "ConstantPotential::operator() - types don't match with ones in arguments.");
    122   const result_t result = params[energy_offset];
    123   return std::vector<result_t>(1, result);
     117  // directly set result as energy constant as independent of arg list
     118  result_t result = params[energy_offset];
     119  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     120      iter != listarguments.end(); ++iter) {
     121    const arguments_t &arguments = *iter;
     122    ASSERT( arguments.size() == 0,
     123        "ConstantPotential::operator() - requires no argument.");
     124    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     125        arguments, getParticleTypes()),
     126        "ConstantPotential::operator() - types don't match with ones in arguments.");
     127    result += 0.;
     128  }
     129  return results_t(1, result);
    124130}
    125131
    126132ConstantPotential::derivative_components_t
    127133ConstantPotential::derivative(
    128     const arguments_t &arguments
     134    const list_of_arguments_t &listarguments
    129135    ) const
    130136{
    131   ASSERT( arguments.size() == 0,
    132       "ConstantPotential::operator() - requires no argument.");
    133   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    134       arguments, getParticleTypes()),
    135       "ConstantPotential::operator() - types don't match with ones in arguments.");
    136   derivative_components_t result(1, 0.);
    137   return result;
     137  result_t result = 0.;
     138  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     139      iter != listarguments.end(); ++iter) {
     140    const arguments_t &arguments = *iter;
     141    ASSERT( arguments.size() == 0,
     142        "ConstantPotential::operator() - requires no argument.");
     143    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     144        arguments, getParticleTypes()),
     145        "ConstantPotential::operator() - types don't match with ones in arguments.");
     146    result += 0.;
     147  }
     148  return derivative_components_t(1, result);
    138149}
    139150
    140151ConstantPotential::results_t
    141152ConstantPotential::parameter_derivative(
    142     const arguments_t &arguments,
     153    const list_of_arguments_t &listarguments,
    143154    const size_t index
    144155    ) const
    145156{
    146   ASSERT( arguments.size() == 0,
    147       "ConstantPotential::parameter_derivative() - requires no argument.");
    148   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    149       arguments, getParticleTypes()),
    150       "ConstantPotential::operator() - types don't match with ones in arguments.");
    151   switch (index) {
    152     case energy_offset:
    153     {
    154       // Maple result: 1
    155       const result_t result = +1.;
    156       return std::vector<result_t>(1, result);
    157       break;
    158     }
    159     default:
    160       break;
     157  // Maple result: 1
     158  result_t result = 1.;   // energy_offset
     159  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     160      iter != listarguments.end(); ++iter) {
     161    const arguments_t &arguments = *iter;
     162    ASSERT( arguments.size() == 0,
     163        "ConstantPotential::parameter_derivative() - requires no argument.");
     164    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     165        arguments, getParticleTypes()),
     166        "ConstantPotential::operator() - types don't match with ones in arguments.");
     167//    switch (index) {
     168//      case energy_offset:
     169//      {
     170//        // Maple result: 1
     171//        result = +1.;
     172//        break;
     173//      }
     174//      default:
     175//        break;
     176//    }
    161177  }
    162   return std::vector<result_t>(1, 0.);
     178  return results_t(1, result);
    163179}
    164180
  • src/Potentials/Specifics/ConstantPotential.hpp

    r16227a re1fe7e  
    3636  friend class PotentialFactory;
    3737  // some repeated typedefs to avoid ambiguities
     38  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    3839  typedef FunctionModel::arguments_t arguments_t;
    3940  typedef FunctionModel::result_t result_t;
     
    9091  /** Evaluates the harmonic potential function for the given arguments.
    9192   *
    92    * @param arguments single distance
     93   * @param listarguments empty list
    9394   * @return value of the potential function
    9495   */
    95   results_t operator()(const arguments_t &arguments) const;
     96  results_t operator()(const list_of_arguments_t &listarguments) const;
    9697
    9798  /** Evaluates the derivative of the potential function.
    9899   *
    99    * @param arguments single distance
     100   * @param listarguments empty list
    100101   * @return vector with derivative with respect to the input degrees of freedom
    101102   */
    102   derivative_components_t derivative(const arguments_t &arguments) const;
     103  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    103104
    104105  /** Evaluates the derivative of the function with the given \a arguments
    105106   * with respect to a specific parameter indicated by \a index.
    106107   *
    107    * \param arguments set of arguments as input variables to the function
     108   * \param listarguments empty list
    108109   * \param index derivative of which parameter
    109110   * \return result vector containing the derivative with respect to the given
    110111   *         input
    111112   */
    112   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     113  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    113114
    114115  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/FourBodyPotential_Torsion.cpp

    r16227a re1fe7e  
    140140FourBodyPotential_Torsion::results_t
    141141FourBodyPotential_Torsion::operator()(
    142     const arguments_t &arguments
     142    const list_of_arguments_t &listarguments
    143143    ) const
    144144{
    145   ASSERT( arguments.size() == getSpecificArgumentCount(),
    146       "FourBodyPotential_Torsion::operator() - requires exactly three arguments.");
    147   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    148       arguments, getParticleTypes()),
    149       "FourBodyPotential_Torsion::operator() - types don't match with ones in arguments.");
    150   const argument_t &r_ij = arguments[0]; // 01
    151   const argument_t &r_ik = arguments[1]; // 02
    152   const argument_t &r_il = arguments[2]; // 03
    153   const argument_t &r_jk = arguments[3]; // 12
    154   const argument_t &r_jl = arguments[4]; // 13
    155   const argument_t &r_kl = arguments[5]; // 23
    156   const result_t result =
    157       params[spring_constant]
    158              * Helpers::pow( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance) - params[equilibrium_distance], 2 );
    159   return std::vector<result_t>(1, result);
     145  result_t result = 0.;
     146  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     147      iter != listarguments.end(); ++iter) {
     148    const arguments_t &arguments = *iter;
     149    ASSERT( arguments.size() == getSpecificArgumentCount(),
     150        "FourBodyPotential_Torsion::operator() - requires exactly three arguments.");
     151    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     152        arguments, getParticleTypes()),
     153        "FourBodyPotential_Torsion::operator() - types don't match with ones in arguments.");
     154    const argument_t &r_ij = arguments[0]; // 01
     155    const argument_t &r_ik = arguments[1]; // 02
     156    const argument_t &r_il = arguments[2]; // 03
     157    const argument_t &r_jk = arguments[3]; // 12
     158    const argument_t &r_jl = arguments[4]; // 13
     159    const argument_t &r_kl = arguments[5]; // 23
     160    result +=
     161        params[spring_constant]
     162               * Helpers::pow( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance)
     163                   - params[equilibrium_distance], 2 );
     164  }
     165  return results_t(1, result);
    160166}
    161167
    162168FourBodyPotential_Torsion::derivative_components_t
    163169FourBodyPotential_Torsion::derivative(
    164     const arguments_t &arguments
     170    const list_of_arguments_t &listarguments
    165171    ) const
    166172{
    167   ASSERT( arguments.size() == getSpecificArgumentCount(),
    168       "FourBodyPotential_Torsion::operator() - requires exactly three arguments.");
    169   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    170       arguments, getParticleTypes()),
    171       "FourBodyPotential_Torsion::operator() - types don't match with ones in arguments.");
    172   derivative_components_t result;
    173   const argument_t &r_ij = arguments[0]; // 01
    174   const argument_t &r_ik = arguments[1]; // 02
    175   const argument_t &r_il = arguments[2]; // 03
    176   const argument_t &r_jk = arguments[3]; // 12
    177   const argument_t &r_jl = arguments[4]; // 13
    178   const argument_t &r_kl = arguments[5]; // 23
    179   result.push_back( 2. * params[spring_constant] * ( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance) - params[equilibrium_distance]) );
    180   ASSERT( result.size() == 1,
    181       "FourBodyPotential_Torsion::operator() - we did not create exactly one component.");
    182   return result;
     173  result_t result = 0.;
     174  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     175      iter != listarguments.end(); ++iter) {
     176    const arguments_t &arguments = *iter;
     177    ASSERT( arguments.size() == getSpecificArgumentCount(),
     178        "FourBodyPotential_Torsion::operator() - requires exactly three arguments.");
     179    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     180        arguments, getParticleTypes()),
     181        "FourBodyPotential_Torsion::operator() - types don't match with ones in arguments.");
     182    const argument_t &r_ij = arguments[0]; // 01
     183    const argument_t &r_ik = arguments[1]; // 02
     184    const argument_t &r_il = arguments[2]; // 03
     185    const argument_t &r_jk = arguments[3]; // 12
     186    const argument_t &r_jl = arguments[4]; // 13
     187    const argument_t &r_kl = arguments[5]; // 23
     188    result +=
     189        2. * params[spring_constant] *
     190        ( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance)
     191            - params[equilibrium_distance]);
     192  }
     193  return derivative_components_t(1, result);
    183194}
    184195
    185196FourBodyPotential_Torsion::results_t
    186197FourBodyPotential_Torsion::parameter_derivative(
    187     const arguments_t &arguments,
     198    const list_of_arguments_t &listarguments,
    188199    const size_t index
    189200    ) const
    190201{
    191   ASSERT( arguments.size() == getSpecificArgumentCount(),
    192       "FourBodyPotential_Torsion::parameter_derivative() - requires exactly three arguments.");
    193   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    194       arguments, getParticleTypes()),
    195       "FourBodyPotential_Torsion::operator() - types don't match with ones in arguments.");
    196   const argument_t &r_ij = arguments[0]; // 01
    197   const argument_t &r_ik = arguments[1]; // 02
    198   const argument_t &r_il = arguments[2]; // 03
    199   const argument_t &r_jk = arguments[3]; // 12
    200   const argument_t &r_jl = arguments[4]; // 13
    201   const argument_t &r_kl = arguments[5]; // 23
    202   switch (index) {
    203     case spring_constant:
    204     {
    205       const result_t result =
    206                  Helpers::pow( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance) - params[equilibrium_distance], 2 );
    207       return std::vector<result_t>(1, result);
    208       break;
     202  result_t result = 0.;
     203  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     204      iter != listarguments.end(); ++iter) {
     205    const arguments_t &arguments = *iter;
     206    ASSERT( arguments.size() == getSpecificArgumentCount(),
     207        "FourBodyPotential_Torsion::parameter_derivative() - requires exactly three arguments.");
     208    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     209        arguments, getParticleTypes()),
     210        "FourBodyPotential_Torsion::operator() - types don't match with ones in arguments.");
     211    const argument_t &r_ij = arguments[0]; // 01
     212    const argument_t &r_ik = arguments[1]; // 02
     213    const argument_t &r_il = arguments[2]; // 03
     214    const argument_t &r_jk = arguments[3]; // 12
     215    const argument_t &r_jl = arguments[4]; // 13
     216    const argument_t &r_kl = arguments[5]; // 23
     217    switch (index) {
     218      case spring_constant:
     219      {
     220        result +=
     221            Helpers::pow( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance) - params[equilibrium_distance], 2 );
     222        break;
     223      }
     224      case equilibrium_distance:
     225      {
     226        result +=
     227            -2. * params[spring_constant]
     228                   * ( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance) - params[equilibrium_distance]);
     229        break;
     230      }
     231      default:
     232        ASSERT(0, "FourBodyPotential_Torsion::parameter_derivative() - derivative to unknown parameter desired.");
     233        break;
    209234    }
    210     case equilibrium_distance:
    211     {
    212       const result_t result =
    213           -2. * params[spring_constant]
    214                  * ( function_theta(r_ij.distance, r_ik.distance, r_il.distance, r_jk.distance, r_jl.distance, r_kl.distance) - params[equilibrium_distance]);
    215       return std::vector<result_t>(1, result);
    216       break;
    217     }
    218     default:
    219       ASSERT(0, "FourBodyPotential_Torsion::parameter_derivative() - derivative to unknown parameter desired.");
    220       break;
    221235  }
    222   return std::vector<result_t>(1);
     236  return results_t(1, result);
    223237}
    224238
  • src/Potentials/Specifics/FourBodyPotential_Torsion.hpp

    r16227a re1fe7e  
    3737  friend class PotentialFactory;
    3838  // some repeated typedefs to avoid ambiguities
     39  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    3940  typedef FunctionModel::arguments_t arguments_t;
    4041  typedef FunctionModel::result_t result_t;
     
    9293  /** Evaluates the harmonic potential function for the given arguments.
    9394   *
    94    * @param arguments single distance
     95   * @param listarguments list of six distances
    9596   * @return value of the potential function
    9697   */
    97   results_t operator()(const arguments_t &arguments) const;
     98  results_t operator()(const list_of_arguments_t &listarguments) const;
    9899
    99100  /** Evaluates the derivative of the potential function.
    100101   *
    101    * @param arguments single distance
     102   * @param listarguments list of six distances
    102103   * @return vector with derivative with respect to the input degrees of freedom
    103104   */
    104   derivative_components_t derivative(const arguments_t &arguments) const;
     105  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    105106
    106107  /** Evaluates the derivative of the function with the given \a arguments
    107108   * with respect to a specific parameter indicated by \a index.
    108109   *
    109    * \param arguments set of arguments as input variables to the function
     110   * \param listarguments list of six distances
    110111   * \param index derivative of which parameter
    111112   * \return result vector containing the derivative with respect to the given
    112113   *         input
    113114   */
    114   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     115  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    115116
    116117  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/ManyBodyPotential_Tersoff.cpp

    r16227a re1fe7e  
    184184ManyBodyPotential_Tersoff::results_t
    185185ManyBodyPotential_Tersoff::operator()(
    186     const arguments_t &arguments
     186    const list_of_arguments_t &listarguments
    187187    ) const
    188188{
    189189//  Info info(__func__);
    190   double result = 0.;
    191   for(arguments_t::const_iterator argiter = arguments.begin();
    192       argiter != arguments.end();
    193       ++argiter) {
    194     const argument_t &r_ij = *argiter;
    195     ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    196         arguments_t(1, r_ij), getParticleTypes()),
    197         "ManyBodyPotential_Tersoff::operator() - types don't match with ones in arguments.");
    198  
    199     const double cutoff = function_cutoff(r_ij.distance);
    200     const double temp = (cutoff == 0.) ?
    201         0. :
    202         cutoff * (
    203             function_prefactor(
    204                 alpha,
    205                 function_eta(r_ij))
     190  result_t result = 0.;
     191  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     192      iter != listarguments.end(); ++iter) {
     193    const arguments_t &arguments = *iter;
     194    for(arguments_t::const_iterator argiter = arguments.begin();
     195        argiter != arguments.end();
     196        ++argiter) {
     197      const argument_t &r_ij = *argiter;
     198      ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     199          arguments_t(1, r_ij), getParticleTypes()),
     200          "ManyBodyPotential_Tersoff::operator() - types don't match with ones in arguments.");
     201
     202      const double cutoff = function_cutoff(r_ij.distance);
     203      const double temp = (cutoff == 0.) ?
     204          0. :
     205          cutoff * (
     206              function_prefactor(
     207                  alpha,
     208                  function_eta(r_ij))
     209              * function_smoother(
     210                  params[A],
     211                  params[lambda],
     212                  r_ij.distance)
     213              +
     214              function_prefactor(
     215                  params[beta],
     216                  function_zeta(r_ij))
     217              * function_smoother(
     218                  -params[B],
     219                  params[mu],
     220                  r_ij.distance)
     221          );
     222      result += temp;
     223    }
     224  }
     225//  LOG(2, "DEBUG: operator()(" << r_ij.distance << ") = " << result);
     226  return results_t(1, result);
     227}
     228
     229ManyBodyPotential_Tersoff::derivative_components_t
     230ManyBodyPotential_Tersoff::derivative(
     231    const list_of_arguments_t &listarguments
     232    ) const
     233{
     234//  Info info(__func__);
     235  return derivative_components_t();
     236}
     237
     238ManyBodyPotential_Tersoff::results_t
     239ManyBodyPotential_Tersoff::parameter_derivative(
     240    const list_of_arguments_t &listarguments,
     241    const size_t index
     242    ) const
     243{
     244//  Info info(__func__);
     245//  ASSERT( arguments.size() == 1,
     246//      "ManyBodyPotential_Tersoff::parameter_derivative() - requires exactly one argument.");
     247
     248  result_t result = 0.;
     249  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     250      iter != listarguments.end(); ++iter) {
     251    const arguments_t &arguments = *iter;
     252    for(arguments_t::const_iterator argiter = arguments.begin();
     253        argiter != arguments.end();
     254        ++argiter) {
     255      const argument_t &r_ij = *argiter;
     256      ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     257          arguments_t(1, r_ij), getParticleTypes()),
     258          "ManyBodyPotential_Tersoff::operator() - types don't match with ones in arguments.");
     259
     260    switch (index) {
     261//      case R:
     262//      {
     263//        result += 0.;
     264//        break;
     265//      }
     266//      case S:
     267//      {
     268//        result += 0.;
     269//        break;
     270//      }
     271      case A:
     272      {
     273        const double cutoff = function_cutoff(r_ij.distance);
     274        result += (cutoff == 0.) ?
     275            0. :
     276            cutoff *
     277                function_prefactor(
     278                    alpha,
     279                    function_eta(r_ij))
     280                * function_smoother(
     281                    1.,
     282                    params[lambda],
     283                    r_ij.distance);
     284//            cutoff * function_prefactor(
     285//                alpha,
     286//                function_eta(r_ij))
     287//            * function_smoother(
     288//                1.,
     289//                params[mu],
     290//                r_ij.distance);
     291        break;
     292      }
     293      case B:
     294      {
     295        const double cutoff = function_cutoff(r_ij.distance);
     296        result += (cutoff == 0.) ?
     297            0. :
     298            cutoff * function_prefactor(
     299                params[beta],
     300                function_zeta(r_ij))
    206301            * function_smoother(
    207                 params[A],
    208                 params[lambda],
    209                 r_ij.distance)
    210             +
     302                -1.,
     303                params[mu],
     304                r_ij.distance);
     305//            cutoff * function_prefactor(
     306//                beta,
     307//                function_zeta(r_ij))
     308//            * function_smoother(
     309//                -params[B],
     310//                params[mu],
     311//                r_ij.distance)/params[B];
     312        break;
     313      }
     314      case lambda:
     315      {
     316        const double cutoff = function_cutoff(r_ij.distance);
     317        result += (cutoff == 0.) ?
     318            0. :
     319            -r_ij.distance * cutoff *
     320                function_prefactor(
     321                    alpha,
     322                    function_eta(r_ij))
     323                * function_smoother(
     324                    params[A],
     325                    params[lambda],
     326                    r_ij.distance);
     327        break;
     328      }
     329      case mu:
     330      {
     331        const double cutoff = function_cutoff(r_ij.distance);
     332        result += (cutoff == 0.) ?
     333            0. :
     334            -r_ij.distance * cutoff *(
    211335            function_prefactor(
    212336                params[beta],
     
    217341                r_ij.distance)
    218342        );
    219     result += temp;
    220   }
    221 //  LOG(2, "DEBUG: operator()(" << r_ij.distance << ") = " << result);
    222   return std::vector<result_t>(1, result);
    223 }
    224 
    225 ManyBodyPotential_Tersoff::derivative_components_t
    226 ManyBodyPotential_Tersoff::derivative(
    227     const arguments_t &arguments
    228     ) const
    229 {
    230 //  Info info(__func__);
    231   return ManyBodyPotential_Tersoff::derivative_components_t();
    232 }
    233 
    234 ManyBodyPotential_Tersoff::results_t
    235 ManyBodyPotential_Tersoff::parameter_derivative(
    236     const arguments_t &arguments,
    237     const size_t index
    238     ) const
    239 {
    240 //  Info info(__func__);
    241 //  ASSERT( arguments.size() == 1,
    242 //      "ManyBodyPotential_Tersoff::parameter_derivative() - requires exactly one argument.");
    243 
    244   double result = 0.;
    245   for(arguments_t::const_iterator argiter = arguments.begin();
    246       argiter != arguments.end();
    247       ++argiter) {
    248     const argument_t &r_ij = *argiter;
    249     ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    250         arguments_t(1, r_ij), getParticleTypes()),
    251         "ManyBodyPotential_Tersoff::operator() - types don't match with ones in arguments.");
    252 
    253   switch (index) {
    254 //    case R:
    255 //    {
    256 //      result += 0.;
    257 //      break;
    258 //    }
    259 //    case S:
    260 //    {
    261 //      result += 0.;
    262 //      break;
    263 //    }
    264     case A:
    265     {
    266       const double cutoff = function_cutoff(r_ij.distance);
    267       result += (cutoff == 0.) ?
    268           0. :
    269           cutoff *
    270               function_prefactor(
    271                   alpha,
    272                   function_eta(r_ij))
    273               * function_smoother(
    274                   1.,
    275                   params[lambda],
    276                   r_ij.distance);
    277 //          cutoff * function_prefactor(
    278 //              alpha,
    279 //              function_eta(r_ij))
    280 //          * function_smoother(
    281 //              1.,
    282 //              params[mu],
    283 //              r_ij.distance);
    284       break;
     343        break;
     344      }
     345//      case lambda3:
     346//      {
     347//        result += 0.;
     348//        break;
     349//      }
     350//      case alpha:
     351//      {
     352//        const double temp =
     353//            pow(alpha*function_eta(r_ij), params[n]);
     354//        const double cutoff = function_cutoff(r_ij.distance);
     355//        result += (cutoff == 0.) || (alpha == 0. )?
     356//            0. :
     357//            function_smoother(
     358//                params[A],
     359//                params[lambda],
     360//                r_ij.distance)
     361//            * (-.5) * alpha * (temp/alpha)
     362//            / (1. + temp)
     363//            ;
     364//        break;
     365//      }
     366//      case chi:
     367//      {
     368//        result += 0.;
     369//        break;
     370//      }
     371//      case omega:
     372//      {
     373//        result += 0.;
     374//        break;
     375//      }
     376      case beta:
     377      {
     378        const double temp =
     379            pow(params[beta]*function_zeta(r_ij), params[n]);
     380        const double cutoff = function_cutoff(r_ij.distance);
     381        result += (cutoff == 0.) || (params[beta] == 0. )?
     382            0. : cutoff *
     383            function_smoother(
     384                -params[B],
     385                params[mu],
     386                r_ij.distance)
     387            * (-.5)
     388            * function_prefactor(
     389                params[beta],
     390                function_zeta(r_ij))
     391            * (temp/params[beta])
     392            / (1. + temp)
     393            ;
     394        break;
     395      }
     396      case n:
     397      {
     398        const double zeta = function_zeta(r_ij);
     399        const double temp = pow( params[beta]*zeta , params[n]);
     400        const double cutoff = function_cutoff(r_ij.distance);
     401        const double tempres = ((cutoff == 0.) || (zeta == 0.)) ? // zeta must be caught if zero due to log
     402            0. : .5 * cutoff *
     403            function_smoother(
     404                -params[B],
     405                params[mu],
     406                r_ij.distance)
     407            * function_prefactor(
     408                params[beta],
     409                function_zeta(r_ij))
     410            * ( log(1.+temp)/(params[n]*params[n]) - temp
     411                * (log(function_zeta(r_ij)) + log(params[beta]))
     412                /(params[n]*(1.+temp)))
     413            ;
     414//        if (tempres != tempres)
     415//    LOG(2, "DEBUG: tempres is NaN.");
     416//        LOG(2, "DEBUG: Adding " << tempres << " for p.d. w.r.t n, temp=" << temp << ", cutoff=" << cutoff);
     417        result += tempres;
     418        break;
     419      }
     420      case c:
     421      {
     422        const double zeta = function_zeta(r_ij);
     423        if (zeta == 0.)
     424          break;
     425        const double temp =
     426            pow(zeta, params[n]-1.) * pow(params[beta],params[n]);
     427        const double cutoff = function_cutoff(r_ij.distance);
     428        const double tempres = (cutoff == 0.) ?
     429            0. : cutoff *
     430            function_smoother(
     431                -params[B],
     432                params[mu],
     433                r_ij.distance)
     434            * function_prefactor(
     435                params[beta],
     436                zeta)
     437             * (-1.) * temp / (1.+temp*zeta);
     438        double factor = function_derivative_c(r_ij);
     439        result += tempres*factor;
     440        if (result != result)
     441          ELOG(1, "result is NaN, zeta=" << zeta << ", temp=" << temp << ", cutoff=" << cutoff << ", tempres=" << tempres << ", factor=" << factor);
     442        break;
     443      }
     444      case d:
     445      {
     446        const double zeta = function_zeta(r_ij);
     447        const double temp =
     448            pow(zeta, params[n]-1.) * pow(params[beta],params[n]);
     449        const double cutoff = function_cutoff(r_ij.distance);
     450        const double tempres = (cutoff == 0.) ?
     451            0. : cutoff *
     452            function_smoother(
     453                -params[B],
     454                params[mu],
     455                r_ij.distance)
     456            * function_prefactor(
     457                params[beta],
     458                zeta)
     459             * (-1.) * temp / (1.+temp*zeta);
     460        double factor = function_derivative_d(r_ij);
     461        result += tempres*factor;
     462        if (result != result)
     463          ELOG(1, "result is NaN, zeta=" << zeta << ", temp=" << temp << ", cutoff=" << cutoff << ", tempres=" << tempres << ", factor=" << factor);
     464        break;
     465      }
     466      case h:
     467      {
     468        const double zeta = function_zeta(r_ij);
     469        const double temp =
     470            pow(zeta, params[n]-1.) * pow(params[beta],params[n]);
     471        const double cutoff = function_cutoff(r_ij.distance);
     472        const double tempres = (cutoff == 0.) ?
     473            0. : cutoff *
     474            function_smoother(
     475                -params[B],
     476                params[mu],
     477                r_ij.distance)
     478            * function_prefactor(
     479                params[beta],
     480                zeta)
     481             * (-1.) * temp / (1.+temp*zeta);
     482        double factor = function_derivative_h(r_ij);
     483        result += tempres*factor;
     484        if (result != result)
     485          ELOG(1, "result is NaN, zeta=" << zeta << ", temp=" << temp << ", cutoff=" << cutoff << ", tempres=" << tempres << ", factor=" << factor);
     486        break;
     487      }
     488      default:
     489        ASSERT(0, "ManyBodyPotential_Tersoff::parameter_derivative() - derivative to unknown parameter desired.");
     490        break;
    285491    }
    286     case B:
    287     {
    288       const double cutoff = function_cutoff(r_ij.distance);
    289       result += (cutoff == 0.) ?
    290           0. :
    291           cutoff * function_prefactor(
    292               params[beta],
    293               function_zeta(r_ij))
    294           * function_smoother(
    295               -1.,
    296               params[mu],
    297               r_ij.distance);
    298 //          cutoff * function_prefactor(
    299 //              beta,
    300 //              function_zeta(r_ij))
    301 //          * function_smoother(
    302 //              -params[B],
    303 //              params[mu],
    304 //              r_ij.distance)/params[B];
    305       break;
     492    if (result != result)
     493      ELOG(1, "result is NaN.");
    306494    }
    307     case lambda:
    308     {
    309       const double cutoff = function_cutoff(r_ij.distance);
    310       result += (cutoff == 0.) ?
    311           0. :
    312           -r_ij.distance * cutoff *
    313               function_prefactor(
    314                   alpha,
    315                   function_eta(r_ij))
    316               * function_smoother(
    317                   params[A],
    318                   params[lambda],
    319                   r_ij.distance);
    320       break;
    321     }
    322     case mu:
    323     {
    324       const double cutoff = function_cutoff(r_ij.distance);
    325       result += (cutoff == 0.) ?
    326           0. :
    327           -r_ij.distance * cutoff *(
    328           function_prefactor(
    329               params[beta],
    330               function_zeta(r_ij))
    331           * function_smoother(
    332               -params[B],
    333               params[mu],
    334               r_ij.distance)
    335       );
    336       break;
    337     }
    338 //    case lambda3:
    339 //    {
    340 //      result += 0.;
    341 //      break;
    342 //    }
    343 //    case alpha:
    344 //    {
    345 //      const double temp =
    346 //          pow(alpha*function_eta(r_ij), params[n]);
    347 //      const double cutoff = function_cutoff(r_ij.distance);
    348 //      result += (cutoff == 0.) || (alpha == 0. )?
    349 //          0. :
    350 //          function_smoother(
    351 //              params[A],
    352 //              params[lambda],
    353 //              r_ij.distance)
    354 //          * (-.5) * alpha * (temp/alpha)
    355 //          / (1. + temp)
    356 //          ;
    357 //      break;
    358 //    }
    359 //    case chi:
    360 //    {
    361 //      result += 0.;
    362 //      break;
    363 //    }
    364 //    case omega:
    365 //    {
    366 //      result += 0.;
    367 //      break;
    368 //    }
    369     case beta:
    370     {
    371       const double temp =
    372           pow(params[beta]*function_zeta(r_ij), params[n]);
    373       const double cutoff = function_cutoff(r_ij.distance);
    374       result += (cutoff == 0.) || (params[beta] == 0. )?
    375           0. : cutoff *
    376           function_smoother(
    377               -params[B],
    378               params[mu],
    379               r_ij.distance)
    380           * (-.5)
    381           * function_prefactor(
    382               params[beta],
    383               function_zeta(r_ij))
    384           * (temp/params[beta])
    385           / (1. + temp)
    386           ;
    387       break;
    388     }
    389     case n:
    390     {
    391       const double zeta = function_zeta(r_ij);
    392       const double temp = pow( params[beta]*zeta , params[n]);
    393       const double cutoff = function_cutoff(r_ij.distance);
    394       const double tempres = ((cutoff == 0.) || (zeta == 0.)) ? // zeta must be caught if zero due to log
    395           0. : .5 * cutoff *
    396           function_smoother(
    397               -params[B],
    398               params[mu],
    399               r_ij.distance)
    400           * function_prefactor(
    401               params[beta],
    402               function_zeta(r_ij))
    403           * ( log(1.+temp)/(params[n]*params[n]) - temp
    404               * (log(function_zeta(r_ij)) + log(params[beta]))
    405               /(params[n]*(1.+temp)))
    406           ;
    407 //      if (tempres != tempres)
    408 //      LOG(2, "DEBUG: tempres is NaN.");
    409 //      LOG(2, "DEBUG: Adding " << tempres << " for p.d. w.r.t n, temp=" << temp << ", cutoff=" << cutoff);
    410       result += tempres;
    411       break;
    412     }
    413     case c:
    414     {
    415       const double zeta = function_zeta(r_ij);
    416       if (zeta == 0.)
    417         break;
    418       const double temp =
    419           pow(zeta, params[n]-1.) * pow(params[beta],params[n]);
    420       const double cutoff = function_cutoff(r_ij.distance);
    421       const double tempres = (cutoff == 0.) ?
    422           0. : cutoff *
    423           function_smoother(
    424               -params[B],
    425               params[mu],
    426               r_ij.distance)
    427           * function_prefactor(
    428               params[beta],
    429               zeta)
    430            * (-1.) * temp / (1.+temp*zeta);
    431       double factor = function_derivative_c(r_ij);
    432       result += tempres*factor;
    433       if (result != result)
    434         ELOG(1, "result is NaN, zeta=" << zeta << ", temp=" << temp << ", cutoff=" << cutoff << ", tempres=" << tempres << ", factor=" << factor);
    435       break;
    436     }
    437     case d:
    438     {
    439       const double zeta = function_zeta(r_ij);
    440       const double temp =
    441           pow(zeta, params[n]-1.) * pow(params[beta],params[n]);
    442       const double cutoff = function_cutoff(r_ij.distance);
    443       const double tempres = (cutoff == 0.) ?
    444           0. : cutoff *
    445           function_smoother(
    446               -params[B],
    447               params[mu],
    448               r_ij.distance)
    449           * function_prefactor(
    450               params[beta],
    451               zeta)
    452            * (-1.) * temp / (1.+temp*zeta);
    453       double factor = function_derivative_d(r_ij);
    454       result += tempres*factor;
    455       if (result != result)
    456         ELOG(1, "result is NaN, zeta=" << zeta << ", temp=" << temp << ", cutoff=" << cutoff << ", tempres=" << tempres << ", factor=" << factor);
    457       break;
    458     }
    459     case h:
    460     {
    461       const double zeta = function_zeta(r_ij);
    462       const double temp =
    463           pow(zeta, params[n]-1.) * pow(params[beta],params[n]);
    464       const double cutoff = function_cutoff(r_ij.distance);
    465       const double tempres = (cutoff == 0.) ?
    466           0. : cutoff *
    467           function_smoother(
    468               -params[B],
    469               params[mu],
    470               r_ij.distance)
    471           * function_prefactor(
    472               params[beta],
    473               zeta)
    474            * (-1.) * temp / (1.+temp*zeta);
    475       double factor = function_derivative_h(r_ij);
    476       result += tempres*factor;
    477       if (result != result)
    478         ELOG(1, "result is NaN, zeta=" << zeta << ", temp=" << temp << ", cutoff=" << cutoff << ", tempres=" << tempres << ", factor=" << factor);
    479       break;
    480     }
    481     default:
    482       ASSERT(0, "ManyBodyPotential_Tersoff::parameter_derivative() - derivative to unknown parameter desired.");
    483       break;
    484   }
    485   if (result != result)
    486     ELOG(1, "result is NaN.");
    487495  }
    488496  return results_t(1,-result);
  • src/Potentials/Specifics/ManyBodyPotential_Tersoff.hpp

    r16227a re1fe7e  
    3838  friend class PotentialFactory;
    3939  // some repeated typedefs to avoid ambiguities
     40  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    4041  typedef FunctionModel::arguments_t arguments_t;
    4142  typedef FunctionModel::result_t result_t;
     
    106107  /** Evaluates the Tersoff potential for the given arguments.
    107108   *
    108    * @param arguments single distance
     109   * @param listarguments list of distances
    109110   * @return value of the potential function
    110111   */
    111   results_t operator()(const arguments_t &arguments) const;
     112  results_t operator()(const list_of_arguments_t &listarguments) const;
    112113
    113114  /** Evaluates the derivative of the Tersoff potential with respect to the
    114115   * input variables.
    115116   *
    116    * @param arguments single distance
     117   * @param listarguments list of distances
    117118   * @return vector with components of the derivative
    118119   */
    119   derivative_components_t derivative(const arguments_t &arguments) const;
     120  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    120121
    121122  /** Evaluates the derivative of the function with the given \a arguments
    122123   * with respect to a specific parameter indicated by \a index.
    123124   *
    124    * \param arguments set of arguments as input variables to the function
     125   * \param listarguments list of distances
    125126   * \param index derivative of which parameter
    126127   * \return result vector containing the derivative with respect to the given
    127128   *         input
    128129   */
    129   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     130  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    130131
    131132  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/PairPotential_Harmonic.cpp

    r16227a re1fe7e  
    117117PairPotential_Harmonic::results_t
    118118PairPotential_Harmonic::operator()(
    119     const arguments_t &arguments
     119    const list_of_arguments_t &listarguments
    120120    ) const
    121121{
    122   ASSERT( arguments.size() == 1,
    123       "PairPotential_Harmonic::operator() - requires exactly one argument.");
    124   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    125       arguments, getParticleTypes()),
    126       "PairPotential_Harmonic::operator() - types don't match with ones in arguments.");
    127   const argument_t &r_ij = arguments[0];
    128   const result_t result =
    129       params[spring_constant]
    130              * Helpers::pow( r_ij.distance - params[equilibrium_distance], 2 );
    131   return std::vector<result_t>(1, result);
     122  result_t result = 0.;
     123  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     124      iter != listarguments.end(); ++iter) {
     125    const arguments_t &arguments = *iter;
     126    ASSERT( arguments.size() == 1,
     127        "PairPotential_Harmonic::operator() - requires exactly one argument.");
     128    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     129        arguments, getParticleTypes()),
     130        "PairPotential_Harmonic::operator() - types don't match with ones in arguments.");
     131    const argument_t &r_ij = arguments[0];
     132    result +=
     133        params[spring_constant]
     134               * Helpers::pow( r_ij.distance - params[equilibrium_distance], 2 );
     135  }
     136  return results_t(1, result);
    132137}
    133138
    134139PairPotential_Harmonic::derivative_components_t
    135140PairPotential_Harmonic::derivative(
    136     const arguments_t &arguments
     141    const list_of_arguments_t &listarguments
    137142    ) const
    138143{
    139   ASSERT( arguments.size() == 1,
    140       "PairPotential_Harmonic::operator() - requires exactly one argument.");
    141   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    142       arguments, getParticleTypes()),
    143       "PairPotential_Harmonic::operator() - types don't match with ones in arguments.");
    144   derivative_components_t result;
    145   const argument_t &r_ij = arguments[0];
    146   result.push_back( 2. * params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]) );
    147   ASSERT( result.size() == 1,
    148       "PairPotential_Harmonic::operator() - we did not create exactly one component.");
    149   return result;
     144  result_t result = 0.;
     145  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     146      iter != listarguments.end(); ++iter) {
     147    const arguments_t &arguments = *iter;
     148      ASSERT( arguments.size() == 1,
     149          "PairPotential_Harmonic::operator() - requires exactly one argument.");
     150      ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     151          arguments, getParticleTypes()),
     152          "PairPotential_Harmonic::operator() - types don't match with ones in arguments.");
     153      const argument_t &r_ij = arguments[0];
     154      result +=
     155          2. * params[spring_constant] *
     156          ( r_ij.distance - params[equilibrium_distance]);
     157  }
     158  return derivative_components_t(1, result);
    150159}
    151160
    152161PairPotential_Harmonic::results_t
    153162PairPotential_Harmonic::parameter_derivative(
    154     const arguments_t &arguments,
     163    const list_of_arguments_t &listarguments,
    155164    const size_t index
    156165    ) const
    157166{
    158   ASSERT( arguments.size() == 1,
    159       "PairPotential_Harmonic::parameter_derivative() - requires exactly one argument.");
    160   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    161       arguments, getParticleTypes()),
    162       "PairPotential_Harmonic::operator() - types don't match with ones in arguments.");
    163   const argument_t &r_ij = arguments[0];
    164   switch (index) {
    165     case spring_constant:
    166     {
    167       const result_t result =
    168                  Helpers::pow( r_ij.distance - params[equilibrium_distance], 2 );
    169       return std::vector<result_t>(1, result);
    170       break;
     167  result_t result = 0.;
     168  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     169      iter != listarguments.end(); ++iter) {
     170    const arguments_t &arguments = *iter;
     171    ASSERT( arguments.size() == 1,
     172        "PairPotential_Harmonic::parameter_derivative() - requires exactly one argument.");
     173    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     174        arguments, getParticleTypes()),
     175        "PairPotential_Harmonic::operator() - types don't match with ones in arguments.");
     176    const argument_t &r_ij = arguments[0];
     177    switch (index) {
     178      case spring_constant:
     179      {
     180        result +=
     181            Helpers::pow( r_ij.distance - params[equilibrium_distance], 2 );
     182        break;
     183      }
     184      case equilibrium_distance:
     185      {
     186        result +=
     187            -2. * params[spring_constant]
     188                   * ( r_ij.distance - params[equilibrium_distance]);
     189        break;
     190      }
     191      default:
     192        ASSERT(0, "PairPotential_Harmonic::parameter_derivative() - derivative to unknown parameter desired.");
     193        break;
    171194    }
    172     case equilibrium_distance:
    173     {
    174       const result_t result =
    175           -2. * params[spring_constant]
    176                  * ( r_ij.distance - params[equilibrium_distance]);
    177       return std::vector<result_t>(1, result);
    178       break;
    179     }
    180     default:
    181       ASSERT(0, "PairPotential_Harmonic::parameter_derivative() - derivative to unknown parameter desired.");
    182       break;
    183195  }
    184 
    185   return PairPotential_Harmonic::results_t(1, 0.);
     196  return results_t(1, result);
    186197}
    187198
  • src/Potentials/Specifics/PairPotential_Harmonic.hpp

    r16227a re1fe7e  
    3535  friend class PotentialFactory;
    3636  // some repeated typedefs to avoid ambiguities
     37  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    3738  typedef FunctionModel::arguments_t arguments_t;
    3839  typedef FunctionModel::result_t result_t;
     
    9091  /** Evaluates the harmonic potential function for the given arguments.
    9192   *
    92    * @param arguments single distance
     93   * @param listarguments list of single distances
    9394   * @return value of the potential function
    9495   */
    95   results_t operator()(const arguments_t &arguments) const;
     96  results_t operator()(const list_of_arguments_t &listarguments) const;
    9697
    9798  /** Evaluates the derivative of the potential function.
    9899   *
    99    * @param arguments single distance
     100   * @param listarguments list of single distances
    100101   * @return vector with derivative with respect to the input degrees of freedom
    101102   */
    102   derivative_components_t derivative(const arguments_t &arguments) const;
     103  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    103104
    104105  /** Evaluates the derivative of the function with the given \a arguments
    105106   * with respect to a specific parameter indicated by \a index.
    106107   *
    107    * \param arguments set of arguments as input variables to the function
     108   * \param listarguments list of single distances
    108109   * \param index derivative of which parameter
    109110   * \return result vector containing the derivative with respect to the given
    110111   *         input
    111112   */
    112   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     113  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    113114
    114115  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/PairPotential_LennardJones.cpp

    r16227a re1fe7e  
    123123PairPotential_LennardJones::results_t
    124124PairPotential_LennardJones::operator()(
    125     const arguments_t &arguments
     125    const list_of_arguments_t &listarguments
    126126    ) const
    127127{
    128   ASSERT( arguments.size() == 1,
    129       "PairPotential_LennardJones::operator() - requires one argument.");
    130   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    131       arguments, getParticleTypes()),
    132       "PairPotential_LennardJones::operator() - types don't match with ones in arguments.");
    133   const double &r = arguments[0].distance;
    134   const double temp = Helpers::pow(params[sigma]/r, 6);
    135   const result_t result = 4.*params[epsilon] * (temp*temp - temp);
    136   return std::vector<result_t>(1, result);
     128  result_t result = 0.;
     129  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     130      iter != listarguments.end(); ++iter) {
     131    const arguments_t &arguments = *iter;
     132    ASSERT( arguments.size() == 1,
     133        "PairPotential_LennardJones::operator() - requires one argument.");
     134    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     135        arguments, getParticleTypes()),
     136        "PairPotential_LennardJones::operator() - types don't match with ones in arguments.");
     137    const double &r = arguments[0].distance;
     138    const double temp = Helpers::pow(params[sigma]/r, 6);
     139    result += 4.*params[epsilon] * (temp*temp - temp);
     140  }
     141  return results_t(1, result);
    137142}
    138143
    139144PairPotential_LennardJones::derivative_components_t
    140145PairPotential_LennardJones::derivative(
    141     const arguments_t &arguments
     146    const list_of_arguments_t &listarguments
    142147    ) const
    143148{
    144   ASSERT( arguments.size() == 1,
    145       "PairPotential_LennardJones::operator() - requires no argument.");
    146   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    147       arguments, getParticleTypes()),
    148       "PairPotential_LennardJones::operator() - types don't match with ones in arguments.");
    149   const double &r = arguments[0].distance;
    150149  const double sigma6 = Helpers::pow(params[sigma], 6);
    151   const result_t result =
    152       4.*params[epsilon] * (
    153           sigma6*sigma6*(-12.) / Helpers::pow(r,13)
    154           - sigma6*(-6.) /Helpers::pow(r,7)
    155       );
    156   derivative_components_t results(1, result);
    157   return results;
     150  result_t result = 0.;
     151  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     152      iter != listarguments.end(); ++iter) {
     153    const arguments_t &arguments = *iter;
     154    ASSERT( arguments.size() == 1,
     155        "PairPotential_LennardJones::operator() - requires no argument.");
     156    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     157        arguments, getParticleTypes()),
     158        "PairPotential_LennardJones::operator() - types don't match with ones in arguments.");
     159    const double &r = arguments[0].distance;
     160    result +=
     161        4.*params[epsilon] * (
     162            sigma6*sigma6*(-12.) / Helpers::pow(r,13)
     163            - sigma6*(-6.) /Helpers::pow(r,7)
     164        );
     165  }
     166  return derivative_components_t(1, result);
    158167}
    159168
    160169PairPotential_LennardJones::results_t
    161170PairPotential_LennardJones::parameter_derivative(
    162     const arguments_t &arguments,
     171    const list_of_arguments_t &listarguments,
    163172    const size_t index
    164173    ) const
    165174{
    166   ASSERT( arguments.size() == 1,
    167       "PairPotential_LennardJones::parameter_derivative() - requires no argument.");
    168   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    169       arguments, getParticleTypes()),
    170       "PairPotential_LennardJones::operator() - types don't match with ones in arguments.");
    171   const double &r = arguments[0].distance;
    172   switch (index) {
    173     case epsilon:
    174     {
    175       const double temp = Helpers::pow(params[sigma]/r, 6);
    176       const result_t result = 4. * (temp*temp - temp);
    177       return std::vector<result_t>(1, result);
    178       break;
     175  result_t result = 0.;
     176  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     177      iter != listarguments.end(); ++iter) {
     178    const arguments_t &arguments = *iter;
     179    ASSERT( arguments.size() == 1,
     180        "PairPotential_LennardJones::parameter_derivative() - requires no argument.");
     181    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     182        arguments, getParticleTypes()),
     183        "PairPotential_LennardJones::operator() - types don't match with ones in arguments.");
     184    const double &r = arguments[0].distance;
     185    switch (index) {
     186      case epsilon:
     187      {
     188        const double temp = Helpers::pow(params[sigma]/r, 6);
     189        result += 4. * (temp*temp - temp);
     190        break;
     191      }
     192      case sigma:
     193      {
     194        const double r6 = Helpers::pow(r, 6);
     195        result +=
     196            4.*params[epsilon] * (
     197                12. * Helpers::pow(params[sigma],11)/(r6*r6)
     198                - 6. * Helpers::pow(params[sigma],5)/r6
     199            );
     200        break;
     201      }
     202      default:
     203        break;
    179204    }
    180     case sigma:
    181     {
    182       const double r6 = Helpers::pow(r, 6);
    183       const result_t result =
    184           4.*params[epsilon] * (
    185               12. * Helpers::pow(params[sigma],11)/(r6*r6)
    186               - 6. * Helpers::pow(params[sigma],5)/r6
    187           );
    188       return std::vector<result_t>(1, result);
    189       break;
    190     }
    191     default:
    192       break;
    193205  }
    194   return std::vector<result_t>(1, 0.);
     206  return results_t(1, result);
    195207}
    196208
  • src/Potentials/Specifics/PairPotential_LennardJones.hpp

    r16227a re1fe7e  
    3636  friend class PotentialFactory;
    3737  // some repeated typedefs to avoid ambiguities
     38  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    3839  typedef FunctionModel::arguments_t arguments_t;
    3940  typedef FunctionModel::result_t result_t;
     
    9293  /** Evaluates the harmonic potential function for the given arguments.
    9394   *
    94    * @param arguments single distance
     95   * @param listarguments list of single distances
    9596   * @return value of the potential function
    9697   */
    97   results_t operator()(const arguments_t &arguments) const;
     98  results_t operator()(const list_of_arguments_t &listarguments) const;
    9899
    99100  /** Evaluates the derivative of the potential function.
    100101   *
    101    * @param arguments single distance
     102   * @param listarguments list of single distances
    102103   * @return vector with derivative with respect to the input degrees of freedom
    103104   */
    104   derivative_components_t derivative(const arguments_t &arguments) const;
     105  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    105106
    106107  /** Evaluates the derivative of the function with the given \a arguments
    107108   * with respect to a specific parameter indicated by \a index.
    108109   *
    109    * \param arguments set of arguments as input variables to the function
     110   * \param listarguments list of single distances
    110111   * \param index derivative of which parameter
    111112   * \return result vector containing the derivative with respect to the given
    112113   *         input
    113114   */
    114   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     115  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    115116
    116117  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/PairPotential_Morse.cpp

    r16227a re1fe7e  
    124124PairPotential_Morse::results_t
    125125PairPotential_Morse::operator()(
    126     const arguments_t &arguments
     126    const list_of_arguments_t &listarguments
    127127    ) const
    128128{
    129   ASSERT( arguments.size() == 1,
    130       "PairPotential_Morse::operator() - requires exactly one argument.");
    131   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    132       arguments, getParticleTypes()),
    133       "PairPotential_Morse::operator() - types don't match with ones in arguments.");
    134   const argument_t &r_ij = arguments[0];
    135   // Maple: f(r,D,k,R,c) := D * (1 - exp(-k*(r-R)))^(2)+c
    136   const result_t result =
    137       params[dissociation_energy] * Helpers::pow( 1.
    138           - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])),2);
     129  result_t result = 0.;
     130  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     131      iter != listarguments.end(); ++iter) {
     132    const arguments_t &arguments = *iter;
     133    ASSERT( arguments.size() == 1,
     134        "PairPotential_Morse::operator() - requires exactly one argument.");
     135    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     136        arguments, getParticleTypes()),
     137        "PairPotential_Morse::operator() - types don't match with ones in arguments.");
     138    const argument_t &r_ij = arguments[0];
     139    // Maple: f(r,D,k,R,c) := D * (1 - exp(-k*(r-R)))^(2)+c
     140    result +=
     141        params[dissociation_energy] * Helpers::pow( 1.
     142            - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])),2);
     143  }
    139144  return std::vector<result_t>(1, result);
    140145}
     
    142147PairPotential_Morse::derivative_components_t
    143148PairPotential_Morse::derivative(
    144     const arguments_t &arguments
     149    const list_of_arguments_t &listarguments
    145150    ) const
    146151{
    147   ASSERT( arguments.size() == 1,
    148       "PairPotential_Morse::operator() - requires exactly one argument.");
    149   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    150       arguments, getParticleTypes()),
    151       "PairPotential_Morse::operator() - types don't match with ones in arguments.");
    152   derivative_components_t result;
    153   const argument_t &r_ij = arguments[0];
    154   // Maple result: 2*D*(1-exp(-k*(r-R)))*k*exp(-k*(r-R))
    155   result.push_back(
    156       2. * params[dissociation_energy]
    157       * ( 1. - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])))
    158       * (- params[spring_constant]) * exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]))
    159   );
    160   ASSERT( result.size() == 1,
    161       "PairPotential_Morse::operator() - we did not create exactly one component.");
    162   return result;
     152  result_t result = 0.;
     153  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     154      iter != listarguments.end(); ++iter) {
     155    const arguments_t &arguments = *iter;
     156    ASSERT( arguments.size() == 1,
     157        "PairPotential_Morse::operator() - requires exactly one argument.");
     158    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     159        arguments, getParticleTypes()),
     160        "PairPotential_Morse::operator() - types don't match with ones in arguments.");
     161    const argument_t &r_ij = arguments[0];
     162    // Maple result: 2*D*(1-exp(-k*(r-R)))*k*exp(-k*(r-R))
     163    result +=
     164        2. * params[dissociation_energy]
     165        * ( 1. - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])))
     166        * (- params[spring_constant]) * exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]));
     167  }
     168  return derivative_components_t(1, result);
    163169}
    164170
    165171PairPotential_Morse::results_t
    166172PairPotential_Morse::parameter_derivative(
    167     const arguments_t &arguments,
     173    const list_of_arguments_t &listarguments,
    168174    const size_t index
    169175    ) const
    170176{
    171   ASSERT( arguments.size() == 1,
    172       "PairPotential_Morse::parameter_derivative() - requires exactly one argument.");
    173   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    174       arguments, getParticleTypes()),
    175       "PairPotential_Morse::operator() - types don't match with ones in arguments.");
    176   const argument_t &r_ij = arguments[0];
    177   switch (index) {
    178     case spring_constant:
    179     {
    180       // Maple result: -2*D*(1-exp(-k*(r-R)))*(-r+R)*exp(-k*(r-R))
    181       const result_t result =
    182           - 2. * params[dissociation_energy]
    183           * ( 1. - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])))
    184           * (- r_ij.distance + params[equilibrium_distance])
    185           * exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]))
    186           ;
    187       return std::vector<result_t>(1, result);
    188       break;
     177  result_t result = 0.;
     178  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     179      iter != listarguments.end(); ++iter) {
     180    const arguments_t &arguments = *iter;
     181    ASSERT( arguments.size() == 1,
     182        "PairPotential_Morse::parameter_derivative() - requires exactly one argument.");
     183    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     184        arguments, getParticleTypes()),
     185        "PairPotential_Morse::operator() - types don't match with ones in arguments.");
     186    const argument_t &r_ij = arguments[0];
     187    switch (index) {
     188      case spring_constant:
     189      {
     190        // Maple result: -2*D*(1-exp(-k*(r-R)))*(-r+R)*exp(-k*(r-R))
     191        result +=
     192            - 2. * params[dissociation_energy]
     193            * ( 1. - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])))
     194            * (- r_ij.distance + params[equilibrium_distance])
     195            * exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]))
     196            ;
     197        break;
     198      }
     199      case equilibrium_distance:
     200      {
     201        // Maple result: -2*D*(1-exp(-k*(r-R)))*k*exp(-k*(r-R))
     202        result +=
     203            - 2. * params[dissociation_energy]
     204            * ( 1. - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])))
     205            * params[spring_constant] * exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]))
     206            ;
     207        break;
     208      }
     209      case dissociation_energy:
     210      {
     211        // Maple result: (1-exp(-k*(r-R)))^2
     212        result +=
     213            Helpers::pow(1.
     214                - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])),2);
     215        break;
     216      }
     217      default:
     218        ASSERT(0, "PairPotential_Morse::parameter_derivative() - derivative to unknown parameter desired.");
     219        break;
    189220    }
    190     case equilibrium_distance:
    191     {
    192       // Maple result: -2*D*(1-exp(-k*(r-R)))*k*exp(-k*(r-R))
    193       const result_t result =
    194           - 2. * params[dissociation_energy]
    195           * ( 1. - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])))
    196           * params[spring_constant] * exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance]))
    197           ;
    198       return std::vector<result_t>(1, result);
    199       break;
    200     }
    201     case dissociation_energy:
    202     {
    203       // Maple result: (1-exp(-k*(r-R)))^2
    204       const result_t result =
    205           Helpers::pow(1.
    206               - exp( - params[spring_constant] * ( r_ij.distance - params[equilibrium_distance])),2);
    207       return std::vector<result_t>(1, result);
    208       break;
    209     }
    210     default:
    211       ASSERT(0, "PairPotential_Morse::parameter_derivative() - derivative to unknown parameter desired.");
    212       break;
    213221  }
    214   return std::vector<result_t>(1, 0.);
     222  return results_t(1, result);
    215223}
    216224
  • src/Potentials/Specifics/PairPotential_Morse.hpp

    r16227a re1fe7e  
    3535  friend class PotentialFactory;
    3636  // some repeated typedefs to avoid ambiguities
     37  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    3738  typedef FunctionModel::arguments_t arguments_t;
    3839  typedef FunctionModel::result_t result_t;
     
    9192  /** Evaluates the harmonic potential function for the given arguments.
    9293   *
    93    * @param arguments single distance
     94   * @param listarguments list of single distances
    9495   * @return value of the potential function
    9596   */
    96   results_t operator()(const arguments_t &arguments) const;
     97  results_t operator()(const list_of_arguments_t &listarguments) const;
    9798
    9899  /** Evaluates the derivative of the potential function.
    99100   *
    100    * @param arguments single distance
     101   * @param listarguments list of single distances
    101102   * @return vector with derivative with respect to the input degrees of freedom
    102103   */
    103   derivative_components_t derivative(const arguments_t &arguments) const;
     104  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    104105
    105106  /** Evaluates the derivative of the function with the given \a arguments
    106107   * with respect to a specific parameter indicated by \a index.
    107108   *
    108    * \param arguments set of arguments as input variables to the function
     109   * \param listarguments list of single distances
    109110   * \param index derivative of which parameter
    110111   * \return result vector containing the derivative with respect to the given
    111112   *         input
    112113   */
    113   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     114  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    114115
    115116  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/ThreeBodyPotential_Angle.cpp

    r16227a re1fe7e  
    136136ThreeBodyPotential_Angle::results_t
    137137ThreeBodyPotential_Angle::operator()(
    138     const arguments_t &arguments
     138    const list_of_arguments_t &listarguments
    139139    ) const
    140140{
    141   ASSERT( arguments.size() == 3,
    142       "ThreeBodyPotential_Angle::operator() - requires exactly three arguments.");
    143   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    144       arguments, getParticleTypes()),
    145       "ThreeBodyPotential_Angle::operator() - types don't match with ones in arguments.");
    146   const argument_t &r_ij = arguments[0]; // 01
    147   const argument_t &r_jk = arguments[2]; // 12
    148   const argument_t &r_ik = arguments[1]; // 02
    149   const result_t result =
    150       params[spring_constant]
    151              * Helpers::pow( function_theta(r_ij.distance, r_jk.distance, r_ik.distance) - params[equilibrium_distance], 2 );
    152   return std::vector<result_t>(1, result);
     141  result_t result = 0.;
     142  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     143      iter != listarguments.end(); ++iter) {
     144    const arguments_t &arguments = *iter;
     145    ASSERT( arguments.size() == 3,
     146        "ThreeBodyPotential_Angle::operator() - requires exactly three arguments.");
     147    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     148        arguments, getParticleTypes()),
     149        "ThreeBodyPotential_Angle::operator() - types don't match with ones in arguments.");
     150    const argument_t &r_ij = arguments[0]; // 01
     151    const argument_t &r_jk = arguments[2]; // 12
     152    const argument_t &r_ik = arguments[1]; // 02
     153    result +=
     154        params[spring_constant]
     155               * Helpers::pow( function_theta(r_ij.distance, r_jk.distance, r_ik.distance)
     156                   - params[equilibrium_distance], 2 );
     157  }
     158  return results_t(1, result);
    153159}
    154160
    155161ThreeBodyPotential_Angle::derivative_components_t
    156162ThreeBodyPotential_Angle::derivative(
    157     const arguments_t &arguments
     163    const list_of_arguments_t &listarguments
    158164    ) const
    159165{
    160   ASSERT( arguments.size() == 3,
    161       "ThreeBodyPotential_Angle::operator() - requires exactly three arguments.");
    162   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    163       arguments, getParticleTypes()),
    164       "ThreeBodyPotential_Angle::operator() - types don't match with ones in arguments.");
    165   derivative_components_t result;
    166   const argument_t &r_ij = arguments[0]; //01
    167   const argument_t &r_jk = arguments[2]; //12
    168   const argument_t &r_ik = arguments[1]; //02
    169   result.push_back( 2. * params[spring_constant] * ( function_theta(r_ij.distance, r_jk.distance, r_ik.distance) - params[equilibrium_distance]) );
    170   ASSERT( result.size() == 1,
    171       "ThreeBodyPotential_Angle::operator() - we did not create exactly one component.");
    172   return result;
     166  result_t result = 0.;
     167  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     168      iter != listarguments.end(); ++iter) {
     169    const arguments_t &arguments = *iter;
     170    ASSERT( arguments.size() == 3,
     171        "ThreeBodyPotential_Angle::operator() - requires exactly three arguments.");
     172    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     173        arguments, getParticleTypes()),
     174        "ThreeBodyPotential_Angle::operator() - types don't match with ones in arguments.");
     175    const argument_t &r_ij = arguments[0]; //01
     176    const argument_t &r_jk = arguments[2]; //12
     177    const argument_t &r_ik = arguments[1]; //02
     178    result +=
     179        2. * params[spring_constant] *
     180        ( function_theta(r_ij.distance, r_jk.distance, r_ik.distance)
     181            - params[equilibrium_distance]);
     182  }
     183  return derivative_components_t(1, result);
    173184}
    174185
    175186ThreeBodyPotential_Angle::results_t
    176187ThreeBodyPotential_Angle::parameter_derivative(
    177     const arguments_t &arguments,
     188    const list_of_arguments_t &listarguments,
    178189    const size_t index
    179190    ) const
    180191{
    181   ASSERT( arguments.size() == 3,
    182       "ThreeBodyPotential_Angle::parameter_derivative() - requires exactly three arguments.");
    183   ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
    184       arguments, getParticleTypes()),
    185       "ThreeBodyPotential_Angle::operator() - types don't match with ones in arguments.");
    186   const argument_t &r_ij = arguments[0]; //01
    187   const argument_t &r_jk = arguments[2]; //12
    188   const argument_t &r_ik = arguments[1]; //02
    189   switch (index) {
    190     case spring_constant:
    191     {
    192       const result_t result =
    193                  Helpers::pow( function_theta(r_ij.distance, r_jk.distance, r_ik.distance) - params[equilibrium_distance], 2 );
    194       return std::vector<result_t>(1, result);
    195       break;
     192  result_t result = 0.;
     193  for(list_of_arguments_t::const_iterator iter = listarguments.begin();
     194      iter != listarguments.end(); ++iter) {
     195    const arguments_t &arguments = *iter;
     196    ASSERT( arguments.size() == 3,
     197        "ThreeBodyPotential_Angle::parameter_derivative() - requires exactly three arguments.");
     198    ASSERT( ParticleTypeChecker::checkArgumentsAgainstParticleTypes(
     199        arguments, getParticleTypes()),
     200        "ThreeBodyPotential_Angle::operator() - types don't match with ones in arguments.");
     201    const argument_t &r_ij = arguments[0]; //01
     202    const argument_t &r_jk = arguments[2]; //12
     203    const argument_t &r_ik = arguments[1]; //02
     204    switch (index) {
     205      case spring_constant:
     206      {
     207        result +=
     208            Helpers::pow( function_theta(r_ij.distance, r_jk.distance, r_ik.distance) - params[equilibrium_distance], 2 );
     209        break;
     210      }
     211      case equilibrium_distance:
     212      {
     213        result +=
     214            -2. * params[spring_constant]
     215                   * ( function_theta(r_ij.distance, r_jk.distance, r_ik.distance) - params[equilibrium_distance]);
     216        break;
     217      }
     218      default:
     219        ASSERT(0, "ThreeBodyPotential_Angle::parameter_derivative() - derivative to unknown parameter desired.");
     220        break;
    196221    }
    197     case equilibrium_distance:
    198     {
    199       const result_t result =
    200           -2. * params[spring_constant]
    201                  * ( function_theta(r_ij.distance, r_jk.distance, r_ik.distance) - params[equilibrium_distance]);
    202       return std::vector<result_t>(1, result);
    203       break;
    204     }
    205     default:
    206       ASSERT(0, "ThreeBodyPotential_Angle::parameter_derivative() - derivative to unknown parameter desired.");
    207       break;
    208222  }
    209   return std::vector<result_t>(1);
     223  return results_t(1, result);
    210224}
    211225
  • src/Potentials/Specifics/ThreeBodyPotential_Angle.hpp

    r16227a re1fe7e  
    3535  friend class PotentialFactory;
    3636  // some repeated typedefs to avoid ambiguities
     37  typedef FunctionModel::list_of_arguments_t list_of_arguments_t;
    3738  typedef FunctionModel::arguments_t arguments_t;
    3839  typedef FunctionModel::result_t result_t;
     
    9091  /** Evaluates the harmonic potential function for the given arguments.
    9192   *
    92    * @param arguments single distance
     93   * @param listarguments list of three distances
    9394   * @return value of the potential function
    9495   */
    95   results_t operator()(const arguments_t &arguments) const;
     96  results_t operator()(const list_of_arguments_t &listarguments) const;
    9697
    9798  /** Evaluates the derivative of the potential function.
    9899   *
    99    * @param arguments single distance
     100   * @param listarguments list of three distances
    100101   * @return vector with derivative with respect to the input degrees of freedom
    101102   */
    102   derivative_components_t derivative(const arguments_t &arguments) const;
     103  derivative_components_t derivative(const list_of_arguments_t &listarguments) const;
    103104
    104105  /** Evaluates the derivative of the function with the given \a arguments
    105106   * with respect to a specific parameter indicated by \a index.
    106107   *
    107    * \param arguments set of arguments as input variables to the function
     108   * \param listarguments list of three distances
    108109   * \param index derivative of which parameter
    109110   * \return result vector containing the derivative with respect to the given
    110111   *         input
    111112   */
    112   results_t parameter_derivative(const arguments_t &arguments, const size_t index) const;
     113  results_t parameter_derivative(const list_of_arguments_t &listarguments, const size_t index) const;
    113114
    114115  /** Returns the functor that converts argument_s into the
  • src/Potentials/Specifics/unittests/ConstantPotentialUnitTest.cpp

    r16227a re1fe7e  
    8282      Helpers::isEqual(
    8383          offset,
    84           constant( FunctionModel::arguments_t() )[0]
     84          constant( FunctionModel::list_of_arguments_t() )[0]
    8585      )
    8686  );
     
    9797          0.,
    9898          constant.derivative(
    99               FunctionModel::arguments_t()
     99              FunctionModel::list_of_arguments_t()
    100100          )[0]
    101101      )
     
    114114          1.,
    115115          constant.parameter_derivative(
    116               FunctionModel::arguments_t(),
     116              FunctionModel::list_of_arguments_t(),
    117117              0
    118118          )[0]
  • src/Potentials/Specifics/unittests/FourBodyPotential_ImproperUnitTest.cpp

    r16227a re1fe7e  
    121121  for (size_t index = 0; index < input.size(); ++index) {
    122122    const FourBodyPotential_Improper::results_t result =
    123         angle( input[index] );
     123        angle( FunctionModel::list_of_arguments_t(1, input[index]) );
    124124    CPPUNIT_ASSERT(
    125125        Helpers::isEqual(
     
    144144          0.,
    145145          angle.derivative(
    146               input[5]
     146              FunctionModel::list_of_arguments_t(1, input[5])
    147147          )[0],
    148148          10.
     
    165165          0.,
    166166          angle.parameter_derivative(
    167               input[5],
     167              FunctionModel::list_of_arguments_t(1, input[5]),
    168168              0
    169169          )[0],
     
    175175          0.,
    176176          angle.parameter_derivative(
    177               input[5],
     177              FunctionModel::list_of_arguments_t(1, input[5]),
    178178              1
    179179          )[0],
  • src/Potentials/Specifics/unittests/FourBodyPotential_TorsionUnitTest.cpp

    r16227a re1fe7e  
    121121  for (size_t index = 0; index < input.size(); ++index) {
    122122    const FourBodyPotential_Torsion::results_t result =
    123         angle( input[index] );
     123        angle( FunctionModel::list_of_arguments_t(1, input[index]) );
    124124    CPPUNIT_ASSERT(
    125125        Helpers::isEqual(
     
    144144          0.,
    145145          angle.derivative(
    146               input[5]
     146              FunctionModel::list_of_arguments_t(1, input[5])
    147147          )[0],
    148148          10.
     
    165165          0.,
    166166          angle.parameter_derivative(
    167               input[5],
     167              FunctionModel::list_of_arguments_t(1, input[5]),
    168168              0
    169169          )[0],
     
    175175          0.,
    176176          angle.parameter_derivative(
    177               input[5],
     177              FunctionModel::list_of_arguments_t(1, input[5]),
    178178              1
    179179          )[0],
  • src/Potentials/Specifics/unittests/ManyBodyPotential_TersoffUnitTest.cpp

    r16227a re1fe7e  
    299299        arg.globalid = index; // this is needed for the triplefunction to the configuration
    300300        FunctionModel::arguments_t args(1,arg);
    301         const ManyBodyPotential_Tersoff::results_t res = tersoff(args);
     301        const ManyBodyPotential_Tersoff::results_t res =
     302            tersoff( FunctionModel::list_of_arguments_t(1, args) );
    302303        temp += res[0];
    303304      }
     
    336337//          0.,
    337338//          tersoff.derivative(
    338 //              FunctionModel::arguments_t(1,argument_t(1.))
     339//              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,argument_t(1.)))
    339340//          )[0]
    340341//      )
     
    361362//          0.,
    362363//          tersoff.parameter_derivative(
    363 //              FunctionModel::arguments_t(1,argument_t(1.)),
     364//              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,argument_t(1.))),
    364365//              0
    365366//          )[0]
     
    370371//          0.,
    371372//          tersoff.parameter_derivative(
    372 //              FunctionModel::arguments_t(1,argument_t(1.)),
     373//              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,argument_t(1.))),
    373374//              1
    374375//          )[0]
     
    379380//          1.,
    380381//          tersoff.parameter_derivative(
    381 //              FunctionModel::arguments_t(1,argument_t(1.)),
     382//              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,argument_t(1.))),
    382383//              2
    383384//          )[0]
  • src/Potentials/Specifics/unittests/PairPotential_HarmonicUnitTest.cpp

    r16227a re1fe7e  
    102102        Helpers::isEqual(
    103103            output[index],
    104             harmonic( FunctionModel::arguments_t(1,arg) )[0]
     104            harmonic( FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)) )[0]
    105105        )
    106106    );
     
    122122          0.,
    123123          harmonic.derivative(
    124               FunctionModel::arguments_t(1,arg)
     124              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg))
    125125          )[0]
    126126      )
     
    143143          0.,
    144144          harmonic.parameter_derivative(
    145               FunctionModel::arguments_t(1,arg),
     145              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    146146              0
    147147          )[0]
     
    152152          0.,
    153153          harmonic.parameter_derivative(
    154               FunctionModel::arguments_t(1,arg),
     154              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    155155              1
    156156          )[0]
  • src/Potentials/Specifics/unittests/PairPotential_LennardJonesUnitTest.cpp

    r16227a re1fe7e  
    118118        Helpers::isEqual(
    119119            output[index],
    120             lj( FunctionModel::arguments_t(1,arg) )[0],
     120            lj( FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)) )[0],
    121121            1.e-4/std::numeric_limits<double>::epsilon() // only compare four digits
    122122        )
     
    140140          0.,
    141141          lj.derivative(
    142               FunctionModel::arguments_t(1,arg)
     142              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg))
    143143          )[0],
    144144          1.e+6
     
    162162          -1.,
    163163          lj.parameter_derivative(
    164               FunctionModel::arguments_t(1,arg),
     164              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    165165              0
    166166          )[0],
     
    172172          0.,
    173173          lj.parameter_derivative(
    174               FunctionModel::arguments_t(1,arg),
     174              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    175175              1
    176176          )[0],
  • src/Potentials/Specifics/unittests/PairPotential_MorseUnitTest.cpp

    r16227a re1fe7e  
    120120        Helpers::isEqual(
    121121            output[index],
    122             Morse( FunctionModel::arguments_t(1,arg) )[0],
     122            Morse( FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)) )[0],
    123123            1.e-4/std::numeric_limits<double>::epsilon() // only compare four digits
    124124        )
     
    141141          0.,
    142142          Morse.derivative(
    143               FunctionModel::arguments_t(1,arg)
     143              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg))
    144144          )[0],
    145145          1.e+6
     
    163163          0.,
    164164          Morse.parameter_derivative(
    165               FunctionModel::arguments_t(1,arg),
     165              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    166166              0
    167167          )[0],
     
    173173          0.,
    174174          Morse.parameter_derivative(
    175               FunctionModel::arguments_t(1,arg),
     175              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    176176              1
    177177          )[0],
     
    183183          0.,
    184184          Morse.parameter_derivative(
    185               FunctionModel::arguments_t(1,arg),
     185              FunctionModel::list_of_arguments_t(1, FunctionModel::arguments_t(1,arg)),
    186186              2
    187187          )[0],
  • src/Potentials/Specifics/unittests/ThreeBodyPotential_AngleUnitTest.cpp

    r16227a re1fe7e  
    112112  for (size_t index = 0; index < input.size(); ++index) {
    113113    const ThreeBodyPotential_Angle::results_t result =
    114         angle( input[index] );
     114        angle( FunctionModel::list_of_arguments_t(1, input[index]) );
    115115    CPPUNIT_ASSERT(
    116116        Helpers::isEqual(
     
    135135          0.,
    136136          angle.derivative(
    137               input[5]
     137              FunctionModel::list_of_arguments_t(1, input[5])
    138138          )[0],
    139139          10.
     
    156156          0.,
    157157          angle.parameter_derivative(
    158               input[5],
     158              FunctionModel::list_of_arguments_t(1, input[5]),
    159159              0
    160160          )[0],
     
    166166          0.,
    167167          angle.parameter_derivative(
    168               input[5],
     168              FunctionModel::list_of_arguments_t(1, input[5]),
    169169              1
    170170          )[0],
  • src/Potentials/helpers.hpp

    r16227a re1fe7e  
    108108  }
    109109
     110  inline FunctionModel::list_of_arguments_t
     111  returnEmptyListArguments()
     112  {
     113    return FunctionModel::list_of_arguments_t();
     114  }
     115
    110116}; /* namespace Helpers */
    111117
  • src/Potentials/unittests/CompoundPotentialUnitTest.cpp

    r16227a re1fe7e  
    159159    argument_t firstarg(argument_t::indices_t(0,1), argument_t::types_t(0,1), input[index]);
    160160    argument_t secondarg(argument_t::indices_t(0,2), argument_t::types_t(0,1), input[index]);
    161     FunctionModel::arguments_t args;
    162     args += firstarg,secondarg;
    163     const double result = compound( args )[0];
     161    FunctionModel::list_of_arguments_t listargs;
     162    listargs += FunctionModel::arguments_t(1,firstarg),FunctionModel::arguments_t(1,secondarg);
     163    const double result = compound( listargs )[0];
    164164    CPPUNIT_ASSERT(
    165165        Helpers::isEqual(
     
    205205  argument_t firstarg(argument_t::indices_t(0,1), argument_t::types_t(0,1), c);
    206206  argument_t secondarg(argument_t::indices_t(0,2), argument_t::types_t(0,1), c);
    207   FunctionModel::arguments_t args;
    208   args += firstarg,secondarg;
    209   {
    210     const double result =
    211         compound.parameter_derivative(
    212             args,
     207  FunctionModel::list_of_arguments_t listargs;
     208  listargs += FunctionModel::arguments_t(1,firstarg),FunctionModel::arguments_t(1,secondarg);
     209  {
     210    const double result =
     211        compound.parameter_derivative(
     212            listargs,
    213213            0)[0];
    214214    CPPUNIT_ASSERT(
     
    223223    const double result =
    224224        compound.parameter_derivative(
    225             args,
     225            listargs,
    226226            1)[0];
    227227    CPPUNIT_ASSERT(
     
    236236    const double result =
    237237        compound.parameter_derivative(
    238             args,
     238            listargs,
    239239            2)[0];
    240240    CPPUNIT_ASSERT(
     
    249249    const double result =
    250250        compound.parameter_derivative(
    251             args,
     251            listargs,
    252252            3)[0];
    253253    CPPUNIT_ASSERT(
Note: See TracChangeset for help on using the changeset viewer.