1
0
mirror of https://github.com/paboyle/Grid.git synced 2025-04-09 21:50:45 +01:00

Cleaning up files and directories

This commit is contained in:
Guido Cossu 2017-01-26 15:24:49 +00:00
parent 7996f06335
commit e863a948e3
15 changed files with 587 additions and 551 deletions

View File

@ -89,7 +89,6 @@ Author: paboyle <paboyle@ph.ed.ac.uk>
#include <Grid/qcd/hmc/HMCModules.h> #include <Grid/qcd/hmc/HMCModules.h>
#include <Grid/qcd/modules/mods.h> #include <Grid/qcd/modules/mods.h>
#include <Grid/qcd/hmc/HMCResourceManager.h> #include <Grid/qcd/hmc/HMCResourceManager.h>
#include <Grid/qcd/hmc/HmcRunner.h>
#include <Grid/qcd/hmc/GenericHMCrunner.h> #include <Grid/qcd/hmc/GenericHMCrunner.h>
#include <Grid/qcd/hmc/HMCRunnerModule.h> #include <Grid/qcd/hmc/HMCRunnerModule.h>

View File

@ -4,10 +4,11 @@ Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/qcd/action/ActionBase.h Source file: ./lib/qcd/action/ActionBase.h
Copyright (C) 2015 Copyright (C) 2015-2016
Author: Peter Boyle <paboyle@ph.ed.ac.uk> Author: Peter Boyle <paboyle@ph.ed.ac.uk>
Author: neo <cossu@post.kek.jp> Author: neo <cossu@post.kek.jp>
Author: Guido Cossu <guido.cossu@ed.ac.uk>
This program is free software; you can redistribute it and/or modify This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by it under the terms of the GNU General Public License as published by
@ -27,97 +28,29 @@ See the full license in the file "LICENSE" in the top level distribution
directory directory
*************************************************************************************/ *************************************************************************************/
/* END LEGAL */ /* END LEGAL */
#ifndef QCD_ACTION_BASE
#define QCD_ACTION_BASE #ifndef ACTION_BASE_H
#define ACTION_BASE_H
namespace Grid { namespace Grid {
namespace QCD { namespace QCD {
template <class GaugeField> template <class GaugeField >
class Action { class Action
{
public: public:
bool is_smeared = false; bool is_smeared = false;
// Heatbath? // Heatbath?
virtual void refresh(const GaugeField& U, GridParallelRNG& pRNG) = 0; // refresh pseudofermions virtual void refresh(const GaugeField& U, GridParallelRNG& pRNG) = 0; // refresh pseudofermions
virtual RealD S(const GaugeField& U) = 0; // evaluate the action virtual RealD S(const GaugeField& U) = 0; // evaluate the action
virtual void deriv(const GaugeField& U, GaugeField& dSdU) = 0; // evaluate the action derivative virtual void deriv(const GaugeField& U, GaugeField& dSdU) = 0; // evaluate the action derivative
virtual std::string action_name() = 0; // return the action name virtual std::string action_name() = 0; // return the action name
virtual std::string LogParameters() = 0; // prints action parameters virtual std::string LogParameters() = 0; // prints action parameters
virtual ~Action(){} virtual ~Action(){}
}; };
// Indexing of tuple types
template <class T, class Tuple>
struct Index;
template <class T, class... Types>
struct Index<T, std::tuple<T, Types...>> {
static const std::size_t value = 0;
};
template <class T, class U, class... Types>
struct Index<T, std::tuple<U, Types...>> {
static const std::size_t value = 1 + Index<T, std::tuple<Types...>>::value;
};
template <class Field, class Repr = NoHirep >
struct ActionLevel {
public:
unsigned int multiplier;
// Fundamental repr actions separated because of the smearing
typedef Action<Field>* ActPtr;
// construct a tuple of vectors of the actions for the corresponding higher
// representation fields
typedef typename AccessTypes<Action, Repr>::VectorCollection action_collection;
action_collection actions_hirep;
typedef typename AccessTypes<Action, Repr>::FieldTypeCollection action_hirep_types;
std::vector<ActPtr>& actions;
explicit ActionLevel(unsigned int mul = 1) : actions(std::get<0>(actions_hirep)), multiplier(mul) {
// initialize the hirep vectors to zero.
// apply(this->resize, actions_hirep, 0); //need a working resize
assert(mul >= 1);
}
// void push_back(ActPtr ptr) { actions.push_back(ptr); }
template < class GenField >
void push_back(Action<GenField>* ptr) {
// insert only in the correct vector
std::get< Index < GenField, action_hirep_types>::value >(actions_hirep).push_back(ptr);
};
template < class ActPtr>
static void resize(ActPtr ap, unsigned int n){
ap->resize(n);
}
//template <std::size_t I>
//auto getRepresentation(Repr& R)->decltype(std::get<I>(R).U) {return std::get<I>(R).U;}
// Loop on tuple for a callable function
template <std::size_t I = 1, typename Callable, typename ...Args>
inline typename std::enable_if<I == std::tuple_size<action_collection>::value, void>::type apply(
Callable, Repr& R,Args&...) const {}
template <std::size_t I = 1, typename Callable, typename ...Args>
inline typename std::enable_if<I < std::tuple_size<action_collection>::value, void>::type apply(
Callable fn, Repr& R, Args&... arguments) const {
fn(std::get<I>(actions_hirep), std::get<I>(R.rep), arguments...);
apply<I + 1>(fn, R, arguments...);
}
};
template <class GaugeField, class R>
using ActionSet = std::vector<ActionLevel<GaugeField, R> >;
} }
} }
#endif
#endif // ACTION_BASE_H

View File

@ -1,68 +1,80 @@
/************************************************************************************* /*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/qcd/action/ActionParams.h Source file: ./lib/qcd/action/ActionParams.h
Copyright (C) 2015 Copyright (C) 2015
Author: Peter Boyle <paboyle@ph.ed.ac.uk> Author: Peter Boyle <paboyle@ph.ed.ac.uk>
Author: paboyle <paboyle@ph.ed.ac.uk> Author: paboyle <paboyle@ph.ed.ac.uk>
Author: Guido Cossu <guido.cossu@ed.ac.uk>
This program is free software; you can redistribute it and/or modify This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or the Free Software Foundation; either version 2 of the License, or
(at your option) any later version. (at your option) any later version.
This program is distributed in the hope that it will be useful, This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details. GNU General Public License for more details.
You should have received a copy of the GNU General Public License along You should have received a copy of the GNU General Public License along
with this program; if not, write to the Free Software Foundation, Inc., with this program; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution
directory
*************************************************************************************/
/* END LEGAL */
See the full license in the file "LICENSE" in the top level distribution directory
*************************************************************************************/
/* END LEGAL */
#ifndef GRID_QCD_ACTION_PARAMS_H #ifndef GRID_QCD_ACTION_PARAMS_H
#define GRID_QCD_ACTION_PARAMS_H #define GRID_QCD_ACTION_PARAMS_H
namespace Grid { namespace Grid {
namespace QCD { namespace QCD {
// These can move into a params header and be given MacroMagic serialisation // These can move into a params header and be given MacroMagic serialisation
struct GparityWilsonImplParams { struct GparityWilsonImplParams {
bool overlapCommsCompute; bool overlapCommsCompute;
std::vector<int> twists; std::vector<int> twists;
GparityWilsonImplParams () : twists(Nd,0), overlapCommsCompute(false) {}; GparityWilsonImplParams() : twists(Nd, 0), overlapCommsCompute(false){};
};
}; struct WilsonImplParams {
bool overlapCommsCompute;
WilsonImplParams() : overlapCommsCompute(false){};
};
struct WilsonImplParams { struct OneFlavourRationalParams : Serializable {
bool overlapCommsCompute; GRID_SERIALIZABLE_CLASS_MEMBERS(OneFlavourRationalParams,
WilsonImplParams() : overlapCommsCompute(false) {}; RealD, lo,
}; RealD, hi,
int, MaxIter,
struct OneFlavourRationalParams : Serializable { RealD, tolerance,
GRID_SERIALIZABLE_CLASS_MEMBERS(OneFlavourRationalParams, int, degree,
RealD, lo, int, precision);
RealD, hi,
int, MaxIter, // MaxIter and tolerance, vectors??
RealD, tolerance,
int, degree, // constructor
int, precision); OneFlavourRationalParams( RealD _lo = 0.0,
RealD _hi = 1.0,
// MaxIter and tolerance, vectors?? int _maxit = 1000,
public: RealD tol = 1.0e-8,
OneFlavourRationalParams (RealD _lo = 0.0, RealD _hi = 0.0 , int _degree = 10,
int _maxit = 1000, RealD tol = 1.0e-8, int _precision = 64)
int _degree = 10, int _precision = 64) : : lo(_lo),
lo(_lo), hi(_hi), MaxIter(_maxit), tolerance(tol), degree(_degree), precision(_precision) hi(_hi),
{}; MaxIter(_maxit),
}; tolerance(tol),
degree(_degree),
}} precision(_precision){};
};
}
}
#endif #endif

116
lib/qcd/action/ActionSet.h Normal file
View File

@ -0,0 +1,116 @@
/*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/qcd/action/ActionSet.h
Copyright (C) 2015
Author: Peter Boyle <paboyle@ph.ed.ac.uk>
Author: neo <cossu@post.kek.jp>
This program 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 2 of the License, or
(at your option) any later version.
This program 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, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution
directory
*************************************************************************************/
/* END LEGAL */
#ifndef ACTION_SET_H
#define ACTION_SET_H
namespace Grid {
// Should drop this namespace here
namespace QCD {
//////////////////////////////////
// Indexing of tuple types
//////////////////////////////////
template <class T, class Tuple>
struct Index;
template <class T, class... Types>
struct Index<T, std::tuple<T, Types...>> {
static const std::size_t value = 0;
};
template <class T, class U, class... Types>
struct Index<T, std::tuple<U, Types...>> {
static const std::size_t value = 1 + Index<T, std::tuple<Types...>>::value;
};
////////////////////////////////////////////
// Action Level
// Action collection
// in a integration level
// (for multilevel integration schemes)
////////////////////////////////////////////
template <class Field, class Repr = NoHirep >
struct ActionLevel {
public:
unsigned int multiplier;
// Fundamental repr actions separated because of the smearing
typedef Action<Field>* ActPtr;
// construct a tuple of vectors of the actions for the corresponding higher
// representation fields
typedef typename AccessTypes<Action, Repr>::VectorCollection action_collection;
typedef typename AccessTypes<Action, Repr>::FieldTypeCollection action_hirep_types;
action_collection actions_hirep;
std::vector<ActPtr>& actions;
explicit ActionLevel(unsigned int mul = 1) :
actions(std::get<0>(actions_hirep)), multiplier(mul) {
// initialize the hirep vectors to zero.
// apply(this->resize, actions_hirep, 0); //need a working resize
assert(mul >= 1);
}
template < class GenField >
void push_back(Action<GenField>* ptr) {
// insert only in the correct vector
std::get< Index < GenField, action_hirep_types>::value >(actions_hirep).push_back(ptr);
};
template <class ActPtr>
static void resize(ActPtr ap, unsigned int n) {
ap->resize(n);
}
// Loop on tuple for a callable function
template <std::size_t I = 1, typename Callable, typename ...Args>
inline typename std::enable_if<I == std::tuple_size<action_collection>::value, void>::type apply(Callable, Repr& R,Args&...) const {}
template <std::size_t I = 1, typename Callable, typename ...Args>
inline typename std::enable_if<I < std::tuple_size<action_collection>::value, void>::type apply(Callable fn, Repr& R, Args&... arguments) const {
fn(std::get<I>(actions_hirep), std::get<I>(R.rep), arguments...);
apply<I + 1>(fn, R, arguments...);
}
};
// Define the ActionSet
template <class GaugeField, class R>
using ActionSet = std::vector<ActionLevel<GaugeField, R> >;
} // QCD
} // Grid
#endif // ACTION_SET_H

View File

@ -1,4 +1,4 @@
/************************************************************************************* /*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid Grid physics library, www.github.com/paboyle/Grid
@ -29,7 +29,8 @@ Author: paboyle <paboyle@ph.ed.ac.uk>
See the full license in the file "LICENSE" in the top level distribution directory See the full license in the file "LICENSE" in the top level distribution directory
*************************************************************************************/ *************************************************************************************/
/* END LEGAL */ /* END LEGAL */
#ifndef GRID_QCD_ACTIONS_H #ifndef GRID_QCD_ACTIONS_H
#define GRID_QCD_ACTIONS_H #define GRID_QCD_ACTIONS_H
@ -41,12 +42,14 @@ Author: paboyle <paboyle@ph.ed.ac.uk>
// Abstract base interface // Abstract base interface
//////////////////////////////////////////// ////////////////////////////////////////////
#include <Grid/qcd/action/ActionBase.h> #include <Grid/qcd/action/ActionBase.h>
#include <Grid/qcd/action/ActionSet.h>
#include <Grid/qcd/action/ActionParams.h> #include <Grid/qcd/action/ActionParams.h>
//////////////////////////////////////////// ////////////////////////////////////////////
// Utility functions // Utility functions
//////////////////////////////////////////// ////////////////////////////////////////////
#include <Grid/qcd/action/gauge/GaugeImpl.h> //#include <Grid/qcd/action/gauge/GaugeImplTypes.h>
#include <Grid/qcd/action/gauge/GaugeImplementations.h>
#include <Grid/qcd/utils/WilsonLoops.h> #include <Grid/qcd/utils/WilsonLoops.h>
#include <Grid/qcd/action/fermion/WilsonCompressor.h> //used by all wilson type fermions #include <Grid/qcd/action/fermion/WilsonCompressor.h> //used by all wilson type fermions
@ -63,7 +66,7 @@ Author: paboyle <paboyle@ph.ed.ac.uk>
//////////////////////////////////////////// ////////////////////////////////////////////
// Scalar Actions // Scalar Actions
//////////////////////////////////////////// ////////////////////////////////////////////
#include <Grid/qcd/action/scalar/scalarImpl.h> #include <Grid/qcd/action/scalar/ScalarImpl.h>
#include <Grid/qcd/action/scalar/ScalarAction.h> #include <Grid/qcd/action/scalar/ScalarAction.h>
namespace Grid { namespace Grid {

View File

@ -0,0 +1,144 @@
/*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/qcd/action/gauge/GaugeImpl.h
Copyright (C) 2015
Author: paboyle <paboyle@ph.ed.ac.uk>
This program 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 2 of the License, or
(at your option) any later version.
This program 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, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution
directory
*************************************************************************************/
/* END LEGAL */
#ifndef GRID_GAUGE_IMPL_TYPES_H
#define GRID_GAUGE_IMPL_TYPES_H
namespace Grid {
namespace QCD {
////////////////////////////////////////////////////////////////////////
// Implementation dependent gauge types
////////////////////////////////////////////////////////////////////////
#define INHERIT_GIMPL_TYPES(GImpl) \
typedef typename GImpl::Simd Simd; \
typedef typename GImpl::LinkField GaugeLinkField; \
typedef typename GImpl::Field GaugeField; \
typedef typename GImpl::SiteField SiteGaugeField; \
typedef typename GImpl::SiteLink SiteGaugeLink;
#define INHERIT_FIELD_TYPES(Impl) \
typedef typename Impl::Simd Simd; \
typedef typename Impl::SiteField SiteField; \
typedef typename Impl::Field Field;
// hardcodes the exponential approximation in the template
template <class S, int Nrepresentation = Nc, int Nexp = 12 > class GaugeImplTypes {
public:
typedef S Simd;
template <typename vtype> using iImplGaugeLink = iScalar<iScalar<iMatrix<vtype, Nrepresentation>>>;
template <typename vtype> using iImplGaugeField = iVector<iScalar<iMatrix<vtype, Nrepresentation>>, Nd>;
typedef iImplGaugeLink<Simd> SiteLink;
typedef iImplGaugeField<Simd> SiteField;
typedef Lattice<SiteLink> LinkField;
typedef Lattice<SiteField> Field;
// Guido: we can probably separate the types from the HMC functions
// this will create 2 kind of implementations
// probably confusing the users
// Now keeping only one class
// Move this elsewhere? FIXME
static inline void AddLink(Field &U, LinkField &W,
int mu) { // U[mu] += W
PARALLEL_FOR_LOOP
for (auto ss = 0; ss < U._grid->oSites(); ss++) {
U._odata[ss]._internal[mu] =
U._odata[ss]._internal[mu] + W._odata[ss]._internal;
}
}
///////////////////////////////////////////////////////////
// Move these to another class
// HMC auxiliary functions
static inline void generate_momenta(Field& P, GridParallelRNG& pRNG){
// specific for SU gauge fields
LinkField Pmu(P._grid);
Pmu = zero;
for (int mu = 0; mu < Nd; mu++) {
SU<Nrepresentation>::GaussianFundamentalLieAlgebraMatrix(pRNG, Pmu);
PokeIndex<LorentzIndex>(P, Pmu, mu);
}
}
static inline Field projectForce(Field& P){
return Ta(P);
}
static inline void update_field(Field& P, Field& U, double ep){
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
auto Pmu = PeekIndex<LorentzIndex>(P, mu);
Umu = expMat(Pmu, ep, Nexp) * Umu;
PokeIndex<LorentzIndex>(U, ProjectOnGroup(Umu), mu);
}
}
static inline RealD FieldSquareNorm(Field& U){
LatticeComplex Hloc(U._grid);
Hloc = zero;
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
Hloc += trace(Umu * Umu);
}
Complex Hsum = sum(Hloc);
return Hsum.real();
}
static inline void HotConfiguration(GridParallelRNG &pRNG, Field &U) {
SU<Nc>::HotConfiguration(pRNG, U);
}
static inline void TepidConfiguration(GridParallelRNG &pRNG, Field &U) {
SU<Nc>::TepidConfiguration(pRNG, U);
}
static inline void ColdConfiguration(GridParallelRNG &pRNG, Field &U) {
SU<Nc>::ColdConfiguration(pRNG, U);
}
};
typedef GaugeImplTypes<vComplex, Nc> GimplTypesR;
typedef GaugeImplTypes<vComplexF, Nc> GimplTypesF;
typedef GaugeImplTypes<vComplexD, Nc> GimplTypesD;
typedef GaugeImplTypes<vComplex, SU<Nc>::AdjointDimension> GimplAdjointTypesR;
typedef GaugeImplTypes<vComplexF, SU<Nc>::AdjointDimension> GimplAdjointTypesF;
typedef GaugeImplTypes<vComplexD, SU<Nc>::AdjointDimension> GimplAdjointTypesD;
} // QCD
} // Grid
#endif // GRID_GAUGE_IMPL_TYPES_H

View File

@ -26,107 +26,14 @@ See the full license in the file "LICENSE" in the top level distribution
directory directory
*************************************************************************************/ *************************************************************************************/
/* END LEGAL */ /* END LEGAL */
#ifndef GRID_QCD_GAUGE_IMPL_H #ifndef GRID_QCD_GAUGE_IMPLEMENTATIONS_H
#define GRID_QCD_GAUGE_IMPL_H #define GRID_QCD_GAUGE_IMPLEMENTATIONS_H
#include "GaugeImplTypes.h"
namespace Grid { namespace Grid {
namespace QCD { namespace QCD {
////////////////////////////////////////////////////////////////////////
// Implementation dependent gauge types
////////////////////////////////////////////////////////////////////////
template <class Gimpl> class WilsonLoops;
//
#define INHERIT_GIMPL_TYPES(GImpl) \
typedef typename GImpl::Simd Simd; \
typedef typename GImpl::LinkField GaugeLinkField; \
typedef typename GImpl::Field GaugeField; \
typedef typename GImpl::SiteField SiteGaugeField; \
typedef typename GImpl::SiteLink SiteGaugeLink;
#define INHERIT_FIELD_TYPES(Impl) \
typedef typename Impl::Simd Simd; \
typedef typename Impl::SiteField SiteField; \
typedef typename Impl::Field Field;
// hard codes the exponential approximation in the template
template <class S, int Nrepresentation = Nc, int Nexp = 12 > class GaugeImplTypes {
public:
typedef S Simd;
template <typename vtype>
using iImplGaugeLink = iScalar<iScalar<iMatrix<vtype, Nrepresentation>>>;
template <typename vtype>
using iImplGaugeField = iVector<iScalar<iMatrix<vtype, Nrepresentation>>, Nd>;
typedef iImplGaugeLink<Simd> SiteLink;
typedef iImplGaugeField<Simd> SiteField;
typedef Lattice<SiteLink> LinkField;
typedef Lattice<SiteField> Field;
// Move this elsewhere? FIXME
static inline void AddLink(Field &U, LinkField &W,
int mu) { // U[mu] += W
PARALLEL_FOR_LOOP
for (auto ss = 0; ss < U._grid->oSites(); ss++) {
U._odata[ss]._internal[mu] =
U._odata[ss]._internal[mu] + W._odata[ss]._internal;
}
}
///////////////////////////////////////////////////////////
// Move these to another class
// HMC auxiliary functions
static inline void generate_momenta(Field& P, GridParallelRNG& pRNG){
// specific for SU gauge fields
LinkField Pmu(P._grid);
Pmu = zero;
for (int mu = 0; mu < Nd; mu++) {
SU<Nrepresentation>::GaussianFundamentalLieAlgebraMatrix(pRNG, Pmu);
PokeIndex<LorentzIndex>(P, Pmu, mu);
}
}
static inline Field projectForce(Field& P){
return Ta(P);
}
static inline void update_field(Field& P, Field& U, double ep){
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
auto Pmu = PeekIndex<LorentzIndex>(P, mu);
Umu = expMat(Pmu, ep, Nexp) * Umu;
PokeIndex<LorentzIndex>(U, ProjectOnGroup(Umu), mu);
}
}
static inline RealD FieldSquareNorm(Field& U){
LatticeComplex Hloc(U._grid);
Hloc = zero;
for (int mu = 0; mu < Nd; mu++) {
auto Umu = PeekIndex<LorentzIndex>(U, mu);
Hloc += trace(Umu * Umu);
}
Complex Hsum = sum(Hloc);
return Hsum.real();
}
static inline void HotConfiguration(GridParallelRNG &pRNG, Field &U) {
SU<Nc>::HotConfiguration(pRNG, U);
}
static inline void TepidConfiguration(GridParallelRNG &pRNG, Field &U) {
SU<Nc>::TepidConfiguration(pRNG, U);
}
static inline void ColdConfiguration(GridParallelRNG &pRNG, Field &U) {
SU<Nc>::ColdConfiguration(pRNG, U);
}
};
// Composition with smeared link, bc's etc.. probably need multiple inheritance // Composition with smeared link, bc's etc.. probably need multiple inheritance
// Variable precision "S" and variable Nc // Variable precision "S" and variable Nc
template <class GimplTypes> class PeriodicGaugeImpl : public GimplTypes { template <class GimplTypes> class PeriodicGaugeImpl : public GimplTypes {
@ -222,14 +129,6 @@ public:
static inline bool isPeriodicGaugeField(void) { return false; } static inline bool isPeriodicGaugeField(void) { return false; }
}; };
typedef GaugeImplTypes<vComplex, Nc> GimplTypesR;
typedef GaugeImplTypes<vComplexF, Nc> GimplTypesF;
typedef GaugeImplTypes<vComplexD, Nc> GimplTypesD;
typedef GaugeImplTypes<vComplex, SU<Nc>::AdjointDimension> GimplAdjointTypesR;
typedef GaugeImplTypes<vComplexF, SU<Nc>::AdjointDimension> GimplAdjointTypesF;
typedef GaugeImplTypes<vComplexD, SU<Nc>::AdjointDimension> GimplAdjointTypesD;
typedef PeriodicGaugeImpl<GimplTypesR> PeriodicGimplR; // Real.. whichever prec typedef PeriodicGaugeImpl<GimplTypesR> PeriodicGimplR; // Real.. whichever prec
typedef PeriodicGaugeImpl<GimplTypesF> PeriodicGimplF; // Float typedef PeriodicGaugeImpl<GimplTypesF> PeriodicGimplF; // Float
typedef PeriodicGaugeImpl<GimplTypesD> PeriodicGimplD; // Double typedef PeriodicGaugeImpl<GimplTypesD> PeriodicGimplD; // Double
@ -241,6 +140,8 @@ typedef PeriodicGaugeImpl<GimplAdjointTypesD> PeriodicGimplAdjD; // Double
typedef ConjugateGaugeImpl<GimplTypesR> ConjugateGimplR; // Real.. whichever prec typedef ConjugateGaugeImpl<GimplTypesR> ConjugateGimplR; // Real.. whichever prec
typedef ConjugateGaugeImpl<GimplTypesF> ConjugateGimplF; // Float typedef ConjugateGaugeImpl<GimplTypesF> ConjugateGimplF; // Float
typedef ConjugateGaugeImpl<GimplTypesD> ConjugateGimplD; // Double typedef ConjugateGaugeImpl<GimplTypesD> ConjugateGimplD; // Double
} }
} }

View File

@ -29,51 +29,48 @@ Author: paboyle <paboyle@ph.ed.ac.uk>
directory directory
*************************************************************************************/ *************************************************************************************/
/* END LEGAL */ /* END LEGAL */
#ifndef SCALAR_ACTION_H #ifndef SCALAR_ACTION_H
#define SCALAR_ACTION_H #define SCALAR_ACTION_H
namespace Grid { namespace Grid {
namespace QCD { // FIXME drop the QCD namespace everywhere here
////////////////////////////////////////////////////////////////////////
// Wilson Gauge Action .. should I template the Nc etc..
////////////////////////////////////////////////////////////////////////
template <class Impl>
class ScalarAction : public Action<typename Impl::Field> {
public:
INHERIT_FIELD_TYPES(Impl);
private:
RealD mass_square;
RealD lambda;
public:
ScalarAction(RealD ms, RealD l) : mass_square(ms), lambda(l){};
virtual std::string action_name(){return "ScalarAction";}
virtual void refresh(const Field &U, template <class Impl>
GridParallelRNG &pRNG){}; // noop as no pseudoferms class ScalarAction : public QCD::Action<typename Impl::Field> {
public:
virtual RealD S(const Field &p) { INHERIT_FIELD_TYPES(Impl);
return (mass_square * 0.5 + Nd) * ScalarObs<Impl>::sumphisquared(p) +
(lambda / 24.) * ScalarObs<Impl>::sumphifourth(p) + private:
ScalarObs<Impl>::sumphider(p); RealD mass_square;
RealD lambda;
public:
ScalarAction(RealD ms, RealD l) : mass_square(ms), lambda(l){};
virtual std::string action_name(){return "ScalarAction";}
virtual void refresh(const Field &U,
GridParallelRNG &pRNG){}; // noop as no pseudoferms
virtual RealD S(const Field &p) {
return (mass_square * 0.5 + QCD::Nd) * ScalarObs<Impl>::sumphisquared(p) +
(lambda / 24.) * ScalarObs<Impl>::sumphifourth(p) +
ScalarObs<Impl>::sumphider(p);
};
virtual void deriv(const Field &p,
Field &force) {
Field tmp(p._grid);
Field p2(p._grid);
ScalarObs<Impl>::phisquared(p2, p);
tmp = -(Cshift(p, 0, -1) + Cshift(p, 0, 1));
for (int mu = 1; mu < QCD::Nd; mu++) tmp -= Cshift(p, mu, -1) + Cshift(p, mu, 1);
force=+(mass_square + 2. * QCD::Nd) * p + (lambda / 6.) * p2 * p + tmp;
};
}; };
} // Grid
virtual void deriv(const Field &p, #endif // SCALAR_ACTION_H
Field &force) {
Field tmp(p._grid);
Field p2(p._grid);
ScalarObs<Impl>::phisquared(p2, p);
tmp = -(Cshift(p, 0, -1) + Cshift(p, 0, 1));
for (int mu = 1; mu < Nd; mu++) tmp -= Cshift(p, mu, -1) + Cshift(p, mu, 1);
force=+(mass_square + 2. * Nd) * p + (lambda / 6.) * p2 * p + tmp;
};
};
}
}
#endif

View File

@ -0,0 +1,57 @@
#ifndef SCALAR_IMPL
#define SCALAR_IMPL
namespace Grid {
//namespace QCD {
template <class S>
class ScalarImplTypes {
public:
typedef S Simd;
template <typename vtype>
using iImplField = iScalar<iScalar<iScalar<vtype> > >;
typedef iImplField<Simd> SiteField;
typedef Lattice<SiteField> Field;
static inline void generate_momenta(Field& P, GridParallelRNG& pRNG){
gaussian(pRNG, P);
}
static inline Field projectForce(Field& P){return P;}
static inline void update_field(Field& P, Field& U, double ep){
U += P*ep;
}
static inline RealD FieldSquareNorm(Field& U){
return (- sum(trace(U*U))/2.0);
}
static inline void HotConfiguration(GridParallelRNG &pRNG, Field &U) {
gaussian(pRNG, U);
}
static inline void TepidConfiguration(GridParallelRNG &pRNG, Field &U) {
gaussian(pRNG, U);
}
static inline void ColdConfiguration(GridParallelRNG &pRNG, Field &U) {
U = 1.0;
}
};
typedef ScalarImplTypes<vReal> ScalarImplR;
typedef ScalarImplTypes<vRealF> ScalarImplF;
typedef ScalarImplTypes<vRealD> ScalarImplD;
//}
}
#endif

View File

@ -1,56 +0,0 @@
#ifndef SCALAR_IMPL
#define SCALAR_IMPL
namespace Grid {
namespace QCD {
template <class S>
class ScalarImplTypes {
public:
typedef S Simd;
template <typename vtype>
using iImplField = iScalar<iScalar<iScalar<vtype> > >;
typedef iImplField<Simd> SiteField;
typedef Lattice<SiteField> Field;
static inline void generate_momenta(Field& P, GridParallelRNG& pRNG){
gaussian(pRNG, P);
}
static inline Field projectForce(Field& P){return P;}
static inline void update_field(Field& P, Field& U, double ep){
U += P*ep;
}
static inline RealD FieldSquareNorm(Field& U){
return (- sum(trace(U*U))/2.0);
}
static inline void HotConfiguration(GridParallelRNG &pRNG, Field &U) {
gaussian(pRNG, U);
}
static inline void TepidConfiguration(GridParallelRNG &pRNG, Field &U) {
gaussian(pRNG, U);
}
static inline void ColdConfiguration(GridParallelRNG &pRNG, Field &U) {
U = 1.0;
}
};
typedef ScalarImplTypes<vReal> ScalarImplR;
typedef ScalarImplTypes<vRealF> ScalarImplF;
typedef ScalarImplTypes<vRealD> ScalarImplD;
}}
#endif

View File

@ -1,191 +0,0 @@
/*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/qcd/hmc/HmcRunner.h
Copyright (C) 2015
Author: paboyle <paboyle@ph.ed.ac.uk>
This program 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 2 of the License, or
(at your option) any later version.
This program 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, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution
directory
*************************************************************************************/
/* END LEGAL */
#ifndef HMC_RUNNER
#define HMC_RUNNER
namespace Grid {
namespace QCD {
// Class for HMC specific for gauge theories
template <class Gimpl, class RepresentationsPolicy = NoHirep >
class NerscHmcRunnerTemplate {
public:
INHERIT_GIMPL_TYPES(Gimpl);
enum StartType_t { ColdStart, HotStart, TepidStart, CheckpointStart };
ActionSet<GaugeField, RepresentationsPolicy> TheAction;
GridCartesian *UGrid;
GridCartesian *FGrid;
GridRedBlackCartesian *UrbGrid;
GridRedBlackCartesian *FrbGrid;
virtual void BuildTheAction(int argc, char **argv) = 0; // necessary?
void Run(int argc, char **argv) {
StartType_t StartType = HotStart;
std::string arg;
if (GridCmdOptionExists(argv, argv + argc, "--StartType")) {
arg = GridCmdOptionPayload(argv, argv + argc, "--StartType");
if (arg == "HotStart") {
StartType = HotStart;
} else if (arg == "ColdStart") {
StartType = ColdStart;
} else if (arg == "TepidStart") {
StartType = TepidStart;
} else if (arg == "CheckpointStart") {
StartType = CheckpointStart;
} else {
std::cout << GridLogError << "Unrecognized option in --StartType\n";
std::cout << GridLogError << "Valid [HotStart, ColdStart, TepidStart, CheckpointStart]\n";
assert(0);
}
}
int StartTraj = 0;
if (GridCmdOptionExists(argv, argv + argc, "--StartTrajectory")) {
arg = GridCmdOptionPayload(argv, argv + argc, "--StartTrajectory");
std::vector<int> ivec(0);
GridCmdOptionIntVector(arg, ivec);
StartTraj = ivec[0];
}
int NumTraj = 1;
if (GridCmdOptionExists(argv, argv + argc, "--Trajectories")) {
arg = GridCmdOptionPayload(argv, argv + argc, "--Trajectories");
std::vector<int> ivec(0);
GridCmdOptionIntVector(arg, ivec);
NumTraj = ivec[0];
}
int NumThermalizations = 10;
if (GridCmdOptionExists(argv, argv + argc, "--Thermalizations")) {
arg = GridCmdOptionPayload(argv, argv + argc, "--Thermalizations");
std::vector<int> ivec(0);
GridCmdOptionIntVector(arg, ivec);
NumThermalizations = ivec[0];
}
GridSerialRNG sRNG;
GridParallelRNG pRNG(UGrid);
LatticeGaugeField U(UGrid); // change this to an extended field (smearing class)?
std::vector<int> SerSeed({1, 2, 3, 4, 5});
std::vector<int> ParSeed({6, 7, 8, 9, 10});
// Create integrator, including the smearing policy
// Smearing policy, only defined for Nc=3
/*
std::cout << GridLogDebug << " Creating the Stout class\n";
double rho = 0.1; // smearing parameter, now hardcoded
int Nsmear = 1; // number of smearing levels
Smear_Stout<Gimpl> Stout(rho);
std::cout << GridLogDebug << " Creating the SmearedConfiguration class\n";
//SmearedConfiguration<Gimpl> SmearingPolicy(UGrid, Nsmear, Stout);
std::cout << GridLogDebug << " done\n";
*/
//////////////
NoSmearing<Gimpl> SmearingPolicy;
typedef MinimumNorm2<Gimpl, NoSmearing<Gimpl>, RepresentationsPolicy >
IntegratorType; // change here to change the algorithm
IntegratorParameters MDpar(40, 1.0);
IntegratorType MDynamics(UGrid, MDpar, TheAction, SmearingPolicy);
// Checkpoint strategy
NerscHmcCheckpointer<Gimpl> Checkpoint(std::string("ckpoint_lat"),
std::string("ckpoint_rng"), 1);
PlaquetteLogger<Gimpl> PlaqLog(std::string("plaq"));
HMCparameters HMCpar;
HMCpar.StartTrajectory = StartTraj;
HMCpar.Trajectories = NumTraj;
HMCpar.NoMetropolisUntil = NumThermalizations;
if (StartType == HotStart) {
// Hot start
HMCpar.MetropolisTest = true;
sRNG.SeedFixedIntegers(SerSeed);
pRNG.SeedFixedIntegers(ParSeed);
SU<Nc>::HotConfiguration(pRNG, U);
} else if (StartType == ColdStart) {
// Cold start
HMCpar.MetropolisTest = true;
sRNG.SeedFixedIntegers(SerSeed);
pRNG.SeedFixedIntegers(ParSeed);
SU<Nc>::ColdConfiguration(pRNG, U);
} else if (StartType == TepidStart) {
// Tepid start
HMCpar.MetropolisTest = true;
sRNG.SeedFixedIntegers(SerSeed);
pRNG.SeedFixedIntegers(ParSeed);
SU<Nc>::TepidConfiguration(pRNG, U);
} else if (StartType == CheckpointStart) {
HMCpar.MetropolisTest = true;
// CheckpointRestart
Checkpoint.CheckpointRestore(StartTraj, U, sRNG, pRNG);
}
// Attach the gauge field to the smearing Policy and create the fill the
// smeared set
// notice that the unit configuration is singular in this procedure
std::cout << GridLogMessage << "Filling the smeared set\n";
SmearingPolicy.set_Field(U);
HybridMonteCarlo<IntegratorType> HMC(HMCpar, MDynamics, sRNG, pRNG, U);
HMC.AddObservable(&Checkpoint);
HMC.AddObservable(&PlaqLog);
// Run it
HMC.evolve();
}
};
typedef NerscHmcRunnerTemplate<PeriodicGimplR> NerscHmcRunner;
typedef NerscHmcRunnerTemplate<PeriodicGimplF> NerscHmcRunnerF;
typedef NerscHmcRunnerTemplate<PeriodicGimplD> NerscHmcRunnerD;
typedef NerscHmcRunnerTemplate<PeriodicGimplR> PeriodicNerscHmcRunner;
typedef NerscHmcRunnerTemplate<PeriodicGimplF> PeriodicNerscHmcRunnerF;
typedef NerscHmcRunnerTemplate<PeriodicGimplD> PeriodicNerscHmcRunnerD;
typedef NerscHmcRunnerTemplate<ConjugateGimplR> ConjugateNerscHmcRunner;
typedef NerscHmcRunnerTemplate<ConjugateGimplF> ConjugateNerscHmcRunnerF;
typedef NerscHmcRunnerTemplate<ConjugateGimplD> ConjugateNerscHmcRunnerD;
template <class RepresentationsPolicy>
using NerscHmcRunnerHirep = NerscHmcRunnerTemplate<PeriodicGimplR, RepresentationsPolicy>;
}
}
#endif

View File

@ -0,0 +1,120 @@
/*************************************************************************************
Grid physics library, www.github.com/paboyle/Grid
Source file: ./lib/qcd/modules/Registration.h
Copyright (C) 2016
Author: Guido Cossu <guido.cossu@ed.ac.uk>
This program 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 2 of the License, or
(at your option) any later version.
This program 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, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
See the full license in the file "LICENSE" in the top level distribution
directory
*************************************************************************************/
/* END LEGAL */
#ifndef MODULES_REGISTRATION_H
#define MODULES_REGISTRATION_H
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Actions
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
static Registrar<QCD::WilsonGMod, HMC_LGTActionModuleFactory<gauge_string, Serialiser > > __WGmodXMLInit("Wilson");
static Registrar<QCD::SymanzikGMod, HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __SymGmodXMLInit("Symanzik");
static Registrar<QCD::IwasakiGMod, HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __IwGmodXMLInit("Iwasaki");
static Registrar<QCD::DBW2GMod, HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __DBW2GmodXMLInit("DBW2");
static Registrar<QCD::RBCGMod, HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __RBCGmodXMLInit("RBC");
static Registrar<QCD::PlaqPlusRectangleGMod, HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __PPRectGmodXMLInit("PlaqPlusRect");
// FIXME more general implementation
static Registrar<QCD::TwoFlavourFModule<FermionImplementationPolicy> ,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __TwoFlavourFmodXMLInit("TwoFlavours");
static Registrar<QCD::TwoFlavourRatioFModule<FermionImplementationPolicy> ,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __TwoFlavourRatioFmodXMLInit("TwoFlavoursRatio");
static Registrar<QCD::TwoFlavourEOFModule<FermionImplementationPolicy> ,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __TwoFlavourEOFmodXMLInit("TwoFlavoursEvenOdd");
static Registrar<QCD::TwoFlavourRatioEOFModule<FermionImplementationPolicy>,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __TwoFlavourRatioEOFmodXMLInit("TwoFlavoursEvenOddRatio");
static Registrar<QCD::OneFlavourFModule<FermionImplementationPolicy> ,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __OneFlavourFmodXMLInit("OneFlavour");
static Registrar<QCD::OneFlavourEOFModule<FermionImplementationPolicy> ,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __OneFlavourEOFmodXMLInit("OneFlavourEvenOdd");
static Registrar<QCD::OneFlavourRatioFModule<FermionImplementationPolicy> ,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __OneFlavourRatioFmodXMLInit("OneFlavourRatio");
static Registrar<QCD::OneFlavourRatioEOFModule<FermionImplementationPolicy>,
HMC_LGTActionModuleFactory<gauge_string, Serialiser> > __OneFlavourRatioEOFmodXMLInit("OneFlavourEvenOddRatio");
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Solvers
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Now a specific registration with a fermion field
// here must instantiate CG and CR for every new fermion field type (macro!!)
static Registrar< ConjugateGradientModule<QCD::WilsonFermionR::FermionField>,
HMC_SolverModuleFactory<solver_string, QCD::WilsonFermionR::FermionField, Serialiser> > __CGWFmodXMLInit("ConjugateGradientWF");
static Registrar< ConjugateResidualModule<QCD::WilsonFermionR::FermionField>,
HMC_SolverModuleFactory<solver_string, QCD::WilsonFermionR::FermionField, Serialiser> > __CRWFmodXMLInit("ConjugateResidualWF");
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Fermion operators
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
static Registrar< QCD::WilsonFermionModule<FermionImplementationPolicy>,
HMC_FermionOperatorModuleFactory<fermionop_string, FermionImplementationPolicy, Serialiser> > __WilsonFOPmodXMLInit("Wilson");
static Registrar< QCD::MobiusFermionModule<FermionImplementationPolicy>,
HMC_FermionOperatorModuleFactory<fermionop_string, FermionImplementationPolicy, Serialiser> > __MobiusFOPmodXMLInit("Mobius");
static Registrar< QCD::DomainWallFermionModule<FermionImplementationPolicy>,
HMC_FermionOperatorModuleFactory<fermionop_string, FermionImplementationPolicy, Serialiser> > __DWFOPmodXMLInit("DomainWall");
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Observables
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
static Registrar<QCD::PlaquetteMod<ImplementationPolicy>, HMC_ObservablesModuleFactory<observable_string, Serialiser> > __OBSPLmodXMLInit("Plaquette");
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Checkpointers
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
static Registrar<QCD::BinaryCPModule<ImplementationPolicy>, HMC_CPModuleFactory<cp_string, ImplementationPolicy, Serialiser> > __CPBinarymodXMLInit("Binary");
static Registrar<QCD::NerscCPModule<ImplementationPolicy> , HMC_CPModuleFactory<cp_string, ImplementationPolicy, Serialiser> > __CPNerscmodXMLInit("Nersc");
#ifdef HAVE_LIME
static Registrar<QCD::ILDGCPModule<ImplementationPolicy> , HMC_CPModuleFactory<cp_string, ImplementationPolicy, Serialiser> > __CPILDGmodXMLInit("ILDG");
#endif
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Integrators
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
static Registrar< HMCLeapFrog<ImplementationPolicy, RepresentationPolicy, Serialiser> , HMCRunnerModuleFactory<hmc_string, Serialiser> > __HMCLFmodXMLInit("LeapFrog");
static Registrar< HMCMinimumNorm2<ImplementationPolicy, RepresentationPolicy, Serialiser> , HMCRunnerModuleFactory<hmc_string, Serialiser> > __HMCMN2modXMLInit("MinimumNorm2");
static Registrar< HMCForceGradient<ImplementationPolicy, RepresentationPolicy, Serialiser> , HMCRunnerModuleFactory<hmc_string, Serialiser> > __HMCFGmodXMLInit("ForceGradient");
typedef HMCRunnerModuleFactory<hmc_string, Serialiser > HMCModuleFactory;
#endif

View File

@ -4,7 +4,7 @@
#include <Grid/qcd/representations/adjoint.h> #include <Grid/qcd/representations/adjoint.h>
#include <Grid/qcd/representations/two_index.h> #include <Grid/qcd/representations/two_index.h>
#include <Grid/qcd/representations/fundamental.h> #include <Grid/qcd/representations/fundamental.h>
#include <Grid/qcd/action/scalar/scalarImpl.h> #include <Grid/qcd/action/scalar/ScalarImpl.h>
#include <tuple> #include <tuple>
#include <utility> #include <utility>

View File

@ -6,10 +6,7 @@
Copyright (C) 2015 Copyright (C) 2015
Author: Azusa Yamaguchi <ayamaguc@staffmail.ed.ac.uk>
Author: Peter Boyle <paboyle@ph.ed.ac.uk>
Author: neo <cossu@post.kek.jp> Author: neo <cossu@post.kek.jp>
Author: paboyle <paboyle@ph.ed.ac.uk>
This program is free software; you can redistribute it and/or modify This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by it under the terms of the GNU General Public License as published by
@ -32,7 +29,9 @@ directory
#ifndef SCALAR_OBJS_H #ifndef SCALAR_OBJS_H
#define SCALAR_OBJS_H #define SCALAR_OBJS_H
namespace Grid { namespace Grid {
namespace QCD {
// FIXME drop the QCD namespace in Nd
// Scalar field obs // Scalar field obs
template <class Impl> template <class Impl>
@ -62,7 +61,7 @@ class ScalarObs {
static void phider(typename Impl::Field &fsq, static void phider(typename Impl::Field &fsq,
const typename Impl::Field &f) { const typename Impl::Field &f) {
fsq = Cshift(f, 0, -1) * f; fsq = Cshift(f, 0, -1) * f;
for (int mu = 1; mu < Nd; mu++) fsq += Cshift(f, mu, -1) * f; for (int mu = 1; mu < QCD::Nd; mu++) fsq += Cshift(f, mu, -1) * f;
} }
////////////////////////////////////////////////// //////////////////////////////////////////////////
@ -72,7 +71,7 @@ class ScalarObs {
static RealD sumphider(const typename Impl::Field &f) { static RealD sumphider(const typename Impl::Field &f) {
typename Impl::Field tmp(f._grid); typename Impl::Field tmp(f._grid);
tmp = Cshift(f, 0, -1) * f; tmp = Cshift(f, 0, -1) * f;
for (int mu = 1; mu < Nd; mu++) { for (int mu = 1; mu < QCD::Nd; mu++) {
tmp += Cshift(f, mu, -1) * f; tmp += Cshift(f, mu, -1) * f;
} }
return -sum(trace(tmp)); return -sum(trace(tmp));
@ -90,7 +89,8 @@ class ScalarObs {
return sum(trace(tmp)); return sum(trace(tmp));
} }
}; };
}
} }
#endif #endif

View File

@ -31,36 +31,35 @@ directory
namespace Grid{ namespace Grid{
// Put this section in a separate header // Put this section in a separate header
// ifdefs ?? Local makefile suggestion , policy as make parameter // ifdefs ?? Local makefile suggestion , policy as make parameter
typedef QCD::PeriodicGimplR ImplementationPolicy; typedef QCD::PeriodicGimplR ImplementationPolicy;
typedef QCD::NoHirep RepresentationPolicy; typedef QCD::WilsonImplR FermionImplementationPolicy;
typedef QCD::WilsonFermionR FermionImplementation; typedef QCD::NoHirep RepresentationPolicy;
typedef Grid::XmlReader Serialiser;
/////////////////////////////////////////////////////////////////////// // Register all object names
// Put all registrations in an header file #include "Grid/qcd/modules/Registration.h"
static Registrar< HMCLeapFrog<ImplementationPolicy, RepresentationPolicy, XmlReader> , HMCRunnerModuleFactory<hmc_string, XmlReader> > __HMCLFmodXMLInit("LeapFrog");
static Registrar< HMCMinimumNorm2<ImplementationPolicy, RepresentationPolicy, XmlReader> , HMCRunnerModuleFactory<hmc_string, XmlReader> > __HMCMN2modXMLInit("MinimumNorm2"); } // Grid
static Registrar< HMCForceGradient<ImplementationPolicy, RepresentationPolicy, XmlReader> , HMCRunnerModuleFactory<hmc_string, XmlReader> > __HMCFGmodXMLInit("ForceGradient");
}
int main(int argc, char **argv) { int main(int argc, char **argv) {
using namespace Grid; using namespace Grid;
using namespace Grid::QCD; using namespace Grid::QCD;
Grid_init(&argc, &argv); Grid_init(&argc, &argv);
int threads = GridThread::GetThreads(); int threads = GridThread::GetThreads();
// here make a routine to print all the relevant information on the run // here make a routine to print all the relevant information on the run
std::cout << GridLogMessage << "Grid is setup to use " << threads << " threads" << std::endl; std::cout << GridLogMessage << "Grid is setup to use " << threads << " threads" << std::endl;
// Typedefs to simplify notation // Typedefs to simplify notation
typedef Grid::XmlReader InputFileReader; //typedef XmlReader InputFileReader;
// Reader, file should come from command line // Reader, file should come from command line
InputFileReader Reader("input.wilson_gauge.params.xml"); Serialiser Reader("input.wilson_gauge.params.xml");
// Test HMC factory (put in an external file) // Test HMC factory (put in an external file)
auto &HMCfactory = HMCRunnerModuleFactory<hmc_string, InputFileReader >::getInstance(); auto &HMCfactory = HMCModuleFactory::getInstance();
// Simplify this step (IntergratorName field?) // Simplify this step (IntergratorName field?)
HMCparameters HMCpar(Reader); HMCparameters HMCpar(Reader);
@ -71,5 +70,7 @@ int main(int argc, char **argv) {
HMCmodule->getPtr()->Run(); HMCmodule->getPtr()->Run();
Grid_finalize(); Grid_finalize();
return 0;
} // main } // main