Skip to content

Commit

Permalink
Merge pull request #1637 from glotzerlab/external_magnetic_field
Browse files Browse the repository at this point in the history
External magnetic field
  • Loading branch information
joaander authored Nov 16, 2023
2 parents da999a3 + 0cb657c commit 3db0b6b
Show file tree
Hide file tree
Showing 12 changed files with 371 additions and 26 deletions.
3 changes: 2 additions & 1 deletion hoomd/md/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -87,6 +87,7 @@ set(_md_headers ActiveForceComputeGPU.h
EvaluatorSpecialPairLJ.h
EvaluatorSpecialPairCoulomb.h
EvaluatorExternalElectricField.h
EvaluatorExternalMagneticField.h
EvaluatorExternalPeriodic.h
EvaluatorPairALJ.h
EvaluatorPairBuckingham.h
Expand Down Expand Up @@ -426,7 +427,7 @@ foreach(_evaluator ${_triplets})
endif()
endforeach()

set(_external_evaluators Periodic ElectricField)
set(_external_evaluators Periodic ElectricField MagneticField)

foreach(_evaluator ${_external_evaluators})
configure_file(export_PotentialExternal.cc.inc
Expand Down
15 changes: 14 additions & 1 deletion hoomd/md/EvaluatorExternalElectricField.h
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@

#include "hoomd/BoxDim.h"
#include "hoomd/HOOMDMath.h"
#include "hoomd/VectorMath.h"
#include <math.h>

/*! \file EvaluatorExternalElectricField.h
Expand Down Expand Up @@ -75,13 +76,19 @@ class EvaluatorExternalElectricField
\param params per-type parameters of external potential
*/
DEVICE EvaluatorExternalElectricField(Scalar3 X,
quat<Scalar> q,
const BoxDim& box,
const param_type& params,
const field_type& field)
: m_pos(X), m_box(box), m_E(params.E)
{
}

DEVICE static bool isAnisotropic()
{
return false;
}

//! ExternalElectricField needs charges
DEVICE static bool needsCharge()
{
Expand All @@ -106,10 +113,12 @@ class EvaluatorExternalElectricField

//! Evaluate the force, energy and virial
/*! \param F force vector
\param T torque vector
\param energy value of the energy
\param virial array of six scalars for the upper triangular virial tensor
*/
DEVICE void evalForceEnergyAndVirial(Scalar3& F, Scalar& energy, Scalar* virial)
DEVICE void
evalForceTorqueEnergyAndVirial(Scalar3& F, Scalar3& T, Scalar& energy, Scalar* virial)
{
F = m_qi * m_E;
energy = -m_qi * dot(m_E, m_pos);
Expand All @@ -120,6 +129,10 @@ class EvaluatorExternalElectricField
virial[3] = F.y * m_pos.y;
virial[4] = F.y * m_pos.z;
virial[5] = F.z * m_pos.z;

T.x = Scalar(0.0);
T.y = Scalar(0.0);
T.z = Scalar(0.0);
}

#ifndef __HIPCC__
Expand Down
164 changes: 164 additions & 0 deletions hoomd/md/EvaluatorExternalMagneticField.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,164 @@
// Copyright (c) 2009-2023 The Regents of the University of Michigan.
// Part of HOOMD-blue, released under the BSD 3-Clause License.

#ifndef __EVALUATOR_EXTERNAL_MAGNETIC_FIELD_H__
#define __EVALUATOR_EXTERNAL_MAGNETIC_FIELD_H__

#ifndef __HIPCC__
#include <string>
#endif

#include "hoomd/BoxDim.h"
#include "hoomd/HOOMDMath.h"
#include "hoomd/VectorMath.h"
#include <math.h>

/*! \file EvaluatorExternalMagneticField.h
\brief Defines the external potential evaluator to induce a magnetic field
*/

// need to declare these class methods with __device__ qualifiers when building in nvcc
// DEVICE is __host__ __device__ when included in nvcc and blank when included into the host
// compiler
#ifdef __HIPCC__
#define DEVICE __device__
#else
#define DEVICE
#endif

namespace hoomd
{
namespace md
{
//! Class for evaluating an magnetic field
/*! <b>General Overview</b>
The external potential \f$V(\theta) \f$ is implemented using the following formula:
\f[
V(\theta}) = - \vec{B} \cdot \vec{n}_i(\theta)
\f]
where \f$B\f$ is the strength of the magnetic field and \f$\vec{n}_i\f$ is the magnetic moment
of particle i.
*/
class EvaluatorExternalMagneticField
{
public:
//! type of parameters this external potential accepts
struct param_type
{
Scalar3 B;
Scalar3 mu;

#ifndef __HIPCC__
param_type() : B(make_scalar3(0, 0, 0)), mu(make_scalar3(0, 0, 0)) { }

param_type(pybind11::dict params)
{
pybind11::tuple py_B = params["B"];
B.x = pybind11::cast<Scalar>(py_B[0]);
B.y = pybind11::cast<Scalar>(py_B[1]);
B.z = pybind11::cast<Scalar>(py_B[2]);
pybind11::tuple py_mu = params["mu"];
mu.x = pybind11::cast<Scalar>(py_mu[0]);
mu.y = pybind11::cast<Scalar>(py_mu[1]);
mu.z = pybind11::cast<Scalar>(py_mu[2]);
}

pybind11::dict toPython()
{
pybind11::dict d;
d["B"] = pybind11::make_tuple(B.x, B.y, B.z);
d["mu"] = pybind11::make_tuple(mu.x, mu.y, mu.z);
return d;
}
#endif // ifndef __HIPCC__
} __attribute__((aligned(16)));

typedef void* field_type;

//! Constructs the external field evaluator
/*! \param X position of particle
\param box box dimensions
\param params per-type parameters of external potential
*/
DEVICE EvaluatorExternalMagneticField(Scalar3 X,
quat<Scalar> q,
const BoxDim& box,
const param_type& params,
const field_type& field)
: m_q(q), m_B(params.B), m_mu(params.mu)
{
}

DEVICE static bool isAnisotropic()
{
return true;
}

//! ExternalMagneticField needs charges
DEVICE static bool needsCharge()
{
return false;
}

//! Accept the optional charge value
/*! \param qi Charge of particle i
*/
DEVICE void setCharge(Scalar qi) { }

//! Declares additional virial contributions are needed for the external field
/*! No contribution
*/
DEVICE static bool requestFieldVirialTerm()
{
return false;
}

//! Evaluate the force, energy and virial
/*! \param F force vector
\param T torque vector
\param energy value of the energy
\param virial array of six scalars for the upper triangular virial tensor
*/
DEVICE void
evalForceTorqueEnergyAndVirial(Scalar3& F, Scalar3& T, Scalar& energy, Scalar* virial)
{
vec3<Scalar> dir = rotate(m_q, m_mu);

vec3<Scalar> T_vec = cross(dir, m_B);

T.x = T_vec.x;
T.y = T_vec.y;
T.z = T_vec.z;

energy = -dot(dir, m_B);

F.x = Scalar(0.0);
F.y = Scalar(0.0);
F.z = Scalar(0.0);

for (unsigned int i = 0; i < 6; i++)
virial[i] = Scalar(0.0);
}

#ifndef __HIPCC__
//! Get the name of this potential
/*! \returns The potential name.
*/
static std::string getName()
{
return std::string("b_field");
}
#endif

protected:
quat<Scalar> m_q; //!< Particle orientation
vec3<Scalar> m_B; //!< Magnetic field vector (box frame).
vec3<Scalar> m_mu; //!< Magnetic dipole moment (particle frame).
};

} // end namespace md
} // end namespace hoomd

#endif // __EVALUATOR_EXTERNAL_LAMELLAR_H__
16 changes: 15 additions & 1 deletion hoomd/md/EvaluatorExternalPeriodic.h
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@

#include "hoomd/BoxDim.h"
#include "hoomd/HOOMDMath.h"
#include "hoomd/VectorMath.h"
#include <math.h>

/*! \file EvaluatorExternalPeriodic.h
Expand Down Expand Up @@ -89,6 +90,7 @@ class EvaluatorExternalPeriodic
\param params per-type parameters of external potential
*/
DEVICE EvaluatorExternalPeriodic(Scalar3 X,
quat<Scalar> q,
const BoxDim& box,
const param_type& params,
const field_type& field)
Expand All @@ -97,6 +99,11 @@ class EvaluatorExternalPeriodic
{
}

DEVICE static bool isAnisotropic()
{
return false;
}

//! External Periodic doesn't need charges
DEVICE static bool needsCharge()
{
Expand All @@ -117,17 +124,24 @@ class EvaluatorExternalPeriodic

//! Evaluate the force, energy and virial
/*! \param F force vector
\param T torque vector
\param energy value of the energy
\param virial array of six scalars for the upper triangular virial tensor
*/
DEVICE void evalForceEnergyAndVirial(Scalar3& F, Scalar& energy, Scalar* virial)
DEVICE void
evalForceTorqueEnergyAndVirial(Scalar3& F, Scalar3& T, Scalar& energy, Scalar* virial)
{
Scalar3 a2 = make_scalar3(0, 0, 0);
Scalar3 a3 = make_scalar3(0, 0, 0);

F.x = Scalar(0.0);
F.y = Scalar(0.0);
F.z = Scalar(0.0);

T.x = Scalar(0.0);
T.y = Scalar(0.0);
T.z = Scalar(0.0);

energy = Scalar(0.0);

// For this potential, since it uses scaled positions, the virial is always zero.
Expand Down
19 changes: 17 additions & 2 deletions hoomd/md/EvaluatorWalls.h
Original file line number Diff line number Diff line change
Expand Up @@ -134,11 +134,20 @@ template<class evaluator> class EvaluatorWalls
typedef wall_type field_type;

//! Constructs the external wall potential evaluator
DEVICE EvaluatorWalls(Scalar3 pos, const BoxDim& box, const param_type& p, const field_type& f)
DEVICE EvaluatorWalls(Scalar3 pos,
quat<Scalar> q,
const BoxDim& box,
const param_type& p,
const field_type& f)
: m_pos(pos), m_field(f), m_params(p)
{
}

DEVICE static bool isAnisotropic()
{
return false;
}

//! Charges not supported by walls evals
DEVICE static bool needsCharge()
{
Expand Down Expand Up @@ -229,11 +238,17 @@ template<class evaluator> class EvaluatorWalls
}

//! Generates force and energy from standard evaluators using wall geometry functions
DEVICE void evalForceEnergyAndVirial(Scalar3& F, Scalar& energy, Scalar* virial)
DEVICE void
evalForceTorqueEnergyAndVirial(Scalar3& F, Scalar3& T, Scalar& energy, Scalar* virial)
{
F.x = Scalar(0.0);
F.y = Scalar(0.0);
F.z = Scalar(0.0);

T.x = Scalar(0.0);
T.y = Scalar(0.0);
T.z = Scalar(0.0);

energy = Scalar(0.0);
// initialize virial
for (unsigned int i = 0; i < 6; i++)
Expand Down
Loading

0 comments on commit 3db0b6b

Please sign in to comment.