Open Qmin
0.8.0
GPU-accelerated Q-tensor-based liquid crystal simulations
|
A landau-de gennes q-tensor framework force computer...currently working with the one-constant approximation for the distortion term. More...
#include <landauDeGennesLC.h>
Public Member Functions | |
landauDeGennesLC (bool _neverGPU=false) | |
landauDeGennesLC (scalar _A, scalar _B, scalar _C, scalar _L1) | |
landauDeGennesLC (scalar _A, scalar _B, scalar _C, scalar _L1, scalar _L2, scalar _L3, scalar _L4, scalar _L6) | |
landauDeGennesLC (scalar _A, scalar _B, scalar _C, scalar _L1, scalar _L2, scalar _L3orWavenumber, distortionEnergyType _type) | |
void | baseInitialization () |
set up a few basic things (common force tuners, number of energy components, etc.) More... | |
virtual void | setModel (shared_ptr< cubicLattice > _model) |
The model setting creates an additional data structure to help with 2- or 3- constant approximation. More... | |
virtual void | computeForces (GPUArray< dVec > &forces, bool zeroOutForce=true, int type=0) |
use the "type" flag to select either bulk or boundary routines More... | |
virtual void | computeForceGPU (GPUArray< dVec > &forces, bool zeroOutForce=true) |
As an example of usage, we'll implement an n-Vector model force w/ nearest-neighbor interactions. More... | |
void | setPhaseConstants (scalar _a=-1, scalar _b=-12.325581395, scalar _c=10.058139535) |
void | setElasticConstants (scalar _l1=2.32, scalar _l2=0, scalar _l3=0, scalar _l4=0, scalar _l6=0) |
void | setNumberOfConstants (distortionEnergyType _type) |
virtual void | computeForceCPU (GPUArray< dVec > &forces, bool zeroOutForce=true, int type=0) |
virtual void | computeObjectForces (int objectIdx) |
compute the forces on the objects in the system More... | |
virtual void | computeFirstDerivatives () |
Precompute the first derivatives at all of the LC Sites. More... | |
virtual void | computeStressTensors (GPUArray< int > &sites, GPUArray< Matrix3x3 > &stress) |
compute the stress tensors at the given set of sites More... | |
virtual void | computeBoundaryForcesCPU (GPUArray< dVec > &forces, bool zeroOutForce) |
virtual void | computeBoundaryForcesGPU (GPUArray< dVec > &forces, bool zeroOutForce) |
virtual void | computeEorHFieldForcesCPU (GPUArray< dVec > &forces, bool zeroOutForce, scalar3 field, scalar anisotropicSusceptibility, scalar vacuumPermeability) |
virtual void | computeEorHFieldForcesGPU (GPUArray< dVec > &forces, bool zeroOutForce, scalar3 field, scalar anisotropicSusceptibility, scalar vacuumPermeability) |
virtual void | computeEnergyCPU (bool verbose=false) |
As an example of usage, we'll implement an n-Vector model force w/ nearest-neighbor interactions. More... | |
virtual void | computeEnergyGPU (bool verbose=false) |
void | printTuners () |
void | setEField (scalar3 field, scalar eps, scalar eps0, scalar deltaEps) |
void | setHField (scalar3 field, scalar chi, scalar _mu0, scalar _deltaChi) |
virtual scalar | getClassSize () |
virtual void | computeForceCPU (GPUArray< dVec > &forces, bool zeroOutForce=true) |
As an example of usage, we'll implement an n-Vector model force w/ nearest-neighbor interactions. More... | |
void | setJ (scalar _j) |
virtual scalar | computeEnergy (bool verbose=false) |
compute the energy associated with this force More... | |
virtual void | setModel (shared_ptr< simpleModel > _model) |
virtual function to allow the model to be a derived class More... | |
virtual void | setForceParameters (vector< scalar > ¶ms) |
some generic function to set parameters More... | |
void | setSimulation (shared_ptr< basicSimulation > _sim) |
set the simulation More... | |
virtual MatrixDxD | computePressureTensor () |
compute the system-averaged pressure tensor; return identity if the force hasn't defined this yet More... | |
virtual void | setGPU (bool _useGPU=true) |
Enforce GPU operation. More... | |
void | setNeighborList (shared_ptr< neighborList > _neighbor) |
tell the force to use a neighbor list More... | |
virtual void | setNThreads (int n) |
allow for setting multiple threads More... | |
Public Attributes | |
vector< scalar > | energyComponents |
A vector storing the components of energy (phase,distortion,anchoring) More... | |
GPUArray< scalar > | energyDensity |
the free energy density at each lattice site More... | |
GPUArray< scalar > | energyDensityReduction |
A helper array for energy reductions. More... | |
GPUArray< scalar3 > | objectForceArray |
the force from stresses at the surface of an object More... | |
shared_ptr< kernelTuner > | forceTuner |
kernelTuner object More... | |
shared_ptr< basicSimulation > | sim |
A pointer to the governing simulation. More... | |
shared_ptr< simpleModel > | model |
A pointer to a simpleModel that the updater acts on. More... | |
bool | useGPU |
whether the updater does its work on the GPU or not More... | |
bool | neverGPU |
whether the updater never does work on the GPU More... | |
scalar | energy |
Forces might update the total energy associated with them. More... | |
GPUArray< scalar > | energyPerParticle |
on the gpu, this is per particle and then a reduction can be called More... | |
bool | useNeighborList |
does the force get an assist from a neighbor list? More... | |
shared_ptr< neighborList > | neighbors |
a pointer to a neighbor list the force might use More... | |
int | nThreads =1 |
number of threads to use if compiled with openmp More... | |
Protected Member Functions | |
virtual void | computeL1BulkCPU (GPUArray< dVec > &forces, bool zeroOutForce) |
Compute L1 distortion terms in the bulk and the phase force. More... | |
virtual void | computeL1BoundaryCPU (GPUArray< dVec > &forces, bool zeroOutForce) |
Compute L1 distortion terms at boundaries and the phase force. More... | |
virtual void | computeAllDistortionTermsBulkCPU (GPUArray< dVec > &forces, bool zeroOutForce) |
Compute all distortion terms in the bulk and the phase force. More... | |
virtual void | computeAllDistortionTermsBoundaryCPU (GPUArray< dVec > &forces, bool zeroOutForce) |
Compute all distortion terms at boundaries and the phase force. More... | |
Protected Attributes | |
scalar | A |
constants, etc. More... | |
scalar | B |
scalar | C |
scalar | L1 |
scalar | L2 |
scalar | L3 |
scalar | L4 |
scalar | L6 |
scalar | q0 |
scalar3 | Efield |
scalar | deltaEpsilon |
scalar | epsilon0 |
scalar | epsilon |
scalar3 | Hfield |
scalar | deltaChi |
scalar | Chi |
scalar | mu0 |
distortionEnergyType | numberOfConstants |
number of elastic constants More... | |
bool | computeEfieldContribution |
switches for extra parts of the energy/force calculations More... | |
bool | computeHfieldContribution |
GPUArray< cubicLatticeDerivativeVector > | forceCalculationAssist |
for 2- and 3- constant approximations, the force calculation is helped by first pre-computing first derivatives More... | |
shared_ptr< kernelTuner > | forceAssistTuner |
performance for the first derivative calculation More... | |
shared_ptr< kernelTuner > | boundaryForceTuner |
performance for the boundary force kernel More... | |
shared_ptr< kernelTuner > | l24ForceTuner |
performance for the l24 force kernel More... | |
shared_ptr< kernelTuner > | fieldForceTuner |
performance for the E/H field force kernel More... | |
shared_ptr< cubicLattice > | lattice |
scalar | J |
if all lattice interactions are uniform More... | |
A landau-de gennes q-tensor framework force computer...currently working with the one-constant approximation for the distortion term.
landauDeGennesLC::landauDeGennesLC | ( | bool | _neverGPU = false | ) |
References baseInitialization(), oneConstant, and setNumberOfConstants().
landauDeGennesLC::landauDeGennesLC | ( | scalar | _A, |
scalar | _B, | ||
scalar | _C, | ||
scalar | _L1, | ||
scalar | _L2, | ||
scalar | _L3, | ||
scalar | _L4, | ||
scalar | _L6 | ||
) |
References baseInitialization(), multiConstant, and setNumberOfConstants().
landauDeGennesLC::landauDeGennesLC | ( | scalar | _A, |
scalar | _B, | ||
scalar | _C, | ||
scalar | _L1, | ||
scalar | _L2, | ||
scalar | _L3orWavenumber, | ||
distortionEnergyType | _type | ||
) |
References baseInitialization(), and setNumberOfConstants().
void landauDeGennesLC::baseInitialization | ( | ) |
set up a few basic things (common force tuners, number of energy components, etc.)
References boundaryForceTuner, computeEfieldContribution, computeHfieldContribution, energyComponents, fieldForceTuner, forceAssistTuner, baseLatticeForce::forceTuner, l24ForceTuner, and force::useNeighborList.
Referenced by landauDeGennesLC().
|
virtual |
The model setting creates an additional data structure to help with 2- or 3- constant approximation.
Reimplemented from baseLatticeForce.
References energyDensity, energyDensityReduction, baseLatticeForce::lattice, force::model, multiConstant, numberOfConstants, and GPUArray< T >::resize().
|
virtual |
use the "type" flag to select either bulk or boundary routines
Reimplemented from baseLatticeForce.
References computeForceCPU(), computeForceGPU(), and force::useGPU.
|
virtual |
As an example of usage, we'll implement an n-Vector model force w/ nearest-neighbor interactions.
Reimplemented from baseLatticeForce.
References A, B, C, computeBoundaryForcesGPU(), computeEfieldContribution, computeEorHFieldForcesGPU(), computeFirstDerivatives(), computeHfieldContribution, ArrayHandle< T >::data, deltaChi, deltaEpsilon, access_location::device, Efield, epsilon0, forceCalculationAssist, baseLatticeForce::forceTuner, gpu_qTensor_multiConstantForce(), gpu_qTensor_oneConstantForce(), Hfield, L1, L2, L3, L4, L6, baseLatticeForce::lattice, mu0, multiConstant, numberOfConstants, oneConstant, access_mode::read, and access_mode::readwrite.
Referenced by computeForces().
void landauDeGennesLC::setNumberOfConstants | ( | distortionEnergyType | _type | ) |
References multiConstant, and numberOfConstants.
Referenced by landauDeGennesLC().
|
virtual |
References computeAllDistortionTermsBoundaryCPU(), computeAllDistortionTermsBulkCPU(), computeBoundaryForcesCPU(), computeEfieldContribution, computeEorHFieldForcesCPU(), computeFirstDerivatives(), computeHfieldContribution, computeL1BoundaryCPU(), computeL1BulkCPU(), deltaChi, deltaEpsilon, Efield, epsilon0, Hfield, baseLatticeForce::lattice, mu0, multiConstant, numberOfConstants, and oneConstant.
Referenced by computeForces().
|
virtual |
compute the forces on the objects in the system
References computeStressTensors(), ArrayHandle< T >::data, access_location::device, GPUArray< T >::getNumElements(), gpu_qTensor_computeObjectForceFromStresses(), access_location::host, baseLatticeForce::lattice, make_scalar3(), objectForceArray, access_mode::overwrite, access_mode::read, GPUArray< T >::resize(), scalar3, and force::useGPU.
|
virtual |
Precompute the first derivatives at all of the LC Sites.
References ArrayHandle< T >::data, access_location::device, forceAssistTuner, forceCalculationAssist, GPUArray< T >::getNumElements(), gpu_qTensor_firstDerivatives(), access_location::host, idx, baseLatticeForce::lattice, force::neighbors, access_mode::read, access_mode::readwrite, GPUArray< T >::resize(), and force::useGPU.
Referenced by computeForceCPU(), computeForceGPU(), and computeStressTensors().
|
virtual |
compute the stress tensors at the given set of sites
expression from "Hierarchical self-assembly of nematic colloidal superstructures" PHYSICAL REVIEW E 77, 061706 (2008)
References computeEnergyCPU(), computeFirstDerivatives(), ArrayHandle< T >::data, energyDensity, forceCalculationAssist, GPUArray< T >::getNumElements(), access_location::host, L1, n, numberOfConstants, oneConstant, access_mode::overwrite, access_mode::read, GPUArray< T >::resize(), Matrix3x3::set(), Matrix3x3::x11, Matrix3x3::x22, and Matrix3x3::x33.
Referenced by computeObjectForces().
|
virtual |
References computeBoundaryForce(), ArrayHandle< T >::data, baseLatticeForce::lattice, make_dVec(), and force::neighbors.
Referenced by computeForceCPU().
|
virtual |
|
virtual |
References ArrayHandle< T >::data, baseLatticeForce::lattice, make_dVec(), force::neighbors, and scalar.
Referenced by computeForceCPU().
|
virtual |
As an example of usage, we'll implement an n-Vector model force w/ nearest-neighbor interactions.
Reimplemented from baseLatticeForce.
References A, B, C, Chi, computeBoundaryEnergy(), computeEfieldContribution, computeHfieldContribution, ArrayHandle< T >::data, deltaChi, deltaEpsilon, Efield, force::energy, energyComponents, energyDensity, epsilon, epsilon0, Hfield, L1, L2, L3, L4, L6, baseLatticeForce::lattice, mu0, force::neighbors, scalar, TrQ2(), TrQ2Squared(), and TrQ3().
Referenced by computeStressTensors().
|
virtual |
Reimplemented from baseLatticeForce.
References A, B, C, Chi, computeEfieldContribution, computeHfieldContribution, ArrayHandle< T >::data, deltaChi, deltaEpsilon, access_location::device, Efield, force::energy, energyDensity, energyDensityReduction, epsilon, epsilon0, getNumBlocksAndThreads(), gpu_computeAllEnergyTerms(), gpuReduction(), Hfield, L1, L2, L3, L4, L6, baseLatticeForce::lattice, mu0, access_mode::read, access_mode::readwrite, and scalar.
|
inline |
References computeEfieldContribution, deltaEpsilon, Efield, epsilon, and epsilon0.
|
inline |
References Chi, computeHfieldContribution, deltaChi, Hfield, and mu0.
|
inlinevirtual |
Reimplemented from baseLatticeForce.
References energyComponents, energyDensity, forceCalculationAssist, baseLatticeForce::getClassSize(), GPUArray< T >::getNumElements(), objectForceArray, and scalar.
|
protectedvirtual |
Compute L1 distortion terms in the bulk and the phase force.
References A, B, lcForce::bulkL1Force(), C, ArrayHandle< T >::data, derivativeTrQ2(), derivativeTrQ2Squared(), derivativeTrQ3(), force::force(), access_location::host, L1, baseLatticeForce::lattice, access_mode::read, and scalar.
Referenced by computeForceCPU().
|
protectedvirtual |
Compute L1 distortion terms at boundaries and the phase force.
References A, B, lcForce::boundaryL1Force(), lcForce::bulkL1Force(), C, ArrayHandle< T >::data, derivativeTrQ2(), derivativeTrQ2Squared(), derivativeTrQ3(), force::force(), access_location::host, L1, baseLatticeForce::lattice, access_mode::read, and scalar.
Referenced by computeForceCPU().
|
protectedvirtual |
Compute all distortion terms in the bulk and the phase force.
References A, B, lcForce::bulkL1Force(), lcForce::bulkL2Force(), lcForce::bulkL3Force(), lcForce::bulkL4Force(), lcForce::bulkL6Force(), C, ArrayHandle< T >::data, derivativeTrQ2(), derivativeTrQ2Squared(), derivativeTrQ3(), force::force(), forceCalculationAssist, access_location::host, L1, L2, L3, L4, L6, baseLatticeForce::lattice, access_mode::read, and scalar.
Referenced by computeForceCPU().
|
protectedvirtual |
Compute all distortion terms at boundaries and the phase force.
References A, B, lcForce::boundaryL1Force(), lcForce::boundaryL2Force(), lcForce::boundaryL3Force(), lcForce::boundaryL4Force(), lcForce::boundaryL6Force(), lcForce::bulkL1Force(), lcForce::bulkL2Force(), lcForce::bulkL3Force(), lcForce::bulkL4Force(), lcForce::bulkL6Force(), C, ArrayHandle< T >::data, derivativeTrQ2(), derivativeTrQ2Squared(), derivativeTrQ3(), force::force(), forceCalculationAssist, lcForce::getBoundaryCase(), access_location::host, L1, L2, L3, L4, L6, baseLatticeForce::lattice, access_mode::read, and scalar.
Referenced by computeForceCPU().
|
virtualinherited |
As an example of usage, we'll implement an n-Vector model force w/ nearest-neighbor interactions.
References ArrayHandle< T >::data, baseLatticeForce::J, baseLatticeForce::lattice, make_dVec(), and force::neighbors.
Referenced by baseLatticeForce::computeForces().
|
inlineinherited |
References baseLatticeForce::J.
|
inlinevirtualinherited |
compute the energy associated with this force
Reimplemented from force.
References baseLatticeForce::computeEnergyCPU(), baseLatticeForce::computeEnergyGPU(), force::energy, and force::useGPU.
|
inlinevirtualinherited |
virtual function to allow the model to be a derived class
References force::model.
|
virtualinherited |
some generic function to set parameters
|
inlineinherited |
set the simulation
References force::sim.
|
inlinevirtualinherited |
compute the system-averaged pressure tensor; return identity if the force hasn't defined this yet
|
inlinevirtualinherited |
Enforce GPU operation.
References force::neighbors, force::useGPU, and force::useNeighborList.
|
inlineinherited |
tell the force to use a neighbor list
References force::neighbors, and force::useNeighborList.
|
inlinevirtualinherited |
allow for setting multiple threads
References n, and force::nThreads.
vector<scalar> landauDeGennesLC::energyComponents |
A vector storing the components of energy (phase,distortion,anchoring)
Referenced by baseInitialization(), computeEnergyCPU(), and getClassSize().
the free energy density at each lattice site
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeStressTensors(), getClassSize(), landauDeGennesLC(), and setModel().
A helper array for energy reductions.
Referenced by computeEnergyGPU(), landauDeGennesLC(), and setModel().
the force from stresses at the surface of an object
Referenced by computeObjectForces(), getClassSize(), and landauDeGennesLC().
|
protected |
constants, etc.
Referenced by computeAllDistortionTermsBoundaryCPU(), computeAllDistortionTermsBulkCPU(), computeEnergyCPU(), computeEnergyGPU(), computeForceGPU(), computeL1BoundaryCPU(), computeL1BulkCPU(), and setPhaseConstants().
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setEField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setEField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setEField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), and setEField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setHField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setHField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), and setHField().
|
protected |
Referenced by computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setHField().
|
protected |
number of elastic constants
Referenced by computeForceCPU(), computeForceGPU(), computeStressTensors(), printTuners(), setModel(), and setNumberOfConstants().
|
protected |
switches for extra parts of the energy/force calculations
Referenced by baseInitialization(), computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setEField().
|
protected |
Referenced by baseInitialization(), computeEnergyCPU(), computeEnergyGPU(), computeForceCPU(), computeForceGPU(), and setHField().
|
protected |
for 2- and 3- constant approximations, the force calculation is helped by first pre-computing first derivatives
Referenced by computeAllDistortionTermsBoundaryCPU(), computeAllDistortionTermsBulkCPU(), computeFirstDerivatives(), computeForceGPU(), computeStressTensors(), getClassSize(), and landauDeGennesLC().
|
protected |
performance for the first derivative calculation
Referenced by baseInitialization(), computeFirstDerivatives(), and printTuners().
|
protected |
performance for the boundary force kernel
Referenced by baseInitialization(), computeBoundaryForcesGPU(), computeEorHFieldForcesGPU(), and printTuners().
|
protected |
performance for the l24 force kernel
Referenced by baseInitialization().
|
protected |
performance for the E/H field force kernel
Referenced by baseInitialization(), and computeEorHFieldForcesGPU().
|
inherited |
kernelTuner object
Referenced by baseInitialization(), baseLatticeForce::baseLatticeForce(), baseLatticeForce::computeForceGPU(), computeForceGPU(), and printTuners().
|
protectedinherited |
Referenced by computeAllDistortionTermsBoundaryCPU(), computeAllDistortionTermsBulkCPU(), computeBoundaryForcesCPU(), computeBoundaryForcesGPU(), baseLatticeForce::computeEnergyCPU(), computeEnergyCPU(), computeEnergyGPU(), computeEorHFieldForcesCPU(), computeEorHFieldForcesGPU(), computeFirstDerivatives(), baseLatticeForce::computeForceCPU(), computeForceCPU(), baseLatticeForce::computeForceGPU(), computeForceGPU(), computeL1BoundaryCPU(), computeL1BulkCPU(), computeObjectForces(), printTuners(), setModel(), and baseLatticeForce::setModel().
|
protectedinherited |
if all lattice interactions are uniform
Referenced by baseLatticeForce::baseLatticeForce(), baseLatticeForce::computeEnergyCPU(), baseLatticeForce::computeForceCPU(), baseLatticeForce::computeForceGPU(), and baseLatticeForce::setJ().
|
inherited |
A pointer to the governing simulation.
Referenced by force::setSimulation().
|
inherited |
A pointer to a simpleModel that the updater acts on.
Referenced by setModel(), force::setModel(), and baseLatticeForce::setModel().
|
inherited |
whether the updater does its work on the GPU or not
Referenced by baseLatticeForce::computeEnergy(), computeFirstDerivatives(), baseLatticeForce::computeForces(), computeForces(), computeObjectForces(), force::force(), and force::setGPU().
|
inherited |
whether the updater never does work on the GPU
Referenced by landauDeGennesLC().
|
inherited |
Forces might update the total energy associated with them.
Referenced by baseLatticeForce::computeEnergy(), baseLatticeForce::computeEnergyCPU(), computeEnergyCPU(), baseLatticeForce::computeEnergyGPU(), and computeEnergyGPU().
on the gpu, this is per particle and then a reduction can be called
Referenced by force::getClassSize(), and landauDeGennesLC().
|
inherited |
does the force get an assist from a neighbor list?
Referenced by baseInitialization(), baseLatticeForce::baseLatticeForce(), force::force(), force::setGPU(), and force::setNeighborList().
|
inherited |
a pointer to a neighbor list the force might use
Referenced by computeBoundaryForcesCPU(), baseLatticeForce::computeEnergyCPU(), computeEnergyCPU(), computeEorHFieldForcesCPU(), computeFirstDerivatives(), baseLatticeForce::computeForceCPU(), force::setGPU(), and force::setNeighborList().
|
inherited |
number of threads to use if compiled with openmp
Referenced by force::setNThreads().