| [201199] | 1 | /* | 
|---|
|  | 2 | * Project: MoleCuilder | 
|---|
|  | 3 | * Description: creates and alters molecular systems | 
|---|
|  | 4 | * Copyright (C)  2010 University of Bonn. All rights reserved. | 
|---|
|  | 5 | * Please see the LICENSE file or "Copyright notice" in builder.cpp for details. | 
|---|
|  | 6 | */ | 
|---|
|  | 7 |  | 
|---|
|  | 8 | /** | 
|---|
|  | 9 | * \file potentials.dox | 
|---|
|  | 10 | * | 
|---|
|  | 11 | * Created on: Nov 28, 2012 | 
|---|
|  | 12 | *    Author: heber | 
|---|
|  | 13 | */ | 
|---|
|  | 14 |  | 
|---|
|  | 15 | /** \page potentials Empirical Potentials and FunctionModels | 
|---|
|  | 16 | * | 
|---|
|  | 17 | * On this page we explain what is meant with the Potentials sub folder. | 
|---|
|  | 18 | * | 
|---|
|  | 19 | * First, we are based on fragmenting a molecular system, i.e. dissecting its | 
|---|
|  | 20 | * bond structure into connected subgraphs, calculating the energies of the | 
|---|
|  | 21 | * fragments (ab-initio) and summing up to a good approximation of the total | 
|---|
|  | 22 | * energy of the whole system. | 
|---|
|  | 23 | * Second, having calculated these energies, there quickly comes up the thought | 
|---|
|  | 24 | * that one actually calculates quite similar systems all time and if one could | 
|---|
|  | 25 | * not cache results in an intelligent (i.e. interpolating) fashion ... | 
|---|
|  | 26 | * | 
|---|
|  | 27 | * That's where so-called empirical potentials come into play. They are | 
|---|
|  | 28 | * functions depending on a number of "fitted" parameters and the variable | 
|---|
|  | 29 | * distances within a molecular fragment (i.e. the bond lengths) in order to | 
|---|
|  | 30 | * give a value for the total energy without the need to solve a complex | 
|---|
|  | 31 | * ab-initio model. | 
|---|
|  | 32 | * | 
|---|
|  | 33 | * Empirical potentials have been thought of by fellows such as Lennard-Jones, | 
|---|
|  | 34 | * Morse, Tersoff, Stillinger and Weber, etc. And in their honor, the | 
|---|
|  | 35 | * potential form is named after its inventor. Hence, we speak e.g. of a | 
|---|
|  | 36 | * Lennard-Jones potential. | 
|---|
|  | 37 | * | 
|---|
|  | 38 | * So, what we have to do in order to cache results is the following procedure: | 
|---|
|  | 39 | * -# gather similar fragments | 
|---|
|  | 40 | * -# perform a fit procedure to obtain the parameters for the empirical | 
|---|
|  | 41 | *    potential | 
|---|
|  | 42 | * -# evaluate the potential instead of an ab-initio calculation | 
|---|
|  | 43 | * | 
|---|
|  | 44 | * The terms we use, model the classes that are implemented: | 
|---|
|  | 45 | * -# EmpiricalPotential: Contains the interface to a function that can be | 
|---|
|  | 46 | *    evaluated given a number of arguments_t, i.e. distances. Also, one might | 
|---|
|  | 47 | *    want to evaluate derivatives. | 
|---|
|  | 48 | * -# FunctionModel: Is a function that can be fitted, i.e. that has internal | 
|---|
|  | 49 | *    parameters to be set and got. | 
|---|
|  | 50 | * -# argument_t: The Argument stores not only the distance but also the index | 
|---|
|  | 51 | *    pair of the associated atoms and also their charges, to let the potential | 
|---|
|  | 52 | *    check on validity. | 
|---|
|  | 53 | * -# SerializablePotential: Eventually, one wants to store to or parse from | 
|---|
|  | 54 | *    a file all potential parameters. This functionality is encoded in this | 
|---|
|  | 55 | *    class. | 
|---|
|  | 56 | * -# HomologyGraph: "Similar" fragments in our case have to have the same bond | 
|---|
|  | 57 | *    graph. It is stored in the HomologyGraph that acts as representative | 
|---|
|  | 58 | * -# HomologyContainer: This container combines, in multimap fashion, all | 
|---|
|  | 59 | *    similar fragments with their energies together, with the HomologyGraph | 
|---|
|  | 60 | *    as their "key". | 
|---|
|  | 61 | * -# TrainingData: Here, we combine InputVector and OutputVector that contain | 
|---|
|  | 62 | *    the set of distances required for the FunctionModel (e.g. only a single | 
|---|
|  | 63 | *    distance/argument for a pair potential, three for an angle potential, | 
|---|
|  | 64 | *    etc.) and also the expected OutputVector. This in combination with the | 
|---|
|  | 65 | *    FunctionModel is the basis for the non-linear regression used for the | 
|---|
|  | 66 | *    fitting procedure. | 
|---|
|  | 67 | * -# Extractors: These set of functions yield the set of distances from a | 
|---|
|  | 68 | *    given fragment that is stored in the HomologyContainer. | 
|---|
|  | 69 | * -# FunctionApproximation: Contains the interface to the levmar package where | 
|---|
|  | 70 | *    the Levenberg-Marquardt (Newton + Trust region) algorithm is used to | 
|---|
|  | 71 | *    perform the fit. | 
|---|
|  | 72 | * | 
|---|
|  | 73 | * \section potentials-howto Howto use the potentials | 
|---|
|  | 74 | * | 
|---|
|  | 75 | *  We just give a brief run-down in terms of code on how to use the potentials. | 
|---|
|  | 76 | *  Here, we just describe what to do in order to perform the fitting. | 
|---|
|  | 77 | * | 
|---|
|  | 78 | *  \code | 
|---|
|  | 79 | *  // we need the homology container and the representative graph we want to | 
|---|
|  | 80 | *  // fit to. | 
|---|
|  | 81 | *  HomologyContainer homologies; | 
|---|
|  | 82 | *  const HomologyGraph graph = getSomeGraph(homologies); | 
|---|
|  | 83 | *  Fragment::charges_t h2o; | 
|---|
|  | 84 | *  h2o += 8,1,1; | 
|---|
|  | 85 | *  // TrainingData needs so called Extractors to get the required distances | 
|---|
|  | 86 | *  // from the stored fragment. These are functions are bound. | 
|---|
|  | 87 | *  TrainingData AngleData( | 
|---|
|  | 88 | *      boost::bind(&Extractors::gatherDistancesFromFragment, | 
|---|
|  | 89 | *          boost::bind(&Fragment::getPositions, _1), | 
|---|
|  | 90 | *          boost::bind(&Fragment::getCharges, _1), | 
|---|
|  | 91 | *          boost::cref(h2o), | 
|---|
|  | 92 | *          _2) | 
|---|
|  | 93 | *      ); | 
|---|
|  | 94 | *  // now we extract the distances and energies and store them | 
|---|
|  | 95 | *  AngleData(homologies.getHomologousGraphs(graph)); | 
|---|
|  | 96 | *  // give ParticleTypes of this potential to make it unique | 
|---|
|  | 97 | *  PairPotential_Angle::ParticleTypes_t types = | 
|---|
|  | 98 | *        boost::assign::list_of<PairPotential_Angle::ParticleType_t> | 
|---|
|  | 99 | *        (8)(1)(1) | 
|---|
|  | 100 | *      ; | 
|---|
|  | 101 | *  PairPotential_Angle angle(types); | 
|---|
|  | 102 | *  // give initial parameter | 
|---|
|  | 103 | *  FunctionModel::parameters_t params(PairPotential_Angle::MAXPARAMS, 0.); | 
|---|
|  | 104 | *  ... set some initial parameters | 
|---|
|  | 105 | *  angle.setParameters(params); | 
|---|
|  | 106 | * | 
|---|
|  | 107 | *  // use the potential as a FunctionModel along with prepared TrainingData | 
|---|
|  | 108 | *  FunctionModel &model = angle; | 
|---|
|  | 109 | *  FunctionApproximation approximator(AngleData, model); | 
|---|
|  | 110 | *  approximator(FunctionApproximation::ParameterDerivative); | 
|---|
|  | 111 | * | 
|---|
|  | 112 | *  // obtain resulting parameters and check remaining L_2 and L_max error | 
|---|
|  | 113 | *  angleparams = model.getParameters(); | 
|---|
|  | 114 | *  LOG(1, "INFO: L2sum = " << AngleData(model) | 
|---|
|  | 115 | *      << ", LMax = " << AngleData(model) << "."); | 
|---|
|  | 116 | *  \endcode | 
|---|
|  | 117 | * | 
|---|
|  | 118 | *  The evaluation of the fitted potential is then trivial, e.g. | 
|---|
|  | 119 | *  \code | 
|---|
|  | 120 | *  // constructed someplace | 
|---|
|  | 121 | *  PairPotential_Angle angle(...); | 
|---|
|  | 122 | * | 
|---|
|  | 123 | *  // evaluate | 
|---|
|  | 124 | *  FunctionModel::arguments_t args; | 
|---|
|  | 125 | *  .. initialise args to the desired distances | 
|---|
|  | 126 | *  const double value = angle(args)[0]; // output is a vector! | 
|---|
|  | 127 | *  \endcode | 
|---|
|  | 128 | * | 
|---|
|  | 129 | * \date 2012-11-28 | 
|---|
|  | 130 | */ | 
|---|