src/gromacs/mdlib/calcmu.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2004, The GROMACS development team.
   6  * Copyright (c) 2012,2014,2015,2017,2018 by the GROMACS development team.
   7  * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
   8  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   9  * and including many others, as listed in the AUTHORS file in the
  10  * top-level source directory and at http://www.gromacs.org.
  11  *
  12  * GROMACS is free software; you can redistribute it and/or
  13  * modify it under the terms of the GNU Lesser General Public License
  14  * as published by the Free Software Foundation; either version 2.1
  15  * of the License, or (at your option) any later version.
  16  *
  17  * GROMACS is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  20  * Lesser General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU Lesser General Public
  23  * License along with GROMACS; if not, see
  24  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  25  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  26  *
  27  * If you want to redistribute modifications to GROMACS, please
  28  * consider that scientific software is very special. Version
  29  * control is crucial - bugs must be traceable. We will be happy to
  30  * consider code for inclusion in the official distribution, but
  31  * derived work must not be called official GROMACS. Details are found
  32  * in the README & COPYING files - if they are missing, get the
  33  * official version at http://www.gromacs.org.
  34  *
  35  * To help us fund GROMACS development, we humbly ask that you cite
  36  * the research papers on the package. Check out http://www.gromacs.org.
  37  */
  38 /* This file is completely threadsafe - keep it that way! */
  39 #include "gmxpre.h"
  40
  41 #include "calcmu.h"
  42
  43 #include <cstdio>
  44 #include <cstdlib>
  45
  46 #include "gromacs/gmxlib/network.h"
  47 #include "gromacs/math/units.h"
  48 #include "gromacs/math/vec.h"
  49 #include "gromacs/mdlib/gmx_omp_nthreads.h"
  50 #include "gromacs/utility/arrayref.h"
  51
  52 void calc_mu(int                            start,
  53              int                            homenr,
  54              gmx::ArrayRef<const gmx::RVec> x,
  55              gmx::ArrayRef<const real>      q,
  56              gmx::ArrayRef<const real>      qB,
  57              bool                           havePerturbedCharges,
  58              dvec                           mu,
  59              dvec                           mu_B)
  60 {
  61     int    end, m;
  62     double mu_x, mu_y, mu_z;
  63
  64     end = start + homenr;
  65
  66     mu_x = mu_y = mu_z = 0.0;
  67 #pragma omp parallel for reduction(+: mu_x, mu_y, mu_z) schedule(static) \
  68     num_threads(gmx_omp_nthreads_get(ModuleMultiThread::Default))
  69     for (int i = start; i < end; i++)
  70     {
  71         // Trivial OpenMP region that cannot throw
  72         mu_x += q[i] * x[i][XX];
  73         mu_y += q[i] * x[i][YY];
  74         mu_z += q[i] * x[i][ZZ];
  75     }
  76     mu[XX] = mu_x;
  77     mu[YY] = mu_y;
  78     mu[ZZ] = mu_z;
  79
  80     for (m = 0; (m < DIM); m++)
  81     {
  82         mu[m] *= gmx::c_enm2Debye;
  83     }
  84
  85     if (havePerturbedCharges)
  86     {
  87         mu_x = mu_y = mu_z = 0.0;
  88 #pragma omp parallel for reduction(+: mu_x, mu_y, mu_z) schedule(static) \
  89         num_threads(gmx_omp_nthreads_get(ModuleMultiThread::Default))
  90         for (int i = start; i < end; i++)
  91         {
  92             // Trivial OpenMP region that cannot throw
  93             mu_x += qB[i] * x[i][XX];
  94             mu_y += qB[i] * x[i][YY];
  95             mu_z += qB[i] * x[i][ZZ];
  96         }
  97         mu_B[XX] = mu_x * gmx::c_enm2Debye;
  98         mu_B[YY] = mu_y * gmx::c_enm2Debye;
  99         mu_B[ZZ] = mu_z * gmx::c_enm2Debye;
 100     }
 101     else
 102     {
 103         copy_dvec(mu, mu_B);
 104     }
 105 }
 106
 107 gmx_bool read_mu(FILE* fp, rvec mu, real* vol)
 108 {
 109     /* For backward compatibility */
 110     real mmm[4];
 111
 112     if (fread(mmm, static_cast<size_t>(4 * sizeof(real)), 1, fp) != 1)
 113     {
 114         return FALSE;
 115     }
 116
 117     copy_rvec(mmm, mu);
 118     *vol = mmm[3];
 119
 120     return TRUE;
 121 }