diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index de7749e2788..21d00c80e5c 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -472,9 +472,9 @@ run_tutorials: run_doxygen: <<: *global_job_definition stage: additional_checks - image: gitlab.icp.uni-stuttgart.de:4567/espressomd/docker/ubuntu-python3:cuda-9.0 + image: gitlab.icp.uni-stuttgart.de:4567/espressomd/docker/ubuntu-python3:cuda-10.1 needs: - - cuda9-maxset + - cuda10-maxset when: on_success only: - python diff --git a/doc/doxygen/CMakeLists.txt b/doc/doxygen/CMakeLists.txt index 044c4052ffd..74d70346e75 100644 --- a/doc/doxygen/CMakeLists.txt +++ b/doc/doxygen/CMakeLists.txt @@ -7,6 +7,17 @@ if(DOXYGEN_FOUND) DEPENDS EspressoConfig ) + # transform BibTeX DOI fields into URL fields (bibliographic + # styles available to Doxygen do not process the DOI field) + add_custom_command( + OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/bibliography.bib + COMMAND sed -r + "'s_^ *doi *= *([^0-9]+)(10\\.[0-9]+)_url=\\1https://dx.doi.org/\\2_'" + ${CMAKE_CURRENT_SOURCE_DIR}/bibliography.bib > + ${CMAKE_CURRENT_BINARY_DIR}/bibliography.bib + DEPENDS ${CMAKE_CURRENT_SOURCE_DIR}/bibliography.bib + ) + set(DOXYGEN_IN ${CMAKE_CURRENT_SOURCE_DIR}/Doxyfile.in) set(DOXYGEN_OUT ${CMAKE_CURRENT_BINARY_DIR}/Doxyfile) @@ -17,6 +28,7 @@ if(DOXYGEN_FOUND) WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} COMMENT "Generating API documentation with Doxygen" DEPENDS ${CMAKE_CURRENT_BINARY_DIR}/doxy-features + ${CMAKE_CURRENT_BINARY_DIR}/bibliography.bib VERBATIM ) diff --git a/doc/doxygen/Doxyfile.in b/doc/doxygen/Doxyfile.in index 92794e4530f..75e06ceb96c 100644 --- a/doc/doxygen/Doxyfile.in +++ b/doc/doxygen/Doxyfile.in @@ -585,7 +585,7 @@ LAYOUT_FILE = # of the bibliography can be controlled using LATEX_BIB_STYLE. To use this # feature you need bibtex and perl available in the search path. -CITE_BIB_FILES = +CITE_BIB_FILES = bibliography.bib #--------------------------------------------------------------------------- # configuration options related to warning and progress messages @@ -700,7 +700,8 @@ EXCLUDE_SYMLINKS = NO # against the file with absolute path, so to exclude all test directories # for example use the pattern */test/* -EXCLUDE_PATTERNS = +EXCLUDE_PATTERNS = @CMAKE_SOURCE_DIR@/*/unit_tests/* \ + @CMAKE_SOURCE_DIR@/*/tests/* # The EXCLUDE_SYMBOLS tag can be used to specify one or more symbol names # (namespaces, classes, functions, etc.) that should be excluded from the @@ -1315,7 +1316,7 @@ LATEX_SOURCE_CODE = NO # bibliography, e.g. plainnat, or ieeetr. The default style is "plain". See # http://en.wikipedia.org/wiki/BibTeX for more info. -LATEX_BIB_STYLE = plain +LATEX_BIB_STYLE = plainnat #--------------------------------------------------------------------------- # configuration options related to the RTF output diff --git a/doc/doxygen/bibliography.bib b/doc/doxygen/bibliography.bib new file mode 100644 index 00000000000..3a215288175 --- /dev/null +++ b/doc/doxygen/bibliography.bib @@ -0,0 +1,494 @@ +% Encoding: US-ASCII + +@Book{abramowitz65a, + title = {{Handbook of mathematical functions: With formulas, graphs, and mathematical tables}}, + author = {Abramowitz, M. and Stegun, I.}, + publisher = {Dover Publications Inc.}, + year = {1965}, + address = {New York}, + edition = {9th}, + isbn = {9780486612720}, + url = {https://store.doverpublications.com/0486612724.html}, +} + +@Article{ahlrichs99a, + title = {{Simulation of a single polymer chain in solution by combining lattice Boltzmann and molecular dynamics}}, + author = {Ahlrichs, P. and D\"{u}nweg, B.}, + journal = {Journal of Chemical Physics}, + year = {1999}, + number = {17}, + pages = {8225--8239}, + volume = {111}, + doi = {10.1063/1.480156}, +} + +@Article{andersen83a, + title = {{Rattle: A "velocity" version of the shake algorithm for molecular dynamics calculations}}, + author = {Andersen, Hans C.}, + journal = {Journal of Computational Physics}, + year = {1983}, + pages = {24--34}, + volume = {52}, + doi = {10.1016/0021-9991(83)90014-1}, +} + +@Article{arnold02a, + author = {Arnold, Axel and Holm, Christian}, + title = {{{MMM2D}: A fast and accurate summation method for electrostatic interactions in {2D} slab geometries}}, + journal = {Computer Physics Communications}, + year = {2002}, + volume = {148}, + number = {3}, + pages = {327--348}, + doi = {10.1016/s0010-4655(02)00586-6}, +} + +@Article{arnold13c, + title = {{Efficient algorithms for electrostatic interactions including dielectric contrasts}}, + author = {Arnold, Axel and Breitsprecher, Konrad and Fahrenberger, Florian and Kesselheim, Stefan and Lenz, Olaf and Holm, Christian}, + journal = {Entropy}, + year = {2013}, + number = {11}, + pages = {4569--4588}, + volume = {15}, + doi = {10.3390/e15114569}, + issn = {1099-4300}, +} + +@Article{brodka04a, + author = {Br\'{o}dka, A.}, + title = {{Ewald summation method with electrostatic layer correction for interactions of point dipoles in slab geometry}}, + journal = {Chemical Physics Letters}, + year = {2004}, + volume = {400}, + number = {1-3}, + pages = {62--67}, + doi = {10.1016/j.cplett.2004.10.086}, +} + +@InCollection{burtscher11a, +author = {Burtscher, Martin and Pingali, Keshav}, +chapter = {6}, +title = {{An efficient CUDA implementation of the tree-based Barnes Hut n-body algorithm}}, +editor = {Hwu, Wen-mei W.}, +booktitle = {{GPU Computing Gems Emerald Edition}}, +publisher = {Morgan Kaufmann}, +address = {Boston}, +pages = {75-92}, +year = {2011}, +series = {Applications of GPU Computing Series}, +isbn = {978-0-12-384988-5}, +doi = {10.1016/B978-0-12-384988-5.00006-1}, +} + +@Article{cerda08d, + author = {Cerd\`{a}, Juan J. and Ballenegger, Vincent and Lenz, Olaf and Holm, Christian}, + title = {{{P3M} algorithm for dipolar interactions}}, + journal = {Journal of Chemical Physics}, + year = {2008}, + volume = {129}, + pages = {234104}, + doi = {10.1063/1.3000389}, +} + +@PhdThesis{deserno00b, + author = {Deserno, Markus}, + title = {Counterion condensation for rigid linear polyelectrolytes}, + year = {2000}, + month = FEB, + url = {http://archimed.uni-mainz.de/pub/2000/0018}, + school = {Universit\"{a}t Mainz}, +} + +@InProceedings{deserno00e, + author = {Deserno, Markus and Holm, Christian and Limbach, Hans J\"{o}rg}, + title = {{How to mesh up Ewald sums}}, + booktitle = {{Molecular Dynamics on Parallel Computers}}, + year = {2000}, + editor = {Esser, R. and Grassberger, P. and Grotendorst, J. and Lewerenz, M.}, + pages = {319}, + doi = {10.1142/9789812793768_0023}, + address = {Singapore}, +} + +@Article{deserno98a, + author = {Deserno, Markus and Holm, Christian}, + title = {{How to mesh up Ewald sums. I. A theoretical and numerical comparison of various particle mesh routines}}, + journal = {Journal of Chemical Physics}, + year = {1998}, + volume = {109}, + pages = {7678}, + doi = {10.1063/1.477414}, +} + +@Article{deserno98b, + author = {Deserno, Markus and Holm, Christian}, + title = {{How to mesh up Ewald sums. II. An accurate error estimate for the Particle-Particle-Particle-Mesh algorithm}}, + journal = {Journal of Chemical Physics}, + year = {1998}, + volume = {109}, + pages = {7694}, + doi = {10.1063/1.477415}, +} + +@Article{dhumieres09a, + title = {{Viscosity independent numerical errors for lattice Boltzmann models: from recurrence equations to ``magic'' collision numbers}}, + author = {d'Humi\`{e}res, Dominique and Ginzburg, Irina}, + journal = {Computers \& Mathematics with Applications}, + year = {2009}, + number = {5}, + pages = {823--840}, + volume = {58}, + doi = {10.1016/j.camwa.2009.02.008} +} + +@InCollection{duenweg09a, + title = {{Lattice Boltzmann simulations of soft matter systems}}, + author = {D\"{u}nweg, B. and Ladd, A. J. C.}, + booktitle = {{Advanced Computer Simulation Approaches for Soft Matter Sciences III}}, + publisher = {Springer-Verlag Berlin}, + year = {2009}, + address = {Berlin, Germany}, + pages = {89--166}, + series = {Advances in Polymer Science}, + type = {Review}, + volume = {221}, + doi = {10.1007/12_2008_4}, + issn = {0065-3195}, +} + +@Article{dunweg07a, + author = {D\"{u}nweg, Burkhard and Schiller, Ulf D. and Ladd, Anthony J. C.}, + title = {{Statistical mechanics of the fluctuating lattice Boltzmann equation}}, + journal = {Physical Review E}, + year = {2007}, + volume = {76}, + number = {3}, + pages = {036704}, + doi = {10.1103/PhysRevE.76.036704}, +} + +@Article{dupin07a, + title = {{Modeling the flow of dense suspensions of deformable particles in three dimensions}}, + author = {Dupin, Michael M. and Halliday, Ian and Care, Chris M. and Alboul, Lyuba and Munn, Lance L.}, + journal = {Physical Review E}, + volume = {75}, + number = {6}, + pages = {066707}, + numpages = {17}, + year = {2007}, + publisher = {American Physical Society}, + doi = {10.1103/PhysRevE.75.066707}, +} + +@Article{dupin08a, +author = {Dupin, M. M. and Halliday, I. and Care, C. M. and Munn, L. L.}, +title = {{Lattice Boltzmann modelling of blood cell dynamics}}, +journal = {International Journal of Computational Fluid Dynamics}, +volume = {22}, +number = {7}, +pages = {481-492}, +year = {2008}, +doi = {10.1080/10618560802238242}, +} + +@Article{essmann95a, + title = {{A smooth Particle Mesh Ewald method}}, + author = {Essmann, U. and Perera, L. and Berkowitz, M. L. and Darden, T. and Lee, H. and Pedersen, L.}, + journal = {Journal of Chemical Physics}, + year = {1995}, + pages = {8577}, + volume = {103}, + doi = {10.1063/1.470117}, +} + +@Article{ewald21a, + title = {{Die Berechnung optischer und elektrostatischer Gitterpotentiale}}, + author = {Ewald, P. P.}, + journal = {Annales de Physique}, + year = {1921}, + number = {3}, + pages = {253--287}, + volume = {369}, + doi = {10.1002/andp.19213690304}, +} + +@Book{frenkel02b, + author = {Frenkel, Daan and Smit, Berend}, + title = {{Understanding molecular simulation: From algorithms to applications}}, + year = {2002}, + edition = {2nd}, + publisher = {Academic Press}, + isbn = {978-0-12-267351-1}, + doi = {10.1016/B978-0-12-267351-1.X5000-7}, + address = {San Diego}, +} + +@Book{goldstein80a, + title = {{Classical Mechanics}}, + author = {Goldstein, Herbert}, + publisher = {Addison-Wesley}, + year = {1980}, + edition = {2nd}, + isbn = {9780201029185}, + url = {https://www.pearson.com/us/higher-education/product/Goldstein-Classical-Mechanics-2nd-Edition/9780201029185.html}, +} + +@Article{gompper96a, +author={Gompper, G. and Kroll, D. M.}, +journal={Journal de Physique I France}, +number={10}, +pages={1305-1320}, +title={{Random surface discretizations and the renormalization of the bending rigidity}}, +volume={6}, +year={1996}, +doi={10.1051/jp1:1996246}, +} + +@Article{guo02a, + title = {{Discrete lattice effects on the forcing term in the lattice Boltzmann method}}, + author = {Guo, Zhaoli and Zheng, Chuguang and Shi, Baochang}, + journal = {Physical Review E}, + year = {2002}, + number = {4}, + pages = {046308}, + volume = {65}, + doi = {10.1103/PhysRevE.65.046308}, +} + +@Book{hockney88a, + title = {{Computer simulation using particles}}, + author = {Hockney, R. W. and Eastwood, J. W.}, + publisher = {CRC Press}, + year = {1988}, + isbn = {9780852743928}, + url = {https://www.crcpress.com/Computer-Simulation-Using-Particles/Hockney-Eastwood/p/book/9780852743928}, +} + +@Article{kolafa92a, + title = {{Cutoff errors in the Ewald summation formulae for point charge systems}}, + author = {Kolafa, Jiri and Perram, John W.}, + journal = {Molecular Simulation}, + year = {1992}, + number = {5}, + pages = {351--368}, + volume = {9}, + doi = {10.1080/08927029208049126}, +} + +@Book{kruger12a, +author={Kr\"{u}ger, Timm}, +title={{Computer simulation study of collective phenomena in dense suspensions of red blood cells under shear}}, +year={2012}, +publisher={Vieweg+Teubner Verlag}, +address={Wiesbaden}, +isbn={978-3-8348-2376-2}, +doi={10.1007/978-3-8348-2376-2}, +} + +@Article{jancigova16a, +author = {Jan\v{c}igov\'{a}, Iveta and Cimr\'{a}k, Ivan}, +title = {{Non-uniform force allocation for area preservation in spring network models}}, +journal = {International Journal for Numerical Methods in Biomedical Engineering}, +volume = {32}, +number = {10}, +pages = {e02757}, +year = {2016}, +doi = {10.1002/cnm.2757}, +} + +@Article{ladd01a, + title = {{Lattice-Boltzmann simulations of particle-fluid suspensions}}, + author = {Ladd, A. J. C. and Verberg, R.}, + journal = {Journal of Statistical Physics}, + year = {2001}, + number = {5/6}, + pages = {1191--1251}, + volume = {104}, + doi = {10.1023/A:1010414013942}, +} + +@Article{ladd94a, + title = {{Numerical simulations of particulate suspensions via a discretized Boltzmann equation. Part 1. Theoretical foundation}}, + author = {Ladd, A. J. C.}, + journal = {Journal of Fluid Mechanics}, + year = {1994}, + pages = {285--309}, + volume = {271}, + doi = {10.1017/S0022112094001771}, +} + +@Article{marsili10a, +author = {Marsili, Simone and Signorini, Giorgio Federico and Chelli, Riccardo and Marchi, Massimo and Procacci, Piero}, +title = {{{ORAC}: A molecular dynamics simulation program to explore free energy surfaces in biomolecular systems at the atomistic level}}, +journal = {Journal of Computational Chemistry}, +volume = {31}, +number = {5}, +pages = {1106-1116}, +doi = {10.1002/jcc.21388}, +year = {2010}, +} + +@Article{moss96a, +author={Moss, G. P.}, +journal={Pure and Applied Chemistry}, +number={12}, +pages={2193-2222}, +title={{Basic terminology of stereochemistry}}, +volume={68}, +year={1996}, +doi={10.1351/pac199668122193}, +} + +@Article{neumann85b, +author = {Neumann, Martin}, +title = {{The dielectric constant of water. Computer simulations with the {MCY} potential}}, +journal = {The Journal of Chemical Physics}, +volume = {82}, +number = {12}, +pages = {5663-5672}, +year = {1985}, +doi = {10.1063/1.448553}, +} + +@Article{reed92a, + title = {{Monte Carlo study of titration of linear polyelectrolytes}}, + author = {Reed, C. E. and Reed, W. F.}, + journal = {Journal of Chemical Physics}, + year = {1992}, + pages = {1609}, + volume = {96}, + doi = {10.1063/1.462145}, +} + +@Article{smith94c, + title={{The reaction ensemble method for the computer simulation of chemical and phase equilibria. I. Theory and basic examples}}, + author={Smith, W. R. and Triska, B.}, + journal={The Journal of Chemical Physics}, + volume={100}, + number={4}, + pages={3019--3027}, + year={1994}, + doi={10.1063/1.466443}, +} + +@Article{soddemann03a, + title = {{Dissipative particle dynamics: A useful thermostat for equilibrium and nonequilibrium molecular dynamics simulations}}, + author = {Soddemann, T. and D\"{u}nweg, B. and Kremer, K.}, + journal = {Physical Review E}, + year = {2003}, + number = {4}, + pages = {46702}, + volume = {68}, + doi = {10.1103/PhysRevE.68.046702}, +} + +@Article{sonnenschein85a, +author = {Sonnenschein, Roland}, +title = {{An improved algorithm for molecular dynamics simulation of rigid molecules}}, +journal = {Journal of Computational Physics}, +volume = {59}, +number = {2}, +pages = {347-350}, +year = {1985}, +doi = {10.1016/0021-9991(85)90151-2}, +} + +@Article{swope92a, +author = {Swope, William C. and Ferguson, David M.}, +title = {{Alternative expressions for energies and forces due to angle bending and torsional energy}}, +journal = {Journal of Computational Chemistry}, +volume = {13}, +number = {5}, +pages = {585-594}, +year = {1992}, +doi = {10.1002/jcc.540130508}, +} + +@Article{thole81a, + author = {Thole, B. T.}, + title = {{Molecular polarizabilities calculated with a modified dipole interaction}}, + journal = {Chemical Physics}, + year = {1981}, + volume = {59}, + number = {3}, + pages = {341--350}, + doi = {10.1016/0301-0104(81)85176-2}, +} + +@Article{troester05a, + title = {{Wang-Landau sampling with self-adaptive range}}, + author = {Tr\"{o}ster, Andreas and Dellago, Christoph}, + journal = {Physical Review E}, + volume = {71}, + number = {6}, + pages = {066705}, + year = {2005}, + doi = {10.1103/PhysRevE.71.066705}, +} + +@Article{tyagi10a, + author = {Tyagi, Sandeep and S\"{u}zen, Mehmet and Sega, Marcello and Barbosa, Marcia C. and Kantorovich, Sofia S. and Holm, Christian}, + title = {{An iterative, fast, linear-scaling method for computing induced charges on arbitrary dielectric boundaries}}, + journal = {The Journal of Chemical Physics}, + year = {2010}, + volume = {132}, + pages = {154112}, + doi = {10.1063/1.3376011}, +} + +@Article{usta05a, + title = {{Lattice-Boltzmann simulations of the dynamics of polymer solutions in periodic and confined geometries}}, + author = {Usta, O. B. and Ladd, A. J. C. and Butler, J. E.}, + journal = {Journal of Chemical Physics}, + year = {2005}, + pages = {094902}, + volume = {122}, + doi = {10.1063/1.1854151}, +} + +@Article{wang01a, + author = {Wang, Zuowei and Holm, Christian}, + title = {{Estimate of the cutoff errors in the Ewald summation for dipolar systems}}, + journal = {The Journal of Chemical Physics}, + year = {2001}, + volume = {115}, + pages = {6351--6359}, + doi = {10.1063/1.1398588}, +} + +@Article{yan02b, +author = {Yan, Qiliang and Faller, Roland and {de Pablo}, Juan J.}, +title = {{Density-of-states Monte Carlo method for simulation of fluids}}, +journal = {The Journal of Chemical Physics}, +volume = {116}, +number = {20}, +pages = {8745-8749}, +year = {2002}, +doi = {10.1063/1.1463055}, +} + +@Article{yeh99a, + title = {{Ewald summation for systems with slab geometry}}, + author = {Yeh, In-Chul and Berkowitz, Max L.}, + journal = {Journal of Chemical Physics}, + year = {1999}, + number = {7}, + pages = {3155--3162}, + volume = {111}, + doi = {10.1063/1.479595}, +} + +@InProceedings{zhang01b, +author={Zhang, Cha and Chen, Tsuhan}, +booktitle={{Proceedings 2001 International Conference on Image Processing (Cat. No.01CH37205)}}, +title={{Efficient feature extraction for {2D/3D} objects in mesh representation}}, +year={2001}, +volume={3}, +pages={935-938}, +doi={10.1109/ICIP.2001.958278}, +} + +@Comment{jabref-meta: databaseType:biblatex;} + +@Comment{jabref-meta: saveOrderConfig:specified;bibtexkey;false;bibtexkey;false;bibtexkey;false;} diff --git a/doc/sphinx/advanced_methods.rst b/doc/sphinx/advanced_methods.rst index df11b2d04a9..58d33e0a439 100644 --- a/doc/sphinx/advanced_methods.rst +++ b/doc/sphinx/advanced_methods.rst @@ -1701,7 +1701,7 @@ An example script can be found here: * `Reaction ensemble / constant pH ensemble `_ -The reaction ensemble :cite:`smith94a,turner2008simulation` allows to simulate +The reaction ensemble :cite:`smith94c,turner2008simulation` allows to simulate chemical reactions which can be represented by the general equation: .. math:: @@ -1766,7 +1766,7 @@ In the *forward* reaction, the appropriate number of reactants (given by products is inserted into the system. In the *backward* reaction, reactants and products exchange their roles. The acceptance probability :math:`P^{\xi}` for move from state :math:`o` to :math:`n` reaction -ensemble is given by the criterion :cite:`smith94a` +ensemble is given by the criterion :cite:`smith94c` .. math:: diff --git a/doc/sphinx/zrefs.bib b/doc/sphinx/zrefs.bib index 00321264741..5750d701fb5 100644 --- a/doc/sphinx/zrefs.bib +++ b/doc/sphinx/zrefs.bib @@ -855,7 +855,7 @@ @ARTICLE{smith81a timestamp = {2007.06.13} } -@article{smith94a, +@article{smith94c, title={The reaction ensemble method for the computer simulation of chemical and phase equilibria. {I.} {T}heory and basic examples}, author={Smith, WR and Triska, B}, journal={J. Chem. Phys.}, diff --git a/maintainer/files_with_header.sh b/maintainer/files_with_header.sh index 5bb674b534a..15d1a088e04 100755 --- a/maintainer/files_with_header.sh +++ b/maintainer/files_with_header.sh @@ -1,4 +1,4 @@ -# Copyright (C) 2012-2018 The ESPResSo project +# Copyright (C) 2012-2019 The ESPResSo project # Copyright (C) 2012 Olaf Lenz # # This file is part of ESPResSo. diff --git a/src/config/config.hpp b/src/config/config.hpp index 130e29d5cbb..87f69f5577b 100644 --- a/src/config/config.hpp +++ b/src/config/config.hpp @@ -22,11 +22,11 @@ #define ESPRESSO_CONFIG_HPP /** \file - - This file contains the defaults for Espresso. To modify them, add - an appropriate line in myconfig.hpp. To find a list of features that - can be compiled into Espresso, refer to myconfig-sample.hpp or to - the documentation of the features. + * + * This file contains the defaults for ESPResSo. To modify them, add + * an appropriate line in myconfig.hpp. To find a list of features that + * can be compiled into ESPResSo, refer to myconfig-sample.hpp or to + * the documentation of the features. */ /* Prevent C++ bindings in MPI (there is a DataType called LB in there) */ @@ -36,44 +36,47 @@ #include "config-features.hpp" /** P3M: Default for number of interpolation points of the charge - assignment function. */ + * assignment function. + */ #ifndef P3M_N_INTERPOL #define P3M_N_INTERPOL 32768 #endif -/** P3M: Default for boundary condition: Epsilon of the surrounding - medium. */ +/** P3M: Default for boundary condition: Epsilon of the surrounding medium. */ #ifndef P3M_EPSILON #define P3M_EPSILON 0.0 #endif -/** P3M: Default for boundary condition in magnetic calculations */ +/** P3M: Default for boundary condition in magnetic calculations. */ #ifndef P3M_EPSILON_MAGNETIC #define P3M_EPSILON_MAGNETIC 0.0 #endif /** P3M: Default for offset of first mesh point from the origin (left - down corner of the simulation box. */ + * down corner of the simulation box). + */ #ifndef P3M_MESHOFF #define P3M_MESHOFF 0.5 #endif /** P3M: Number of Brillouin zones taken into account - in the calculation of the optimal influence function (aliasing - sums). */ + * in the calculation of the optimal influence function (aliasing sums). + */ #ifndef P3M_BRILLOUIN #define P3M_BRILLOUIN 0 #endif /** P3M: Maximal mesh size that will be checked. The current setting - limits the memory consumption to below 1GB, which is probably - reasonable for a while. */ + * limits the memory consumption to below 1GB, which is probably + * reasonable for a while. + */ #ifndef P3M_MAX_MESH #define P3M_MAX_MESH 128 #endif -/** Whether to use the approximation of Abramowitz/Stegun - AS_erfc_part() for \f$\exp(d^2) erfc(d)\f$, or the C function erfc - in P3M and Ewald summation. */ +/** Whether to use the approximation of Abramowitz/Stegun @cite abramowitz65a + * AS_erfc_part() for \f$\exp(d^2) erfc(d)\f$, or the C function erfc + * in P3M and Ewald summation. + */ #ifndef USE_ERFC_APPROXIMATION #define USE_ERFC_APPROXIMATION 1 #endif @@ -83,38 +86,39 @@ #define ROUND_ERROR_PREC 1.0e-14 #endif -/** Tiny angle cutoff for sinus calculations */ +/** Tiny angle cutoff for sinus calculations. */ #ifndef TINY_SIN_VALUE #define TINY_SIN_VALUE 1e-10 #endif -/** Tiny angle cutoff for cosine calculations */ +/** Tiny angle cutoff for cosine calculations. */ #ifndef TINY_COS_VALUE #define TINY_COS_VALUE 0.9999999999 #endif -/** Tiny length cutoff */ +/** Tiny length cutoff. */ #ifndef TINY_LENGTH_VALUE #define TINY_LENGTH_VALUE 0.0001 #endif -/** Tiny oif elasticity cutoff */ +/** Tiny oif elasticity cutoff. */ #ifndef TINY_OIF_ELASTICITY_COEFFICIENT #define TINY_OIF_ELASTICITY_COEFFICIENT 1e-10 #endif -/** Small oif membrane collision cutoff */ +/** Small oif membrane collision cutoff. */ #ifndef SMALL_OIF_MEMBRANE_CUTOFF #define SMALL_OIF_MEMBRANE_CUTOFF 0.05 #endif -/** maximal number of iterations in the RATTLE algorithm before it bails out. */ +/** Maximal number of iterations in the RATTLE algorithm before it bails out. */ #ifndef SHAKE_MAX_ITERATIONS #define SHAKE_MAX_ITERATIONS 1000 #endif -/** maximal number of objects in the object-in-fluid framework. */ +/** Maximal number of objects in the object-in-fluid framework. */ #ifndef MAX_OBJECTS_IN_FLUID #define MAX_OBJECTS_IN_FLUID 10000 #endif -/* Mathematical constants, from gcc's math.hpp */ +/** @name Mathematical constants, from gcc's math.hpp */ +/*@{*/ #ifndef M_PI #define M_E 2.7182818284590452353602874713526625L /* e */ #define M_LOG2E 1.4426950408889634073599246810018921L /* log_2 e */ @@ -130,5 +134,6 @@ #define M_SQRT2 1.4142135623730950488016887242096981L /* sqrt(2) */ #define M_SQRT1_2 0.7071067811865475244008443621048490L /* 1/sqrt(2) */ #endif +/*@}*/ #endif diff --git a/src/core/Observable_stat.hpp b/src/core/Observable_stat.hpp index 86af2b71198..a15453957e4 100644 --- a/src/core/Observable_stat.hpp +++ b/src/core/Observable_stat.hpp @@ -22,10 +22,11 @@ #include struct Observable_stat { - /** Status flag for observable calculation. For 'analyze energy': 0 - re-initialize observable struct, else everything is fine, - calculation can start. For 'analyze pressure' and 'analyze - p_inst': 0 or !(1+v_comp) re-initialize, else all OK. */ + /** Status flag for observable calculation. For 'analyze energy': 0 + * re-initialize observable struct, else everything is fine, + * calculation can start. For 'analyze pressure' and 'analyze p_inst': + * 0 or !(1+v_comp) re-initialize, else all OK. + */ int init_status; /** Array for observables on each node. */ diff --git a/src/core/accumulators/Correlator.hpp b/src/core/accumulators/Correlator.hpp index e0cfd370736..7b18362ba51 100644 --- a/src/core/accumulators/Correlator.hpp +++ b/src/core/accumulators/Correlator.hpp @@ -16,21 +16,20 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see . */ -/* * @File Header file for the correlation class +/** @file * - * This module allow to compute correlations (and other two time averages) on + * This module computes correlations (and other two time averages) on * the fly and from files. * - * The basic idea is that the user can write arbitrary function A and B that + * The basic idea is that the user can write arbitrary functions A and B that * can depend on e.g. particle coordinates or whatever state of the MD box. - * These function can be vector valued, always indicated by dim_A and dim_B. + * These functions can be vector valued, always indicated by dim_A and dim_B. * * The way they can be correlated is can be formulated as a (vector valued) * operation on A and B. One example would be to calculate the product * component by component, and if A and B are both the particle velocities, - * then one would obtain - * { ... } + * then one would obtain { + * ... }. * * The idea of the algorithm is not to keep all As and Bs in memory, but * feed the algorithm (or the container object) successively by new As and Bs @@ -39,21 +38,21 @@ * discarded. * * To save memory, increase statistics and make the calculation possible for - * many orders of magnitude in time, the blocking algorithm in Frenkels book - * (p.91) + * many orders of magnitude in time, the blocking algorithm in + * @cite frenkel02b (algorithm 8, chapter 4.4.2) * is applied. Thus not all As and Bs of the whole "past" are stored but - * some of them are blocked. In this implementation always a blocking based - * on 2 is applied: All As and Bs not older than a certain tau_lin are stored - * as they were, those which are older are not entire stored, but only + * some of them are blocked. In this implementation, a blocking based on 2 is + * always applied: all As and Bs not older than a certain tau_lin are stored + * as they were, those which are older are not entirely stored, but only * their compressed average value. All As and Bs older than 2*tau_lin - * are compressed in blocks of four etc. + * are compressed in blocks of four, etc. * * This leads to a hierarchical "history": on level 0 the last tau_lin values - * are stored. This is done in a cyclic array: The newest is always appended + * are stored. This is done in a cyclic array: the newest is always appended * at the end of the array, and if the array length is reached, values * are appended at the beginning, overwriting older values. We therefore - * have to carry the index newest[0] which tells, what is the index of the - * newest value of A and B. + * have to carry the index newest[0] which tells, what is the index + * of the newest value of A and B. * * As soon as the first row of As and Bs is full, two of them are * compressed by calculating the arithmetic mean, and stored in the @@ -68,29 +67,29 @@ * This allows to have a "history" over many orders of magnitude * in time, without the full memory effort. * - * Correlations are only calculated on each level: For + * Correlations are only calculated on each level: for * tau=1,2,..,tau_lin the values are taken from level 1. * For tau=tau_lin, tau_lin+2, .., 2*tau_lin we take the values * from level 2. On level 2 we also have values for 0,..tau_lin-2, * but these are discarded as we have "better" estimates on level 1. * * The functions A and B can get extra arguments. This can e.g. be an - * integer describing to the "type" of the particles to be considered, + * integer describing the "type" of the particles to be considered, * or in case of files, it is a file_data_source pointer, which tells * the function from where to read. These arguments have to be a single * pointer to an object that carries all relevant information * to obtain A and B (from the MD Box or somewhere else). * * The correlation has to be initialized with all necessary information, - * i.e. all function pointers, the dimensions of A and B and their dimensions - * etc. - * When new A and B values should be processed, a single call of - * double_correlation_get_data(..) with the correlation object as an + * i.e. all function pointers, the dimensions of A and B and their dimensions, + * etc. When new A and B values should be processed, a single call of + * double_correlation_get_data() with the correlation object as an * argument is enough to update A and B and the correlation estimate. * * Eventually the correlation can be printed. - * - * + */ + +/* * There is a lot of stuff to do: * ============================== * -> Write input functions that take TCL arrays for A and B to @@ -110,8 +109,6 @@ * -> Write some nice samples * -> Write a destructor * -> Finally: write the users guide - * - * */ #ifndef _STATISTICS_CORRELATION_H #define _STATISTICS_CORRELATION_H @@ -132,41 +129,37 @@ namespace Accumulators { -/** The struct that is used to calculate correlations * +/** The main correlator class * - * Data organization: - * We use a ring-like way to manage the data: at the beginning we have a linear - * array - * which we fill from index 0 to tau_lin. The index newest[i] always indicates - * the latest - * entry of the hierarchic "past" For every new entry in is incremented and if - * tau_lin is reached, - * it starts again from the beginning. + * Data organization: + * We use a ring-like way to manage the data: at the beginning we have a + * linear array, which we fill from index 0 to @c tau_lin. The index + * newest[i] always indicates the latest entry of the hierarchic + * "past" For every new entry in is incremented and if @c tau_lin is reached, + * it starts again from the beginning. */ class Correlator : public AccumulatorBase { using obs_ptr = std::shared_ptr; public: - /** - * The initialization procedure for the correlation object. All important - * parameters have to be specified - * at the same time. They can not be change later, so every instance of the - * correlation class - * has to be fed with correct data from the very beginning. + /** The initialization procedure for the correlation object. All important + * parameters have to be specified at the same time. They cannot be changed + * later, so every instance of the correlation class has to be fed with + * correct data from the very beginning. * - * @param delta_N The number of time steps between subsequent updates - * @param tau_lin The linear part of the correlation function. - * @param tau_max maximal time delay tau to sample - * @param obs1 First observable to correlate - * @param obs2 Second observable to correlate - * @param corr_operation how to correlate the two observables A and B - * (this has no default) - * @param compress1_ how the A values should be compressed (usually - * the linear compression method) - * @param compress2_ how the B values should be compressed (usually - * the linear compression method) - * @param correlation_args_ optional arguments for the correlation function - * (currently only used when @p corr_operation is "fcs_acf") + * @param delta_N The number of time steps between subsequent updates + * @param tau_lin The linear part of the correlation function. + * @param tau_max maximal time delay tau to sample + * @param obs1 First observable to correlate + * @param obs2 Second observable to correlate + * @param corr_operation how to correlate the two observables A and B + * (this has no default) + * @param compress1_ how the A values should be compressed (usually + * the linear compression method) + * @param compress2_ how the B values should be compressed (usually + * the linear compression method) + * @param correlation_args_ optional arguments for the correlation function + * (currently only used when @p corr_operation is "fcs_acf") * */ Correlator(int tau_lin, double tau_max, int delta_N, std::string compress1_, @@ -188,37 +181,29 @@ class Correlator : public AccumulatorBase { public: /** The function to process a new datapoint of A and B * - * First the function finds out if it necessary to make some space for the new - * entries of A and B. - * Then, if necessary, it compresses old Values of A and B to make for the new - * value. Finally - * The new values of A and B are stored in A[newest[0]] and B[newest[0]], - * where the newest indices - * have been increased before. Finally the correlation estimate is updated. - * TODO: Not all - * the correlation estimates have to be updated. - * + * First the function finds out if it is necessary to make some space for + * the new entries of A and B. Then, if necessary, it compresses old values + * of A and B to make room for the new value. Finally, the new values of A + * and B are stored in A[newest[0]] and B[newest[0]], + * where the newest indices have been increased before. Finally, + * the correlation estimate is updated. + * TODO: Not all the correlation estimates have to be updated. */ void update() override; /** At the end of data collection, go through the whole hierarchy and - * correlate data left there - * - * This works pretty much the same as get_data, but does not feed on new data, - * just uses what - * is already available. + * correlate data left there. * + * This works pretty much the same as get_data, but does not feed on new + * data, just uses what is already available. */ int finalize(); /** Return an estimate of the integrated correlation time * * We calculate the correlation time for each dim_corr by normalizing the - * correlation, - * integrating it and finding out where C(tau)=tau; - * + * correlation, integrating it and finding out where C(tau)=tau; */ - int get_correlation_time(double *correlation_time); /** Return correlation result */ @@ -242,57 +227,56 @@ class Correlator : public AccumulatorBase { return corr_operation_name; } - /* Partial serialization of state that is not accessible - via the interface. */ + /** Partial serialization of state that is not accessible via the interface. + */ std::string get_internal_state() const; void set_internal_state(std::string const &); private: - unsigned int finalized; // non-zero of correlation is finalized - unsigned int t; // global time in number of frames + unsigned int finalized; ///< flag whether non-zero of correlation is finalized + unsigned int t; ///< global time in number of frames - Utils::Vector3d - m_correlation_args; // additional arguments, which the correlation - // may need (currently only used by fcs_acf) + Utils::Vector3d m_correlation_args; ///< additional arguments, which the + ///< correlation may need (currently + ///< only used by fcs_acf) - int hierarchy_depth; // maximum level of data compression - int m_tau_lin; // number of frames in the linear correlation + int hierarchy_depth; ///< maximum level of data compression + int m_tau_lin; ///< number of frames in the linear correlation int m_dim_corr; - double m_dt; // time interval at which samples arrive - double - m_tau_max; // maximum time, for which the correlation should be calculated + double m_dt; ///< time interval at which samples arrive + double m_tau_max; ///< maximum time for which the correlation should be + ///< calculated std::string compressA_name; std::string compressB_name; - std::string corr_operation_name; + std::string corr_operation_name; ///< Name of the correlation operator - int m_n_result; // the total number of result values + int m_n_result; ///< the total number of result values std::shared_ptr A_obs; std::shared_ptr B_obs; - std::vector tau; // time differences + std::vector tau; ///< time differences boost::multi_array, 2> A; boost::multi_array, 2> B; - boost::multi_array result; // output quantity + boost::multi_array result; ///< output quantity - // The actual allocated storage space - std::vector - n_sweeps; // number of correlation sweeps at a particular value of tau - std::vector n_vals; // number of data values already present at - // a particular value of tau - std::vector - newest; // index of the newest entry in each hierarchy level + /// number of correlation sweeps at a particular value of tau + std::vector n_sweeps; + /// number of data values already present at a particular value of tau + std::vector n_vals; + /// index of the newest entry in each hierarchy level + std::vector newest; - std::vector A_accumulated_average; // all A values are added up here - std::vector B_accumulated_average; // all B values are added up here - unsigned int n_data; // a counter to calculated averages and variances + std::vector A_accumulated_average; ///< all A values are added up here + std::vector B_accumulated_average; ///< all B values are added up here + unsigned int n_data; ///< a counter for calculated averages and variances double m_last_update; - unsigned int dim_A; // dimensionality of A - unsigned int dim_B; + unsigned int dim_A; ///< dimensionality of A + unsigned int dim_B; ///< dimensionality of B using correlation_operation_type = std::vector (*)(std::vector const &, diff --git a/src/core/actor/DipolarBarnesHut_cuda.cu b/src/core/actor/DipolarBarnesHut_cuda.cu index 575b3af54c2..a7859f2d732 100644 --- a/src/core/actor/DipolarBarnesHut_cuda.cu +++ b/src/core/actor/DipolarBarnesHut_cuda.cu @@ -17,6 +17,9 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see . */ +/** @file + * The method concept is revealed within @cite burtscher11a + */ #include "cuda_wrapper.hpp" @@ -30,13 +33,6 @@ typedef float dds_float; #ifdef DIPOLAR_BARNES_HUT -//////////////////////////////////////////////////////////////////////////////// -// The method concept is revealed within: M. Burtscher, K. Pingali, in: GPU -// Gems’11: GPU Computing Gems Emerald Edition, 2011. An Efficient CUDA -// Implementation of the Tree-Based Barnes Hut n-Body Algorithm -// http://iss.ices.utexas.edu/Publications/Papers/burtscher11.pdf -//////////////////////////////////////////////////////////////////////////////// - #include "DipolarBarnesHut_cuda.cuh" #define IND (blockDim.x * blockIdx.x + threadIdx.x) diff --git a/src/core/actor/Mmm1dgpuForce_cuda.cu b/src/core/actor/Mmm1dgpuForce_cuda.cu index c2deaafabe0..b4ce6b99f6c 100644 --- a/src/core/actor/Mmm1dgpuForce_cuda.cu +++ b/src/core/actor/Mmm1dgpuForce_cuda.cu @@ -28,7 +28,7 @@ #ifdef MMM1D_GPU -// the code is mostly multi-GPU capable, but Espresso is not yet +// the code is mostly multi-GPU capable, but ESPResSo is not yet const int deviceCount = 1; float multigpu_factors[] = {1.0}; #define cudaSetDevice(d) diff --git a/src/core/bonded_interactions/bonded_interactions.dox b/src/core/bonded_interactions/bonded_interactions.dox index 9fe42a0e1d0..c769a7b8d39 100644 --- a/src/core/bonded_interactions/bonded_interactions.dox +++ b/src/core/bonded_interactions/bonded_interactions.dox @@ -210,9 +210,7 @@ * * @subsection bondedIA_angle_force General expressions for the forces * - * This section uses the particle force expressions derived in - * Swope, Ferguson *J. Comput. Chem.* **1992**, *13*(5): 585-594, - * doi:[10.1002/jcc.540130508](https://doi.org/10.1002/jcc.540130508). + * This section uses the particle force expressions derived in @cite swope92a. * * The gradient of the potential at particle @f$ i @f$ is given by the chain * rule in equation 6: diff --git a/src/core/bonded_interactions/dihedral.hpp b/src/core/bonded_interactions/dihedral.hpp index 47c9b9c2333..656691e2d8e 100644 --- a/src/core/bonded_interactions/dihedral.hpp +++ b/src/core/bonded_interactions/dihedral.hpp @@ -22,7 +22,7 @@ #define DIHEDRAL_H /** \file * Routines to calculate the dihedral energy or/and - * force for a particle quadruple. Note that usage of dihedrals + * force for a particle quadruple. Note that usage of dihedrals * increases the interaction range of bonded interactions to 2 times * the maximal bond length! * diff --git a/src/core/cells.cpp b/src/core/cells.cpp index 19f2eebf344..bc6c7293b3b 100644 --- a/src/core/cells.cpp +++ b/src/core/cells.cpp @@ -457,7 +457,7 @@ void cells_update_ghosts() { ? CELL_GLOBAL_EXCHANGE : CELL_NEIGHBOR_EXCHANGE; - /* Communication step: number of ghosts and ghost information */ + /* Communication step: number of ghosts and ghost information */ cells_resort_particles(global); } else diff --git a/src/core/cells.hpp b/src/core/cells.hpp index 05f4b56bef5..af8bfdba5f6 100644 --- a/src/core/cells.hpp +++ b/src/core/cells.hpp @@ -24,37 +24,33 @@ * This file contains everything related to the cell structure / cell * system. * - * The cell system (\ref Cell Structure) describes how particles are + * The cell system (\ref CellStructure) describes how particles are * distributed on the cells and how particles of different cells * (regardless if they reside on the same or different nodes) * interact with each other. The following cell systems are implemented: * - *
    - *
  • domain decomposition: The simulation box is divided spatially + * - domain decomposition: The simulation box is divided spatially * into cells (see \ref domain_decomposition.hpp). This is suitable for * short range interactions. - *
  • nsquare: The particles are distributed equally on all nodes + * - nsquare: The particles are distributed equally on all nodes * regardless their spatial position (see \ref nsquare.hpp). This is - * suitable for long range interactions that can not be treated by a + * suitable for long range interactions that cannot be treated by a * special method like P3M (see \ref p3m.hpp). - *
  • layered: in x and y directions, it uses a nsquared type of + * - layered: in x and y directions, it uses a nsquared type of * interaction calculation, but in z it has a domain decomposition * into layers. - *
* * Some structures are common to all cell systems: * - *
    - *
  • All cells, real cells as well as ghost cells, are stored in the + * - All cells, real cells as well as ghost cells, are stored in the * vector \ref cells::cells whose size has to be changed with * \ref realloc_cells. - *
  • There are two lists of cell pointers to access particles and + * - There are two lists of cell pointers to access particles and * ghost particles on a node: \ref local_cells contains pointers to * all cells containing the particles physically residing on that * node. \ref ghost_cells contains pointers to all cells containing * the ghost particles of that node. The size of these lists has to be * changed with \ref realloc_cellplist - *
*/ #include diff --git a/src/core/cuda_common_cuda.cu b/src/core/cuda_common_cuda.cu index 3a61f07c742..9d83d8327f8 100644 --- a/src/core/cuda_common_cuda.cu +++ b/src/core/cuda_common_cuda.cu @@ -263,7 +263,7 @@ void gpu_init_particle_comm() { if (cuda_check_gpu(0) != ES_OK) { runtimeWarningMsg() << "CUDA device 0 is not capable of running ESPResSo but is used " - "by default. Espresso has detected a CUDA capable card but it " + "by default. ESPResSo has detected a CUDA capable card but it " "is not the one used by ESPResSo by default. Please set the " "GPU to use by setting System.cuda_init_handle.device. A list " "of available GPUs is available through " @@ -291,7 +291,7 @@ void copy_part_data_to_gpu(ParticleRange particles) { global_part_vars_host.number_of_particles) { cuda_mpi_get_particles(particles, particle_data_host); - /* get espresso md particle values*/ + /* get espressomd particle values */ if (this_node == 0) cudaMemcpyAsync(particle_data_device, particle_data_host, global_part_vars_host.number_of_particles * diff --git a/src/core/cuda_init.cpp b/src/core/cuda_init.cpp index 2a21af76567..990828b7f1c 100644 --- a/src/core/cuda_init.cpp +++ b/src/core/cuda_init.cpp @@ -28,9 +28,8 @@ #include #include -/** Helper class force device set +/** Helper class force device set. */ - struct CompareDevices { bool operator()(const EspressoGpuDevice &a, const EspressoGpuDevice &b) const { @@ -43,12 +42,11 @@ struct CompareDevices { } }; -/** Gather list of CUDA devices on all nodes on the master node - It relies on MPI_Get_processor_name() to get a unique identifier of - the physical node, as opposed to the logical rank of which there can - be more than one on one node. +/** Gather list of CUDA devices on all nodes on the master node. + * It relies on MPI_Get_processor_name() to get a unique identifier of + * the physical node, as opposed to the logical rank of which there can + * be more than one on one node. */ - std::vector cuda_gather_gpus() { int n_gpus = cuda_get_n_gpus(); char proc_name[MPI_MAX_PROCESSOR_NAME]; @@ -76,7 +74,7 @@ std::vector cuda_gather_gpus() { } } - /** Update n_gpus to number of usable devices */ + /* Update n_gpus to number of usable devices */ n_gpus = devices.size(); if (this_node == 0) { diff --git a/src/core/cuda_init.hpp b/src/core/cuda_init.hpp index c26d9b6bd6c..1ab5a2ea553 100644 --- a/src/core/cuda_init.hpp +++ b/src/core/cuda_init.hpp @@ -25,10 +25,10 @@ #include -/** Struct to hold information relevant to Espresso - about GPUs. Should contain only fixed length plain - old datatypes, as it is intended for MPI communication */ - +/** Struct to hold information relevant to ESPResSo + * about GPUs. Should contain only fixed-length plain + * old datatypes, as it is intended for MPI communication. + */ struct EspressoGpuDevice { /* Local CUDA device id */ int id; @@ -51,65 +51,63 @@ struct EspressoGpuDevice { */ void cuda_init(); -/** get the number of CUDA devices. - - @return the number of GPUs, or -1 if CUDA could not be - initialized. The error message from CUDA can be found in \ref - cuda_error. -*/ +/** Get the number of CUDA devices. + * + * @return the number of GPUs, or -1 if CUDA could not be + * initialized. The error message from CUDA can be found in \ref + * cuda_error. + */ int cuda_get_n_gpus(); -/** check that a given GPU is capable of what we need, that is, at - least compute capability 1.1. - - @param dev CUDA device number - @return \ref ES_OK if and only if the GPU with the given id is - usable for CUDA computations. Only devices with compute - capability of 1.1 or higher are ok, since atomic operations are - required for CUDA-LB. -*/ +/** Check that a given GPU is capable of what we need, that is, at + * least compute capability 1.1. + * + * @param dev CUDA device number + * @return \ref ES_OK if and only if the GPU with the given id is + * usable for CUDA computations. Only devices with compute + * capability of 1.1 or higher are ok, since atomic operations are + * required for CUDA-LB. + */ int cuda_check_gpu(int dev); -/** get the name of a CUDA device. - - @param dev the CUDA device number to ask the name for - @param name a buffer to write the name to, at least 64 characters -*/ +/** Get the name of a CUDA device. + * + * @param dev the CUDA device number to ask the name for + * @param name a buffer to write the name to, at least 64 characters + */ void cuda_get_gpu_name(int dev, char name[64]); -/** choose a device for future CUDA computations. - - @param dev the device to use - @return \ref ES_OK on success, \ref ES_ERROR else. The error - message from CUDA can be found in \ref cuda_error. -*/ +/** Choose a device for future CUDA computations. + * + * @param dev the device to use + * @return \ref ES_OK on success, \ref ES_ERROR else. The error + * message from CUDA can be found in \ref cuda_error. + */ int cuda_set_device(int dev); -/** get the current CUDA device. - - @return the current device's number or -1 if an error occurred. The error - message from CUDA can be found in \ref cuda_error. -*/ +/** Get the current CUDA device. + * + * @return the current device's number or -1 if an error occurred. The error + * message from CUDA can be found in \ref cuda_error. + */ int cuda_get_device(); /** Test if actual CUDA device works. - @return \ref ES_OK on success, \ref ES_ERROR else. - The error message from CUDA can be found in \ref cuda_error. -*/ - + * @return \ref ES_OK on success, \ref ES_ERROR else. + * The error message from CUDA can be found in \ref cuda_error. + */ int cuda_test_device_access(); /** Gather unique list of CUDA devices on all nodes - @return vector of device on master, empty vector on other nodes. -*/ - + * @return vector of device on master, empty vector on other nodes. + */ std::vector cuda_gather_gpus(); /** Get properties of a CUDA device */ int cuda_get_device_props(int dev, EspressoGpuDevice &d); -/** current error message of CUDA. */ +/** Current error message of CUDA. */ extern const char *cuda_error; #endif diff --git a/src/core/domain_decomposition.hpp b/src/core/domain_decomposition.hpp index 8c21ddf1081..caa2ee6971c 100644 --- a/src/core/domain_decomposition.hpp +++ b/src/core/domain_decomposition.hpp @@ -98,7 +98,7 @@ extern int max_num_cells; /** Minimal number of cells per node. This is mainly to avoid excessively large * numbers of particles per cell, which will result in really large Verlet - * lists and eventually crash Espresso. + * lists and eventually crash ESPResSo. */ extern int min_num_cells; @@ -123,7 +123,7 @@ void dd_on_geometry_change(int flags, const Utils::Vector3i &grid, /** Initialize the topology. The argument is a list of cell pointers, * containing particles that have to be sorted into new cells. The - * particles might not belong to this node. This procedure is used + * particles might not belong to this node. This procedure is used * when particle data or cell structure has changed and the cell * structure has to be reinitialized. This also includes setting up * the cell_structure array. diff --git a/src/core/dpd.hpp b/src/core/dpd.hpp index 5da75ace6e2..e129588f77c 100644 --- a/src/core/dpd.hpp +++ b/src/core/dpd.hpp @@ -21,10 +21,9 @@ #ifndef CORE_DPD_HPP #define CORE_DPD_HPP /** \file - * Routines to use dpd as thermostat or pair force - * T. Soddemann, B. Duenweg and K. Kremer, Phys. Rev. E 68, 046702 (2003) + * Routines to use dpd as thermostat or pair force @cite soddemann03a * - * Implementation in forces.cpp. + * Implementation in @ref dpd.cpp. */ #include "config.hpp" diff --git a/src/core/electrostatics_magnetostatics/elc.cpp b/src/core/electrostatics_magnetostatics/elc.cpp index 9258dff0004..cd768a9e2f9 100644 --- a/src/core/electrostatics_magnetostatics/elc.cpp +++ b/src/core/electrostatics_magnetostatics/elc.cpp @@ -229,6 +229,9 @@ void distribute(int size) { /* dipole terms */ /*****************************************************************/ +/** Calculate the dipole force. + * See @cite yeh99a. + */ static void add_dipole_force(const ParticleRange &particles) { double pref = coulomb.prefactor * 4 * M_PI * ux * uy * uz; int size = 3; @@ -270,7 +273,7 @@ static void add_dipole_force(const ParticleRange &particles) { distribute(size); - // Yeh+Berkowitz dipole term + // Yeh + Berkowitz dipole term @cite yeh99a field_tot = gblcblk[0]; // Const. potential contribution @@ -290,6 +293,9 @@ static void add_dipole_force(const ParticleRange &particles) { } } +/** Calculate the dipole energy. + * See @cite yeh99a. + */ static double dipole_energy(const ParticleRange &particles) { double pref = coulomb.prefactor * 2 * M_PI * ux * uy * uz; int size = 7; @@ -332,7 +338,7 @@ static double dipole_energy(const ParticleRange &particles) { distribute(size); - // Yeh + Berkowitz term + // Yeh + Berkowitz term @cite yeh99a double eng = 2 * pref * (Utils::sqr(gblcblk[2]) + gblcblk[2] * gblcblk[3]); if (!elc_params.neutralize) { @@ -351,7 +357,7 @@ static double dipole_energy(const ParticleRange &particles) { } /* counter the P3M homogeneous background contribution to the - boundaries. We never need that, since a homogeneous background + boundaries. We never need that, since a homogeneous background spanning the artificial boundary layers is aphysical. */ eng += pref * (-(gblcblk[1] * gblcblk[4] + gblcblk[0] * gblcblk[5]) - (1. - 2. / 3.) * gblcblk[0] * gblcblk[1] * diff --git a/src/core/electrostatics_magnetostatics/fft.cpp b/src/core/electrostatics_magnetostatics/fft.cpp index 95eb36cb14f..064f7ebd0be 100644 --- a/src/core/electrostatics_magnetostatics/fft.cpp +++ b/src/core/electrostatics_magnetostatics/fft.cpp @@ -21,7 +21,7 @@ /** \file * * Routines, row decomposition, data structures and communication for the - * 3D-FFT. + * 3D-FFT. * */ @@ -169,12 +169,12 @@ find_comm_groups(Utils::Vector3i const &grid1, Utils::Vector3i const &grid2, return group; } -/** Calculate the local fft mesh. Calculate the local mesh (loc_mesh) - * of a node at position (n_pos) in a node grid (n_grid) for a global - * mesh of size (mesh) and a mesh offset (mesh_off (in mesh units)) - * and store also the first point (start) of the local mesh. +/** Calculate the local fft mesh. Calculate the local mesh (@p loc_mesh) + * of a node at position (@p n_pos) in a node grid (@p n_grid) for a global + * mesh of size (@p mesh) and a mesh offset (@p mesh_off (in mesh units)) + * and store also the first point (@p start) of the local mesh. * - * \param[in] n_pos Position of the node in n_grid. + * \param[in] n_pos Position of the node in @p n_grid. * \param[in] n_grid node grid. * \param[in] mesh global mesh dimensions. * \param[in] mesh_off global mesh offset (see \ref p3m_data_struct). @@ -206,23 +206,24 @@ int calc_local_mesh(const int *n_pos, const int *n_grid, const int *mesh, } /** Calculate a send (or recv.) block for grid communication during a - * decomposition change. Calculate the send block specification + * decomposition change. Calculate the send block specification * (block = lower left corner and upper right corner) which a node at - * position (pos1) in the actual node grid (grid1) has to send to - * another node at position (pos2) in the desired node grid - * (grid2). The global mesh, subject to communication, is specified - * via its size (mesh) and its mesh offset (mesh_off (in mesh - * units)). + * position (@p pos1) in the actual node grid (@p grid1) has to send to + * another node at position (@p pos2) in the desired node grid (@p grid2). + * The global mesh, subject to communication, is specified via its size + * (@p mesh) and its mesh offset (@p mesh_off (in mesh units)). * * For the calculation of a receive block you have to change the arguments in - * the following way:
pos1 - position of receiving node in the desired - * node grid.
grid1 - desired node grid.
pos2 - position of the node - * you intend to receive the data from in the actual node grid.
grid2 - - * actual node grid.
+ * the following way: + * - @p pos1: position of receiving node in the desired node grid. + * - @p grid1: desired node grid. + * - @p pos2: position of the node you intend to receive the data from in the + * actual node grid. + * - @p grid2: actual node grid. * - * \param[in] pos1 Position of send node in grid1. + * \param[in] pos1 Position of send node in @p grid1. * \param[in] grid1 node grid 1. - * \param[in] pos2 Position of recv node in grid2. + * \param[in] pos2 Position of recv node in @p grid2. * \param[in] grid2 node grid 2. * \param[in] mesh global mesh dimensions. * \param[in] mesh_off global mesh offset (see \ref p3m_data_struct). @@ -249,24 +250,22 @@ int calc_send_block(const int *pos1, const int *grid1, const int *pos2, return size; } -/** pack a block with dimensions (size[0] * size[1] * size[2]) starting - * at start[3] of an input 3d-grid with dimension dim[3] into an - * output 3d-grid with dimensions (size[2] * size[0] * size[1]) with - * a simultaneous one-fold permutation of the indices. +/** Pack a block with dimensions size[0] * size[1] * size[2] starting + * at @p start of an input 3D-grid with dimension @p dim into an output + * 3D-grid with dimensions size[2] * size[0] * size[1] with + * a simultaneous one-fold permutation of the indices. The permutation is + * defined as: slow_in -> fast_out, mid_in ->slow_out, fast_in -> mid_out. * - * The permutation is defined as: - * slow_in -> fast_out, mid_in ->slow_out, fast_in -> mid_out + * An element (i0_in, i1_in, i2_in) is then + * (i0_out = i1_in-start[1], i1_out = i2_in-start[2], + * i2_out = i0_in-start[0]) and for the linear indices we have: + * - li_in = i2_in + size[2] * (i1_in + (size[1]*i0_in)) + * - li_out = i2_out + size[0] * (i1_out + (size[2]*i0_out)) * - * An element (i0_in , i1_in , i2_in ) is then - * (i0_out = i1_in-start[1], i1_out = i2_in-start[2], i2_out = i0_in-start[0]) - * and for the linear indices we have:
li_in = - * i2_in + size[2] * (i1_in + (size[1]*i0_in))
li_out = i2_out + - * size[0] * (i1_out + (size[2]*i0_out)) + * For index definition see \ref fft_pack_block. * - * For index definition see \ref fft_pack_block. - * - * \param[in] in pointer to input 3d-grid. - * \param[out] out pointer to output 3d-grid (block). + * \param[in] in input 3D-grid. + * \param[out] out output 3D-grid (block). * \param[in] start start index of the block in the in-grid. * \param[in] size size of the block (=dimension of the out-grid). * \param[in] dim size of the in-grid. @@ -275,7 +274,7 @@ int calc_send_block(const int *pos1, const int *grid1, const int *pos2, void pack_block_permute1(double const *const in, double *const out, const int *start, const int *size, const int *dim, int element) { - /* slow,mid and fast changing indices for input grid */ + /* slow, mid and fast changing indices for input grid */ int s, m, f, e; /* linear index of in grid, linear index of out grid */ int li_in, li_out = 0; @@ -303,24 +302,22 @@ void pack_block_permute1(double const *const in, double *const out, } } -/** pack a block with dimensions (size[0] * size[1] * size[2]) starting - * at start[3] of an input 3d-grid with dimension dim[3] into an - * output 3d-grid with dimensions (size[2] * size[0] * size[1]), this - * is a simultaneous two-fold permutation of the indices. - * - * The permutation is defined as: - * slow_in -> mid_out, mid_in ->fast_out, fast_in -> slow_out +/** Pack a block with dimensions size[0] * size[1] * size[2] starting + * at @p start of an input 3D-grid with dimension @p dim into an output + * 3D-grid with dimensions size[2] * size[0] * size[1] with + * a simultaneous two-fold permutation of the indices. The permutation is + * defined as: slow_in -> mid_out, mid_in ->fast_out, fast_in -> slow_out. * - * An element (i0_in , i1_in , i2_in ) is then - * (i0_out = i2_in-start[2], i1_out = i0_in-start[0], i2_out = i1_in-start[1]) - * and for the linear indices we have:
li_in = - * i2_in + size[2] * (i1_in + (size[1]*i0_in))
li_out = i2_out + - * size[0] * (i1_out + (size[2]*i0_out)) + * An element (i0_in, i1_in, i2_in) is then + * (i0_out = i2_in-start[2], i1_out = i0_in-start[0], + * i2_out = i1_in-start[1]) and for the linear indices we have: + * - li_in = i2_in + size[2] * (i1_in + (size[1]*i0_in)) + * - li_out = i2_out + size[0] * (i1_out + (size[2]*i0_out)) * - * For index definition see \ref fft_pack_block. + * For index definition see \ref fft_pack_block. * - * \param[in] in pointer to input 3d-grid. - * \param[out] out pointer to output 3d-grid (block). + * \param[in] in input 3D-grid. + * \param[out] out output 3D-grid (block). * \param[in] start start index of the block in the in-grid. * \param[in] size size of the block (=dimension of the out-grid). * \param[in] dim size of the in-grid. @@ -329,7 +326,7 @@ void pack_block_permute1(double const *const in, double *const out, void pack_block_permute2(double const *const in, double *const out, const int *start, const int *size, const int *dim, int element) { - /* slow,mid and fast changing indices for input grid */ + /* slow, mid and fast changing indices for input grid */ int s, m, f, e; /* linear index of in grid, linear index of out grid */ int li_in, li_out = 0; @@ -423,14 +420,14 @@ void back_grid_comm(fft_forw_plan plan_f, fft_back_plan plan_b, } } -/** calculate 'best' mapping between a 2d and 3d grid. - * This we need for the communication from 3d domain decomposition - * to 2d row decomposition. - * The dimensions of the 2d grid are resorted, if necessary, in a way - * that they are multiples of the 3d grid dimensions. - * \param g3d 3d grid. - * \param g2d 2d grid. - * \param mult factors between 3d and 2d grid dimensions +/** Calculate 'best' mapping between a 2D and 3D grid. + * Required for the communication from 3D domain decomposition + * to 2D row decomposition. + * The dimensions of the 2D grid are resorted, if necessary, in a way + * that they are multiples of the 3D grid dimensions. + * \param g3d 3D grid. + * \param g2d 2D grid. + * \param mult factors between 3D and 2D grid dimensions * \return index of the row direction [0,1,2]. */ int map_3don2d_grid(int const g3d[3], int g2d[3], int mult[3]) { @@ -478,7 +475,7 @@ int map_3don2d_grid(int const g3d[3], int g2d[3], int mult[3]) { return row_dir; } -/** calculate most square 2d grid. */ +/** Calculate most square 2D grid. */ void calc_2d_grid(int n, int grid[3]) { for (auto i = static_cast(std::sqrt(n)); i >= 1; i--) { if (n % i == 0) { diff --git a/src/core/electrostatics_magnetostatics/icc.cpp b/src/core/electrostatics_magnetostatics/icc.cpp index 295d8c7b22c..4bbeb488acf 100644 --- a/src/core/electrostatics_magnetostatics/icc.cpp +++ b/src/core/electrostatics_magnetostatics/icc.cpp @@ -20,9 +20,11 @@ */ /** \file - Detailed Information about the method is included in the corresponding header - file \ref icc.hpp. -*/ + * Functions to compute the electric field acting on the induced charges, + * excluding forces other than the electrostatic ones. Detailed information + * about the ICC* method is included in the corresponding header file + * \ref icc.hpp. + */ #include "icc.hpp" @@ -51,22 +53,21 @@ iccp3m_struct iccp3m_cfg; -/* functions that are used in icc* to compute the electric field acting on the - * induced charges, excluding forces other than the electrostatic ones */ void init_forces_iccp3m(const ParticleRange &particles, const ParticleRange &ghosts_particles); -/** Calculation of the electrostatic forces between source charges (= real - * charges) and wall charges. For each electrostatic method the proper functions - * for short and long range parts are called. Long Range Parts are calculated - * directly, short range parts need helper functions according to the particle - * data organisation. A modified version of \ref force_calc in \ref forces.hpp. +/** Calculate the electrostatic forces between source charges (= real charges) + * and wall charges. For each electrostatic method, the proper functions + * for short- and long-range parts are called. Long-range parts are calculated + * directly, short-range parts need helper functions according to the particle + * data organisation. This is a modified version of \ref force_calc. */ void force_calc_iccp3m(const ParticleRange &particles, const ParticleRange &ghost_particles); -/** Variant of add_non_bonded_pair_force where only Coulomb - * contributions are calculated */ +/** Variant of @ref add_non_bonded_pair_force where only Coulomb + * contributions are calculated + */ inline void add_non_bonded_pair_force_iccp3m(Particle &p1, Particle &p2, Utils::Vector3d const &d, double dist, double dist2) { diff --git a/src/core/electrostatics_magnetostatics/icc.hpp b/src/core/electrostatics_magnetostatics/icc.hpp index 97365a7e231..b44616239e0 100644 --- a/src/core/electrostatics_magnetostatics/icc.hpp +++ b/src/core/electrostatics_magnetostatics/icc.hpp @@ -21,30 +21,27 @@ /** \file * * ICCP3M is a method that allows to take into account the influence - * of arbitrarily shaped dielectric interfaces. The dielectric + * of arbitrarily shaped dielectric interfaces. The dielectric * properties of a dielectric medium in the bulk of the simulation * box are taken into account by reproducing the jump in the electric - * field at the inface with charge surface segments. The charge + * field at the interface with charge surface segments. The charge * density of the surface segments have to be determined - * self-consistently using an iterative scheme. It can at presently - * - despite its name - be used with P3M, ELCP3M, MMM2D and MMM1D. For - * details see: S. Tyagi, M. Suzen, M. Sega, M. Barbosa, S. S. Kantorovich, - * C. Holm: An iterative, fast, linear-scaling method for computing induced - * charges on arbitrary dielectric boundaries, J. Chem. Phys. 2010, 132, - * p. 154112, doi:10.1063/1.3376011 + * self-consistently using an iterative scheme. It can at present - + * despite its name - be used with P3M, ELCP3M, MMM2D and MMM1D. For + * details see: @cite tyagi10a * - * To set up ICCP3M first the dielectric boundary has to be modelled - * by espresso particles 0..n where n has to be passed as a parameter + * To set up ICCP3M, first the dielectric boundary has to be modeled + * by ESPResSo particles 0..n where n has to be passed as a parameter * to ICCP3M. This is still a bit inconvenient, as it forces the user * to reserve the first n particle ids to wall charges, but as the - * other parts of espresso do not suffer from a limitation like this, + * other parts of ESPResSo do not suffer from a limitation like this, * it can be tolerated. * * For the determination of the induced charges only the forces * acting on the induced charges has to be determined. As P3M and the * other Coulomb solvers calculate all mutual forces, the force * calculation was modified to avoid the calculation of the short - * range part of the source-source force calculation. For different + * range part of the source-source force calculation. For different * particle data organisation schemes this is performed differently. */ @@ -58,21 +55,21 @@ #include #include -/* iccp3m data structures*/ +/** ICCP3M data structure */ struct iccp3m_struct { - int n_ic; /* Last induced id (can not be smaller then 2) */ - int num_iteration = 30; /* Number of max iterations */ - double eout = 1; /* Dielectric constant of the bulk */ - std::vector areas; /* Array of area of the grid elements */ - std::vector - ein; /* Array of dielectric constants at each surface element */ - std::vector sigma; /* Surface Charge density */ - double convergence = 1e-2; /* Convergence criterion */ - std::vector normals; /* Surface normal vectors */ - Utils::Vector3d ext_field = {0, 0, 0}; /* External field */ - double relax = 0.7; /* relaxation parameter for iterative */ - int citeration = 0; /* current number of iterations*/ - int first_id = 0; + int n_ic; /**< Last induced id (cannot be smaller than 2) */ + int num_iteration = 30; /**< Number of max iterations */ + double eout = 1; /**< Dielectric constant of the bulk */ + std::vector areas; /**< Array of area of the grid elements */ + std::vector ein; /**< Array of dielectric constants at each surface + element */ + std::vector sigma; /**< Surface charge density */ + double convergence = 1e-2; /**< Convergence criterion */ + std::vector normals; /**< Surface normal vectors */ + Utils::Vector3d ext_field = {0, 0, 0}; /**< External field */ + double relax = 0.7; /**< relaxation parameter for iteration */ + int citeration = 0; /**< current number of iterations */ + int first_id = 0; /**< id of the first particle in the dielectric boundary */ template void serialize(Archive &ar, long int /* version */) { @@ -90,10 +87,10 @@ struct iccp3m_struct { ar &citeration; } }; -extern iccp3m_struct iccp3m_cfg; /* global variable with ICCP3M configuration */ +extern iccp3m_struct iccp3m_cfg; /**< Global state of the ICCP3M solver */ /** The main iterative scheme, where the surface element charges are calculated - * self-consistently. + * self-consistently. */ int iccp3m_iteration(const ParticleRange &particles, const ParticleRange &ghost_particles); diff --git a/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp b/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp index 7cfe5c4597f..85f823ee5d9 100644 --- a/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp +++ b/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp @@ -301,7 +301,7 @@ magnetic_dipolar_direct_sum_calculations(bool force_flag, bool energy_flag, } /* of j and i */ } /* end of the area of calculation */ - /* set the forces, and torques of the particles within Espresso */ + /* set the forces, and torques of the particles within ESPResSo */ if (force_flag) { int dip_particles2 = 0; diff --git a/src/core/electrostatics_magnetostatics/mdlc_correction.cpp b/src/core/electrostatics_magnetostatics/mdlc_correction.cpp index b1d9c919c23..4977284c31c 100644 --- a/src/core/electrostatics_magnetostatics/mdlc_correction.cpp +++ b/src/core/electrostatics_magnetostatics/mdlc_correction.cpp @@ -19,18 +19,6 @@ * along with this program. If not, see . */ /** \file - * code for calculating the MDLC (magnetic dipolar layer correction). - * Purpose: get the corrections for dipolar 3D algorithms - * when applied to a slab geometry and dipolar - * particles. DLC & co - * Article: A. Brodka, Chemical Physics Letters 400, 62-67 (2004). - * - * We also include a tuning function that returns the - * cut-off necessary to attend a certain accuracy. - * - * Restrictions: the slab must be such that the z is the short - * direction. Otherwise we get trash. - * */ #include "electrostatics_magnetostatics/mdlc_correction.hpp" @@ -76,7 +64,7 @@ inline double g2_DLC_dip(double g, double x) { return a; } -/* Compute Mx, My, Mz and Mtotal */ +/** Compute Mx, My, Mz and Mtotal */ double slab_dip_count_mu(double *mt, double *mx, double *my, const ParticleRange &particles) { Utils::Vector3d node_sums{}; @@ -103,10 +91,8 @@ double slab_dip_count_mu(double *mt, double *mx, double *my, return Mz; } -/** - Compute the dipolar DLC corrections for forces and torques. - Algorithm implemented accordingly to the paper of A. Brodka, Chem. Phys. - Lett. 400, 62-67, (2004). +/** Compute the dipolar DLC corrections for forces and torques. + * Algorithm implemented accordingly to @cite brodka04a. */ double get_DLC_dipolar(int kcut, std::vector &fs, std::vector &ts, @@ -264,10 +250,8 @@ double get_DLC_dipolar(int kcut, std::vector &fs, return energy; } -/** - Compute the dipolar DLC corrections - Algorithm implemented accordingly to the paper of A. Brodka, Chem. Phys. - Lett. 400, 62-67, (2004). +/** Compute the dipolar DLC corrections + * Algorithm implemented accordingly to @cite brodka04a. */ double get_DLC_energy_dipolar(int kcut, const ParticleRange &particles) { @@ -369,13 +353,11 @@ void add_mdlc_force_corrections(const ParticleRange &particles) { get_DLC_dipolar(dip_DLC_kcut, dip_DLC_f, dip_DLC_t, particles); // Now we compute the correction like Yeh and Klapp to take into account - // the fact that you are using a - // 3D PBC method which uses spherical summation instead of slab-wise - // summation. - // Slab-wise summation is the one - // required to apply DLC correction. This correction is often called SDC = - // Shape Dependent Correction. - // See Brodka, Chem. Phys. Lett. 400, 62, (2004). + // the fact that you are using a 3D PBC method which uses spherical + // summation instead of slab-wise summation. + // Slab-wise summation is the one required to apply DLC correction. + // This correction is often called SDC = Shape Dependent Correction. + // See @cite brodka04a. mz = slab_dip_count_mu(&mtot, &mx, &my, particles); @@ -431,13 +413,11 @@ double add_mdlc_energy_corrections(const ParticleRange &particles) { // \n",dip_DLC_energy); // Now we compute the correction like Yeh and Klapp to take into account - // the fact that you are using a - // 3D PBC method which uses spherical summation instead of slab-wise - // summation. - // Slab-wise summation is the one - // required to apply DLC correction. This correction is often called SDC = - // Shape Dependent Correction. - // See Brodka, Chem. Phys. Lett. 400, 62, (2004). + // the fact that you are using a 3D PBC method which uses spherical + // summation instead of slab-wise summation. + // Slab-wise summation is the one required to apply DLC correction. + // This correction is often called SDC = Shape Dependent Correction. + // See @cite brodka04a. mz = slab_dip_count_mu(&mtot, &mx, &my, particles); @@ -465,21 +445,16 @@ double add_mdlc_energy_corrections(const ParticleRange &particles) { return 0.0; } -/** - Subroutine to compute the cut-off (NCUT) necessary in the DLC dipolar part - to get a certain accuracy (acc). We assume particles to have all them a - same - value of the dipolar momentum modulus (mu_max). mu_max is taken as the - largest value of - mu inside the system. If we assume the gap has a width gap_size (within which - there is no particles) - - Lz=h+gap_size - - BE CAREFUL: (1) We assume the short distance for the slab to be in the Z - direction - (2) You must also tune the other 3D method to the same accuracy, otherwise - it has no sense to have a good accurate result for DLC-dipolar. +/** Compute the cut-off in the DLC dipolar part to get a certain accuracy. + * We assume particles to have all the same value of the dipolar momentum + * modulus (@ref mu_max). @ref mu_max is taken as the largest value of mu + * inside the system. If we assume the gap has a width @c gap_size (within + * which there is no particles): Lz = h + gap_size + * + * BE CAREFUL: + * 1. We assume the short distance for the slab to be in the *z*-direction + * 2. You must also tune the other 3D method to the same accuracy, otherwise + * it makes no sense to have an accurate result for DLC-dipolar. */ int mdlc_tune(double error) { double de, n, gc, lz, lx, a, fa1, fa2, fa0, h; diff --git a/src/core/electrostatics_magnetostatics/mdlc_correction.hpp b/src/core/electrostatics_magnetostatics/mdlc_correction.hpp index 768899d5ea0..9db866c0c51 100644 --- a/src/core/electrostatics_magnetostatics/mdlc_correction.hpp +++ b/src/core/electrostatics_magnetostatics/mdlc_correction.hpp @@ -22,13 +22,15 @@ * main header-file for MDLC (magnetic dipolar layer correction). * * Developer: Joan J. Cerda. + * * Purpose: get the corrections for dipolar 3D algorithms * when applied to a slab geometry and dipolar * particles. DLC & co - * Article: A. Brodka, Chemical Physics Letters 400, 62-67 (2004). * - * We also include a tuning function that returns the - * cut-off necessary to attend a certain accuracy. + * Article: @cite brodka04a + * + * We also include a tuning function that returns the + * cut-off necessary to attend a certain accuracy. * * Restrictions: the slab must be such that the z is the short * direction. Otherwise we get trash. @@ -60,7 +62,7 @@ struct DLC_struct { */ double gap_size; - /** Flag whether #far_cut was set by the user, or calculated by Espresso. + /** Flag whether #far_cut was set by the user, or calculated by ESPResSo. * In the latter case, the cutoff will be adapted if important parameters, * such as the box dimensions, change. */ diff --git a/src/core/electrostatics_magnetostatics/mmm-common.cpp b/src/core/electrostatics_magnetostatics/mmm-common.cpp index e459fc1137c..d2c95392e88 100644 --- a/src/core/electrostatics_magnetostatics/mmm-common.cpp +++ b/src/core/electrostatics_magnetostatics/mmm-common.cpp @@ -18,17 +18,6 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see . */ -/** \file - Common parts of the MMM family of methods for the - electrostatic interaction, MMM1D, MMM2D and ELC. This file contains the - code for the polygamma expansions used for the near formulas of MMM1D and - MMM2D. - - The expansion of the polygamma functions is fairly easy and follows directly - from Abramowitz and Stegun. For details, see Axel Arnold and Christian Holm, - "MMM2D: A fast and accurate summation method for electrostatic interactions - in 2D slab geometries", Comp. Phys. Comm., 148/3(2002),327-348. -*/ #include "config.hpp" diff --git a/src/core/electrostatics_magnetostatics/mmm-common.hpp b/src/core/electrostatics_magnetostatics/mmm-common.hpp index 8ab21a0c0a8..8672f1bae6f 100644 --- a/src/core/electrostatics_magnetostatics/mmm-common.hpp +++ b/src/core/electrostatics_magnetostatics/mmm-common.hpp @@ -19,8 +19,14 @@ * along with this program. If not, see . */ /** \file - modified polygamma functions. See Arnold,Holm 2002 -*/ + * Common parts of the MMM family of methods for the electrostatic + * interaction, MMM1D, MMM2D and ELC. This file contains the code for the + * polygamma expansions used for the near formulas of MMM1D and MMM2D. + * + * The expansion of the polygamma functions is fairly easy and follows + * directly from @cite abramowitz65a. For details, see @cite arnold02a. + */ + #ifndef MMM_COMMON_H #define MMM_COMMON_H @@ -29,7 +35,7 @@ #include -/** \name Math Constants */ +/** \name Math constants */ /*@{*/ #define C_2PI (2 * M_PI) #define C_GAMMA (0.57721566490153286060651209008) @@ -37,21 +43,21 @@ #define C_2PISQR (C_2PI * C_2PI) /*@}*/ -/** table of the Taylor expansions of the modified polygamma functions */ +/** Table of the Taylor expansions of the modified polygamma functions */ extern std::vector modPsi; extern int n_modPsi; -/** modified polygamma for even order 2*n, n >= 0 */ +/** Modified polygamma for even order 2*n, n >= 0 */ inline double mod_psi_even(int n, double x) { return evaluateAsTaylorSeriesAt(&modPsi[2 * n], x * x); } -/** modified polygamma for odd order 2*n+1, n>= 0 */ +/** Modified polygamma for odd order 2*n+1, n>= 0 */ inline double mod_psi_odd(int n, double x) { return x * evaluateAsTaylorSeriesAt(&modPsi[2 * n + 1], x * x); } -/** create the both the even and odd polygamma functions up to order 2*n */ +/** Create the both the even and odd polygamma functions up to order 2*n */ void create_mod_psi_up_to(int n); #endif diff --git a/src/core/electrostatics_magnetostatics/mmm2d.cpp b/src/core/electrostatics_magnetostatics/mmm2d.cpp index 64e4a332d77..3c3f472e98f 100644 --- a/src/core/electrostatics_magnetostatics/mmm2d.cpp +++ b/src/core/electrostatics_magnetostatics/mmm2d.cpp @@ -527,9 +527,8 @@ static void setup_z_force() { const double pref = coulomb.prefactor * C_2PI * ux * uy; constexpr int e_size = 1, size = 2; - /* there is NO contribution from images here, unlike claimed in Tyagi et al. - Please refer to the Entropy - article of Arnold, Kesselheim, Breitsprecher et al, 2013, for details. */ + /** there is NO contribution from images here, unlike claimed in + * @cite tyagi10a. Please refer to @cite arnold13c for details. */ if (this_node == 0) { clear_vec(blwentry(lclcblk, 0, e_size), e_size); diff --git a/src/core/electrostatics_magnetostatics/mmm2d.hpp b/src/core/electrostatics_magnetostatics/mmm2d.hpp index d194a951984..a8777d289d7 100644 --- a/src/core/electrostatics_magnetostatics/mmm2d.hpp +++ b/src/core/electrostatics_magnetostatics/mmm2d.hpp @@ -19,14 +19,14 @@ * along with this program. If not, see . */ /** \file - MMM2D algorithm for long range Coulomb interaction - in 2d+h geometries. Implementation of the MMM2D method for the calculation - of the electrostatic interaction for two dimensionally periodic systems. - For details on the method see MMM general. The MMM2D method works only with - the layered or nsquared "cell system". The tuning is not automated, since - the only tunable parameter is the cell size, which can be changed easily in - the script interface. Moreover, only a few values make sense to be tested, - since in general the number of cells will be between 5 and 20. + * MMM2D algorithm for long range Coulomb interaction in 2d+h geometries. + * Implementation of the MMM2D method for the calculation + * of the electrostatic interaction for two dimensionally periodic systems. + * For details on the method see MMM general. The MMM2D method works only with + * the layered or nsquared "cell system". The tuning is not automated, since + * the only tunable parameter is the cell size, which can be changed easily in + * the script interface. Moreover, only a few values make sense to be tested, + * since in general the number of cells will be between 5 and 20. */ #ifndef MMM2D_H #define MMM2D_H @@ -52,7 +52,7 @@ typedef struct { double far_cut; /** squared value of #far_cut. */ double far_cut2; - /** Flag whether #far_cut was set by the user, or calculated by Espresso. + /** Flag whether #far_cut was set by the user, or calculated by ESPResSo. * In the latter case, the cutoff will be adapted if important parameters, * such as the box dimensions, change. */ @@ -82,7 +82,7 @@ extern MMM2D_struct mmm2d_params; * Manual setting is probably only good for testing. * @param delta_top @copybrief MMM2D_struct::delta_mid_top * @param delta_bot @copybrief MMM2D_struct::delta_mid_bot - * @param const_pot @copybrief MMM2D_struct::const_pot + * @param const_pot @copybrief MMM2D_struct::const_pot * @param pot_diff @copybrief MMM2D_struct::pot_diff */ int MMM2D_set_params(double maxPWerror, double far_cut, double delta_top, diff --git a/src/core/electrostatics_magnetostatics/p3m-common.hpp b/src/core/electrostatics_magnetostatics/p3m-common.hpp index 33742bd6ee5..60c1f18879f 100644 --- a/src/core/electrostatics_magnetostatics/p3m-common.hpp +++ b/src/core/electrostatics_magnetostatics/p3m-common.hpp @@ -25,29 +25,11 @@ * * We use here a P3M (Particle-Particle Particle-Mesh) method based * on the Ewald summation. Details of the used method can be found in - * Hockney/Eastwood and Deserno/Holm. The file p3m contains only the - * Particle-Mesh part. + * @cite hockney88a and @cite deserno98a @cite deserno98b. The file p3m + * contains only the Particle-Mesh part. * - * Further reading: - * - P. P. Ewald, - * *Die Berechnung optischer und elektrostatischer Gitterpotentiale*, - * Ann. Phys. (64) 253-287, 1921 - * - R. W. Hockney and J. W. Eastwood, - * *Computer simulation using particles*, - * IOP, London, 1988 - * - M. Deserno and C. Holm, - * *How to mesh up Ewald sums I + II*, - * J. Chem. Phys. (109) 7678, 1998; (109) 7694, 1998 - * - M. Deserno, C. Holm and H. J. Limbach, - * *How to mesh up Ewald sums*, - * in Molecular Dynamics on Parallel Computers, - * Ed. R. Esser et al., World Scientific, Singapore, 2000 - * - M. Deserno, - * *Counterion condensation for rigid linear polyelectrolytes*, - * PhD Thesis, Universität Mainz, 2000 - * - J. J. Cerda, - * *P3M for dipolar interactions*, - * J. Chem. Phys (129) 234104, 2008 + * Further reading: @cite ewald21a, @cite hockney88a, @cite deserno98a, + * @cite deserno98b, @cite deserno00e, @cite deserno00b, @cite cerda08d * */ #include "config.hpp" @@ -206,8 +188,8 @@ void p3m_add_block(double const *in, double *out, int const start[3], * most slowly, since it is not damped exponentially) can be * calculated analytically. The result (which depends on the order of * the spline interpolation) can be written as an even trigonometric - * polynomial. The results are tabulated here (The employed formula - * is Eqn. 7.66 in the book of Hockney and Eastwood). + * polynomial. The results are tabulated here (the employed formula + * is eq. (7.66) in @cite hockney88a). */ double p3m_analytic_cotangent_sum(int n, double mesh_i, int cao); diff --git a/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp b/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp index ec8947771a1..7c53d1c6a4f 100644 --- a/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp +++ b/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp @@ -27,7 +27,7 @@ * "charge" by "dipole". In this way one can recognize the similarity of * the functions while avoiding nasty confusions in their use. * - * By default the magnetic epsilon is metallic = 0. + * By default the magnetic epsilon is metallic = 0. * * The corresponding header file is p3m-dipolar.hpp. */ @@ -133,8 +133,8 @@ static bool dp3m_sanity_checks_boxl(); static void dp3m_calc_local_ca_mesh(); /** Interpolate the P-th order charge assignment function from - * Hockney/Eastwood 5-189 (or 8-61). The following charge fractions - * are also tabulated in Deserno/Holm. + * @cite hockney88a 5-189 (or 8-61). The following charge fractions + * are also tabulated in @cite deserno98a @cite deserno98b. */ static void dp3m_interpolate_dipole_assignment_function(); @@ -164,7 +164,7 @@ static void dp3m_compute_constants_energy_dipolar(); * * Calculates the aliasing sums in the nominator and denominator of * the expression for the optimal influence function (see - * Hockney/Eastwood: 8-22, p. 275). + * @cite hockney88a : 8-22, p. 275). * * \param n n-vector for which the aliasing sum is to be performed * \param nominator aliasing sums in the nominator @@ -198,8 +198,9 @@ static void dp3m_tune_aliasing_sums(int nx, int ny, int nz, int mesh, double mesh_i, int cao, double alpha_L_i, double *alias1, double *alias2); -// To compute the value of alpha through a bibisection method from the formula -// 33 of JCP115,6351,(2001). +/** Compute the value of alpha through a bisection method. + * Based on eq. (33) @cite wang01a. + */ double dp3m_rtbisection(double box_size, double prefac, double r_cut_iL, int n_c_part, double sum_q2, double x1, double x2, double xacc, double tuned_accuracy); @@ -531,10 +532,10 @@ int dp3m_set_mesh_offset(double x, double y, double z) { return ES_OK; } -/* We left the handling of the epsilon, due to portability reasons in -the future for the electrical dipoles, or if people wants to do -electrical dipoles alone using the magnetic code .. */ - +/** We left the handling of the epsilon, due to portability reasons in + * the future for the electrical dipoles, or if people want to do + * electrical dipoles alone using the magnetic code. Currently unused. + */ int dp3m_set_eps(double eps) { dp3m.params.epsilon = eps; @@ -1131,6 +1132,7 @@ double calc_surface_term(bool force_flag, bool energy_flag, } /************************************************************/ + void dp3m_gather_fft_grid(double *themesh) { int s_dir, r_dir, evenodd; MPI_Status status; @@ -1435,16 +1437,8 @@ double dp3m_perform_aliasing_sums_energy(const int n[3], double nominator[1]) { /*****************************************************************************/ -/************************************************ - * Functions for dipolar P3M Parameter tuning - * This tuning is based on the P3M tuning of the charges - * which in turn is based on the P3M_tune by M. Deserno - ************************************************/ - #define P3M_TUNE_MAX_CUTS 50 -/*****************************************************************************/ - /** @copybrief p3m_get_accuracy * * The real space error is tuned such that it contributes half of the @@ -1500,8 +1494,6 @@ double dp3m_get_accuracy(int mesh, int cao, double r_cut_iL, double *_alpha_L, return sqrt(Utils::sqr(rs_err) + Utils::sqr(ks_err)); } -/*****************************************************************************/ - /** @copybrief p3m_mcr_time * * @param[in] mesh @copybrief P3MParameters::mesh @@ -1535,8 +1527,6 @@ static double dp3m_mcr_time(int mesh, int cao, double r_cut_iL, return int_time; } -/*****************************************************************************/ - /** @copybrief p3m_mc_time * * The @p _r_cut_iL is determined via a simple bisection. @@ -1645,8 +1635,6 @@ static double dp3m_mc_time(char **log, int mesh, int cao, double r_cut_iL_min, return int_time; } -/*****************************************************************************/ - /** @copybrief p3m_m_time * * @p _cao should contain an initial guess, which is then adapted by stepping @@ -1947,8 +1935,6 @@ int dp3m_adaptive_tune(char **logger) { return ES_OK; } -/*****************************************************************************/ - void dp3m_count_magnetic_particles() { double node_sums[2], tot_sums[2]; @@ -1971,18 +1957,7 @@ void dp3m_count_magnetic_particles() { REGISTER_CALLBACK(dp3m_count_magnetic_particles) -/*****************************************************************************/ - -/* Following are the two functions for computing the error in dipolar P3M and - tune the parameters to minimize the time with the desired accuracy. - - - This functions are called by the functions: dp3m_get_accuracy(). - -*/ - -/*****************************************************************************/ - +/** Calculate the k-space error of dipolar-P3M */ static double dp3m_k_space_error(double box_size, double prefac, int mesh, int cao, int n_c_part, double sum_q2, double alpha_L) { @@ -2012,6 +1987,7 @@ static double dp3m_k_space_error(double box_size, double prefac, int mesh, (box_size * box_size * box_size * box_size); } +/** Tuning dipolar-P3M */ void dp3m_tune_aliasing_sums(int nx, int ny, int nz, int mesh, double mesh_i, int cao, double alpha_L_i, double *alias1, double *alias2) { @@ -2043,17 +2019,12 @@ void dp3m_tune_aliasing_sums(int nx, int ny, int nz, int mesh, double mesh_i, } } -//---------------------------------------------------------- -// Function used to calculate the value of the errors -// for the REAL part of the force in terms of the splitting parameter alpha of -// Ewald -// based on the formulas 33, the paper of Zuowei-HolmJCP, 115,6351,(2001). - -// Please, notice than in this more refined approach we don't use -// formulas 37, but 33 which maintains all the powers in alpha - -//------------------------------------------------------------ - +/** Calculate the value of the errors for the REAL part of the force in terms + * of the splitting parameter alpha of Ewald. Based on eq. (33) @cite wang01a. + * + * Please note that in this more refined approach we don't use + * eq. (37), but eq. (33) which maintains all the powers in alpha. + */ double P3M_DIPOLAR_real_space_error(double box_size, double prefac, double r_cut_iL, int n_c_part, double sum_q2, double alpha_L) { @@ -2082,12 +2053,10 @@ double P3M_DIPOLAR_real_space_error(double box_size, double prefac, return d_error_f; } -/*****************************************************************************/ - -// Using bisection find the root of a function "func-tuned_accuracy/sqrt(2.)" -// known to lie between x1 and x2. The root, returned as rtbis, will be refined -// until its accuracy is +-xacc. - +/** Using bisection, find the root of a function "func-tuned_accuracy/sqrt(2.)" + * known to lie between x1 and x2. The root, returned as rtbis, will be + * refined until its accuracy is \f$\pm\f$ @p xacc. + */ double dp3m_rtbisection(double box_size, double prefac, double r_cut_iL, int n_c_part, double sum_q2, double x1, double x2, double xacc, double tuned_accuracy) { @@ -2413,9 +2382,7 @@ void dp3m_scaleby_box_l() { dp3m_calc_influence_function_energy(); } -/*****************************************************************************/ - -/** function to give the dipolar-P3M energy correction */ +/** Calculate the dipolar-P3M energy correction */ void dp3m_compute_constants_energy_dipolar() { double Eself, Ukp3m; @@ -2434,8 +2401,4 @@ void dp3m_compute_constants_energy_dipolar() { -dp3m.sum_mu2 * (Ukp3m + Eself + 2. * Utils::pi() / 3.); } -/*****************************************************************************/ - -/*****************************************************************************/ - #endif /* DP3M */ diff --git a/src/core/electrostatics_magnetostatics/p3m-dipolar.hpp b/src/core/electrostatics_magnetostatics/p3m-dipolar.hpp index 360ee91739a..c4e0320170f 100644 --- a/src/core/electrostatics_magnetostatics/p3m-dipolar.hpp +++ b/src/core/electrostatics_magnetostatics/p3m-dipolar.hpp @@ -25,13 +25,10 @@ * * We use here a P3M (Particle-Particle Particle-Mesh) method based * on the dipolar Ewald summation. Details of the used method can be found in - * Hockney/Eastwood and Deserno/Holm. The file p3m contains only the - * Particle-Mesh part. + * @cite hockney88a and @cite deserno98a @cite deserno98b. The file p3m + * contains only the Particle-Mesh part. * - * Further reading: - * - J. J. Cerda, - * *P3M for dipolar interactions*, - * J. Chem. Phys (129) 234104, 2008 + * Further reading: @cite cerda08d * * Implementation in p3m-dipolar.cpp. */ @@ -163,10 +160,10 @@ void dp3m_deactivate(); * These parameters are stored in the @ref dp3m object. * * The function utilizes the analytic expression of the error estimate - * for the dipolar P3M method in the paper of J. J. Cerda et al., JCP 2008 in + * for the dipolar P3M method in the paper of @cite cerda08d in * order to obtain the rms error in the force for a system of N randomly - * distributed particles in a cubic box. - * For the real space error the estimate of Kolafa/Perram is used. + * distributed particles in a cubic box. For the real space error, the + * estimate in @cite kolafa92a is used. * * Parameter ranges if not given explicit values via dp3m_set_tune_params(): * - @p mesh is set up such that the number of mesh points is equal to the diff --git a/src/core/electrostatics_magnetostatics/p3m.cpp b/src/core/electrostatics_magnetostatics/p3m.cpp index 378dba79569..2d7211fc976 100644 --- a/src/core/electrostatics_magnetostatics/p3m.cpp +++ b/src/core/electrostatics_magnetostatics/p3m.cpp @@ -20,7 +20,7 @@ */ /** @file * - * The corresponding header file is p3m.hpp. + * The corresponding header file is @ref p3m.hpp. */ #include "p3m.hpp" @@ -149,8 +149,8 @@ static bool p3m_sanity_checks_boxl(); static void p3m_calc_local_ca_mesh(); /** Interpolate the P-th order charge assignment function from - * Hockney/Eastwood 5-189 (or 8-61). The following charge fractions - * are also tabulated in Deserno/Holm. + * @cite hockney88a 5-189 (or 8-61). The following charge fractions + * are also tabulated in @cite deserno98a @cite deserno98b. */ static void p3m_interpolate_charge_assignment_function(); @@ -163,15 +163,15 @@ static void p3m_calc_meshift(); */ static void p3m_calc_differential_operator(); -/** Calculate the optimal influence function of Hockney and Eastwood. +/** Calculate the optimal influence function of @cite hockney88a. * (optimised for force calculations) * * Each node calculates only the values for its domain in k-space * (see fft.plan[3].mesh and fft.plan[3].start). * - * See also: Hockney/Eastwood 8-22 (p275). Note the somewhat + * See also: @cite hockney88a 8-22 (p275). Note the somewhat * different convention for the prefactors, which is described in - * Deserno/Holm. + * @cite deserno98a @cite deserno98b. */ static void p3m_calc_influence_function_force(); @@ -184,7 +184,7 @@ static void p3m_calc_influence_function_energy(); * * Calculate the aliasing sums in the nominator and denominator of * the expression for the optimal influence function (see - * Hockney/Eastwood: 8-22, p. 275). + * @cite hockney88a : 8-22, p. 275). * * \param n n-vector for which the aliasing sum is to be performed. * \param nominator aliasing sums in the nominator. @@ -211,7 +211,7 @@ static double p3m_real_space_error(double prefac, double r_cut_iL, int n_c_part, double sum_q2, double alpha_L); /** Calculate the analytic expression of the error estimate for the - * P3M method in the book of Hockney and Eastwood (Eqn. 8.23) in + * P3M method in @cite hockney88a (eq. (8.23)) in * order to obtain the rms error in the force for a system of N * randomly distributed particles in a cubic box (k-space part). * \param prefac Prefactor of Coulomb interaction. @@ -2236,14 +2236,12 @@ void p3m_scaleby_box_l() { p3m_calc_influence_function_energy(); } +/** @details Calculate the long range electrostatics part of the stress + * tensor. This is part \f$\Pi_{\textrm{dir}, \alpha, \beta}\f$ eq. (2.6) + * in @cite essmann95a. The part \f$\Pi_{\textrm{corr}, \alpha, \beta}\f$ + * eq. (2.8) is not present here since M is the empty set in our simulations. + */ void p3m_calc_kspace_stress(double *stress) { - /** - * Calculates the long range electrostatics part of the stress tensor. This is - * part Pi_{dir, alpha,beta} in the paper by Essmann et al "A smooth particle - * mesh Ewald method", The Journal of Chemical Physics 103, 8577 (1995); - * doi: 10.1063/1.470117. The part Pi_{corr, alpha, beta} in the Essmann paper - * is not present here since M is the empty set in our simulations. - */ if (p3m.sum_q2 > 0) { std::vector node_k_space_stress; std::vector k_space_stress; diff --git a/src/core/electrostatics_magnetostatics/p3m.hpp b/src/core/electrostatics_magnetostatics/p3m.hpp index 1c3bfe6ee21..ae63a559725 100644 --- a/src/core/electrostatics_magnetostatics/p3m.hpp +++ b/src/core/electrostatics_magnetostatics/p3m.hpp @@ -25,28 +25,10 @@ * * We use a P3M (Particle-Particle Particle-Mesh) method based on the * Ewald summation. Details of the used method can be found in - * Hockney/Eastwood and Deserno/Holm. + * @cite hockney88a and @cite deserno98a @cite deserno98b. * - * Further reading: - * - P. P. Ewald, - * *Die Berechnung optischer und elektrostatischer Gitterpotentiale*, - * Ann. Phys. (64) 253-287, 1921 - * - R. W. Hockney and J. W. Eastwood, - * *Computer simulation using particles*, - * IOP, London, 1988 - * - M. Deserno and C. Holm, - * *How to mesh up Ewald sums I + II*, - * J. Chem. Phys. (109) 7678, 1998; (109) 7694, 1998 - * - M. Deserno, C. Holm and H. J. Limbach, - * *How to mesh up Ewald sums*, - * in Molecular Dynamics on Parallel Computers, - * Ed. R. Esser et al., World Scientific, Singapore, 2000 - * - M. Deserno, - * *Counterion condensation for rigid linear polyelectrolytes*, - * PhD Thesis, Universität Mainz, 2000 - * - J. J. Cerda, - * *P3M for dipolar interactions*, - * J. Chem. Phys (129) 234104, 2008 + * Further reading: @cite ewald21a, @cite hockney88a, @cite deserno98a, + * @cite deserno98b, @cite deserno00e, @cite deserno00b, @cite cerda08d. * * Implementation in p3m.cpp. */ @@ -141,7 +123,7 @@ extern p3m_data_struct p3m; * These parameters are stored in the @ref p3m object. * * The function utilizes the analytic expression of the error estimate - * for the P3M method in the book of Hockney and Eastwood (Eqn. 8.23) in + * for the P3M method in @cite hockney88a (eq. (8.23)) in * order to obtain the rms error in the force for a system of N randomly * distributed particles in a cubic box. * For the real space error the estimate of Kolafa/Perram is used. diff --git a/src/core/electrostatics_magnetostatics/p3m_gpu_error_cuda.cu b/src/core/electrostatics_magnetostatics/p3m_gpu_error_cuda.cu index 7a3ebb73a63..cf477771bdb 100644 --- a/src/core/electrostatics_magnetostatics/p3m_gpu_error_cuda.cu +++ b/src/core/electrostatics_magnetostatics/p3m_gpu_error_cuda.cu @@ -42,7 +42,7 @@ using Utils::sqr; #endif /** @todo Extend to higher order. This comes from some 1/sin expansion in - * Hockney/Eastwood + * @cite hockney88a */ template diff --git a/src/core/electrostatics_magnetostatics/reaction_field.hpp b/src/core/electrostatics_magnetostatics/reaction_field.hpp index e4b844f0050..983b7a25999 100644 --- a/src/core/electrostatics_magnetostatics/reaction_field.hpp +++ b/src/core/electrostatics_magnetostatics/reaction_field.hpp @@ -21,11 +21,10 @@ #ifndef REACTION_FIELD_H #define REACTION_FIELD_H /** \file - * Routines to calculate the Reaction Field Energy or/and force - * for a particle pair. - * M. Neumann, J. Chem. Phys 82, 5663 (1985) - * \ref forces.cpp + * Routines to calculate the Reaction Field energy or/and force + * for a particle pair @cite neumann85b. * + * Implementation in \ref reaction_field.cpp */ #include "config.hpp" @@ -35,25 +34,21 @@ /** Structure to hold Reaction Field Parameters. */ typedef struct { - /** ionic strength . */ + /** ionic strength. */ double kappa; - /** epsilon1 (continuum dielectric constant inside) . */ + /** epsilon1 (continuum dielectric constant inside). */ double epsilon1; - /** epsilon2 (continuum dielectric constant outside) . */ + /** epsilon2 (continuum dielectric constant outside). */ double epsilon2; /** Cutoff for Reaction Field interaction. */ double r_cut; - /** B important prefactor . */ + /** B important prefactor. */ double B; } Reaction_field_params; /** Structure containing the Reaction Field parameters. */ extern Reaction_field_params rf_params; -/** \name Functions */ -/************************************************************/ -/*@{*/ - int rf_set_params(double kappa, double epsilon1, double epsilon2, double r_cut); inline void add_rf_coulomb_pair_force_no_cutoff(double const q1q2, @@ -67,8 +62,7 @@ inline void add_rf_coulomb_pair_force_no_cutoff(double const q1q2, force += fac * d; } -/** Computes the Reaction Field pair force and adds this - * force to the particle forces. +/** Compute the Reaction Field pair force. * @param q1q2 Product of the charges on p1 and p2. * @param d Vector pointing from p1 to p2. * @param dist Distance between p1 and p2. @@ -83,6 +77,10 @@ inline void add_rf_coulomb_pair_force(double const q1q2, } } +/** Compute the Reaction Field pair energy. + * @param q1q2 Product of the charges on p1 and p2. + * @param dist Distance between p1 and p2. + */ inline double rf_coulomb_pair_energy_no_cutoff(double const q1q2, double const dist) { double fac; @@ -102,7 +100,6 @@ inline double rf_coulomb_pair_energy(double const q1q2, double const dist) { return 0.0; } -/*@}*/ #endif #endif diff --git a/src/core/electrostatics_magnetostatics/scafacos.cpp b/src/core/electrostatics_magnetostatics/scafacos.cpp index be47ea402df..35befa8632e 100644 --- a/src/core/electrostatics_magnetostatics/scafacos.cpp +++ b/src/core/electrostatics_magnetostatics/scafacos.cpp @@ -314,7 +314,7 @@ void tune() { return; // Tune on the master node will issue mpi calls } } else { - // Espresso is not affected by a short range cutoff. Tune in parallel + // ESPResSo is not affected by a short range cutoff. Tune in parallel scafacos->tune(particles.charges, particles.positions); } } diff --git a/src/core/event.cpp b/src/core/event.cpp index 6594462e25c..33349e4f16e 100644 --- a/src/core/event.cpp +++ b/src/core/event.cpp @@ -324,7 +324,7 @@ void on_boxl_change() { void on_cell_structure_change() { /* Now give methods a chance to react to the change in cell - structure. Most ES methods need to reinitialize, as they depend + structure. Most ES methods need to reinitialize, as they depend on skin, node grid and so on. Only for a change in box length we have separate, faster methods, as this might happen frequently in a NpT simulation. */ diff --git a/src/core/event.hpp b/src/core/event.hpp index 240f8eecdf9..3032f515408 100644 --- a/src/core/event.hpp +++ b/src/core/event.hpp @@ -22,7 +22,7 @@ #define CORE_EVENT_HPP /** \file * This file contains the hook procedures. These are the ones with names - * on_* and are called whenever something is changed in Espresso which + * on_* and are called whenever something is changed in ESPResSo which * might influence other parts. For example, the P3M code has to be * reinitialized whenever the box size changes. The hooking mechanism * allows to keep track of such changes. diff --git a/src/core/ghosts.hpp b/src/core/ghosts.hpp index f97f6d78128..02ee9fd26fb 100644 --- a/src/core/ghosts.hpp +++ b/src/core/ghosts.hpp @@ -22,80 +22,87 @@ #define _GHOSTS_H /** \file * Ghost particles and particle exchange. - -In this file you find everything concerning the exchange of -particle data (particles, ghosts, positions and forces) for short -range interactions between the spacial domains of neighbouring -nodes. - -

How does this work

-The ghost communication transfers data from cells on one node to cells on -another node during the integration process. Note that data can only be -transferred from one cell to one other, and the contents of the other cell will -be overwritten. This communication is invoked by the integrator, at four -different times in an integration step. These steps are reflected in \ref -cell_structure structure, since they are treated differently by different cell -systems:
  • ghost_cells are used to transfer the cell sizes, i.e., make -sure that for all later transfers the target cell has the same size as the -source cell.
  • exchange_ghosts sets up all information on the ghosts that is -necessary. Normally transfers the (shifted) position and particle properties. -
  • update_ghosts is used to update the particle properties if no particles -have been moved between cells. Therefore only the positions are transferred, -otherwise this normally looks pretty much like the exchange_ghosts.
  • -collect_ghost_forces finally is used to transfer back forces that were exerted -on a ghost particle. They are simply added again to the force of the real -particle. The communication process is therefore inverted with respect to -exchange_ghosts and update_ghosts, i.e., sending is replaced by receiving and -the other way round. -
-The particle data that has to be transferred, and especially from where to -where, heavily depends on the cell system. In Espresso, this is abstracted in -form of ghost communicators (\ref GhostCommunicator) and ghost communications -(\ref GhostCommunication). The ghost communicators represent the four -communications above and consist of the data to transfer (which is determined by -their type) and a list of ghost communications. The data types are described by -the particle data classes
  • GHOSTTRANS_PROPRTS transfers the \ref -ParticleProperties
  • GHOSTTRANS_POSITION transfers the \ref ParticlePosition -
  • GHOSTTRANS_POSSHFTD is no transfer class, but marks that an additional -vector, the shift, has to be added to the positions. Can only be used together -with GHOSTTRANS_POSITION.
  • GHOSTTRANS_MOMENTUM transfers the \ref -ParticleMomentum
  • GHOSTTRANS_FORCE transfers the \ref ParticleForce
  • -GHOSTTRANS_PARTNUM transfers the cell sizes -
-Each ghost communication describes a single communication of the local with -another node (or all other nodes). The data transferred can be any number of -cells, there are five communication types:
  • GHOST_SEND sends data to -one other node
  • GHOST_RECV recvs data from one other node. In the case of -forces, they are added up, everything else is overwritten
  • GHOST_BCST sends -data to all nodes
  • GHOST_RDCE recvs data from all nodes. In the case of -forces, they are added up, everything else is overwritten
  • GHOST_LOCL -transfer data from a local cell to another local cell. In this case, the first -half of the cells are the sending cells, the other half are the corresponding -receivers. -
-Note that for the first four communications you have to make sure that when one -node sends to another, that the sender has GHOST_SEND and the receiver -GHOST_RECV. In the case of GHOST_BCST rsp. GHOST_RDCE, all nodes have to have -the same communication type and the same master sender/receiver (just like the -MPI commands). - -A special topic are GHOST_PREFETCH and GHOST_PSTSTORE. For example, if all nodes -broadcast to the other, the naive implementation will be that n_nodes times a -GHOST_BCST is done with different master nodes. But this means that each time -n_nodes-1 nodes wait for the master to construct its send buffer. Therefore -there is the prefetch flag which can be set on a pair of recv/send operations. -If the ghost communication reaches a recv operation with prefetch, the next send -operation (which must have the prefetch set!!) is searched and the send buffer -already created. When sending, this precreated send buffer is used. In the -scenario above, all nodes create the send buffers simultaneously in the first -communication step, thereby reducing the latency a little bit. The pststore is -similar and postpones the write back of received data until a send operation -(with a precreated send buffer) is finished. - -The ghost communicators are created in the init routines of the cell systems, -therefore have a look at \ref dd_topology_init or \ref nsq_topology_init for -further details. -*/ + * + * In this file you find everything concerning the exchange of + * particle data (particles, ghosts, positions and forces) for short + * range interactions between the spacial domains of neighbouring + * nodes. + * + *

How does this work

+ * The ghost communication transfers data from cells on one node to cells on + * another node during the integration process. Note that data can only be + * transferred from one cell to one other, and the contents of the other cell + * will be overwritten. This communication is invoked by the integrator, at + * four different times in an integration step. These steps are reflected in + * @ref cell_structure structure, since they are treated differently by + * different cell systems: + * - @ref ghost_cells are used to transfer the cell sizes, i.e., make sure + * that for all later transfers the target cell has the same size as the + * source cell. + * - @ref CellStructure::exchange_ghosts_comm sets up all information on the + * ghosts that is necessary. Normally transfers the (shifted) position and + * particle properties. + * - @ref cells_update_ghosts is used to update the particle properties if no + * particles have been moved between cells. Therefore only the positions are + * transferred, otherwise this normally looks pretty much like the + * @ref CellStructure::exchange_ghosts_comm. + * - @ref CellStructure::collect_ghost_force_comm finally is used to transfer + * back forces that were exerted on a ghost particle. They are simply added + * again to the force of the real particle. The communication process is + * therefore inverted with respect to + * @ref CellStructure::exchange_ghosts_comm and @ref cells_update_ghosts, + * i.e., sending is replaced by receiving and the other way round. + * + * The particle data that has to be transferred, and especially from where to + * where, heavily depends on the cell system. In ESPResSo, this is abstracted + * in form of ghost communicators (\ref GhostCommunicator) and ghost + * communications (\ref GhostCommunication). The ghost communicators represent + * the four communications above and consist of the data to transfer (which is + * determined by their type) and a list of ghost communications. The data + * types are described by the particle data classes: + * - @ref GHOSTTRANS_PROPRTS transfers the @ref ParticleProperties + * - @ref GHOSTTRANS_POSITION transfers the @ref ParticlePosition + * - @ref GHOSTTRANS_MOMENTUM transfers the @ref ParticleMomentum + * - @ref GHOSTTRANS_FORCE transfers the @ref ParticleForce + * - @ref GHOSTTRANS_PARTNUM transfers the cell sizes + * + * Each ghost communication describes a single communication of the local with + * another node (or all other nodes). The data transferred can be any number + * of cells, there are five communication types: + * - @ref GHOST_SEND sends data to one other node + * - @ref GHOST_RECV recvs data from one other node. In the case of + * forces, they are added up, everything else is overwritten + * - @ref GHOST_BCST sends data to all nodes + * - @ref GHOST_RDCE recvs data from all nodes. In the case of + * forces, they are added up, everything else is overwritten + * - @ref GHOST_LOCL transfer data from a local cell to another local cell. + * In this case, the first half of the cells are the sending cells, the + * other half are the corresponding receivers. + * + * Note that for the first four communications you have to make sure that + * when one node sends to another, that the sender has @ref GHOST_SEND + * and the receiver @ref GHOST_RECV. In the case of @ref GHOST_BCST resp. + * @ref GHOST_RDCE, all nodes have to have the same communication type + * and the same master sender/receiver (just like the MPI commands). + * + * A special topic are @ref GHOST_PREFETCH and @ref GHOST_PSTSTORE. For + * example, if all nodes broadcast to the other, the naive implementation will + * be that @c n_nodes times a @ref GHOST_BCST is done with different master + * nodes. But this means that each time n_nodes - 1 nodes wait for + * the master to construct its send buffer. Therefore there is the prefetch + * flag which can be set on a pair of recv/send operations. If the ghost + * communication reaches a recv operation with prefetch, the next send + * operation (which must have the prefetch set!!) is searched and the send + * buffer already created. When sending, this precreated send buffer is used. + * In the scenario above, all nodes create the send buffers simultaneously in + * the first communication step, thereby reducing the latency a little bit. + * The pststore is similar and postpones the write back of received data + * until a send operation (with a precreated send buffer) is finished. + * + * The ghost communicators are created in the init routines of the cell + * systems, therefore have a look at @ref dd_topology_init or + * @ref nsq_topology_init for further details. + */ #include "Cell.hpp" #include diff --git a/src/core/grid_based_algorithms/electrokinetics_cuda.cu b/src/core/grid_based_algorithms/electrokinetics_cuda.cu index 279a9fcb43a..3a24f714687 100644 --- a/src/core/grid_based_algorithms/electrokinetics_cuda.cu +++ b/src/core/grid_based_algorithms/electrokinetics_cuda.cu @@ -3609,7 +3609,7 @@ LOOKUP_TABLE default\n", int ek_print_vtk_lbforce_density(char *filename) { #if !defined(VIRTUAL_SITES_INERTIALESS_TRACERS) && !defined(EK_DEBUG) - throw std::runtime_error("Please rebuild Espresso with EK_DEBUG"); + throw std::runtime_error("Please rebuild ESPResSo with EK_DEBUG"); #else FILE *fp = fopen(filename, "w"); diff --git a/src/core/grid_based_algorithms/lb.cpp b/src/core/grid_based_algorithms/lb.cpp index b547eb75834..c240bdcfa35 100644 --- a/src/core/grid_based_algorithms/lb.cpp +++ b/src/core/grid_based_algorithms/lb.cpp @@ -63,7 +63,7 @@ using Utils::get_linear_index; #include namespace { -/** Basis of the mode space as described in [Duenweg, Schiller, Ladd] */ +/** Basis of the mode space as described in @cite dunweg07a */ extern constexpr const std::array, 19> e_ki = { {{{1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}}, {{0, 1, -1, 0, 0, 0, 0, 1, -1, 1, -1, 1, -1, 1, -1, 0, 0, 0, 0}}, @@ -85,6 +85,7 @@ extern constexpr const std::array, 19> e_ki = { {{0, -1, -1, 1, 1, -0, -0, 0, 0, 0, 0, 1, 1, 1, 1, -1, -1, -1, -1}}, {{0, -1, -1, -1, -1, 2, 2, 2, 2, 2, 2, -1, -1, -1, -1, -1, -1, -1, -1}}}}; +/** Transposed version of @ref e_ki */ extern constexpr const std::array, 19> e_ki_transposed = { {{{1, 0, 0, 0, -1, 0, -0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0}}, {{1, 1, 0, 0, 0, 1, 1, 0, 0, 0, -2, 0, 0, -0, 0, 0, -2, -1, -1}}, @@ -275,12 +276,12 @@ void lb_reinit_fluid(std::vector &lb_fields, void lb_reinit_parameters(LB_Parameters &lb_parameters) { if (lb_parameters.viscosity > 0.0) { - /* Eq. (80) Duenweg, Schiller, Ladd, PRE 76(3):036704 (2007). */ + /* Eq. (80) @cite dunweg07a. */ lb_parameters.gamma_shear = 1. - 2. / (6. * lb_parameters.viscosity + 1.); } if (lb_parameters.bulk_viscosity > 0.0) { - /* Eq. (81) Duenweg, Schiller, Ladd, PRE 76(3):036704 (2007). */ + /* Eq. (81) @cite dunweg07a. */ lb_parameters.gamma_bulk = 1. - 2. / (9. * lb_parameters.bulk_viscosity + 1.); } @@ -299,7 +300,7 @@ void lb_reinit_parameters(LB_Parameters &lb_parameters) { // gamma_even = 0.0; if (lb_parameters.kT > 0.0) { - /* Eq. (51) Duenweg, Schiller, Ladd, PRE 76(3):036704 (2007). + /* Eq. (51) @cite dunweg07a. * Note that the modes are not normalized as in the paper here! */ double mu = lb_parameters.kT / D3Q19::c_sound_sq * lb_parameters.tau * lb_parameters.tau / @@ -1311,8 +1312,8 @@ Utils::Vector6d lb_calc_stress(std::array const &modes, // Transform the stress tensor components according to the modes that // correspond to those used by U. Schiller. In terms of populations this - // expression then corresponds exactly to those in Eqs. 116 - 121 in the - // Duenweg and Ladd paper, when these are written out in populations. + // expression then corresponds exactly to those in eq. (116)-(121) in + // @cite dunweg07a, when these are written out in populations. // But to ensure this, the expression in Schiller's modes has to be different! Utils::Vector6d stress; @@ -1330,13 +1331,6 @@ Utils::Vector6d lb_calc_stress(std::array const &modes, } #ifdef LB_BOUNDARIES -/** Bounce back boundary conditions. - * The populations that have propagated into a boundary node - * are bounced back to the node they came from. This results - * in no slip boundary conditions. - * - * [cf. Ladd and Verberg, J. Stat. Phys. 104(5/6):1191-1251, 2001] - */ void lb_bounce_back(LB_Fluid &lbfluid, const LB_Parameters &lb_parameters, const std::vector &lb_fields) { int k, i, l; diff --git a/src/core/grid_based_algorithms/lb.hpp b/src/core/grid_based_algorithms/lb.hpp index 27fcbe6c384..29efcd4a244 100644 --- a/src/core/grid_based_algorithms/lb.hpp +++ b/src/core/grid_based_algorithms/lb.hpp @@ -70,7 +70,7 @@ struct LB_FluidNode { Utils::Vector3d force_density; #ifdef VIRTUAL_SITES_INERTIALESS_TRACERS // For particle update, we need the force on the nodes in LBM - // Yet, Espresso resets the force immediately after the LBM update + // Yet, ESPResSo resets the force immediately after the LBM update // Therefore we save it here Utils::Vector3d force_density_buf; #endif @@ -253,9 +253,7 @@ void lb_prepare_communication(HaloCommunicator &halo_comm, /** Bounce back boundary conditions. * The populations that have propagated into a boundary node * are bounced back to the node they came from. This results - * in no slip boundary conditions. - * - * [cf. Ladd and Verberg, J. Stat. Phys. 104(5/6):1191-1251, 2001] + * in no slip boundary conditions, cf. @cite ladd01a. */ void lb_bounce_back(LB_Fluid &lbfluid, const LB_Parameters &lb_parameters, const std::vector &lb_fields); diff --git a/src/core/grid_based_algorithms/lb_boundaries.hpp b/src/core/grid_based_algorithms/lb_boundaries.hpp index 084d2d20a5e..fc2a01ea982 100644 --- a/src/core/grid_based_algorithms/lb_boundaries.hpp +++ b/src/core/grid_based_algorithms/lb_boundaries.hpp @@ -27,7 +27,7 @@ * after the streaming step, in all wall nodes all populations are bounced * back from where they came from. Ulf Schiller spent a lot of time * working on more powerful alternatives, they are to be found in the - * lb_testing branch of espresso until the end of 2010. Now we stripped + * lb_testing branch of ESPResSo until the end of 2010. Now we stripped * down the code to a minimum, as most of it was not sufficiently * understandable. * diff --git a/src/core/grid_based_algorithms/lb_interpolation.cpp b/src/core/grid_based_algorithms/lb_interpolation.cpp index 7b59a8200b0..982a5219662 100644 --- a/src/core/grid_based_algorithms/lb_interpolation.cpp +++ b/src/core/grid_based_algorithms/lb_interpolation.cpp @@ -87,9 +87,8 @@ const Utils::Vector3d lb_lbinterpolation_get_interpolated_velocity(const Utils::Vector3d &pos) { Utils::Vector3d interpolated_u{}; - /* calculate fluid velocity at particle's position - this is done by linear interpolation - (Eq. (11) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) */ + /* Calculate fluid velocity at particle's position. + This is done by linear interpolation (eq. (11) @cite ahlrichs99a) */ lattice_interpolation(lblattice, pos, [&interpolated_u](Lattice::index_t index, double w) { interpolated_u += w * node_u(index); diff --git a/src/core/grid_based_algorithms/lb_particle_coupling.cpp b/src/core/grid_based_algorithms/lb_particle_coupling.cpp index c97e4e39eb0..0cfe9824de8 100644 --- a/src/core/grid_based_algorithms/lb_particle_coupling.cpp +++ b/src/core/grid_based_algorithms/lb_particle_coupling.cpp @@ -120,7 +120,7 @@ namespace { */ void add_md_force(Utils::Vector3d const &pos, Utils::Vector3d const &force) { /* transform momentum transfer to lattice units - (Eq. (12) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) */ + (eq. (12) @cite ahlrichs99a) */ auto const delta_j = -(time_step / lb_lbfluid_get_lattice_speed()) * force; lb_lbinterpolation_add_force_density(pos, delta_j); } @@ -128,18 +128,17 @@ void add_md_force(Utils::Vector3d const &pos, Utils::Vector3d const &force) { /** Coupling of a single particle to viscous fluid with Stokesian friction. * - * Section II.C. Ahlrichs and Duenweg, JCP 111(17):8225 (1999) + * Section II.C. @cite ahlrichs99a * - * @param[in] p The coupled particle. - * @param[in] f_random Additional force to be included. + * @param[in] p The coupled particle. + * @param[in] f_random Additional force to be included. * - * @return The viscous coupling force plus f_random. + * @return The viscous coupling force plus f_random. */ Utils::Vector3d lb_viscous_coupling(Particle const &p, Utils::Vector3d const &f_random) { /* calculate fluid velocity at particle's position - this is done by linear interpolation - (Eq. (11) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) */ + this is done by linear interpolation (eq. (11) @cite ahlrichs99a) */ auto const interpolated_u = lb_lbinterpolation_get_interpolated_velocity(p.r.p) * lb_lbfluid_get_lattice_speed(); @@ -155,9 +154,7 @@ Utils::Vector3d lb_viscous_coupling(Particle const &p, v_drift += p.p.mu_E; #endif - /* calculate viscous force - * (Eq. (9) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) - * */ + /* calculate viscous force (eq. (9) @cite ahlrichs99a) */ auto const force = -lb_lbcoupling_get_gamma() * (p.m.v - v_drift) + f_random; add_md_force(p.r.p, force); @@ -284,7 +281,7 @@ void lb_lbcoupling_calc_particle_lattice_ia( c = ctr_type{{0, 0}}; } - /* Eq. (16) Ahlrichs and Duenweg, JCP 111(17):8225 (1999). + /* Eq. (16) @cite ahlrichs99a. * The factor 12 comes from the fact that we use random numbers * from -0.5 to 0.5 (equally distributed) which have variance 1/12. * time_step comes from the discretization. diff --git a/src/core/grid_based_algorithms/lbgpu.cpp b/src/core/grid_based_algorithms/lbgpu.cpp index 32d03be1c05..7b729aa753b 100644 --- a/src/core/grid_based_algorithms/lbgpu.cpp +++ b/src/core/grid_based_algorithms/lbgpu.cpp @@ -140,28 +140,29 @@ void lb_reinit_fluid_gpu() { } } -/** (Re-)initialize the fluid. */ +/** (Re-)initialize the fluid. + * See @cite dunweg07a and @cite dhumieres09a. + */ void lb_reinit_parameters_gpu() { lbpar_gpu.time_step = (float)time_step; lbpar_gpu.mu = 0.0; if (lbpar_gpu.viscosity > 0.0 && lbpar_gpu.agrid > 0.0 && lbpar_gpu.tau > 0.0) { - /* Eq. (80) Duenweg, Schiller, Ladd, PRE 76(3):036704 (2007). */ + /* Eq. (80) @cite dunweg07a. */ lbpar_gpu.gamma_shear = 1. - 2. / (6. * lbpar_gpu.viscosity + 1.); } if (lbpar_gpu.bulk_viscosity > 0.0) { - /* Eq. (81) Duenweg, Schiller, Ladd, PRE 76(3):036704 (2007). */ + /* Eq. (81) @cite dunweg07a. */ lbpar_gpu.gamma_bulk = 1. - 2. / (9. * lbpar_gpu.bulk_viscosity + 1.); } // By default, gamma_even and gamma_odd are chosen such that the MRT becomes // a TRT with ghost mode relaxation factors that minimize unphysical wall // slip at bounce-back boundaries. For the relation between the gammas - // achieving this, consult - // D. d’Humières, I. Ginzburg, Comp. & Math. w. App. 58(5):823–840 (2009) - // Note that the relaxation operator in Espresso is defined as + // achieving this, consult @cite dhumieres09a. + // Note that the relaxation operator in ESPResSo is defined as // m* = m_eq + gamma * (m - m_eq) // as opposed to this reference, where // m* = m + lambda * (m - m_eq) @@ -175,7 +176,7 @@ void lb_reinit_parameters_gpu() { if (lbpar_gpu.kT > 0.0) { /* fluctuating hydrodynamics ? */ - /* Eq. (51) Duenweg, Schiller, Ladd, PRE 76(3):036704 (2007).*/ + /* Eq. (51) @cite dunweg07a.*/ /* Note that the modes are not normalized as in the paper here! */ lbpar_gpu.mu = lbpar_gpu.kT * lbpar_gpu.tau * lbpar_gpu.tau / D3Q19::c_sound_sq / diff --git a/src/core/grid_based_algorithms/lbgpu_cuda.cu b/src/core/grid_based_algorithms/lbgpu_cuda.cu index a091a0331c0..af5af91432a 100644 --- a/src/core/grid_based_algorithms/lbgpu_cuda.cu +++ b/src/core/grid_based_algorithms/lbgpu_cuda.cu @@ -275,100 +275,108 @@ __device__ __inline__ float calc_mode_x_from_n(LB_nodes_gpu n_a, return 0.0; } -/** Calculate modes from the velocity densities (space-transform) +/** Calculate modes from the velocity densities (space-transform). * @param[in] n_a Local node residing in array a * @param[in] index Node index / thread index * @param[out] mode Local register values mode */ __device__ void calc_m_from_n(LB_nodes_gpu n_a, unsigned int index, Utils::Array &mode) { - // The following convention is used: - // The $\hat{c}_i$ form B. Duenweg's paper are given by: - - /* c_0 = { 0, 0, 0} - c_1 = { 1, 0, 0} - c_2 = {-1, 0, 0} - c_3 = { 0, 1, 0} - c_4 = { 0,-1, 0} - c_5 = { 0, 0, 1} - c_6 = { 0, 0,-1} - c_7 = { 1, 1, 0} - c_8 = {-1,-1, 0} - c_9 = { 1,-1, 0} - c_10 = {-1, 1, 0} - c_11 = { 1, 0, 1} - c_12 = {-1, 0,-1} - c_13 = { 1, 0,-1} - c_14 = {-1, 0, 1} - c_15 = { 0, 1, 1} - c_16 = { 0,-1,-1} - c_17 = { 0, 1,-1} - c_18 = { 0,-1, 1} */ - - // The basis vectors (modes) are constructed as follows - // $m_k = \sum_{i} e_{ki} n_{i}$, where the $e_{ki}$ form a - // linear transformation (matrix) that is given by - - /* $e{ 0,i} = 1$ - $e{ 1,i} = c_{i,x}$ - $e{ 2,i} = c_{i,y}$ - $e{ 3,i} = c_{i,z}$ - $e{ 4,i} = c_{i}^2 - 1$ - $e{ 5,i} = c_{i,x}^2 - c_{i,y}^2$ - $e{ 6,i} = c_{i}^2 - 3*c_{i,z}^2$ - $e{ 7,i} = c_{i,x}*c_{i,y}$ - $e{ 8,i} = c_{i,x}*c_{i,z}$ - $e{ 9,i} = c_{i,y}*c_{i,z}$ - $e{10,i} = (3*c_{i}^2 - 5)*c_{i,x}$ - $e{11,i} = (3*c_{i}^2 - 5)*c_{i,y}$ - $e{12,i} = (3*c_{i}^2 - 5)*c_{i,z}$ - $e{13,i} = (c_{i,y}^2 - c_{i,z}^2)*c_{i,x}$ - $e{14,i} = (c_{i,x}^2 - c_{i,z}^2)*c_{i,y}$ - $e{15,i} = (c_{i,x}^2 - c_{i,y}^2)*c_{i,z}$ - $e{16,i} = 3*c_{i}^2^2 - 6*c_{i}^2 + 1$ - $e{17,i} = (2*c_{i}^2 - 3)*(c_{i,x}^2 - c_{i,y}^2)$ - $e{18,i} = (2*c_{i}^2 - 3)*(c_{i}^2 - 3*c_{i,z}^2)$ */ - - // Such that the transformation matrix is given by - - /* {{ 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, - { 0, 1,-1, 0, 0, 0, 0, 1,-1, 1,-1, 1,-1, 1,-1, 0, 0, 0, 0}, - { 0, 0, 0, 1,-1, 0, 0, 1,-1,-1, 1, 0, 0, 0, 0, 1,-1, 1,-1}, - { 0, 0, 0, 0, 0, 1,-1, 0, 0, 0, 0, 1,-1,-1, 1, 1,-1,-1, 1}, - {-1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, - { 0, 1, 1,-1,-1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,-1,-1,-1,-1}, - { 0, 1, 1, 1, 1,-2,-2, 2, 2, 2, 2,-1,-1,-1,-1,-1,-1,-1,-1}, - { 0, 0, 0, 0, 0, 0, 0, 1, 1,-1,-1, 0, 0, 0, 0, 0, 0, 0, 0}, - { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1,-1,-1, 0, 0, 0, 0}, - { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1,-1,-1}, - { 0,-2, 2, 0, 0, 0, 0, 1,-1, 1,-1, 1,-1, 1,-1, 0, 0, 0, 0}, - { 0, 0, 0,-2, 2, 0, 0, 1,-1,-1, 1, 0, 0, 0, 0, 1,-1, 1,-1}, - { 0, 0, 0, 0, 0,-2, 2, 0, 0, 0, 0, 1,-1,-1, 1, 1,-1,-1, 1}, - { 0, 0, 0, 0, 0, 0, 0, 1,-1, 1,-1,-1, 1,-1, 1, 0, 0, 0, 0}, - { 0, 0, 0, 0, 0, 0, 0, 1,-1,-1, 1, 0, 0, 0, 0,-1, 1,-1, 1}, - { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1,-1,-1, 1,-1, 1, 1,-1}, - { 1,-2,-2,-2,-2,-2,-2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, - { 0,-1,-1, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,-1,-1,-1,-1}, - { 0,-1,-1,-1,-1, 2, 2, 2, 2, 2, 2,-1,-1,-1,-1,-1,-1,-1,-1}} */ - - // With weights - - /* q^{c_{i}} = { 1/3, 1/18, 1/18, 1/18, - 1/18, 1/18, 1/18, 1/36, - 1/36, 1/36, 1/36, 1/36, - 1/36, 1/36, 1/36, 1/36, - 1/36, 1/36, 1/36 } */ - - // Which makes the transformation satisfy the following - // orthogonality condition: - // \sum_{i} q^{c_{i}} e_{ki} e_{li} = w_{k} \delta_{kl}, - // where the weights are: - - /* w_{i} = { 1, 1/3, 1/3, 1/3, - 2/3, 4/9, 4/3, 1/9, - 1/9, 1/9, 2/3, 2/3, - 2/3, 2/9, 2/9, 2/9, - 2, 4/9, 4/3 } */ + /** + * The following convention and equations from @cite duenweg09a are used: + * The \f$\hat{c}_i\f$ are given by: + * + * \f{align*}{ + * c_{ 0} &= ( 0, 0, 0) \\ + * c_{ 1} &= ( 1, 0, 0) \\ + * c_{ 2} &= (-1, 0, 0) \\ + * c_{ 3} &= ( 0, 1, 0) \\ + * c_{ 4} &= ( 0,-1, 0) \\ + * c_{ 5} &= ( 0, 0, 1) \\ + * c_{ 6} &= ( 0, 0,-1) \\ + * c_{ 7} &= ( 1, 1, 0) \\ + * c_{ 8} &= (-1,-1, 0) \\ + * c_{ 9} &= ( 1,-1, 0) \\ + * c_{10} &= (-1, 1, 0) \\ + * c_{11} &= ( 1, 0, 1) \\ + * c_{12} &= (-1, 0,-1) \\ + * c_{13} &= ( 1, 0,-1) \\ + * c_{14} &= (-1, 0, 1) \\ + * c_{15} &= ( 0, 1, 1) \\ + * c_{16} &= ( 0,-1,-1) \\ + * c_{17} &= ( 0, 1,-1) \\ + * c_{18} &= ( 0,-1, 1) + * \f} + * + * The basis vectors (modes) are constructed as follows (eq. (111)): + * \f[m_k = \sum_{i} e_{ki} n_{i}\f] where the \f$e_{ki}\f$ form a + * linear transformation (matrix) that is given by (modified from table 1): + * + * \f{align*}{ + * e_{ 0,i} &= 1 \\ + * e_{ 1,i} &= \hat{c}_{i,x} \\ + * e_{ 2,i} &= \hat{c}_{i,y} \\ + * e_{ 3,i} &= \hat{c}_{i,z} \\ + * e_{ 4,i} &= \hat{c}_{i}^2 - 1 \\ + * e_{ 5,i} &= \hat{c}_{i,x}^2 - \hat{c}_{i,y}^2 \\ + * e_{ 6,i} &= \hat{c}_{i}^2 - 3 \hat{c}_{i,z}^2 \\ + * e_{ 7,i} &= \hat{c}_{i,x} \hat{c}_{i,y} \\ + * e_{ 8,i} &= \hat{c}_{i,x} \hat{c}_{i,z} \\ + * e_{ 9,i} &= \hat{c}_{i,y} \hat{c}_{i,z} \\ + * e_{10,i} &= (3 \hat{c}_{i}^2 - 5) \hat{c}_{i,x} \\ + * e_{11,i} &= (3 \hat{c}_{i}^2 - 5) \hat{c}_{i,y} \\ + * e_{12,i} &= (3 \hat{c}_{i}^2 - 5) \hat{c}_{i,z} \\ + * e_{13,i} &= (\hat{c}_{i,y}^2 - \hat{c}_{i,z}^2) \hat{c}_{i,x} \\ + * e_{14,i} &= (\hat{c}_{i,x}^2 - \hat{c}_{i,z}^2) \hat{c}_{i,y} \\ + * e_{15,i} &= (\hat{c}_{i,x}^2 - \hat{c}_{i,y}^2) \hat{c}_{i,z} \\ + * e_{16,i} &= 3 \hat{c}_{i}^4 - 6 \hat{c}_{i}^2 + 1 \\ + * e_{17,i} &= (2 \hat{c}_{i}^2 - 3) (\hat{c}_{i,x}^2 - \hat{c}_{i,y}^2) \\ + * e_{18,i} &= (2 \hat{c}_{i}^2 - 3) (\hat{c}_{i}^2 - 3 \hat{c}_{i,z}^2) + * \f} + * + * Such that the transformation matrix is given by: + * + * \code{.cpp} + * {{ 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, + * { 0, 1,-1, 0, 0, 0, 0, 1,-1, 1,-1, 1,-1, 1,-1, 0, 0, 0, 0}, + * { 0, 0, 0, 1,-1, 0, 0, 1,-1,-1, 1, 0, 0, 0, 0, 1,-1, 1,-1}, + * { 0, 0, 0, 0, 0, 1,-1, 0, 0, 0, 0, 1,-1,-1, 1, 1,-1,-1, 1}, + * {-1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, + * { 0, 1, 1,-1,-1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,-1,-1,-1,-1}, + * { 0, 1, 1, 1, 1,-2,-2, 2, 2, 2, 2,-1,-1,-1,-1,-1,-1,-1,-1}, + * { 0, 0, 0, 0, 0, 0, 0, 1, 1,-1,-1, 0, 0, 0, 0, 0, 0, 0, 0}, + * { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1,-1,-1, 0, 0, 0, 0}, + * { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1,-1,-1}, + * { 0,-2, 2, 0, 0, 0, 0, 1,-1, 1,-1, 1,-1, 1,-1, 0, 0, 0, 0}, + * { 0, 0, 0,-2, 2, 0, 0, 1,-1,-1, 1, 0, 0, 0, 0, 1,-1, 1,-1}, + * { 0, 0, 0, 0, 0,-2, 2, 0, 0, 0, 0, 1,-1,-1, 1, 1,-1,-1, 1}, + * { 0, 0, 0, 0, 0, 0, 0, 1,-1, 1,-1,-1, 1,-1, 1, 0, 0, 0, 0}, + * { 0, 0, 0, 0, 0, 0, 0, 1,-1,-1, 1, 0, 0, 0, 0,-1, 1,-1, 1}, + * { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1,-1,-1, 1,-1, 1, 1,-1}, + * { 1,-2,-2,-2,-2,-2,-2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, + * { 0,-1,-1, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,-1,-1,-1,-1}, + * { 0,-1,-1,-1,-1, 2, 2, 2, 2, 2, 2,-1,-1,-1,-1,-1,-1,-1,-1}} + * \endcode + * + * with weights + * + * \f[q^{c_{i}} = ( 1/3, 1/18, 1/18, 1/18, + * 1/18, 1/18, 1/18, 1/36, + * 1/36, 1/36, 1/36, 1/36, + * 1/36, 1/36, 1/36, 1/36, + * 1/36, 1/36, 1/36 )\f] + * + * which makes the transformation satisfy the following + * orthogonality condition (eq. (109)): + * \f[\sum_{i} q^{c_{i}} e_{ki} e_{li} = w_{k} \delta_{kl}\f] + * where the weights are: + * + * \f[w_{i} = ( 1, 1/3, 1/3, 1/3, + * 2/3, 4/9, 4/3, 1/9, + * 1/9, 1/9, 2/3, 2/3, + * 2/3, 2/9, 2/9, 2/9, + * 2, 4/9, 4/3 )\f] + */ for (int i = 0; i < 19; ++i) { mode[i] = calc_mode_x_from_n(n_a, index, i); } @@ -448,9 +456,10 @@ __device__ void update_rho_v(Utils::Array const &mode, u_tot[1] += mode[2]; u_tot[2] += mode[3]; - /** if forces are present, the momentum density is redefined to - * include one half-step of the force action. See the - * Chapman-Enskog expansion in [Ladd & Verberg]. */ + /** If forces are present, the momentum density is redefined to + * include one half-step of the force action. See the + * Chapman-Enskog expansion in @cite ladd01a. + */ u_tot[0] += 0.5f * node_f.force_density[0 * para->number_of_nodes + index]; u_tot[1] += 0.5f * node_f.force_density[1 * para->number_of_nodes + index]; @@ -781,9 +790,8 @@ __device__ void calc_n_from_modes_push(LB_nodes_gpu n_b, * * The populations that have propagated into a boundary node * are bounced back to the node they came from. This results - * in no slip boundary conditions. + * in no slip boundary conditions, cf. @cite ladd01a. * - * [cf. Ladd and Verberg, J. Stat. Phys. 104(5/6):1191-1251, 2001] * @param[in] index Node index / thread index * @param[in] n_curr Local node receiving the current node field * @param[in] boundaries Constant velocity at the boundary, set by the user @@ -1143,8 +1151,8 @@ calc_values_in_LB_units(LB_nodes_gpu n_a, Utils::Array &mode, // Transform the stress tensor components according to the modes that // correspond to those used by U. Schiller. In terms of populations this - // expression then corresponds exactly to those in Eqs. 116 - 121 in the - // Duenweg and Ladd paper, when these are written out in populations. + // expression then corresponds exactly to those in eq. (116)-(121) in + // @cite dunweg07a, when these are written out in populations. // But to ensure this, the expression in Schiller's modes has to be // different! @@ -1335,16 +1343,21 @@ __global__ void temperature(LB_nodes_gpu n_a, float *cpu_jsquared, } } -/** - * @param u Distance to grid point in units of agrid - * @retval Value for the interpolation function. - * see Duenweg and Ladd http://arxiv.org/abs/0803.2826 +/** Interpolation kernel. + * See @cite duenweg09a + * @param u Distance to grid point in units of agrid + * @retval Value for the interpolation function. */ __device__ __inline__ float three_point_polynomial_smallerequal_than_half(float u) { return 1.f / 3.f * (1.f + sqrtf(1.f - 3.f * u * u)); } +/** Interpolation kernel. + * See @cite duenweg09a + * @param u Distance to grid point in units of agrid + * @retval Value for the interpolation function. + */ __device__ __inline__ float three_point_polynomial_larger_than_half(float u) { return 1.f / 6.f * (5.f + -3 * fabsf(u) - sqrtf(-2.f + 6.f * fabsf(u) - 3.f * u * u)); @@ -1442,8 +1455,8 @@ velocity_interpolation(LB_nodes_gpu n_a, float *particle_position, return interpolated_u; } -/** - * (Eq. (12) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) +/** Velocity interpolation. + * Eq. (12) @cite ahlrichs99a. * @param[in] n_a Local node residing in array a * @param[in] particle_position Particle position * @param[out] node_index Node index around (8) particle @@ -1456,7 +1469,7 @@ velocity_interpolation(LB_nodes_gpu n_a, float *particle_position, Utils::Array &delta) { Utils::Array left_node_index; Utils::Array temp_delta; - // see ahlrichs + duenweg page 8227 equ (10) and (11) + // Eq. (10) and (11) in @cite ahlrichs99a page 8227 #pragma unroll for (int i = 0; i < 3; ++i) { auto const scaledpos = particle_position[i] / para->agrid - 0.5f; @@ -1507,8 +1520,8 @@ velocity_interpolation(LB_nodes_gpu n_a, float *particle_position, return interpolated_u; } -/** - * (Eq. (12) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) +/** Calculate viscous force. + * Eq. (12) @cite ahlrichs99a. * @param[in] n_a Local node residing in array a * @param[out] delta Weighting of particle position * @param[out] delta_j Weighting of particle momentum @@ -1521,7 +1534,7 @@ velocity_interpolation(LB_nodes_gpu n_a, float *particle_position, * typical) or at the source (1, swimmer only) * @param[in] friction Friction constant for the particle coupling * @tparam no_of_neighbours The number of neighbours to consider for - * interpolation + * interpolation */ template __device__ void calc_viscous_force( @@ -1591,9 +1604,8 @@ __device__ void calc_viscous_force( interpolated_u.z * para->agrid / para->tau; #endif - /** calculate viscous force - * take care to rescale velocities with time_step and transform to MD units - * (Eq. (9) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) */ + /* take care to rescale velocities with time_step and transform to MD units + * (eq. (9) @cite ahlrichs99a) */ /* Viscous force */ float3 viscforce_density{0.0f, 0.0f, 0.0f}; @@ -1611,11 +1623,11 @@ __device__ void calc_viscous_force( #endif if (para->kT > 0.0) { - /** add stochastic force of zero mean (Ahlrichs, Duenweg equ. 15)*/ + /* add stochastic force of zero mean (eq. (15) @cite ahlrichs99a) */ float4 random_floats = random_wrapper_philox( particle_data[part_index].identity, LBQ * 32, philox_counter); - /* lb_coupl_pref is stored in MD units (force) - * Eq. (16) Ahlrichs and Duenweg, JCP 111(17):8225 (1999). + /* lb_coupl_pref is stored in MD units (force). + * Eq. (16) @cite ahlrichs99a. * The factor 12 comes from the fact that we use random numbers * from -0.5 to 0.5 (equally distributed) which have variance 1/12. * time_step comes from the discretization. @@ -1626,9 +1638,8 @@ __device__ void calc_viscous_force( viscforce_density.y += lb_coupl_pref * (random_floats.x - 0.5f); viscforce_density.z += lb_coupl_pref * (random_floats.y - 0.5f); } - /** delta_j for transform momentum transfer to lattice units which is done - in calc_node_force (Eq. (12) Ahlrichs and Duenweg, JCP 111(17):8225 - (1999)) */ + /* delta_j for transform momentum transfer to lattice units which is done + in calc_node_force (eq. (12) @cite ahlrichs99a) */ // only add to particle_force for particle centre <=> (1-flag_cs) = 1 particle_force[3 * part_index + 0] += (1 - flag_cs) * viscforce_density.x; @@ -1658,14 +1669,14 @@ __device__ void calc_viscous_force( } /** Calculate the node force caused by the particles, with atomicAdd due to - * avoiding race conditions - * (Eq. (14) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) + * avoiding race conditions. + * Eq. (14) @cite ahlrichs99a. * @param[in] delta Weighting of particle position * @param[in] delta_j Weighting of particle momentum * @param[in] node_index Node index around (8) particle * @param[out] node_f Node force * @tparam no_of_neighbours The number of neighbours to consider for - * interpolation + * interpolation */ template __device__ void @@ -1689,8 +1700,8 @@ calc_node_force(Utils::Array const &delta, /** Kernel to calculate local populations from hydrodynamic fields. * The mapping is given in terms of the equilibrium distribution. * - * Eq. (2.15) Ladd, J. Fluid Mech. 271, 295-309 (1994) - * Eq. (4) in Berk Usta, Ladd and Butler, JCP 122, 094902 (2005) + * Eq. (2.15) @cite ladd94a. + * Eq. (4) in @cite usta05a. * * @param[out] n_a %Lattice site * @param[out] gpu_check Additional check if GPU kernel are executed @@ -1840,8 +1851,8 @@ __global__ void set_force_density(int single_nodeindex, float *force_density, * from given flow field velocities. The mapping is given in terms of * the equilibrium distribution. * - * Eq. (2.15) Ladd, J. Fluid Mech. 271, 295-309 (1994) - * Eq. (4) in Berk Usta, Ladd and Butler, JCP 122, 094902 (2005) + * Eq. (2.15) @cite ladd94a. + * Eq. (4) in @cite usta05a. * * @param[out] n_a Current nodes array (double buffering!) * @param[in] single_nodeindex Single node index diff --git a/src/core/immersed_boundary/ImmersedBoundaries.cpp b/src/core/immersed_boundary/ImmersedBoundaries.cpp index f9d58a8b89f..0fcf7d39406 100644 --- a/src/core/immersed_boundary/ImmersedBoundaries.cpp +++ b/src/core/immersed_boundary/ImmersedBoundaries.cpp @@ -31,7 +31,7 @@ /** Volume conservation. * Calculate volumes, volume force and add it to each virtual particle. - * This function is called from integrate_vv + * This function is called from integrate_vv. */ void ImmersedBoundaries::volume_conservation() { if (VolumeInitDone && !BoundariesFound) { @@ -89,8 +89,8 @@ int ImmersedBoundaries::volume_conservation_set_params(const int bond_type, // General bond parameters bonded_ia_params[bond_type].type = BONDED_IA_IBM_VOLUME_CONSERVATION; bonded_ia_params[bond_type].num = - 0; // This means that Espresso requires one bond partner. Here we simply - // ignore it, but Espresso cannot handle 0. + 0; // This means that ESPResSo requires one bond partner. Here we simply + // ignore it, but ESPResSo cannot handle 0. // Specific stuff if (softID > MaxNumIBM) { @@ -116,7 +116,9 @@ int ImmersedBoundaries::volume_conservation_set_params(const int bond_type, return ES_OK; } -/** Calculate partial volumes on all compute nodes and call MPI to sum up */ +/** Calculate partial volumes on all compute nodes and call MPI to sum up. + * See @cite zhang01b, @cite dupin08a, @cite kruger12a. + */ void ImmersedBoundaries::calc_volumes() { // Partial volumes for each soft particle, to be summed up @@ -183,23 +185,22 @@ void ImmersedBoundaries::calc_volumes() { return; } - // Unfold position of first node - // this is to get a continuous trajectory with no jumps when box - // boundaries are crossed + // Unfold position of first node. + // This is to get a continuous trajectory with no jumps when box + // boundaries are crossed. auto const x1 = unfolded_position(p1.r.p, p1.l.i, box_geo.length()); auto const x2 = x1 + get_mi_vector(p2->r.p, x1, box_geo); auto const x3 = x1 + get_mi_vector(p3->r.p, x1, box_geo); // Volume of this tetrahedron - // See Cha Zhang et.al. 2001, doi:10.1109/ICIP.2001.958278 - // http://research.microsoft.com/en-us/um/people/chazhang/publications/icip01_ChaZhang.pdf + // See @cite zhang01b // The volume can be negative, but it is not necessarily the "signed // volume" in the above paper (the sign of the real "signed volume" // must be calculated using the normal vector; the result of the // calculation here is simply a term in the sum required to // calculate the volume of a particle). Again, see the paper. This // should be equivalent to the formulation using vector identities - // in Krüger thesis + // in @cite kruger12a const double v321 = x3[0] * x2[1] * x1[2]; const double v231 = x2[0] * x3[1] * x1[2]; @@ -237,11 +238,11 @@ void ImmersedBoundaries::calc_volume_force() { Particle &p1 = cell->part[i]; // Check if particle has a BONDED_IA_IBM_TRIEL and a - // BONDED_IA_IBM_VOLUME_CONSERVATION Basically this loops over all - // triangles, not all particles First round to check for volume - // conservation and virtual Loop over all bonds of this particle Actually - // j loops over the bond-list, i.e. the bond partners (see - // Particle.hpp) + // BONDED_IA_IBM_VOLUME_CONSERVATION. Basically this loops over all + // triangles, not all particles. First round to check for volume + // conservation and virtual. Loop over all bonds of this particle. + // Actually j loops over the bond-list, i.e. the bond partners (see + // Particle.hpp). int softID = -1; double volRef = 0.; double kappaV = 0.; @@ -279,9 +280,9 @@ void ImmersedBoundaries::calc_volume_force() { Particle &p2 = *local_particles[p1.bl.e[j + 1]]; Particle &p3 = *local_particles[p1.bl.e[j + 2]]; - // Unfold position of first node - // this is to get a continuous trajectory with no jumps when box - // boundaries are crossed + // Unfold position of first node. + // This is to get a continuous trajectory with no jumps when box + // boundaries are crossed. auto const x1 = unfolded_position(p1.r.p, p1.l.i, box_geo.length()); // Unfolding seems to work only for the first particle of a triel @@ -290,7 +291,7 @@ void ImmersedBoundaries::calc_volume_force() { auto const a13 = get_mi_vector(p3.r.p, x1, box_geo); // Now we have the true and good coordinates - // Compute force according to eq. C.46 Krüger thesis + // Compute force according to eq. (C.46) in @cite kruger12a. // It is the same as deriving Achim's equation w.r.t x /* const double fact = kappaV * 1/6. * (IBMVolumesCurrent[softID] - volRef) / IBMVolumesCurrent[softID]; @@ -312,8 +313,8 @@ void ImmersedBoundaries::calc_volume_force() { vector_product(x2, x1, n); for (int k=0; k < 3; k++) p3.f.f[k] += fact*n[k];*/ - // This is Dupin 2008. I guess the result will be very similar as - // the code above + // This is eq. (9) in @cite dupin08a. I guess the result will be + // very similar as the code above. auto const n = vector_product(a12, a13); const double ln = n.norm(); const double A = 0.5 * ln; @@ -328,7 +329,7 @@ void ImmersedBoundaries::calc_volume_force() { p3.f.f += force; } // Iterate, increase by the number of partners of this bond + 1 for - // bond type + // bond type. j += iaparams.num + 1; } } diff --git a/src/core/immersed_boundary/ibm_tribend.cpp b/src/core/immersed_boundary/ibm_tribend.cpp index 22d409066ce..4f5632d2009 100644 --- a/src/core/immersed_boundary/ibm_tribend.cpp +++ b/src/core/immersed_boundary/ibm_tribend.cpp @@ -87,6 +87,7 @@ IBM_Tribend_CalcForce(Particle const &p1, Particle const &p2, return std::make_tuple(force1, force2, force3, force4); } +/** @details See @cite gompper96a and @cite kruger12a. */ int IBM_Tribend_SetParams(const int bond_type, const int ind1, const int ind2, const int ind3, const int ind4, const double kb, const bool flat) { @@ -143,8 +144,8 @@ int IBM_Tribend_SetParams(const int bond_type, const int ind1, const int ind2, bonded_ia_params[bond_type].p.ibm_tribend.theta0 = theta0; // NOTE: This is the bare bending modulus used by the program. // If triangle pairs appear only once, the total bending force should get a - // factor 2 For the numerical model, a factor sqrt(3) should be added, see - // Gompper&Kroll J. Phys. 1996 and Krüger thesis This is an approximation, + // factor 2. For the numerical model, a factor sqrt(3) should be added, see + // @cite gompper96a and @cite kruger12a. This is an approximation, // it holds strictly only for a sphere bonded_ia_params[bond_type].p.ibm_tribend.kb = kb; } diff --git a/src/core/immersed_boundary/ibm_tribend.hpp b/src/core/immersed_boundary/ibm_tribend.hpp index 6661d1bdab4..493c3308593 100644 --- a/src/core/immersed_boundary/ibm_tribend.hpp +++ b/src/core/immersed_boundary/ibm_tribend.hpp @@ -24,8 +24,9 @@ #include "bonded_interactions/bonded_interaction_data.hpp" -// This function is used to set the parameters -// Also calculates and stores the reference state +/** Set the IBM Tribend parameters. + * Also calculate and store the reference state. + */ int IBM_Tribend_SetParams(int bond_type, int ind1, int ind2, int ind3, int ind4, double kb, bool flat); diff --git a/src/core/immersed_boundary/ibm_triel.cpp b/src/core/immersed_boundary/ibm_triel.cpp index 4725fba1827..cbbe2365b7f 100644 --- a/src/core/immersed_boundary/ibm_triel.cpp +++ b/src/core/immersed_boundary/ibm_triel.cpp @@ -74,7 +74,7 @@ IBM_Triel_CalcForce(Particle const &p1, Particle const &p2, Particle const &p3, // Calculate the current shape of the triangle (l,lp,cos(phi),sin(phi)); // l = length between 1 and 3 - // get_mi_vector is an Espresso function which considers PBC + // get_mi_vector is an ESPResSo function which considers PBC auto const vec2 = get_mi_vector(p3.r.p, p1.r.p, box_geo); auto const l = vec2.norm(); @@ -104,7 +104,7 @@ IBM_Triel_CalcForce(Particle const &p1, Particle const &p2, Particle const &p3, const double b2 = iaparams.p.ibm_triel.b2; const double A0 = iaparams.p.ibm_triel.area0; - // Displacement gradient tensor D: Eq. (C.9) Krüger thesis + // Displacement gradient tensor D: eq. (C.9) in @cite kruger12a const double Dxx = lp / lp0; const double Dxy = ((l / l0 * cosPhi) - (lp / lp0 * cosPhi0)) / sinPhi0; const double Dyx = 0.0; @@ -120,7 +120,7 @@ IBM_Triel_CalcForce(Particle const &p1, Particle const &p2, Particle const &p3, const double i1 = (Gxx + Gyy) - 2; const double i2 = ((Gxx * Gyy) - (Gxy * Gyx)) - 1; - // Derivatives of energy density E used in chain rule below: Eq. (C.14) + // Derivatives of energy density E used in chain rule below: eq. (C.14) double dEdI1; double dEdI2; if (iaparams.p.ibm_triel.elasticLaw == tElasticLaw::NeoHookean) { diff --git a/src/core/immersed_boundary/ibm_triel.hpp b/src/core/immersed_boundary/ibm_triel.hpp index 66503daa91a..3577727c96f 100644 --- a/src/core/immersed_boundary/ibm_triel.hpp +++ b/src/core/immersed_boundary/ibm_triel.hpp @@ -23,18 +23,22 @@ #include "bonded_interactions/bonded_interaction_data.hpp" #include "config.hpp" -// This function is used to set the parameters -// Also calculates and stores the reference state +/** Set the IBM Triel parameters. + * Also calculate and store the reference state. + */ int IBM_Triel_SetParams(int bond_type, int ind1, int ind2, int ind3, double maxDist, tElasticLaw elasticLaw, double k1, double k2); -// For reading checkpoints. -// Idea: * parameters are set in the run-continue script -// * also reference shape is recomputed there -// * only pass two values here to check consistency + +/** Update the IBM Triel parameters from a checkpoint. + * Ideas: + * - parameters are set in the run-continue script + * - also reference shape is recomputed there + * - only pass two values here to check consistency + */ int IBM_Triel_ResetParams(int bond_type, double k1, double l0); -/** Calculate the forces +/** Calculate the forces. * @return the forces on @p p1, @p p2, @p p3 */ boost::optional> diff --git a/src/core/io/writer/h5md_core.cpp b/src/core/io/writer/h5md_core.cpp index 97ef8a2ae8a..590ac16cb5b 100644 --- a/src/core/io/writer/h5md_core.cpp +++ b/src/core/io/writer/h5md_core.cpp @@ -117,7 +117,7 @@ void File::InitFile() { if (check_for_H5MD_structure(m_filename)) { /* * If the file exists and has a valid H5MD structure, lets create a - * backup of it. This has the advantage, that the new file can + * backup of it. This has the advantage, that the new file can * just be deleted if the simulation crashes at some point and we * still have a valid trajectory, we can start from. */ diff --git a/src/core/metadynamics.cpp b/src/core/metadynamics.cpp index 5723d8e8c0d..cee4ff007ae 100644 --- a/src/core/metadynamics.cpp +++ b/src/core/metadynamics.cpp @@ -29,7 +29,7 @@ /** \file * - * This file contains routines to perform metadynamics. Right now, the + * This file contains routines to perform metadynamics. Right now, the * reaction coordinate is defined between two particles (either distance * or z-projected distance). Note that these * particles can be virtual sites, in order to handle molecules. @@ -99,7 +99,7 @@ void meta_init() { /** Metadynamics main function: * - Calculate reaction coordinate - * - Update profile and biased force + * - Update profile and biased force @cite marsili10a * - apply external force */ void meta_perform(const ParticleRange &particles) { @@ -143,8 +143,7 @@ void meta_perform(const ParticleRange &particles) { } /* Now update free energy profile - * Here, we're following the functional form of - * Marsili et al., J Comp Chem, 31 (2009). + * Here, we're following the functional form of @cite marsili10a. * Instead of Gaussians, we use so-called Lucy's functions */ for (int i = 0; i < meta_xi_num_bins; ++i) { diff --git a/src/core/metadynamics.hpp b/src/core/metadynamics.hpp index 97d21868e47..2346b4afcf5 100644 --- a/src/core/metadynamics.hpp +++ b/src/core/metadynamics.hpp @@ -27,7 +27,7 @@ /** \file * - * This file contains routines to perform metadynamics. Right now, the + * This file contains routines to perform metadynamics. Right now, the * reaction coordinate is defined between two particles. Note that these * particles can be virtual sites, in order to handle molecules. * diff --git a/src/core/nonbonded_interactions/thole.hpp b/src/core/nonbonded_interactions/thole.hpp index ceb7f41fc0f..0caa98d5fd3 100644 --- a/src/core/nonbonded_interactions/thole.hpp +++ b/src/core/nonbonded_interactions/thole.hpp @@ -22,6 +22,7 @@ #define _THOLE_H /** \file * Routines to calculate the Thole damping potential between particle pairs. + * See @cite thole81a. * * Implementation in \ref thole.cpp. */ @@ -47,7 +48,7 @@ inline Utils::Vector3d thole_pair_force(Particle const &p1, Particle const &p2, if (thole_s != 0 && thole_q1q2 != 0 && !(pair_bond_enum_exists_between(p1, p2, BONDED_IA_THERMALIZED_DIST))) { - // Calc damping function (see doi.org/10.1016/0301-0104(81)85176-2) + // Calc damping function (see @cite thole81a) // S(r) = 1.0 - (1.0 + thole_s*r/2.0) * exp(-thole_s*r); // Calc F = - d/dr ( S(r)*q1q2/r) = // -(1/2)*(-2+(r^2*s^2+2*r*s+2)*exp(-s*r))*q1q2/r^2 Everything before diff --git a/src/core/npt.hpp b/src/core/npt.hpp index 745e494bfab..2dd2587d3c9 100644 --- a/src/core/npt.hpp +++ b/src/core/npt.hpp @@ -33,7 +33,7 @@ typedef struct { double piston; /** inverse of \ref piston */ double inv_piston; - /** isotropic volume. Note that we use the term volume throughout, + /** isotropic volume. Note that we use the term volume throughout, * although for a 2d or 1d system we mean Area and Length respectively */ double volume; diff --git a/src/core/nsquare.hpp b/src/core/nsquare.hpp index cc8cd0ca229..ef10bd7a4b6 100644 --- a/src/core/nsquare.hpp +++ b/src/core/nsquare.hpp @@ -25,7 +25,7 @@ * This file contains the code for a simple n^2 particle loop. * * The nsquare cell system performs a full n^2 particle interaction - * calculation over the simulation box. Therefore every node just + * calculation over the simulation box. Therefore every node just * has a single cell containing all local particles plus one ghost * cell per other node. The communication is done via broadcasts * (exchange_ghosts and update_ghosts) and reduce operations diff --git a/src/core/object-in-fluid/oif_global_forces.hpp b/src/core/object-in-fluid/oif_global_forces.hpp index 4132d985ced..abed668f9d3 100644 --- a/src/core/object-in-fluid/oif_global_forces.hpp +++ b/src/core/object-in-fluid/oif_global_forces.hpp @@ -20,7 +20,7 @@ #define _OBJECT_IN_FLUID_OIF_GLOBAL_FORCES_H /** \file * Routines to calculate the OIF global forces energy or/and and force - * for a particle triple (triangle from mesh). (Dupin2007) + * for a particle triple (triangle from mesh). See @cite dupin07a. */ #include diff --git a/src/core/object-in-fluid/oif_local_forces.hpp b/src/core/object-in-fluid/oif_local_forces.hpp index 73534bbacf7..d92d7aef5c6 100644 --- a/src/core/object-in-fluid/oif_local_forces.hpp +++ b/src/core/object-in-fluid/oif_local_forces.hpp @@ -20,9 +20,11 @@ #define _OBJECT_IN_FLUID_OIF_LOCAL_FORCES_H /** \file - * Routines to calculate the OIF local forces - * for a particle quadruple (two neighboring triangles with common edge). - * (Dupin2007) \ref forces.cpp + * Routines to calculate the OIF local forces for a particle quadruple + * (two neighboring triangles with common edge). + * See @cite dupin07a. + * + * Implementation in \ref oif_local_forces.cpp */ #include "Particle.hpp" @@ -37,13 +39,15 @@ int oif_local_forces_set_params(int bond_type, double r0, double ks, double A01, double A02, double kal, double kvisc); -inline double KS(double lambda) { // Defined by (19) from Dupin2007 +/** @details see eq. (19) in @cite dupin07a */ +inline double KS(double lambda) { double res; res = (pow(lambda, 0.5) + pow(lambda, -2.5)) / (lambda + pow(lambda, -3.)); return res; } -/** Compute the OIF local forces (Dupin2007). +/** Compute the OIF local forces. + * See @cite dupin07a, @cite jancigova16a. * @param p2 %Particle of triangle 1. * @param p1 , p3 Particles common to triangle 1 and triangle 2. * @param p4 %Particle of triangle 2. @@ -147,16 +151,11 @@ calc_oif_local(Particle const &p2, Particle const &p1, Particle const &p3, } /* local area - for both triangles - only 1/3 of calculated forces are added, because each triangle will enter - this calculation 3 times (one time per edge) - - Proportional distribution of forces, implemented according to the - article I.Jancigova, I.Cimrak, Non-uniform force allocation for area - preservation in spring network models, International Journal for Numerical - Methods in Biomedical Engineering, DOI: 10.1002/cnm.2757 - - */ + * for both triangles, only 1/3 of calculated forces are added, because each + * triangle will enter this calculation 3 times (one time per edge). + * Proportional distribution of forces, implemented according to + * @cite jancigova16a. + */ if (iaparams.p.oif_local_forces.kal > TINY_OIF_ELASTICITY_COEFFICIENT) { auto handle_triangle = [](double kal, double A0, Utils::Vector3d const &fp1, diff --git a/src/core/observables/BondDihedrals.hpp b/src/core/observables/BondDihedrals.hpp index b49b83733f2..3fb1cee0e12 100644 --- a/src/core/observables/BondDihedrals.hpp +++ b/src/core/observables/BondDihedrals.hpp @@ -32,10 +32,8 @@ namespace Observables { * along the chain, in radians. * * The sign of the dihedral angles follow the IUPAC nomenclature for the - * Newman projection, specifically section "Torsion Angle" pages 2220-2221 in - * G. P. Moss, *Basic Terminology of Stereochemistry* (IUPAC Recommendations - * 1996), *Pure and Applied Chemistry* **1996**, *68*(12): 2193-2222, - * doi:[10.1351/pac199668122193](https://doi.org/10.1351%2Fpac199668122193) + * Newman projection, specifically section "Torsion Angle" pages 2220-2221 + * in @cite moss96a. * */ class BondDihedrals : public PidObservable { diff --git a/src/core/particle_data.hpp b/src/core/particle_data.hpp index 3cb85f60aab..a221849a703 100644 --- a/src/core/particle_data.hpp +++ b/src/core/particle_data.hpp @@ -60,7 +60,7 @@ enum { #ifdef EXTERNAL_FORCES /** - * \ref ParticleProperties::ext_flag "ext_flag" value for fixed coordinate + * \ref ParticleProperties::ext_flag "ext_flag" value for fixed coordinate * coord. */ #define COORD_FIXED(coord) (2u << (coord)) @@ -74,11 +74,11 @@ enum { ************************************************/ /** Highest particle number seen so far. If you leave out some - particle numbers, this number might be higher than the - true number of particles. On the other hand, if you start - your particle numbers at 0, the total number of particles - is larger by 1. -*/ + * particle numbers, this number might be higher than the + * true number of particles. On the other hand, if you start + * your particle numbers at 0, the total number of particles + * is larger by 1. + */ extern int max_seen_particle; /** total number of particles on all nodes. */ extern int n_part; @@ -86,7 +86,8 @@ extern int n_part; extern bool swimming_particles_exist; /** id->particle mapping on all nodes. This is used to find partners - of bonded interactions. */ + * of bonded interactions. + */ extern std::vector local_particles; /************************************************ @@ -102,42 +103,43 @@ void free_particle(Particle *part); /* Functions acting on Particle Lists */ /************************************************/ -/** Append a particle at the end of a particle List. - reallocates particles if necessary! - This procedure does not care for \ref local_particles. - \param l List to append the particle to. - \param part Particle to append. */ +/** Append a particle at the end of a particle list. + * Reallocate particles if necessary! + * This procedure does not care for \ref local_particles. + * \param l List to append the particle to. + * \param part Particle to append. + */ void append_unindexed_particle(ParticleList *l, Particle &&part); -/** Append a particle at the end of a particle List. - reallocates particles if necessary! - This procedure cares for \ref local_particles. - \param plist List to append the particle to. - \param part Particle to append. - \return Pointer to new location of the particle. */ +/** Append a particle at the end of a particle list. + * Reallocate particles if necessary! + * This procedure cares for \ref local_particles. + * \param plist List to append the particle to. + * \param part Particle to append. + * \return Pointer to new location of the particle. + */ Particle *append_indexed_particle(ParticleList *plist, Particle &&part); -/** Remove a particle from one particle List and append it to another. - Refill the sourceList with last particle and update its entry in - local_particles. reallocates particles if necessary. This - procedure does not care for \ref local_particles. - NOT IN USE AT THE MOMENT. - \param destList List where the particle is appended. - \param sourceList List where the particle will be removed. - \param ind Index of the particle in the sourceList. - \return Pointer to new location of the particle. +/** Remove a particle from one particle list and append it to another. + * Refill the @p sourceList with last particle and update its entry in + * local_particles. Reallocate particles if necessary. This + * procedure does not care for \ref local_particles. + * \param destList List where the particle is appended. + * \param sourceList List where the particle will be removed. + * \param ind Index of the particle in the @p sourceList. + * \return Pointer to new location of the particle. */ Particle *move_unindexed_particle(ParticleList *destList, ParticleList *sourceList, int ind); -/** Remove a particle from one particle List and append it to another. - Refill the sourceList with last particle and update its entry in - local_particles. Reallocates particles if necessary. This - procedure cares for \ref local_particles. - \param destList List where the particle is appended. - \param sourceList List where the particle will be removed. - \param ind Index of the particle in the sourceList. - \return Pointer to new location of the particle. +/** Remove a particle from one particle list and append it to another. + * Refill the @p sourceList with last particle and update its entry in + * local_particles. Reallocate particles if necessary. This + * procedure cares for \ref local_particles. + * \param destList List where the particle is appended. + * \param sourceList List where the particle will be removed. + * \param ind Index of the particle in the @p sourceList. + * \return Pointer to new location of the particle. */ Particle *move_indexed_particle(ParticleList *destList, ParticleList *sourceList, int ind); @@ -165,7 +167,7 @@ void realloc_local_particles(int part); * * @param part the identity of the particle to fetch * @return Pointer to copy of particle if it exists, - * nullptr otherwise; + * nullptr otherwise; */ const Particle &get_particle_data(int part); @@ -175,8 +177,7 @@ const Particle &get_particle_data(int part); * If the range is larger than the cache size, only * the particle that fit into the cache are fetched. * - * @param ids Ids of the particles that should be - * fetched. + * @param ids Ids of the particles that should be fetched. */ void prefetch_particle_data(std::vector ids); @@ -435,8 +436,8 @@ void remove_all_bonds_to(int part); * Move a particle to a new position. If it does not exist, it is created. * The position must be on the local node! * - * @param id the identity of the particle to move - * @param pos its new position + * @param id the identity of the particle to move + * @param pos its new position * @param _new if true, the particle is allocated, else has to exists already * * @return Pointer to the particle. @@ -501,7 +502,7 @@ void recv_particles(ParticleList *particles, int node); */ inline bool do_nonbonded(Particle const &p1, Particle const &p2) { /* check for particle 2 in particle 1's exclusion list. The exclusion list is - symmetric, so this is sufficient. */ + * symmetric, so this is sufficient. */ return std::none_of(p1.el.begin(), p1.el.end(), [&p2](int id) { return p2.p.identity == id; }); } diff --git a/src/core/rattle.hpp b/src/core/rattle.hpp index 76827ab733f..05af070ad7e 100644 --- a/src/core/rattle.hpp +++ b/src/core/rattle.hpp @@ -22,14 +22,12 @@ #define RATTLE_H /** \file - * RATTLE Algorithm (Rattle: A "Velocity" Version of the - * Shake Algorithm for Molecular Dynamics Calculations, H.C Andersen, J Comp - * Phys, 52, 24-34, 1983) + * RATTLE algorithm (@cite andersen83a). * - * For more information see \ref rattle.cpp "rattle.cpp". + * For more information see \ref rattle.cpp. */ -/** number of rigid bonds */ +/** Number of rigid bonds. */ extern int n_rigidbonds; #include "cells.hpp" @@ -43,13 +41,14 @@ void save_old_pos(const ParticleRange &particles, const ParticleRange &ghost_particles); /** Propagate velocity and position while using SHAKE algorithm for bond - * constraint.*/ + * constraint. + */ void correct_pos_shake(const ParticleRange &particles); -/** Correction of current velocities using RATTLE algorithm*/ +/** Correction of current velocities using RATTLE algorithm. */ void correct_vel_shake(); -/** set the parameter for a rigid, aka RATTLE bond */ +/** Set the parameter for a rigid, aka RATTLE, bond. */ int rigid_bond_set_params(int bond_type, double d, double p_tol, double v_tol); #endif diff --git a/src/core/reaction_ensemble.cpp b/src/core/reaction_ensemble.cpp index 507bafd92c8..ae98f5e5f59 100644 --- a/src/core/reaction_ensemble.cpp +++ b/src/core/reaction_ensemble.cpp @@ -141,7 +141,7 @@ void ReactionAlgorithm::add_reaction( new_reaction.nu_bar = calculate_nu_bar(new_reaction.reactant_coefficients, new_reaction.product_coefficients); - // make espresso count the particle numbers which take part in the reactions + // make ESPResSo count the particle numbers which take part in the reactions for (int reactant_type : new_reaction.reactant_types) init_type_map(reactant_type); for (int product_type : new_reaction.product_types) @@ -1172,9 +1172,9 @@ int WangLandauReactionEnsemble::initialize_wang_landau() { return ES_OK; } -/** - * Calculates the expression in the acceptance probability of the Wang-Landau - * reaction ensemble +/** Calculate the expression in the acceptance probability of the Wang-Landau + * reaction ensemble. + * Modify Boltzmann factor according to Wang-Landau algorithm in @cite yan02b. */ double WangLandauReactionEnsemble::calculate_acceptance_probability( SingleReaction ¤t_reaction, double E_pot_old, double E_pot_new, @@ -1197,26 +1197,19 @@ double WangLandauReactionEnsemble::calculate_acceptance_probability( } else { // pass } - // look whether the proposed state lies in the reaction coordinate space gamma - // and add the Wang-Landau modification factor, this is a bit nasty due to the - // energy collective variable case (memory layout of storage array of the - // histogram and the wang_landau_potential values is "cuboid") + // Check whether the proposed state lies in the reaction coordinate space + // gamma and add the Wang-Landau modification factor, this is a bit nasty + // due to the energy collective variable case (memory layout of storage + // array of the histogram and the wang_landau_potential values is "cuboid"). if (old_state_index >= 0 && new_state_index >= 0) { if (histogram[new_state_index] >= 0 && histogram[old_state_index] >= 0) { - bf = std::min(1.0, - bf * exp(wang_landau_potential[old_state_index] - - wang_landau_potential[new_state_index])); // modify - // Boltzmann - // factor - // according to Wang-Landau - // algorithm, according to - // grand canonical simulation - // paper "Density-of-states - // Monte Carlo method for - // simulation of fluids" - // this makes the new state being accepted with the conditional + // Modify Boltzmann factor (bf) according to Wang-Landau algorithm + // @cite yan02b. + // This makes the new state being accepted with the conditional // probability bf (bf is a transition probability = conditional - // probability from the old state to move to the new state) + // probability from the old state to move to the new state). + bf = std::min(1.0, bf * exp(wang_landau_potential[old_state_index] - + wang_landau_potential[new_state_index])); } else { if (histogram[new_state_index] >= 0 && histogram[old_state_index] < 0) bf = 10; // this makes the reaction get accepted, since we found a state @@ -1242,21 +1235,16 @@ double WangLandauReactionEnsemble::calculate_acceptance_probability( return bf; } -/** Performs a randomly selected reaction using the Wang-Landau algorithm. - * - * make sure to perform additional configuration changing steps, after the - * reaction step! like in Density-of-states Monte Carlo method for simulation - * of fluids Yan, De Pablo. this can be done with MD in the case of the - * no-energy-reweighting case, or with the functions - * do_global_mc_move_for_particles_of_type +/** Perform a randomly selected reaction using the Wang-Landau algorithm. * - * perform additional Monte Carlo moves to sample configurational - * partition function according to "Density-of-states Monte Carlo method - * for simulation of fluids" + * Make sure to perform additional configuration changing steps, after the + * reaction step! Like in @cite yan02b. This can be done with MD in the case + * of the no-energy-reweighting case, or with the function + * @ref ReactionAlgorithm::do_global_mc_move_for_particles_of_type. * - * do as many steps as needed to get to a new conformation (compare - * Density-of-states Monte Carlo method for simulation of fluids Yan, - * De Pablo) + * Perform additional Monte Carlo moves to sample configurational + * partition function according to @cite yan02b. Do as many steps + * as needed to get to a new conformation. */ int WangLandauReactionEnsemble::do_reaction(int reaction_steps) { m_WL_tries += reaction_steps; @@ -1536,13 +1524,12 @@ int WangLandauReactionEnsemble:: return index; } -/** remove bins from the range of to be sampled values if they have not been +/** Remove bins from the range of to be sampled values if they have not been * sampled. - * use with caution otherwise you produce unphysical results, do only use + * Use with caution otherwise you produce unphysical results, do only use * when you know what you want to do. This can make Wang-Landau converge on a - * reduced set gamma. use this function e.g. in do_reaction_wang_landau() for - * the diprotonic acid compare "Wang-Landau sampling with self-adaptive range" - * by Troester and Dellago + * reduced set gamma. Use this function e.g. in do_reaction_wang_landau(). For + * the diprotonic acid compare with @cite troester05a. */ void WangLandauReactionEnsemble::remove_bins_that_have_not_been_sampled() { int removed_bins = 0; @@ -1779,7 +1766,9 @@ WidomInsertion::measure_excess_chemical_potential(int reaction_id) { /** * Calculates the whole product of factorial expressions which occur in the - * reaction ensemble acceptance probability + * reaction ensemble acceptance probability. + * + * See @cite smith94c. */ double calculate_factorial_expression(SingleReaction ¤t_reaction, @@ -1791,9 +1780,9 @@ calculate_factorial_expression(SingleReaction ¤t_reaction, int N_i0 = old_particle_numbers[current_reaction.reactant_types[i]]; factorial_expr = factorial_expr * factorial_Ni0_divided_by_factorial_Ni0_plus_nu_i( - N_i0, nu_i); // zeta = 1 (see smith paper) since - // we only perform one reaction at - // one call of the function + N_i0, nu_i); // zeta = 1 (see @cite smith94c) + // since we only perform one reaction + // at one call of the function } // factorial contribution of products for (int i = 0; i < current_reaction.product_types.size(); i++) { @@ -1801,9 +1790,9 @@ calculate_factorial_expression(SingleReaction ¤t_reaction, int N_i0 = old_particle_numbers[current_reaction.product_types[i]]; factorial_expr = factorial_expr * factorial_Ni0_divided_by_factorial_Ni0_plus_nu_i( - N_i0, nu_i); // zeta = 1 (see smith paper) since - // we only perform one reaction at - // one call of the function + N_i0, nu_i); // zeta = 1 (see @cite smith94c) + // since we only perform one reaction + // at one call of the function } return factorial_expr; } diff --git a/src/core/reaction_ensemble.hpp b/src/core/reaction_ensemble.hpp index b3fc94c9a85..f73605cbb12 100644 --- a/src/core/reaction_ensemble.hpp +++ b/src/core/reaction_ensemble.hpp @@ -232,14 +232,14 @@ class ReactionAlgorithm { ///////////////////////////// actual declaration of specific reaction algorithms -/** Reaction ensemble method according to smith94x. +/** Reaction ensemble method. * Works for the reaction ensemble at constant volume and temperature. For the - * reaction ensemble at constant pressure additionally employ a barostat! + * reaction ensemble at constant pressure, additionally employ a barostat! * NOTE: a chemical reaction consists of a forward and backward reaction. * Here both reactions have to be defined separately. The extent of the * reaction is here chosen to be +1. If the reaction trial move for a * dissociation of HA is accepted then there is one more dissociated ion - * pair H+ and A- + * pair H+ and A-. Implementation of @cite smith94c. */ class ReactionEnsemble : public ReactionAlgorithm { public: @@ -358,7 +358,7 @@ class WangLandauReactionEnsemble : public ReactionAlgorithm { }; /** - * Constant-pH Ensemble, for derivation see Reed and Reed 1992. + * Constant-pH Ensemble, for derivation see @cite reed92a. * For the constant pH reactions you need to provide the deprotonation and * afterwards the corresponding protonation reaction (in this order). If you * want to deal with multiple reactions do it multiple times. Note that there is diff --git a/src/core/rotation.cpp b/src/core/rotation.cpp index 73fad9d5158..46ea60a879c 100644 --- a/src/core/rotation.cpp +++ b/src/core/rotation.cpp @@ -43,7 +43,8 @@ #include /** Calculate the derivatives of the quaternion and angular acceleration - * for a given particle + * for a given particle. + * See @cite sonnenschein85a * @param[in] p %Particle * @param[out] Qd First derivative of the particle quaternion * @param[out] Qdd Second derivative of the particle quaternion @@ -54,8 +55,7 @@ static void define_Qdd(Particle const &p, double Qd[4], double Qdd[4], double S[3], double Wd[3]) { /* calculate the first derivative of the quaternion */ - /* Taken from "An improved algorithm for molecular dynamics simulation of - * rigid molecules", Sonnenschein, Roland (1985), Eq. 4.*/ + /* Eq. (4) @cite sonnenschein85a */ Qd[0] = 0.5 * (-p.r.quat[1] * p.m.omega[0] - p.r.quat[2] * p.m.omega[1] - p.r.quat[3] * p.m.omega[2]); @@ -69,8 +69,7 @@ static void define_Qdd(Particle const &p, double Qd[4], double Qdd[4], p.r.quat[0] * p.m.omega[2]); /* Calculate the angular acceleration. */ - /* Taken from "An improved algorithm for molecular dynamics simulation of - * rigid molecules", Sonnenschein, Roland (1985), Eq. 5.*/ + /* Eq. (5) @cite sonnenschein85a */ if (p.p.rotation & ROTATION_X) Wd[0] = (p.f.torque[0] + p.m.omega[1] * p.m.omega[2] * (p.p.rinertia[1] - p.p.rinertia[2])) / @@ -93,8 +92,7 @@ static void define_Qdd(Particle const &p, double Qd[4], double Qdd[4], auto const S1 = Qd[0] * Qd[0] + Qd[1] * Qd[1] + Qd[2] * Qd[2] + Qd[3] * Qd[3]; /* Calculate the second derivative of the quaternion. */ - /* Taken from "An improved algorithm for molecular dynamics simulation of - * rigid molecules", Sonnenschein, Roland (1985), Eq. 8.*/ + /* Eq. (8) @cite sonnenschein85a */ Qdd[0] = 0.5 * (-p.r.quat[1] * Wd[0] - p.r.quat[2] * Wd[1] - p.r.quat[3] * Wd[2]) - p.r.quat[0] * S1; @@ -133,8 +131,7 @@ void propagate_omega_quat_particle(Particle &p) { define_Qdd(p, Qd.data(), Qdd.data(), S.data(), Wd.data()); - /* Taken from "On the numerical integration of motion for rigid polyatomics: - * The modified quaternion approach", Omeylan, Igor (1998), Eq. 12.*/ + /* Eq. (12) @cite sonnenschein85a. */ auto const lambda = 1 - S[0] * time_step_squared_half - sqrt(1 - time_step_squared * diff --git a/src/core/specfunc.cpp b/src/core/specfunc.cpp index 4a1702d36bb..691e0d14771 100644 --- a/src/core/specfunc.cpp +++ b/src/core/specfunc.cpp @@ -18,12 +18,6 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see . */ -/** \file - Special functions, see \ref specfunc.hpp "specfunc.hpp" -*/ -#include "specfunc.hpp" -#include "polynom.hpp" -#include /* Original gsl header * specfunc/bessel_K0.cpp @@ -45,36 +39,31 @@ * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */ -/* Original Author: G. Jungman */ +/* Original Author: G. Jungman */ + +/** \file + * Implementation of \ref specfunc.hpp. + */ +#include "specfunc.hpp" +#include "polynom.hpp" +#include /************************************************ * chebychev expansions ************************************************/ /* Note that the first coefficient already includes the constant offsets */ -/* based on SLATEC bk0(), bk0e() */ - -/* chebyshev expansions - - series for bk0 on the interval 0. to 4.00000d+00 - with weighted error 3.57e-19 - log weighted error 18.45 - significant figures required 17.99 - decimal places required 18.97 - - series for ak0 on the interval 1.25000d-01 to 5.00000d-01 - with weighted error 5.34e-17 - log weighted error 16.27 - significant figures required 14.92 - decimal places required 16.89 - - series for ak02 on the interval 0. to 1.25000d-01 - with weighted error 2.34e-17 - log weighted error 16.63 - significant figures required 14.67 - decimal places required 17.20 -*/ - +/** @name Chebyshev expansions based on SLATEC bk0(), bk0e() */ +/*@{*/ +/** Series for @c bk0. + * On the interval 0. to 4.00000d+00 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 3.57e-19 | + * | log weighted error | 18.45 | + * | significant figures required | 17.99 | + * | decimal places required | 18.97 | + */ static double bk0_data[11] = { -.5 - 0.03532739323390276872, 0.3442898999246284869, 0.03597993651536150163, 0.00126461541144692592, @@ -84,6 +73,15 @@ static double bk0_data[11] = { 0.00000000000000000035}; static Polynom bk0_cs{bk0_data}; +/** Series for @c ak0. + * On the interval 1.25000d-01 to 5.00000d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 5.34e-17 | + * | log weighted error | 16.27 | + * | significant figures required | 14.92 | + * | decimal places required | 16.89 | + */ static double ak0_data[17] = { 2.5 - 0.07643947903327941, -0.02235652605699819, 0.00077341811546938, -0.00004281006688886, 0.00000308170017386, -0.00000026393672220, @@ -93,6 +91,15 @@ static double ak0_data[17] = { -0.00000000000000033, 0.00000000000000005}; static Polynom ak0_cs{ak0_data}; +/** Series for @c ak02. + * On the interval 0. to 1.25000d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 2.34e-17 | + * | log weighted error | 16.63 | + * | significant figures required | 14.67 | + * | decimal places required | 17.20 | + */ static double ak02_data[14] = { 2.5 - 0.01201869826307592, -0.00917485269102569, 0.00014445509317750, -0.00000401361417543, 0.00000015678318108, -0.00000000777011043, @@ -100,31 +107,19 @@ static double ak02_data[14] = { -0.00000000000020743, 0.00000000000001925, -0.00000000000000192, 0.00000000000000020, -0.00000000000000002}; static Polynom ak02_cs{ak02_data}; - -/* based on SLATEC besi0 */ - -/* chebyshev expansions - - series for bi0 on the interval 0. to 9.00000d+00 - with weighted error 2.46e-18 - log weighted error 17.61 - significant figures required 17.90 - decimal places required 18.15 - - series for ai0 on the interval 1.25000d-01 to 3.33333d-01 - with weighted error 7.87e-17 - log weighted error 16.10 - significant figures required 14.69 - decimal places required 16.76 - - - series for ai02 on the interval 0. to 1.25000d-01 - with weighted error 3.79e-17 - log weighted error 16.42 - significant figures required 14.86 - decimal places required 17.09 -*/ - +/*@}*/ + +/** @name Chebyshev expansions based on SLATEC besi0() */ +/*@{*/ +/** Series for @c bi0. + * On the interval 0. to 9.00000d+00 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 2.46e-18 | + * | log weighted error | 17.61 | + * | significant figures required | 17.90 | + * | decimal places required | 18.15 | + */ static double bi0_data[12] = { 5.5 - .07660547252839144951, 1.92733795399380827000, .22826445869203013390, .01304891466707290428, .00043442709008164874, .00000942265768600193, @@ -132,6 +127,15 @@ static double bi0_data[12] = { .00000000000009579451, .00000000000000053339, .00000000000000000245}; static Polynom bi0_cs{bi0_data}; +/** Series for @c ai0. + * On the interval 1.25000d-01 to 3.33333d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 7.87e-17 | + * | log weighted error | 16.10 | + * | significant figures required | 14.69 | + * | decimal places required | 16.76 | + */ static double ai0_data[21] = { .75 + .07575994494023796, .00759138081082334, .00041531313389237, .00001070076463439, -.00000790117997921, -.00000078261435014, @@ -142,6 +146,15 @@ static double ai0_data[21] = { .00000000000000314, -.00000000000000071, .00000000000000007}; static Polynom ai0_cs{ai0_data}; +/** Series for @c ai02. + * On the interval 0. to 1.25000d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 3.79e-17 | + * | log weighted error | 16.42 | + * | significant figures required | 14.86 | + * | decimal places required | 17.09 | + */ static double ai02_data[22] = { .75 + .05449041101410882, .00336911647825569, .00006889758346918, .00000289137052082, .00000020489185893, .00000002266668991, @@ -152,30 +165,19 @@ static double ai02_data[22] = { .00000000000000176, -.00000000000000034, -.00000000000000027, .00000000000000003}; static Polynom ai02_cs{ai02_data}; - -/* based on SLATEC besk1(), besk1e() */ - -/* chebyshev expansions - - series for bk1 on the interval 0. to 4.00000d+00 - with weighted error 7.02e-18 - log weighted error 17.15 - significant figures required 16.73 - decimal places required 17.67 - - series for ak1 on the interval 1.25000d-01 to 5.00000d-01 - with weighted error 6.06e-17 - log weighted error 16.22 - significant figures required 15.41 - decimal places required 16.83 - - series for ak12 on the interval 0. to 1.25000d-01 - with weighted error 2.58e-17 - log weighted error 16.59 - significant figures required 15.22 - decimal places required 17.16 -*/ - +/*@}*/ + +/** @name Chebyshev expansions based on SLATEC besk1(), besk1e() */ +/*@{*/ +/** Series for @c bk1. + * On the interval 0. to 4.00000d+00 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 7.02e-18 | + * | log weighted error | 17.15 | + * | significant figures required | 16.73 | + * | decimal places required | 17.67 | + */ static double bk1_data[11] = { 1.5 + 0.0253002273389477705, -0.3531559607765448760, -0.1226111808226571480, -0.0069757238596398643, -0.0001730288957513052, -0.0000024334061415659, @@ -184,6 +186,15 @@ static double bk1_data[11] = { static Polynom bk1_cs{bk1_data}; +/** Series for @c ak1. + * On the interval 1.25000d-01 to 5.00000d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 6.06e-17 | + * | log weighted error | 16.22 | + * | significant figures required | 15.41 | + * | decimal places required | 16.83 | + */ static double ak1_data[17] = { 2.5 + 0.27443134069738830, 0.07571989953199368, -0.00144105155647540, 0.00006650116955125, -0.00000436998470952, 0.00000035402774997, @@ -193,6 +204,15 @@ static double ak1_data[17] = { 0.00000000000000038, -0.00000000000000006}; static Polynom ak1_cs{ak1_data}; +/** Series for @c ak12. + * On the interval 0. to 1.25000d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 2.58e-17 | + * | log weighted error | 16.59 | + * | significant figures required | 15.22 | + * | decimal places required | 17.16 | + */ static double ak12_data[14] = { 2.5 + 0.06379308343739001, 0.02832887813049721, -0.00024753706739052, 0.00000577197245160, -0.00000020689392195, 0.00000000973998344, @@ -200,30 +220,19 @@ static double ak12_data[14] = { 0.00000000000023720, -0.00000000000002176, 0.00000000000000215, -0.00000000000000022, 0.00000000000000002}; static Polynom ak12_cs{ak12_data}; - -/* based on SLATEC besi1(), besi1e() */ - -/* chebyshev expansions - - series for bi1 on the interval 0. to 9.00000d+00 - with weighted error 2.40e-17 - log weighted error 16.62 - significant figures required 16.23 - decimal places required 17.14 - - series for ai1 on the interval 1.25000d-01 to 3.33333d-01 - with weighted error 6.98e-17 - log weighted error 16.16 - significant figures required 14.53 - decimal places required 16.82 - - series for ai12 on the interval 0. to 1.25000d-01 - with weighted error 3.55e-17 - log weighted error 16.45 - significant figures required 14.69 - decimal places required 17.12 -*/ - +/*@}*/ + +/** @name Chebyshev expansions based on SLATEC besi1(), besi1e() */ +/*@{*/ +/** Series for @c bi1. + * On the interval 0. to 9.00000d+00 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 2.40e-17 | + * | log weighted error | 16.62 | + * | significant figures required | 16.23 | + * | decimal places required | 17.14 | + */ static double bi1_data[11] = { 1.75 - 0.001971713261099859, 0.407348876675464810, 0.034838994299959456, 0.001545394556300123, 0.000041888521098377, 0.000000764902676483, @@ -231,6 +240,15 @@ static double bi1_data[11] = { 0.000000000000004741, 0.000000000000000024}; static Polynom bi1_cs{bi1_data}; +/** Series for @c ai1. + * On the interval 1.25000d-01 to 3.33333d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 6.98e-17 | + * | log weighted error | 16.16 | + * | significant figures required | 14.53 | + * | decimal places required | 16.82 | + */ static double ai1_data[21] = { .75 - 0.02846744181881479, -0.01922953231443221, -0.00061151858579437, -0.00002069971253350, 0.00000858561914581, 0.00000104949824671, @@ -241,6 +259,15 @@ static double ai1_data[21] = { -0.00000000000000333, 0.00000000000000071, -0.00000000000000006}; static Polynom ai1_cs{ai1_data}; +/** Series for @c ai12. + * On the interval 0. to 1.25000d-01 + * | Label | Value | + * | ---------------------------: | :------- | + * | with weighted error | 3.55e-17 | + * | log weighted error | 16.45 | + * | significant figures required | 14.69 | + * | decimal places required | 17.12 | + */ static double ai12_data[22] = { .75 + 0.02857623501828014, -0.00976109749136147, -0.00011058893876263, -0.00000388256480887, -0.00000025122362377, -0.00000002631468847, @@ -251,13 +278,10 @@ static double ai12_data[22] = { -0.00000000000000186, 0.00000000000000033, 0.00000000000000028, -0.00000000000000003}; static Polynom ai12_cs{ai12_data}; +/*@}*/ -/************************************************ - * chebychev expansions - ************************************************/ - -/* coefficients for Maclaurin summation in hzeta() - * B_{2j}/(2j)! +/** Coefficients for Maclaurin summation in hzeta(). + * B_{2j}/(2j)! */ static double hzeta_c[15] = { 1.00000000000000000000000000000, 0.083333333333333333333333333333, @@ -269,10 +293,6 @@ static double hzeta_c[15] = { -1.3954464685812523340707686264e-19, 3.5347070396294674716932299778e-21, -8.9535174270375468504026113181e-23}; -/************************************************ - * functions - ************************************************/ - double hzeta(double s, double q) { double max_bits = 54.0; int jmax = 12, kmax = 10; @@ -360,9 +380,10 @@ double K1(double x) { * optimized K0/1 implementations for 10^(-14) precision ***********************************************************/ -/** necessary orders for K0/1 from 2 up to 22 for 10^-14 precision. Note that at - 8 the expansion changes. From 23 to 26 order 2 is used, above order 1. For - the latter cases separate implementations are necessary. */ +/** necessary orders for K0/1 from 2 up to 22 for 10^-14 precision. Note that + * at 8 the expansion changes. From 23 to 26 order 2 is used, above order 1. + * For the latter cases, separate implementations are necessary. + */ static int ak01_orders[] = { /* 2 - 8 */ 11, 11, 10, 10, 9, 9, diff --git a/src/core/specfunc.hpp b/src/core/specfunc.hpp index c6ed021b576..5d9841a34ce 100644 --- a/src/core/specfunc.hpp +++ b/src/core/specfunc.hpp @@ -22,8 +22,8 @@ * This file contains implementations for some special functions which are * needed by the MMM family of algorithms. This are the modified Hurwitz zeta * function and the modified Bessel functions of first and second kind. The - * implementations are based on the GSL code (see \ref specfunc.cpp - * "specfunc.cpp" for the original GSL header). + * implementations are based on the GSL code (see \ref specfunc.cpp for the + * original GSL header). * * The Hurwitz zeta function is evaluated using the Euler-MacLaurin summation * formula, the Bessel functions are evaluated using several different diff --git a/src/core/statistics.cpp b/src/core/statistics.cpp index e022a364aa6..95fa55a0ad0 100644 --- a/src/core/statistics.cpp +++ b/src/core/statistics.cpp @@ -48,7 +48,8 @@ #include /** Previous particle configurations (needed for offline analysis and - correlation analysis) */ + * correlation analysis) + */ std::vector> configs; int n_configs = 0; int n_part_conf = 0; @@ -64,9 +65,7 @@ double mindist(PartCfg &partCfg, IntList const &set1, IntList const &set2) { auto mindist2 = std::numeric_limits::infinity(); for (auto jt = partCfg.begin(); jt != partCfg.end(); ++jt) { - /* check which sets particle j belongs to - bit 0: set1, bit1: set2 - */ + /* check which sets particle j belongs to (bit 0: set1, bit1: set2) */ auto in_set = 0; if (set1.empty() || contains(set1, jt->p.type)) in_set = 1; @@ -228,12 +227,12 @@ void calc_part_distribution(PartCfg &partCfg, int const *p1_types, int n_p1, else inv_bin_width = (double)r_bins / (r_max - r_min); - /* particle loop: p1_types*/ + /* particle loop: p1_types */ for (auto const &p1 : partCfg) { for (t1 = 0; t1 < n_p1; t1++) { if (p1.p.type == p1_types[t1]) { min_dist2 = start_dist2; - /* particle loop: p2_types*/ + /* particle loop: p2_types */ for (auto const &p2 : partCfg) { if (p1 != p2) { for (t2 = 0; t2 < n_p2; t2++) { @@ -302,14 +301,14 @@ void calc_rdf(PartCfg &partCfg, int const *p1_types, int n_p1, inv_bin_width = 1.0 / bin_width; for (i = 0; i < r_bins; i++) rdf[i] = 0.0; - /* particle loop: p1_types*/ + /* particle loop: p1_types */ for (auto it = partCfg.begin(); it != partCfg.end(); ++it) { for (t1 = 0; t1 < n_p1; t1++) { if (it->p.type == p1_types[t1]) { /* distinguish mixed and identical rdf's */ auto jt = mixed_flag ? partCfg.begin() : std::next(it); - /* particle loop: p2_types*/ + /* particle loop: p2_types */ for (; jt != partCfg.end(); ++jt) { for (t2 = 0; t2 < n_p2; t2++) { if (jt->p.type == p2_types[t2]) { @@ -597,7 +596,7 @@ int calc_cylindrical_average( } // Now we turn the counts into densities by dividing by one radial - // bin (binvolume). We also divide the velocities by the counts. + // bin (binvolume). We also divide the velocities by the counts. double binvolume; for (unsigned int type_id = 0; type_id < types.size(); type_id++) { for (int index_radial = 0; index_radial < bins_radial; index_radial++) { diff --git a/src/core/unit_tests/LocalBox_test.cpp b/src/core/unit_tests/LocalBox_test.cpp index 26d899220d4..9bb940c788e 100644 --- a/src/core/unit_tests/LocalBox_test.cpp +++ b/src/core/unit_tests/LocalBox_test.cpp @@ -26,10 +26,7 @@ #include -/** - * @brief Check that the box corners and side length - * agree. - */ +/* Check that the box corners and side length agree. */ template void check_length(LocalBox box) { auto const expected = box.my_right() - box.my_left(); auto const result = box.length(); @@ -60,4 +57,4 @@ BOOST_AUTO_TEST_CASE(constructors) { BOOST_CHECK(boost::equal(boundaries, box.boundary())); check_length(box); } -} \ No newline at end of file +} diff --git a/src/core/unit_tests/MpiCallbacks_test.cpp b/src/core/unit_tests/MpiCallbacks_test.cpp index 6f29fb93005..67c2102ca68 100644 --- a/src/core/unit_tests/MpiCallbacks_test.cpp +++ b/src/core/unit_tests/MpiCallbacks_test.cpp @@ -18,10 +18,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the MpiCallbacks class. - * - */ +/* Unit tests for the MpiCallbacks class. */ #define BOOST_TEST_NO_MAIN #define BOOST_TEST_MODULE MpiCallbacks test @@ -55,7 +52,7 @@ BOOST_AUTO_TEST_CASE(invoke_test) { BOOST_CHECK_EQUAL(f(i, d), (invoke(f, ia))); } -/** +/* * Test that the implementation of callback_model_t * correctly deserialize the parameters and call * the callback with them. diff --git a/src/core/unit_tests/None_test.cpp b/src/core/unit_tests/None_test.cpp index 8ab2ff44e52..ee96084691d 100644 --- a/src/core/unit_tests/None_test.cpp +++ b/src/core/unit_tests/None_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the ScriptInterface::None class. - * - */ +/* Unit tests for the ScriptInterface::None class. */ #define BOOST_TEST_MODULE None test #define BOOST_TEST_DYN_LINK diff --git a/src/core/unit_tests/ParallelScriptInterface_test.cpp b/src/core/unit_tests/ParallelScriptInterface_test.cpp index c227f0cb7bf..4e56ee399a9 100644 --- a/src/core/unit_tests/ParallelScriptInterface_test.cpp +++ b/src/core/unit_tests/ParallelScriptInterface_test.cpp @@ -74,7 +74,7 @@ bool TestClass::destructed = false; std::pair TestClass::last_method_parameters; std::pair TestClass::last_parameter; -/** +/* * Check that instances are created and correctly destroyed on * the slave nodes. */ @@ -99,7 +99,7 @@ BOOST_AUTO_TEST_CASE(ctor_dtor) { BOOST_CHECK(TestClass::destructed); } -/** +/* * Check that parameters are forwarded correctly. */ BOOST_AUTO_TEST_CASE(set_parameter) { diff --git a/src/core/unit_tests/ParticleCache_test.cpp b/src/core/unit_tests/ParticleCache_test.cpp index 5d714e092cd..a951cae5e1c 100644 --- a/src/core/unit_tests/ParticleCache_test.cpp +++ b/src/core/unit_tests/ParticleCache_test.cpp @@ -18,10 +18,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the MpiCallbacks class. - * - */ +/* Unit tests for the MpiCallbacks class. */ #include #include diff --git a/src/core/unit_tests/Particle_test.cpp b/src/core/unit_tests/Particle_test.cpp index 59cbd0cc3b7..93e066031de 100644 --- a/src/core/unit_tests/Particle_test.cpp +++ b/src/core/unit_tests/Particle_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Particle struct. - * - */ +/* Unit tests for the Particle struct. */ #define BOOST_TEST_MODULE Particle test #define BOOST_TEST_DYN_LINK @@ -116,4 +113,4 @@ BOOST_AUTO_TEST_CASE(properties_serialization) { BOOST_CHECK_EQUAL(ia.bytes_read(), expected_size); BOOST_CHECK_EQUAL(out.identity, prop.identity); } -} \ No newline at end of file +} diff --git a/src/core/unit_tests/RuntimeErrorCollector_test.cpp b/src/core/unit_tests/RuntimeErrorCollector_test.cpp index b6530a9fa0a..87c6af8cb9e 100644 --- a/src/core/unit_tests/RuntimeErrorCollector_test.cpp +++ b/src/core/unit_tests/RuntimeErrorCollector_test.cpp @@ -19,10 +19,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the ErrorHandling::RuntimeErrorCollector class. - * - */ +/* Unit tests for the ErrorHandling::RuntimeErrorCollector class. */ #include #include @@ -67,7 +64,7 @@ BOOST_AUTO_TEST_CASE(count) { Testing::reduce_and_check(world, rec.count() == 0); - /** MPI guarantees that size >= 1 and rank 0 exists. */ + /* MPI guarantees that size >= 1 and rank 0 exists. */ if (world.rank() == (world.size() - 1)) { rec.error("Test_error", "Test_functions", "Test_file", 42); } @@ -78,10 +75,10 @@ BOOST_AUTO_TEST_CASE(count) { Testing::reduce_and_check(world, rec.count() == world.size() + 1); - /** There should now be one error and world.size() warnings */ + /* There should now be one error and world.size() warnings */ Testing::reduce_and_check(world, rec.count(RuntimeError::ErrorLevel::ERROR) == 1); - /** All messages are at least WARNING or higher. */ + /* All messages are at least WARNING or higher. */ { /* Beware of the execution order */ int total = rec.count(); @@ -90,7 +87,7 @@ BOOST_AUTO_TEST_CASE(count) { } } -/** +/* * Check the message gathering. Every node generates a runtime error * and a warning. Then we gather the messages * on the master and check if we got the correct messages. Then we @@ -105,24 +102,24 @@ BOOST_AUTO_TEST_CASE(gather) { rec.warning("Test_error", "Test_functions", "Test_file", world.rank()); if (world.rank() == 0) { - /** Gathered error messages */ + /* Gathered error messages */ auto results = rec.gather(); - /** Track how many messages we have seen from which node. */ + /* Track how many messages we have seen from which node. */ std::vector present(world.size()); for (const auto &err : results) { present[err.who()]++; } - /** Check if we got 2 messages from every node. */ + /* Check if we got 2 messages from every node. */ BOOST_CHECK(std::all_of(present.begin(), present.end(), [](int i) { return i == 2; })); - /** Count warnings, should be world.rank() many */ + /* Count warnings, should be world.rank() many */ BOOST_CHECK(std::count_if( results.begin(), results.end(), [](const RuntimeError &e) { return e.level() == RuntimeError::ErrorLevel::WARNING; }) == world.size()); - /** Count errors, should be world.rank() many */ + /* Count errors, should be world.rank() many */ BOOST_CHECK(std::count_if( results.begin(), results.end(), [](const RuntimeError &e) { return e.level() == RuntimeError::ErrorLevel::ERROR; diff --git a/src/core/unit_tests/RuntimeError_test.cpp b/src/core/unit_tests/RuntimeError_test.cpp index 9e1cd76f1e6..8d0aa9ea453 100644 --- a/src/core/unit_tests/RuntimeError_test.cpp +++ b/src/core/unit_tests/RuntimeError_test.cpp @@ -19,10 +19,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the ErrorHandling::RuntimeError class. - * - */ +/* Unit tests for the ErrorHandling::RuntimeError class. */ #include #include @@ -39,7 +36,7 @@ using ErrorHandling::RuntimeError; using std::string; -/** Check constructor and getters */ +/* Check constructor and getters */ BOOST_AUTO_TEST_CASE(values) { RuntimeError::ErrorLevel level = RuntimeError::ErrorLevel::WARNING; string what("Test error"); @@ -57,7 +54,7 @@ BOOST_AUTO_TEST_CASE(values) { BOOST_CHECK(line == err.line()); } -/** Check copy ctor */ +/* Check copy ctor */ BOOST_AUTO_TEST_CASE(def_ctor_and_assignment) { RuntimeError::ErrorLevel level = RuntimeError::ErrorLevel::WARNING; string what("Test error"); @@ -67,7 +64,7 @@ BOOST_AUTO_TEST_CASE(def_ctor_and_assignment) { int line(42); RuntimeError err(level, who, what, function, file, line); - /** Copy ctor */ + /* Copy ctor */ RuntimeError err2(err); // NOLINT (local copy 'err2' of the variable 'err' is // never modified; consider avoiding the copy) @@ -79,7 +76,7 @@ BOOST_AUTO_TEST_CASE(def_ctor_and_assignment) { BOOST_CHECK(line == err2.line()); } -/** Check the serialization */ +/* Check the serialization */ BOOST_AUTO_TEST_CASE(serialization) { std::stringstream ss; boost::archive::text_oarchive oa(ss); @@ -92,16 +89,16 @@ BOOST_AUTO_TEST_CASE(serialization) { int line(21); RuntimeError err(level, who, what, function, file, line); - /** Serialize to string stream */ + /* Serialize to string stream */ oa << err; boost::archive::text_iarchive ia(ss); RuntimeError err2; - /** Deserialize into empty instance */ + /* Deserialize into empty instance */ ia >> err2; - /** Check that the result is equal to the original instance */ + /* Check that the result is equal to the original instance */ BOOST_CHECK((err.level() == err2.level()) && (err.who() == err2.who()) && (err.what() == err2.what()) && (err.function() == err2.function()) && diff --git a/src/core/unit_tests/ScriptInterface_test.cpp b/src/core/unit_tests/ScriptInterface_test.cpp index 2429364a098..a59d8c80b58 100644 --- a/src/core/unit_tests/ScriptInterface_test.cpp +++ b/src/core/unit_tests/ScriptInterface_test.cpp @@ -35,9 +35,7 @@ using namespace ScriptInterface; namespace Testing { -/** - * @brief Mock to test ScriptInterface. - */ +/* Mock to test ScriptInterface. */ struct ScriptInterfaceTest : public ScriptInterface::ScriptInterfaceBase { VariantMap get_parameters() const override { VariantMap ret; @@ -101,7 +99,7 @@ BOOST_AUTO_TEST_CASE(non_copyable) { static_assert(!std::is_copy_assignable::value, ""); } -/** +/* * We check the default implementations of set_parameters * and get_parameter of ScriptInterface (this is the only * logic in the class). diff --git a/src/core/unit_tests/u32_to_u64_test.cpp b/src/core/unit_tests/u32_to_u64_test.cpp index a71a1a26db9..dccc57bff3c 100644 --- a/src/core/unit_tests/u32_to_u64_test.cpp +++ b/src/core/unit_tests/u32_to_u64_test.cpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2018-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2018-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #define BOOST_TEST_MODULE Utils::sgn test #define BOOST_TEST_DYN_LINK diff --git a/src/core/virtual_sites/lb_inertialess_tracers.cpp b/src/core/virtual_sites/lb_inertialess_tracers.cpp index 466b7aa868d..ca95eb24aee 100644 --- a/src/core/virtual_sites/lb_inertialess_tracers.cpp +++ b/src/core/virtual_sites/lb_inertialess_tracers.cpp @@ -60,15 +60,11 @@ bool in_local_domain(Utils::Vector3d const &pos) { } } // namespace -/**************** - IBM_ForcesIntoFluid_CPU - - Puts the calculated force stored on the ibm particles into the fluid by -updating the lbfields structure - Calls couple_trace_to_fluid for each node - Called from the integrate loop right after the forces have been calculated -*****************/ - +/** Put the calculated force stored on the ibm particles into the fluid by + * updating the @ref lbfields structure. + * Called from the integration loop right after the forces have been + * calculated. + */ void IBM_ForcesIntoFluid_CPU() { // Update the forces on the ghost particles ghost_communicator(&cell_structure.exchange_ghosts_comm, GHOSTTRANS_FORCE); @@ -101,12 +97,10 @@ void IBM_ForcesIntoFluid_CPU() { } } -/************* - IBM_UpdateParticlePositions -This function is called from the integrate right after the LB update -Interpolates LB velocity at the particle positions and propagates the particles -**************/ - +/** Interpolate LB velocity at the particle positions and propagate the + * particles. + * Called from the integration loop right after the LB update. + */ void IBM_UpdateParticlePositions(ParticleRange particles) { // Get velocities if (lattice_switch == ActiveLB::CPU) @@ -148,12 +142,7 @@ void IBM_UpdateParticlePositions(ParticleRange particles) { } } -/************* - CoupleIBMParticleToFluid -This function puts the momentum of a given particle into the LB fluid - only for -CPU -**************/ - +/** Put the momentum of a given particle into the LB fluid. */ void CoupleIBMParticleToFluid(Particle *p) { // Convert units from MD to LB double delta_j[3]; @@ -188,23 +177,20 @@ void CoupleIBMParticleToFluid(Particle *p) { } } -/****************** - GetIBMInterpolatedVelocity -Very similar to the velocity interpolation done in standard Espresso, except -that we add the f/2 contribution - only for CPU -*******************/ - +/** Calculate the LB fluid velocity at a particle position. + * Very similar to the velocity interpolation done in standard ESPResSo, + * except that we add the f/2 contribution, cf. @cite guo02a. + * The fluid velocity is obtained by linear interpolation, + * cf. eq. (11) in @cite ahlrichs99a. + */ void GetIBMInterpolatedVelocity(const Utils::Vector3d &pos, double *v, double *forceAdded) { /* determine elementary lattice cell surrounding the particle - and the relative position of the particle in this cell */ + and the relative position of the particle in this cell */ Utils::Vector node_index{}; Utils::Vector6d delta{}; lblattice.map_position_to_lattice(pos, node_index, delta); - /* calculate fluid velocity at particle's position - this is done by linear interpolation - (Eq. (11) Ahlrichs and Duenweg, JCP 111(17):8225 (1999)) */ Utils::Vector3d interpolated_u = {}; // This for the f/2 contribution to the velocity forceAdded[0] = forceAdded[1] = forceAdded[2] = 0; @@ -218,8 +204,8 @@ void GetIBMInterpolatedVelocity(const Utils::Vector3d &pos, double *v, double local_density; Utils::Vector3d local_j; -// This can be done easier without copying the code twice -// We probably can even set the boundary velocity directly + // This can be done more easily without copying the code twice. + // We probably can even set the boundary velocity directly. #ifdef LB_BOUNDARIES if (lbfields[index].boundary) { local_density = lbpar.density; @@ -233,15 +219,14 @@ void GetIBMInterpolatedVelocity(const Utils::Vector3d &pos, double *v, local_density = lbpar.density + modes[0]; // Add the +f/2 contribution!! - // Guo et al. PRE 2002 local_j[0] = modes[1] + f[0] / 2; local_j[1] = modes[2] + f[1] / 2; local_j[2] = modes[3] + f[2] / 2; - // Keep track of the forces that we added to the fluid + // Keep track of the forces that we added to the fluid. // This is necessary for communication because this part is executed - // for real and ghost particles - // Later on we sum the real and ghost contributions + // for real and ghost particles. + // Later on we sum the real and ghost contributions. const double fExt[3] = { lbpar.ext_force_density[0] * pow(lbpar.agrid, 2) * lbpar.tau * lbpar.tau, @@ -281,13 +266,9 @@ void GetIBMInterpolatedVelocity(const Utils::Vector3d &pos, double *v, v[2] *= lbpar.agrid / lbpar.tau; } -/************ - IsHalo -Builds a cache structure which contains a flag for each LB node whether that -node is a halo node or not -Checks for halo - only for CPU -*************/ - +/** Build a cache structure which contains a flag for each LB node whether that + * node is a halo node or not. + */ bool IsHalo(const int indexCheck) { // First call --> build cache if (isHaloCache == nullptr) { @@ -313,16 +294,13 @@ bool IsHalo(const int indexCheck) { return isHaloCache[indexCheck]; } -/**************** - ParticleVelocitiesFromLB_CPU -Get particle velocities from LB and set the velocity field in the particles data -structure -*****************/ - +/** Get particle velocities from LB and set the velocity field in the particles + * data structure. + */ void ParticleVelocitiesFromLB_CPU() { - // Loop over particles in local cells + // Loop over particles in local cells. // Here all contributions are included: velocity, external force and particle - // force + // force. for (int c = 0; c < local_cells.n; c++) { const Cell *const cell = local_cells.cell[c]; Particle *const p = cell->part; @@ -368,14 +346,14 @@ void ParticleVelocitiesFromLB_CPU() { // Now the local particles contain a velocity (stored in the force field) and // the ghosts contain the rest of the velocity in their respective force - // fields + // fields. // We need to add these. Since we have stored them in the force, not the // velocity fields, we can use the standard force communicator and then - // transfer to the velocity afterwards + // transfer to the velocity afterwards. // Note that this overwrites the actual force which would be a problem for - // real particles + // real particles. // This could be solved by keeping a backup of the local forces before this - // operation is attempted + // operation is attempted. ghost_communicator(&cell_structure.collect_ghost_force_comm, GHOSTTRANS_FORCE); diff --git a/src/core/virtual_sites/lb_inertialess_tracers_cuda.cu b/src/core/virtual_sites/lb_inertialess_tracers_cuda.cu index 8fcb9fb69f1..c2b61abaf72 100644 --- a/src/core/virtual_sites/lb_inertialess_tracers_cuda.cu +++ b/src/core/virtual_sites/lb_inertialess_tracers_cuda.cu @@ -20,8 +20,8 @@ // ******* // This is an internal file of the IMMERSED BOUNDARY implementation -// It should not be included by any main Espresso routines -// Functions to be exported for Espresso are in ibm_main.hpp +// It should not be included by any main ESPResSo routines +// Functions to be exported for ESPResSo are in ibm_main.hpp #include "config.hpp" @@ -68,16 +68,12 @@ extern LB_node_force_density_gpu node_f; extern LB_nodes_gpu *current_nodes; // ** These variables are static in lbgpu_cuda.cu, so we need to duplicate them -// here They are initialized in ForcesIntoFluid The pointers are on the host, -// but point into device memory +// here. They are initialized in ForcesIntoFluid. The pointers are on the host, +// but point into device memory. LB_parameters_gpu *para_gpu = nullptr; float *lb_boundary_velocity_IBM = nullptr; -/**************** - IBM_ResetLBForces_GPU -Calls a kernel to reset the forces on the LB nodes to the external force -*****************/ - +/** Call a kernel to reset the forces on the LB nodes to the external force. */ void IBM_ResetLBForces_GPU() { if (this_node == 0) { // Setup for kernel call @@ -93,13 +89,11 @@ void IBM_ResetLBForces_GPU() { } } -/****************** - IBM_ForcesIntoFluid_GPU -Called from integrate_vv to put the forces into the fluid -This must be the first CUDA-IBM function to be called because it also does some -initialization -*******************/ - +/** Transfer particle forces into the LB fluid. + * Called from @ref integrate_vv. + * This must be the first CUDA-IBM function to be called because it also does + * some initialization. + */ void IBM_ForcesIntoFluid_GPU(ParticleRange particles) { // This function does // (1) Gather forces from all particles via MPI @@ -143,10 +137,6 @@ void IBM_ForcesIntoFluid_GPU(ParticleRange particles) { } } -/*************** - InitCUDA_IBM -***************/ - void InitCUDA_IBM(const int numParticles) { if (this_node == 0) // GPU only on master @@ -208,12 +198,11 @@ void InitCUDA_IBM(const int numParticles) { } } -/************** - Calc_m_from_n_IBM -This is our own version of the calc_m_from_n function in lbgpu_cuda.cu -It does exactly the same, but calculates only the first four modes -***************/ - +/** @copybrief calc_m_from_n + * + * This is a re-implementation of @ref calc_m_from_n. It does exactly the + * same, but calculates only the first four modes. + */ __device__ void Calc_m_from_n_IBM(const LB_nodes_gpu n_a, const unsigned int index, float *mode, const LB_parameters_gpu *const paraP) { @@ -275,12 +264,9 @@ __device__ void Calc_m_from_n_IBM(const LB_nodes_gpu n_a, n_a.vd[18 * para.number_of_nodes + index]); } -/************** - ParticleVelocitiesFromLB_GPU -Calls a kernel function to interpolate the velocity at each IBM particle's -position Store velocity in the particle data structure -**************/ - +/** Call a kernel function to interpolate the velocity at each IBM particle's + * position. Store velocity in the particle data structure. + */ void ParticleVelocitiesFromLB_GPU(ParticleRange particles) { // This function performs three steps: // (1) interpolate velocities on GPU @@ -318,10 +304,6 @@ void ParticleVelocitiesFromLB_GPU(ParticleRange particles) { IBM_cuda_mpi_send_velocities(particles); } -/*************** - ForcesIntoFluid_Kernel -****************/ - __global__ void ForcesIntoFluid_Kernel(const IBM_CUDA_ParticleDataInput *const particle_input, LB_node_force_density_gpu node_f, @@ -407,10 +389,6 @@ ForcesIntoFluid_Kernel(const IBM_CUDA_ParticleDataInput *const particle_input, } } -/************** - ParticleVelocitiesFromLB_Kernel -**************/ - __global__ void ParticleVelocitiesFromLB_Kernel( LB_nodes_gpu n_curr, const IBM_CUDA_ParticleDataInput *const particles_input, @@ -532,10 +510,6 @@ __global__ void ParticleVelocitiesFromLB_Kernel( } } -/**************** - ResetLBForces_Kernel -*****************/ - __global__ void ResetLBForces_Kernel(LB_node_force_density_gpu node_f, const LB_parameters_gpu *const paraP) { diff --git a/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.cpp b/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.cpp index 31b6fb2ec99..b8a4c70c9b1 100644 --- a/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.cpp +++ b/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.cpp @@ -17,10 +17,9 @@ * along with this program. If not, see . */ -// ******* // This is an internal file of the IMMERSED BOUNDARY implementation -// It should not be included by any main Espresso routines -// Functions to be exported for Espresso are in ibm_main.hpp +// It should not be included by any main ESPResSo routines +// Functions to be exported for ESPResSo are in ibm_main.hpp #include "config.hpp" @@ -40,14 +39,6 @@ IBM_CUDA_ParticleDataInput *IBM_ParticleDataInput_host = nullptr; IBM_CUDA_ParticleDataOutput *IBM_ParticleDataOutput_host = nullptr; -/***************** - IBM_cuda_mpi_get_particles -Gather particle positions on the master node in order to communicate them to GPU -We transfer all particles (real and virtual), but actually we would only need -the virtual ones -Room for improvement... - *****************/ - namespace { void pack_particles(ParticleRange particles, IBM_CUDA_ParticleDataInput *buffer) { @@ -72,7 +63,11 @@ void pack_particles(ParticleRange particles, } } // namespace -// Analogous to the usual cuda_mpi_get_particles function +/** Gather particle positions on the master node in order to communicate them + * to GPU. We transfer all particles (real and virtual), but actually we would + * only need the virtual ones. Room for improvement... + * Analogous to @ref cuda_mpi_get_particles. + */ void IBM_cuda_mpi_get_particles(ParticleRange particles) { auto const n_part = particles.size(); @@ -91,12 +86,6 @@ void IBM_cuda_mpi_get_particles(ParticleRange particles) { } } -/***************** - IBM_cuda_mpi_send_velocities -Particle velocities have been communicated from GPU, now transmit to all nodes - ******************/ -// Analogous to cuda_mpi_send_forces - namespace { void set_velocities(ParticleRange particles, IBM_CUDA_ParticleDataOutput *buffer) { @@ -111,6 +100,9 @@ void set_velocities(ParticleRange particles, } } // namespace +/** Particle velocities have been communicated from GPU, now transmit to all + * nodes. Analogous to @ref cuda_mpi_send_forces. + */ void IBM_cuda_mpi_send_velocities(ParticleRange particles) { auto const n_part = particles.size(); diff --git a/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.hpp b/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.hpp index 17025592787..80e02c6709a 100644 --- a/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.hpp +++ b/src/core/virtual_sites/lb_inertialess_tracers_cuda_interface.hpp @@ -19,8 +19,8 @@ // ******* // This is an internal file of the IMMERSED BOUNDARY implementation -// It should not be included by any main Espresso routines -// Functions to be exported for Espresso are in ibm_main.hpp +// It should not be included by any main ESPResSo routines +// Functions to be exported for ESPResSo are in ibm_main.hpp #ifndef IBM_CUDA_INTERFACE_HPP #define IBM_CUDA_INTERFACE_HPP diff --git a/src/python/CMakeLists.txt b/src/python/CMakeLists.txt index 5fd523251b1..8a7d5018751 100644 --- a/src/python/CMakeLists.txt +++ b/src/python/CMakeLists.txt @@ -21,7 +21,7 @@ if( IPYTHON_EXECUTABLE ) endif( ) -option(INSTALL_PYPRESSO "Install pypresso script, not needed when Espresso is installed in /usr" ON) +option(INSTALL_PYPRESSO "Install pypresso script, not needed when ESPResSo is installed in /usr" ON) if(INSTALL_PYPRESSO) install(PROGRAMS ${CMAKE_CURRENT_BINARY_DIR}/pypresso DESTINATION ${CMAKE_INSTALL_BINDIR}) diff --git a/src/python/espressomd/MDA_ESP/__init__.py b/src/python/espressomd/MDA_ESP/__init__.py index e23a70d2155..05959da45d9 100644 --- a/src/python/espressomd/MDA_ESP/__init__.py +++ b/src/python/espressomd/MDA_ESP/__init__.py @@ -15,8 +15,8 @@ # You should have received a copy of the GNU General Public License # along with this program. If not, see . """ -This modules allows to expose ESPREsSo's coordinates and particle attributes -to MDAnalysis without need to save information to files. +This modules exposes ESPResSo's coordinates and particle attributes +to MDAnalysis without the need to save any information to files. The main class is :class:`Stream`, which is used to initialize the stream of data to MDAnalysis' readers. These are the topology reader :class:`ESPParser` @@ -28,14 +28,12 @@ >>> import espressomd >>> from espressomd import MDA_ESP >>> import MDAnalysis as mda - >>> # system setup >>> system = espressomd.System() >>> system.time_step = 1. >>> system.cell_system.skin = 1. >>> system.box_l = [10.,10.,10.] >>> system.part.add(id=0,pos=[1.,2.,3.]) - >>> # set up the stream >>> eos = MDA_ESP.Stream(system) >>> # feed Universe with a topology and with coordinates @@ -121,7 +119,7 @@ def trajectory(self): class ESPParser(TopologyReaderBase): """ - An MDAnalysis reader of espresso's topology + An MDAnalysis reader of ESPResSo's topology. """ format = 'ESP' @@ -131,7 +129,7 @@ def __init__(self, filename, **kwargs): # pylint: disable=unused-argument def parse(self): """ - Access ESPResSo data and return the topology object + Access ESPResSo data and return the topology object. Returns ------- @@ -200,7 +198,7 @@ def dimensions(self, box): class ESPReader(SingleFrameReaderBase): """ - An MDAnalysis single frame reader for the stream provided by Stream() + An MDAnalysis single frame reader for the stream provided by :class:`Stream`. """ format = 'ESP' diff --git a/src/python/espressomd/analyze.pyx b/src/python/espressomd/analyze.pyx index a42f0eb8594..3a6e30b2208 100644 --- a/src/python/espressomd/analyze.pyx +++ b/src/python/espressomd/analyze.pyx @@ -412,7 +412,7 @@ class Analysis: # Dict to store the results p = OrderedDict() - # Update in espresso core if necessary + # Update in ESPResSo core if necessary if (analyze.total_pressure.init_status != 1 + v_comp): analyze.update_pressure(v_comp) @@ -524,7 +524,7 @@ class Analysis: # Dict to store the results p = OrderedDict() - # Update in espresso core if necessary + # Update in ESPResSo core if necessary if (analyze.total_p_tensor.init_status != 1 + v_comp): analyze.update_pressure(v_comp) diff --git a/src/python/espressomd/diamond.pyx b/src/python/espressomd/diamond.pyx index 492e9f05ed8..a0a136e4575 100644 --- a/src/python/espressomd/diamond.pyx +++ b/src/python/espressomd/diamond.pyx @@ -115,6 +115,6 @@ cdef class Diamond: pass elif error_code == -3: raise Exception( - "Failed upon creating one of the monomers in Espresso!\nAborting...\n") + "Failed upon creating one of the monomers in ESPResSo!\nAborting...\n") else: raise Exception("Unknown error code: {}".format(error_code)) diff --git a/src/python/espressomd/gen_code_info.py b/src/python/espressomd/gen_code_info.py index 956d324dc78..be7ee393683 100644 --- a/src/python/espressomd/gen_code_info.py +++ b/src/python/espressomd/gen_code_info.py @@ -47,7 +47,7 @@ include "myconfig.pxi" def features(): - \"\"\"Returns list of features compiled into Espresso core\"\"\" + \"\"\"Returns list of features compiled into ESPResSo core\"\"\" f=[] """) diff --git a/src/python/espressomd/interactions.pyx b/src/python/espressomd/interactions.pyx index fc449528d84..007cd243a98 100644 --- a/src/python/espressomd/interactions.pyx +++ b/src/python/espressomd/interactions.pyx @@ -30,7 +30,7 @@ cdef class NonBondedInteraction: """ Represents an instance of a non-bonded interaction, such as Lennard-Jones. Either called with two particle type id, in which case, the interaction - will represent the bonded interaction as it is defined in Espresso core, + will represent the bonded interaction as it is defined in ESPResSo core, or called with keyword arguments describing a new interaction. """ @@ -50,7 +50,7 @@ cdef class NonBondedInteraction: args[0], int) and is_valid_type(args[1], int): self._part_types = args - # Load the parameters currently set in the Espresso core + # Load the parameters currently set in the ESPResSo core self._params = self._get_params_from_es_core() # Or have we been called with keyword args describing the interaction @@ -72,7 +72,7 @@ cdef class NonBondedInteraction: "The constructor has to be called either with two particle type ids (as integer), or with a set of keyword arguments describing a new interaction") def is_valid(self): - """Check, if the data stored in the instance still matches what is in Espresso. + """Check, if the data stored in the instance still matches what is in ESPResSo. """ temp_params = self._get_params_from_es_core() @@ -123,7 +123,7 @@ cdef class NonBondedInteraction: raise ValueError( "At least the following keys have to be given as keyword arguments: " + self.required_keys().__str__()) - # If this instance refers to an interaction defined in the espresso core, + # If this instance refers to an interaction defined in the ESPResSo core, # load the parameters from there if self._part_types[0] >= 0 and self._part_types[1] >= 0: @@ -1701,7 +1701,7 @@ cdef class BondedInteraction: Base class for bonded interactions. Either called with an interaction id, in which case the interaction - will represent the bonded interaction as it is defined in Espresso core, + will represent the bonded interaction as it is defined in ESPResSo core, or called with keyword arguments describing a new interaction. """ @@ -1713,14 +1713,14 @@ cdef class BondedInteraction: # Interaction id as argument if len(args) == 1 and is_valid_type(args[0], int): bond_id = args[0] - # Check if the bond type in Espresso core matches this class + # Check if the bond type in ESPResSo core matches this class if bonded_ia_params[bond_id].type != self.type_number(): raise Exception( - "The bond with this id is not defined as a " + self.type_name() + " bond in the Espresso core.") + "The bond with this id is not defined as a " + self.type_name() + " bond in the ESPResSo core.") self._bond_id = bond_id - # Load the parameters currently set in the Espresso core + # Load the parameters currently set in the ESPResSo core self._params = self._get_params_from_es_core() self._bond_id = bond_id @@ -1745,16 +1745,16 @@ cdef class BondedInteraction: return (self.__class__, (self._bond_id,)) def is_valid(self): - """Check, if the data stored in the instance still matches what is in Espresso. + """Check, if the data stored in the instance still matches what is in ESPResSo. """ - # Check if the bond type in Espresso still matches the bond type saved + # Check if the bond type in ESPResSo still matches the bond type saved # in this class if bonded_ia_params[self._bond_id].type != self.type_number(): return False # check, if the bond parameters saved in the class still match those - # saved in Espresso + # saved in ESPResSo temp_params = self._get_params_from_es_core() if self._params != temp_params: return False @@ -1866,7 +1866,7 @@ class BondedInteractionNotDefined: def __init__(self, *args, **kwargs): raise Exception( - self.__class__.__name__ + " not compiled into Espresso core") + self.__class__.__name__ + " not compiled into ESPResSo core") def type_number(self): raise Exception(("%s has to be defined in myconfig.hpp.") % self.name) @@ -3434,13 +3434,13 @@ class BondedInteractions: raise ValueError( "Index to BondedInteractions[] has to be an integer referring to a bond id") - # Find out the type of the interaction from Espresso + # Find out the type of the interaction from ESPResSo if key >= bonded_ia_params.size(): raise IndexError( "Index to BondedInteractions[] out of range") bond_type = bonded_ia_params[key].type - # Check if the bonded interaction exists in Espresso core + # Check if the bonded interaction exists in ESPResSo core if bond_type == -1: raise ValueError( "The bonded interaction with the id " + str(key) + " is not yet defined.") @@ -3449,7 +3449,7 @@ class BondedInteractions: bond_class = bonded_interaction_classes[bond_type] # And return an instance of it, which refers to the bonded interaction - # id in Espresso + # id in ESPResSo return bond_class(key) def __setitem__(self, key, value): diff --git a/src/python/espressomd/particle_data.pxd b/src/python/espressomd/particle_data.pxd index 48f052b54e6..978308046a8 100644 --- a/src/python/espressomd/particle_data.pxd +++ b/src/python/espressomd/particle_data.pxd @@ -36,7 +36,7 @@ cdef extern from "particle_data.hpp": stdint.uint8_t ROTATION_Z # Note: Conditional compilation is not possible within ctypedef blocks. - # Therefore, only member variables are imported here, which are always compiled into Espresso. + # Therefore, only member variables are imported here, which are always compiled into ESPResSo. # For all other properties, getter-functions have to be used on the c # level. ctypedef struct particle_properties "ParticleProperties": diff --git a/src/python/espressomd/particle_data.pyx b/src/python/espressomd/particle_data.pyx index 7fdfc3db255..c2e9a535138 100644 --- a/src/python/espressomd/particle_data.pyx +++ b/src/python/espressomd/particle_data.pyx @@ -441,7 +441,7 @@ cdef class ParticleHandle: def __set__(self, _q): raise AttributeError( - "Setting the director is not implemented in the C++-core of Espresso.") + "Setting the director is not implemented in the C++-core of ESPResSo.") def __get__(self): self.update_particle_data() @@ -563,7 +563,7 @@ cdef class ParticleHandle: return np.array( [out_direction[0], out_direction[1], out_direction[2]]) -# Charge + # Charge property q: """ Particle charge. @@ -1293,7 +1293,7 @@ cdef class ParticleHandle: See Also -------- add_bond : Delete an unverified bond held by the ``Particle``. - bonds : ``Particle`` property containing a list of all current bonds help by ``Particle``. + bonds : ``Particle`` property containing a list of all current bonds held by ``Particle``. """ @@ -1324,7 +1324,7 @@ cdef class ParticleHandle: See Also -------- delete_bond : Delete an unverified bond held by the ``Particle``. - bonds : ``Particle`` property containing a list of all current bonds help by ``Particle``. + bonds : ``Particle`` property containing a list of all current bonds held by ``Particle``. """ @@ -1366,7 +1366,7 @@ cdef class ParticleHandle: # interactions if bond[0]._bond_id == -1: raise Exception( - "The bonded interaction has not yet been added to the list of active bonds in Espresso.") + "The bonded interaction has not yet been added to the list of active bonds in ESPResSo.") # Validity of the numeric id if bond[0]._bond_id >= bonded_ia_params.size(): @@ -1403,7 +1403,7 @@ cdef class ParticleHandle: See Also -------- - bonds : ``Particle`` property containing a list of all current bonds help by ``Particle``. + bonds : ``Particle`` property containing a list of all current bonds held by ``Particle``. Examples -------- @@ -1492,7 +1492,7 @@ cdef class ParticleHandle: See Also ---------- delete_bond : Delete an unverified bond held by the particle. - bonds : Particle property containing a list of all current bonds help by particle. + bonds : Particle property containing a list of all current bonds held by particle. """ diff --git a/src/python/espressomd/reaction_ensemble.pyx b/src/python/espressomd/reaction_ensemble.pyx index fb78dc3f060..238460757e2 100644 --- a/src/python/espressomd/reaction_ensemble.pyx +++ b/src/python/espressomd/reaction_ensemble.pyx @@ -413,7 +413,7 @@ cdef class ConstantpHEnsemble(ReactionAlgorithm): "The constant pH method is only implemented for reactions with two product types and one adduct type.") if(kwargs["reactant_coefficients"][0] != 1 or kwargs["product_coefficients"][0] != 1 or kwargs["product_coefficients"][1] != 1): raise ValueError( - "All product and reactant coefficients must equal one in the constant pH method as implemented in Espresso.") + "All product and reactant coefficients must equal one in the constant pH method as implemented in ESPResSo.") super().add_reaction(*args, **kwargs) property constant_pH: @@ -665,7 +665,7 @@ cdef class WangLandauReactionEnsemble(ReactionAlgorithm): polymer configuration changing moves need to be implemented in the case of using Wang-Landau with energy reweighting and a polymer in the system. Polymer configuration changing moves had been implemented - before but were removed from espresso. + before but were removed from ESPResSo. """ use_wang_landau = True diff --git a/src/python/espressomd/system.pyx b/src/python/espressomd/system.pyx index 702b6d60aa1..cbcd31dd35c 100644 --- a/src/python/espressomd/system.pyx +++ b/src/python/espressomd/system.pyx @@ -75,7 +75,7 @@ if OIF_GLOBAL_FORCES: cdef bool _system_created = False cdef class System: - """The espresso system class. + """The ESPResSo system class. .. note:: every attribute has to be declared at the class level. This means that methods cannot define an attribute by using @@ -493,13 +493,13 @@ cdef class System: def setup_type_map(self, type_list=None): """ - For using Espresso conveniently for simulations in the grand canonical + For using ESPResSo conveniently for simulations in the grand canonical ensemble, or other purposes, when particles of certain types are created and deleted frequently. Particle ids can be stored in lists for each individual type and so random ids of particles of a certain type can be - drawn. If you want Espresso to keep track of particle ids of a certain type + drawn. If you want ESPResSo to keep track of particle ids of a certain type you have to initialize the method by calling the setup function. After that - Espresso will keep track of particle ids of that type. + ESPResSo will keep track of particle ids of that type. """ if not hasattr(type_list, "__iter__"): diff --git a/src/python/espressomd/version.pyx b/src/python/espressomd/version.pyx index 77d2af070b7..5b2fb368caf 100644 --- a/src/python/espressomd/version.pyx +++ b/src/python/espressomd/version.pyx @@ -19,13 +19,13 @@ from espressomd.utils import to_str def major(): - """Prints the major version of Espresso. + """Prints the major version of ESPResSo. """ return ESPRESSO_VERSION_MAJOR def minor(): - """Prints the minor version of Espresso. + """Prints the minor version of ESPResSo. """ return ESPRESSO_VERSION_MINOR @@ -53,7 +53,7 @@ def git_commit(): def git_state(): """Git state of the build if known, otherwise empty. State is "CLEAN" if the repository - was not changed from git_commit(), "DIRTY" - otherwise. + was not changed from :meth:`git_commit()`, + "DIRTY" otherwise. """ return to_str(GIT_STATE) diff --git a/src/python/espressomd/visualization_mayavi.pyx b/src/python/espressomd/visualization_mayavi.pyx index e66e4c53047..bd3db9b5907 100644 --- a/src/python/espressomd/visualization_mayavi.pyx +++ b/src/python/espressomd/visualization_mayavi.pyx @@ -180,7 +180,7 @@ cdef class mayaviLive: scalars=bonds[:, 6]) def update(self): - """Pull the latest particle information from Espresso. + """Pull the latest particle information from ESPResSo. This is the only function that should be called from the computation thread. It does not call any Mayavi functions unless it is being @@ -292,7 +292,7 @@ cdef class mayaviLive: """Start the GUI event loop. This function blocks until the Mayavi window is closed. So you should - only use it if your Espresso simulation's integrate loop is running in + only use it if your ESPResSo simulation's integrate loop is running in a secondary thread. """ diff --git a/src/python/espressomd/visualization_opengl.pyx b/src/python/espressomd/visualization_opengl.pyx index 41c5933cfcd..d52c6a6c993 100644 --- a/src/python/espressomd/visualization_opengl.pyx +++ b/src/python/espressomd/visualization_opengl.pyx @@ -226,7 +226,7 @@ class openGLLive: # DEFAULT PROPERTIES self.specs = { 'window_size': [800, 800], - 'name': 'Espresso Visualization', + 'name': 'ESPResSo Visualization', 'background_color': [0, 0, 0], @@ -613,7 +613,7 @@ class openGLLive: def update(self): """Update method to be called after integration. - Changes of espresso system can only happen here. + Changes of ESPResSo system can only happen here. """ if self.glutMainLoop_started: diff --git a/src/script_interface/ScriptInterface.dox b/src/script_interface/ScriptInterface.dox index a63d67dd1d0..72b940cd420 100644 --- a/src/script_interface/ScriptInterface.dox +++ b/src/script_interface/ScriptInterface.dox @@ -76,7 +76,7 @@ /// @subsection script_interface_example Example class /// /// As a first example we implement a hello world script object, that does not -/// interact with the Espresso core. This class has one parameter, a single +/// interact with the ESPResSo core. This class has one parameter, a single /// string @c m_name, and has one callable method: @c greet(). /// /// @code{.cpp} diff --git a/src/script_interface/ScriptInterfaceBase.cpp b/src/script_interface/ScriptInterfaceBase.cpp index 10f54cd9907..59f685b4547 100644 --- a/src/script_interface/ScriptInterfaceBase.cpp +++ b/src/script_interface/ScriptInterfaceBase.cpp @@ -57,8 +57,7 @@ ScriptInterfaceBase::make_shared(std::string const &name, const auto id = sp->id(); /* Now get a reference to the corresponding weak_ptr in ObjectId and update - it with our shared ptr, so that everybody uses the same ref count. - */ + * it with our shared ptr, so that everybody uses the same ref count. */ sp->get_instance(id) = sp; return sp; @@ -75,7 +74,7 @@ ScriptInterfaceBase::get_instance(ObjectId id) { * @brief Return a Variant representation of the state of the object. * * This should return the internal state of the instance, so that - * the instance can be restored from this information. The default + * the instance can be restored from this information. The default * implementation stores all the public parameters, including object * parameters that are captured by calling get_state on them. */ diff --git a/src/utils/include/utils/math/AS_erfc_part.hpp b/src/utils/include/utils/math/AS_erfc_part.hpp index 2ab452d63f4..bb8ff55a5e5 100644 --- a/src/utils/include/utils/math/AS_erfc_part.hpp +++ b/src/utils/include/utils/math/AS_erfc_part.hpp @@ -23,9 +23,9 @@ #define UTILS_MATH_AC_ERF_PART_HPP namespace Utils { -/** approximates \f$ \exp(d^2) \mathrm{erfc}(d)\f$ by applying a formula from: - Abramowitz/Stegun: Handbook of Mathematical Functions, Dover - (9. ed.), chapter 7 */ +/** Approximate \f$ \exp(d^2) \mathrm{erfc}(d)\f$ by applying a formula from + * @cite abramowitz65a chapter 7. + */ template constexpr T AS_erfc_part(T d) { T const constexpr a1 = 0.254829592; T const constexpr a2 = -0.284496736; diff --git a/src/utils/include/utils/math/quaternion.hpp b/src/utils/include/utils/math/quaternion.hpp index c3c9e21a0bb..b9ea1884d06 100644 --- a/src/utils/include/utils/math/quaternion.hpp +++ b/src/utils/include/utils/math/quaternion.hpp @@ -1,23 +1,23 @@ /* - Copyright (C) 2010-2018 The ESPResSo project - Copyright (C) 2002,2003,2004,2005,2006,2007,2008,2009,2010 - Max-Planck-Institute for Polymer Research, Theory Group - - This file is part of ESPResSo. - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2010-2019 The ESPResSo project + * Copyright (C) 2002,2003,2004,2005,2006,2007,2008,2009,2010 + * Max-Planck-Institute for Polymer Research, Theory Group + * + * This file is part of ESPResSo. + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #ifndef UTILS_MATH_QUATERNION_HPP #define UTILS_MATH_QUATERNION_HPP diff --git a/src/utils/include/utils/math/rotation_matrix.hpp b/src/utils/include/utils/math/rotation_matrix.hpp index b35d0d4c1b0..8eecc0b650d 100644 --- a/src/utils/include/utils/math/rotation_matrix.hpp +++ b/src/utils/include/utils/math/rotation_matrix.hpp @@ -24,13 +24,12 @@ namespace Utils { /** - * @brief Rotation matrix for quaternion - Taken from "Goldstein - Classical Mechanics" (Chapter 4.6 Eq. 4.47). - - @param q Quaternion - - @return Corresponding rotation matrix -*/ + * @brief Rotation matrix for quaternion. + * Taken from @cite goldstein80a (Chapter 4.6 eq. (4.47)). + * + * @param q Quaternion + * @return Corresponding rotation matrix + */ template Matrix rotation_matrix(Vector const &q) { auto const q0q0 = q[0] * q[0]; auto const q1q1 = q[1] * q[1]; diff --git a/src/utils/include/utils/math/sinc.hpp b/src/utils/include/utils/math/sinc.hpp index 5c1df95ef92..f82f3f553f1 100644 --- a/src/utils/include/utils/math/sinc.hpp +++ b/src/utils/include/utils/math/sinc.hpp @@ -32,13 +32,13 @@ namespace Utils { /** * @brief Calculates the sinc-function as sin(PI*x)/(PI*x). * - * (same convention as in Hockney/Eastwood). In order to avoid + * (same convention as in @cite hockney88a). In order to avoid * divisions by 0, arguments, whose modulus is smaller than epsi, will * be evaluated by an 8th order Taylor expansion of the sinc * function. Note that the difference between sinc(x) and this * expansion is smaller than 0.235e-12, if x is smaller than 0.1. (The * next term in the expansion is the 10th order contribution - * PI^10/39916800 * x^10 = 0.2346...*x^12). This expansion should + * PI^10/39916800 * x^10 = 0.2346...*x^12). This expansion should * also save time, since it reduces the number of function calls to * sin(). */ diff --git a/src/utils/include/utils/math/triangle_functions.hpp b/src/utils/include/utils/math/triangle_functions.hpp index 489de206c6d..30979fbf25a 100644 --- a/src/utils/include/utils/math/triangle_functions.hpp +++ b/src/utils/include/utils/math/triangle_functions.hpp @@ -33,28 +33,30 @@ inline Vector3d get_n_triangle(const Vector3d &P1, const Vector3d &P2, return vector_product(u, v); } -/** Computes the area of triangle between vectors P1,P2,P3, - * by computing the crossproduct P1P2 x P1P3 and taking the half of its norm */ +/** Computes the area of triangle between vectors P1,P2,P3, by computing + * the cross product P1P2 x P1P3 and taking the half of its norm. + */ inline double area_triangle(const Vector3d &P1, const Vector3d &P2, const Vector3d &P3) { return 0.5 * get_n_triangle(P1, P2, P3).norm(); } -/** This function returns the angle btw the triangle p1,p2,p3 and p2,p3,p4. Be - * careful, the angle depends on the orientation of the triangles! You need to - * be sure that the orientation (direction of normal vector) of p1p2p3 is given - * by the cross product p2p1 x p2p3. The orientation of p2p3p4 must be given - * by p2p3 x p2p4. +/** This function returns the angle between the triangle p1,p2,p3 and p2,p3,p4. + * Be careful, the angle depends on the orientation of the triangles! You need + * to be sure that the orientation (direction of normal vector) of p1p2p3 + * is given by the cross product p2p1 x p2p3. The orientation of p2p3p4 must + * be given by p2p3 x p2p4. * - * Example: p1 = (0,0,1), p2 = (0,0,0), p3=(1,0,0), p4=(0,1,0). The + * Example: p1 = (0,0,1), p2 = (0,0,0), p3=(1,0,0), p4=(0,1,0). The * orientation of p1p2p3 should be in the direction (0,1,0) and indeed: p2p1 x * p2p3 = (0,0,1)x(1,0,0) = (0,1,0) This function is called in the beginning - * of the simulation when creating bonds depending on the angle btw the - * triangles, the bending_force. Here, we determine the orientations by - * looping over the triangles and checking the correct orientation. So if you + * of the simulation when creating bonds depending on the angle between the + * triangles, the bending_force. Here, we determine the orientations by + * looping over the triangles and checking the correct orientation. So if you * have the access to the order of particles, you are safe to call this * function with exactly this order. Otherwise you need to check the - * orientations. */ + * orientations. + */ template double angle_btw_triangles(const T1 &P1, const T2 &P2, const T3 &P3, const T4 &P4) { @@ -85,8 +87,8 @@ double angle_btw_triangles(const T1 &P1, const T2 &P2, const T3 &P3, // the orientation, always less or equal to Pi) // is equal to Pi minus angle between the normals - // Now we need to determine, if the angle btw two triangles is less than Pi or - // more than Pi. To do this we check, + // Now we need to determine, if the angle between two triangles is less than + // Pi or more than Pi. To do this we check, // if the point P4 lies in the halfspace given by triangle P1P2P3 and the // normal to this triangle. If yes, we have // angle less than Pi, if not, we have angle more than Pi. diff --git a/src/utils/include/utils/uniform.hpp b/src/utils/include/utils/uniform.hpp index bd48ec94d96..9a508673b52 100644 --- a/src/utils/include/utils/uniform.hpp +++ b/src/utils/include/utils/uniform.hpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2018-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2018-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #ifndef UTILS_UNIFORM_HPP #define UTILS_UNIFORM_HPP diff --git a/src/utils/tests/AS_erfc_part_test.cpp b/src/utils/tests/AS_erfc_part_test.cpp index 0a61fe869be..4fa53581177 100644 --- a/src/utils/tests/AS_erfc_part_test.cpp +++ b/src/utils/tests/AS_erfc_part_test.cpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2017-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2017-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #define BOOST_TEST_MODULE Utils::AS_erfc_part test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/Factory_test.cpp b/src/utils/tests/Factory_test.cpp index e0f93d4c297..ae627a35259 100644 --- a/src/utils/tests/Factory_test.cpp +++ b/src/utils/tests/Factory_test.cpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Utils::Factory class. +/* Unit tests for the Utils::Factory class. * The factory is tested by registering different types of classes * with it (first test), and then checking if instances of those classes can be * made via the Factory (second test). @@ -45,7 +44,7 @@ struct OtherDerivedTestClass : public TestClass {}; } /* namespace Testing */ -/** Check registration of construction functions */ +/* Check registration of construction functions */ BOOST_AUTO_TEST_CASE(regiser_class) { using namespace Testing; typedef Utils::Factory Factory; @@ -66,7 +65,7 @@ BOOST_AUTO_TEST_CASE(regiser_class) { BOOST_CHECK(Factory::has_builder("other_derived_class")); } -/** Check object construction. */ +/* Check object construction. */ BOOST_AUTO_TEST_CASE(make) { using namespace Testing; /* Make a base object */ diff --git a/src/utils/tests/List_test.cpp b/src/utils/tests/List_test.cpp index 04d18b159b0..d3c17292ced 100644 --- a/src/utils/tests/List_test.cpp +++ b/src/utils/tests/List_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Utils::List class. - * - */ +/* Unit tests for the Utils::List class. */ #define BOOST_TEST_MODULE List test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/NumeratedContainer_test.cpp b/src/utils/tests/NumeratedContainer_test.cpp index d5bc0beaf5c..5a3bbecfbdc 100644 --- a/src/utils/tests/NumeratedContainer_test.cpp +++ b/src/utils/tests/NumeratedContainer_test.cpp @@ -19,10 +19,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Utils::NumeratedContainer class. - * - */ +/* Unit tests for the Utils::NumeratedContainer class. */ #define BOOST_TEST_MODULE NumeratedContainerTest #define BOOST_TEST_DYN_LINK @@ -61,14 +58,14 @@ BOOST_AUTO_TEST_CASE(index_handling) { int a = 18, b = 7, c = 3, d = 35325, e = 588; - /** Check adding elements */ + /* Check adding elements */ BOOST_CHECK(container.add(a) == 0); BOOST_CHECK(container.add(b) == 1); BOOST_CHECK(container.add(c) == 2); BOOST_CHECK(container.add(d) == 3); BOOST_CHECK(container.add(e) == 4); - /** Check removing and id reuse */ + /* Check removing and id reuse */ container.remove(1); BOOST_CHECK(container.add(b) == 1); @@ -80,7 +77,7 @@ BOOST_AUTO_TEST_CASE(index_handling) { BOOST_CHECK(container.add(c) == 2); BOOST_CHECK(container.add(d) == 3); - /** Check out-of-order remove */ + /* Check out-of-order remove */ container.remove(3); container.remove(4); container.remove(2); @@ -93,10 +90,10 @@ BOOST_AUTO_TEST_CASE(index_handling) { BOOST_CHECK(container.add(d) == 3); BOOST_CHECK(container.add(e) == 4); - /** Check values */ + /* Check values */ BOOST_CHECK(container[0] == a); BOOST_CHECK(container[4] == e); - /** Check that out-of-bounds check works */ + /* Check that out-of-bounds check works */ BOOST_CHECK_THROW(container[5], std::out_of_range); } diff --git a/src/utils/tests/RunningAverage_test.cpp b/src/utils/tests/RunningAverage_test.cpp index 70cd15a9ee4..b8694905e9e 100644 --- a/src/utils/tests/RunningAverage_test.cpp +++ b/src/utils/tests/RunningAverage_test.cpp @@ -30,7 +30,7 @@ #include "utils/statistics/RunningAverage.hpp" -/** random sequence */ +/* random sequence */ #include "random_sequence.hpp" using namespace Testing; @@ -67,14 +67,14 @@ BOOST_AUTO_TEST_CASE(simple_variance_check) { avg.add_sample(0.0); avg.add_sample(5.0); - /** Var should be \f$ - ^2 = (0^2 + 5.0^2) / 2. - 2.5^2 = 12.5 + /* Var should be \f$ - ^2 = (0^2 + 5.0^2) / 2. - 2.5^2 = 12.5 * - 6.25 = 6.25\f$ */ BOOST_CHECK(std::fabs(avg.avg() - 2.5) <= std::numeric_limits::epsilon()); BOOST_CHECK(std::fabs(avg.var() - 6.25) <= std::numeric_limits::epsilon()); - /** Standard deviation should be sqrt(var()) */ + /* Standard deviation should be sqrt(var()) */ BOOST_CHECK(std::fabs(avg.sig() - std::sqrt(avg.var())) <= std::numeric_limits::epsilon()); } @@ -89,14 +89,14 @@ BOOST_AUTO_TEST_CASE(mean_and_variance) { BOOST_CHECK(running_average.n() == sample_size); - /** Directly calculate the mean from the data */ + /* Directly calculate the mean from the data */ const double m_mean = std::accumulate(std::begin(RandomSequence::values), std::end(RandomSequence::values), 0.0) / sample_size; BOOST_CHECK(std::fabs(running_average.avg() - m_mean) <= 1e-12); - /** Directly calculate the variance from the data */ + /* Directly calculate the variance from the data */ double m_var = 0.0; for (auto const &val : RandomSequence::values) { m_var += (val - m_mean) * (val - m_mean); diff --git a/src/utils/tests/Vector_test.cpp b/src/utils/tests/Vector_test.cpp index f027a6ae916..279ed542fd6 100644 --- a/src/utils/tests/Vector_test.cpp +++ b/src/utils/tests/Vector_test.cpp @@ -19,10 +19,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Utils::Vector class. - * - */ +/* Unit tests for the Utils::Vector class. */ #define BOOST_TEST_MODULE Vector test #define BOOST_TEST_DYN_LINK @@ -39,7 +36,7 @@ #include "utils/Vector.hpp" using Utils::Vector; -/** Number of nontrivial Baxter permutations of length 2n-1. (A001185) */ +/* Number of nontrivial Baxter permutations of length 2n-1. (A001185) */ #define TEST_NUMBERS \ { 0, 1, 1, 7, 21, 112, 456, 2603, 13203 } #define TEST_NUMBERS_PARTIAL_NORM2 \ diff --git a/src/utils/tests/as_const_test.cpp b/src/utils/tests/as_const_test.cpp index 15aa07065fe..34a7d4c9852 100644 --- a/src/utils/tests/as_const_test.cpp +++ b/src/utils/tests/as_const_test.cpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2018-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2018-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #define BOOST_TEST_MODULE Utils::as_const test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/coordinate_transformation.cpp b/src/utils/tests/coordinate_transformation.cpp index 25201d57289..bf4345d5b42 100644 --- a/src/utils/tests/coordinate_transformation.cpp +++ b/src/utils/tests/coordinate_transformation.cpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2017-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2017-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #define BOOST_TEST_MODULE coordinate_transformation test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/for_each_pair_test.cpp b/src/utils/tests/for_each_pair_test.cpp index c0f4596c4bd..977ef18f865 100644 --- a/src/utils/tests/for_each_pair_test.cpp +++ b/src/utils/tests/for_each_pair_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Utils::for_each_pair. - * - */ +/* Unit tests for the Utils::for_each_pair. */ #define BOOST_TEST_MODULE for_each_pair test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/index_test.cpp b/src/utils/tests/index_test.cpp index 2a0eb70620c..a5e3ae590a5 100644 --- a/src/utils/tests/index_test.cpp +++ b/src/utils/tests/index_test.cpp @@ -1,24 +1,22 @@ /* - Copyright (C) 2017-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ - -/** \file int_pow_test.cpp Unit tests for the - * Utils::ravel_index and Utils::unravel_index functions. + * Copyright (C) 2017-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . */ +/* Unit tests for the Utils::ravel_index and Utils::unravel_index functions. */ + #define BOOST_TEST_MODULE Utils::ravel_index test #define BOOST_TEST_DYN_LINK @@ -90,4 +88,4 @@ BOOST_AUTO_TEST_CASE(upper_triangular_test) { BOOST_CHECK(expected == result); } } -} \ No newline at end of file +} diff --git a/src/utils/tests/int_pow_test.cpp b/src/utils/tests/int_pow_test.cpp index c5a5306ca58..072aabfff3b 100644 --- a/src/utils/tests/int_pow_test.cpp +++ b/src/utils/tests/int_pow_test.cpp @@ -1,24 +1,22 @@ /* - Copyright (C) 2017-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ - -/** \file int_pow_test.cpp Unit tests for the - * Utils::int_pow function. + * Copyright (C) 2017-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . */ +/* Unit tests for the Utils::int_pow function. */ + #define BOOST_TEST_MODULE int_pow test #define BOOST_TEST_DYN_LINK #include diff --git a/src/utils/tests/permute_ifield_test.cpp b/src/utils/tests/permute_ifield_test.cpp index db85678047a..123dd5ff68f 100644 --- a/src/utils/tests/permute_ifield_test.cpp +++ b/src/utils/tests/permute_ifield_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file NumeratedContainer_test.cpp Unit tests for the - * Utils::NumeratedContainer class. - * - */ +/* Unit tests for the Utils::NumeratedContainer class. */ #define BOOST_TEST_MODULE Utils::permute test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/sgn_test.cpp b/src/utils/tests/sgn_test.cpp index 96185dc1a4f..a2c80a2bb3d 100644 --- a/src/utils/tests/sgn_test.cpp +++ b/src/utils/tests/sgn_test.cpp @@ -1,24 +1,22 @@ /* - Copyright (C) 2017-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ - -/** \file int_pow_test.cpp Unit tests for the - * Utils::int_pow function. + * Copyright (C) 2017-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . */ +/* Utils::int_pow function. */ + #define BOOST_TEST_MODULE Utils::sgn test #define BOOST_TEST_DYN_LINK #include diff --git a/src/utils/tests/sinc_test.cpp b/src/utils/tests/sinc_test.cpp index 8de12a07828..6fb4c199c68 100644 --- a/src/utils/tests/sinc_test.cpp +++ b/src/utils/tests/sinc_test.cpp @@ -1,24 +1,22 @@ /* - Copyright (C) 2017-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ - -/** \file int_pow_test.cpp Unit tests for the - * Utils::int_pow function. + * Copyright (C) 2017-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . */ +/* Unit tests for the Utils::int_pow function. */ + #define BOOST_TEST_MODULE Utils::sgn test #define BOOST_TEST_DYN_LINK #include diff --git a/src/utils/tests/tensor_product_test.cpp b/src/utils/tests/tensor_product_test.cpp index 25e8dad195f..a86007942d8 100644 --- a/src/utils/tests/tensor_product_test.cpp +++ b/src/utils/tests/tensor_product_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit tests for the Utils::tensor_product class. - * - */ +/* Unit tests for the Utils::tensor_product class. */ #define BOOST_TEST_MODULE Utils::tensor_product #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/tuple_test.cpp b/src/utils/tests/tuple_test.cpp index 3954ebffd9c..dfcb59076f2 100644 --- a/src/utils/tests/tuple_test.cpp +++ b/src/utils/tests/tuple_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit test for Utils tuple algorithms. - * - */ +/* Unit test for Utils tuple algorithms. */ #define BOOST_TEST_MODULE Utils::tuple_test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/type_traits_test.cpp b/src/utils/tests/type_traits_test.cpp index 30a8c5e37cf..fde913f5378 100644 --- a/src/utils/tests/type_traits_test.cpp +++ b/src/utils/tests/type_traits_test.cpp @@ -17,10 +17,7 @@ * along with this program. If not, see . */ -/** \file - * Unit test for Utils Type Traits. - * - */ +/* Unit test for Utils Type Traits. */ #define BOOST_TEST_MODULE type traits tests #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/uniform_test.cpp b/src/utils/tests/uniform_test.cpp index cea9f467fc9..829722b3347 100644 --- a/src/utils/tests/uniform_test.cpp +++ b/src/utils/tests/uniform_test.cpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2018-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2018-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #define BOOST_TEST_MODULE Utils::uniform test #define BOOST_TEST_DYN_LINK diff --git a/src/utils/tests/vec_rotate_test.cpp b/src/utils/tests/vec_rotate_test.cpp index 6c221191f92..db7ea65c807 100644 --- a/src/utils/tests/vec_rotate_test.cpp +++ b/src/utils/tests/vec_rotate_test.cpp @@ -1,19 +1,19 @@ /* - Copyright (C) 2018-2019 The ESPResSo project - - ESPResSo is free software: you can redistribute it and/or modify - it under the terms of the GNU General Public License as published by - the Free Software Foundation, either version 3 of the License, or - (at your option) any later version. - - ESPResSo is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - GNU General Public License for more details. - - You should have received a copy of the GNU General Public License - along with this program. If not, see . -*/ + * Copyright (C) 2018-2019 The ESPResSo project + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ #define BOOST_TEST_MODULE Utils::vec_rotate test #define BOOST_TEST_DYN_LINK diff --git a/testsuite/python/lb_buoyancy_force.py b/testsuite/python/lb_buoyancy_force.py index d65c4d768a2..91a9d55c9b0 100644 --- a/testsuite/python/lb_buoyancy_force.py +++ b/testsuite/python/lb_buoyancy_force.py @@ -1,4 +1,4 @@ -# Copyright (C) 2010-2018 The ESPResSo project +# Copyright (C) 2010-2019 The ESPResSo project # # This file is part of ESPResSo. # diff --git a/testsuite/python/lb_momentum_conservation.py b/testsuite/python/lb_momentum_conservation.py index 906c0f99ae1..87fef085918 100644 --- a/testsuite/python/lb_momentum_conservation.py +++ b/testsuite/python/lb_momentum_conservation.py @@ -1,4 +1,4 @@ -# Copyright (C) 2010-2018 The ESPResSo project +# Copyright (C) 2010-2019 The ESPResSo project # # This file is part of ESPResSo. #