Skip to content
KERNEL_ARCHS = #KERNEL_ARCHS#
qcd-bench: qcd-bench.o $(KERNEL_ARCHS)
#LD# #LDFLAGS# -o #EXECNAME# qcd-bench.o $(KERNEL_ARCHS) #LDLIBS#
kernel_A.a:
cd kernel_A && gmake prep-#KA_PLATFORM# && gmake kernel
kernel_B.a:
cd kernel_B && gmake kernel
kernel_C.a:
cd kernel_C && make kernel
kernel_D.a:
cd kernel_D && make kernel
kernel_E.a:
cd kernel_E && gmake kernel
qcd-bench.o: qcd-bench.c
#MPI_CC# #CFLAGS# -c -o qcd-bench.o qcd-bench.c
clean:
cd kernel_A && gmake clobber
cd kernel_B && make clean
cd kernel_C && make clean
cd kernel_D && make clean
cd kernel_E && gmake clean
#RM# -f qcd-bench.o qcd-bench
#RM# -f $(KERNEL_ARCHS)
#nm -X64 *.a | grep ' T ' | cut -f 1 -d ' ' | sort > all.sort.dat
#nm -X64 *.a | grep ' T ' | cut -f 1 -d ' ' | sort | sort -u > all.sort.unique.dat
#diff all.sort.dat all.sort.unique.dat | grep '<' > doublicates.dat
\ No newline at end of file
replaceInFile ()
{
sed "s/${1}(/${2}(/g" $3 > tmp.dat; cp tmp.dat $3
}
case $1 in
"doublicates") echo "generate doublicate list";
echo "parse files:\n`ls *.a`";
nm -X64 *.a | grep ' T ' | cut -f 1 -d ' ' | sort > all.sort.dat
nm -X64 *.a | grep ' T ' | cut -f 1 -d ' ' | sort | sort -u > all.sort.unique.dat
diff all.sort.dat all.sort.unique.dat | grep '<' | cut -f 2- -d '.'> doublicates.dat
echo "found `wc -l doublicates.dat` doublicates" ;;
"find") echo "looking for $2";
grep -r $2 `find . -name '*.[c|h|f|f90|F90]'` | cut -f 2- -d '<' > find.dat;
cat find.dat;;
"replace") echo "replace $2 by $3 in $4";
replaceInFile $2 $3 $4;;
"replaceAllFiles") echo "replace all $2 by $3";
grep -r $2 `find . -name '*.[c|h|f|f90|F90]'` | cut -f 2- -d '<' | cut -f 1 -d ':' | sort -u > find.dat;
for i in `cat find.dat`
do
echo "replacing $2 by $3 in $i";
replaceInFile $2 $3 $i
done
;;
"replaceAll")
echo "replace all doublicates in this directory, using doublicate list $2 and postfix $3";
for dn in `cat $2`
do
echo "replace $dn"
grep -r $dn `find . -name '*.[c|h|f|f90|F90]'` | cut -f 2- -d '<' | cut -f 1 -d ':' | sort -u > find.dat;
for fn in `cat find.dat`
do
echo "replacing $dn by ${dn}_$3 in $fn";
replaceInFile ${dn} ${dn}_${3} ${fn}
done
done
;;
esac
\ No newline at end of file
#===============================================================================
#
# BQCD -- Berlin Quantum ChromoDynamics programme
#
# Author: Hinnerk Stueben <stueben@zib.de>
#
# Copyright (C) 1998-2006, Hinnerk Stueben, Zuse-Institut Berlin
#
#-------------------------------------------------------------------------------
#
# Makefile
#
#===============================================================================
include Makefile.defs
MODULES_DIR = modules
.SUFFIXES:
.SUFFIXES: .o .F90 .c
.F90.f90:
$(FPP) $(FPPFLAGS) $< > $@
.F90.o:
$(FPP) $(FPPFLAGS) $< > $*.f90
$(F90) -c $(FFLAGS) -I$(MODULES_DIR) $*.f90
MODULES = $(MODULES_DIR)/*.o
OBJS = \
action.o \
cg.o \
checks.o \
$(CKSUM_O) \
conf.o \
conf_info.o \
cooling.o \
dsd.o \
dsg.o \
dsf.o \
dsf1.o \
dsf2.o \
files.o \
flip_bc.o \
hmc.o \
hmc_init_p.o \
hmc_init_phi.o \
hmc_integrator.o \
hmc_forces.o \
hmc_leap_frog.o \
hmc_test.o \
hmc_u.o \
h_mult.o \
index.o \
index2.o \
init_common.o \
init_modules.o \
iteration_count.o \
mc.o \
misc.o \
mre.o \
mtdagmt.o \
m_tilde.o \
polyakov_loop.o \
$(RANDOM_O) \
sc.o \
service.o \
staple.o \
su3.o \
swap.o \
timing.o \
traces.o \
$(UUU_O) \
w_mult.o \
xyzt2i.o
LIBS = d/$(LIBD) comm/$(LIBCOMM) clover/$(LIBCLOVER)
kernel:
cd modules && $(MAKE)
cd d && $(MAKE) fast
cd comm && $(MAKE) fast
cd clover && $(MAKE) fast
$(MAKE) ../kernel_A.a
../kernel_A.a: bqcd.o $(MODULES) $(OBJS) $(LIBS)
$(AR) $(ARFLAGS) ../kernel_A.a *.o d/*.o modules/*.o comm/*.o clover/*.o
bqcd: bqcd.o $(MODULES) $(OBJS) $(LIBS)
$(F90) -o $@ $(LDFLAGS) bqcd.o $(MODULES) $(OBJS) $(LIBS) $(SYSLIBS)
fast:
cd modules && $(MAKE)
cd d && $(MAKE) fast
cd comm && $(MAKE) fast
cd clover && $(MAKE) fast
$(FAST_MAKE) bqcd
clean:
rm -f bqcd.[0-9][0-9][0-9].* diag.[0-9][0-9] core app.rif
rm -f random_test random_test.dump random_test.out
rm -f test_echo
rm -f a.out out out1 out2
rm -f ../kernel-bqcd.a
tidy: clean
rm -f *.[Toid] *.f90 *.mod work.pc work.pcl
clobber: tidy
rm -f bqcd
$(MAKE) clobber_libs
cd modules && $(MAKE) clean
Modules:
cd modules && $(MAKE)
libd:
cd d && $(MAKE)
libclover: $(MODULES)
cd clover && $(MAKE)
libs:
cd d && $(MAKE)
cd comm && $(MAKE)
cd clover && $(MAKE)
clean_libs:
cd d && $(MAKE) clean
cd comm && $(MAKE) clean
cd clover && $(MAKE) clean
clobber_libs:
cd d && $(MAKE) clobber
cd comm && $(MAKE) clobber
cd clover && $(MAKE) clobber
the_ranf_test: ranf.o
$(FPP) $(FPPFLAGS) ranf_test.F90 ranf_test.f90
$(F90) ranf_test.f90 ranf.o
./a.out | diff - ranf_test.reference
test_echo: test_echo.o service.o
$(F90) -o $@ $(LDFLAGS) test_echo.o service.o
prep:
rm -f Makefile.var service.F90
ln -s platform/Makefile-$(PLATFORM).var Makefile.var
ln -s platform/service-$(PLATFORM).F90 service.F90
prep-altix:
$(MAKE) prep PLATFORM=altix
prep-bgl:
$(MAKE) prep PLATFORM=bgl
prep-cray:
$(MAKE) prep PLATFORM=cray
prep-hitachi:
$(MAKE) prep PLATFORM=hitachi
prep-hitachi-omp:
$(MAKE) prep PLATFORM=hitachi-omp
prep-ibm:
$(MAKE) prep PLATFORM=ibm
prep-hp:
$(MAKE) prep PLATFORM=hp
prep-intel:
$(MAKE) prep PLATFORM=intel
prep-nec:
$(MAKE) prep PLATFORM=nec
prep-sun:
$(MAKE) prep PLATFORM=sun
#===============================================================================
d3_buffer_vol = #d3_buffer_vol#
SHELL = #SHELL#
FPP = #FPP#
FPPFLAGS = #FPPFLAGS#
F90 = #MPI_F90#
FFLAGS = #F90FLAGS#
CC = #MPI_CC#
CFLAGS = #CFLAGS#
AR = #AR#
ARFLAGS = #ARFLAGS#
RANLIB = #RANLIB#
LDFLAGS = #LDFLAGS#
SYSLIBS = #SYSLIBS#
FAST_MAKE = #FAST_MAKE#
RM = #RM#
CKSUM_O = #CKSUM_O#
RANDOM_O = #RANDOM_O#
UUU_O = #UUU_O#
LIBD = #KA_LIBD#
LIBCOMM = #KA_LIBCOMM#
LIBCLOVER = #KA_LIBCLOVER#
\ No newline at end of file
#===============================================================================
#
# BQCD -- Berlin Quantum ChromoDynamics programme
#
# Author: Hinnerk Stueben <stueben@zib.de>
#
# Copyright (C) 1998-2006, Hinnerk Stueben, Zuse-Institut Berlin
#
#-------------------------------------------------------------------------------
#
# Makefile
#
#===============================================================================
#include Makefile.var
include Makefile.defs
.SUFFIXES:
.SUFFIXES: .o .F90 .c
.F90.f90:
$(FPP) $(FPPFLAGS) $< > $@
.F90.o:
$(FPP) $(FPPFLAGS) $< > $*.f90
$(F90) -c $(FFLAGS) $*.f90
MODULES_DIR = modules
MODULES = modules/*.o
OBJS = \
action.o \
cg.o \
checks.o \
$(CKSUM_O) \
conf.o \
conf_info.o \
cooling.o \
dsd.o \
dsg.o \
dsf.o \
dsf1.o \
dsf2.o \
files.o \
flip_bc.o \
hmc.o \
hmc_init_p.o \
hmc_init_phi.o \
hmc_integrator.o \
hmc_forces.o \
hmc_leap_frog.o \
hmc_test.o \
hmc_u.o \
h_mult.o \
index.o \
index2.o \
init_common.o \
init_modules.o \
iteration_count.o \
mc.o \
misc.o \
mre.o \
mtdagmt.o \
m_tilde.o \
polyakov_loop.o \
$(RANDOM_O) \
sc.o \
service.o \
staple.o \
su3.o \
swap.o \
timing.o \
traces.o \
$(UUU_O) \
w_mult.o \
xyzt2i.o
LIBS = d/$(LIBD) comm/$(LIBCOMM) clover/$(LIBCLOVER)
bqcd: bqcd.o $(MODULES) $(OBJS) $(LIBS)
$(F90) -o $@ $(LDFLAGS) bqcd.o $(MODULES) $(OBJS) $(LIBS) $(SYSLIBS)
fast:
cd modules && $(MAKE)
cd d && $(MAKE) fast
cd comm && $(MAKE) fast
cd clover && $(MAKE) fast
$(FAST_MAKE) bqcd
mv bqcd #EXECNAME#
# mv bqcd ../BQCD_SGI-ALTIX_cname_SGI-ALTIX.exe
clean:
rm -f bqcd.[0-9][0-9][0-9].* diag.[0-9][0-9] core app.rif
rm -f random_test random_test.dump random_test.out
rm -f test_echo
rm -f a.out out out1 out2
tidy: clean
rm -f *.[Toid] *.f90 *.mod work.pc work.pcl
clobber: tidy
rm -f bqcd
$(MAKE) clobber_libs
cd modules && $(MAKE) clean
Modules:
cd modules && $(MAKE)
libd:
cd d && $(MAKE)
libclover: $(MODULES)
cd clover && $(MAKE)
libs:
cd d && $(MAKE)
cd comm && $(MAKE)
cd clover && $(MAKE)
clean_libs:
cd d && $(MAKE) clean
cd comm && $(MAKE) clean
cd clover && $(MAKE) clean
clobber_libs:
cd d && $(MAKE) clobber
cd comm && $(MAKE) clobber
cd clover && $(MAKE) clobber
the_ranf_test: ranf.o
$(FPP) $(FPPFLAGS) ranf_test.F90 ranf_test.f90
$(F90) ranf_test.f90 ranf.o
./a.out | diff - ranf_test.reference
test_echo: test_echo.o service.o
$(F90) -o $@ $(LDFLAGS) test_echo.o service.o
prep:
rm -f Makefile.var service.F90
ln -s platform/Makefile-$(PLATFORM).var Makefile.var
ln -s platform/service-$(PLATFORM).F90 service.F90
prep-altix:
$(MAKE) prep PLATFORM=altix
prep-bgl:
$(MAKE) prep PLATFORM=bgl
prep-cray:
$(MAKE) prep PLATFORM=cray
prep-hitachi:
$(MAKE) prep PLATFORM=hitachi
prep-hitachi-omp:
$(MAKE) prep PLATFORM=hitachi-omp
prep-ibm:
$(MAKE) prep PLATFORM=ibm
prep-hp:
$(MAKE) prep PLATFORM=hp
prep-intel:
$(MAKE) prep PLATFORM=intel
prep-nec:
$(MAKE) prep PLATFORM=nec
prep-sun:
$(MAKE) prep PLATFORM=sun
#===============================================================================
#===============================================================================
#
# BQCD -- Berlin Quantum ChromoDynamics programme
#
# Author: Hinnerk Stueben <stueben@zib.de>
#
# Copyright (C) 2005, Hinnerk Stueben, Zuse-Institut Berlin
#
#-------------------------------------------------------------------------------
#
# Makefile-altix.var - settings on SGI-Altix
#
#-------------------------------------------------------------------------------
timing = 1
mpi = 1
omp = 1
shmem =
shmempi =
debug =
libd = 2
d3_buffer_vol = 32*32*16*16
#-------------------------------------------------------------------------------
SHELL = /bin/ksh
FPP = mpif90 -g -E
FPP2 = icc -E -C -P
F90 = mpif90
CC = mpicc
AR = ar
RANLIB = echo
MODULES_FLAG = -I$(MODULES_DIR)
MYFLAGS = -DINTEL -DALTIX
FFLAGS_STD= $(MODULES_FLAG)
CFLAGS_STD= -DNamesToLower_
ARFLAGS = rv
LDFLAGS = -Vaxlib
SYSLIBS =
FAST_MAKE = gmake -j 8
CKSUM_O = cksum.o
RANDOM_O = ran.o ranf.o
UUU_O = uuu_f90.o
LIBD =
LIBCOMM = lib_single_pe.a
LIBCLOVER = libclover.a
#-------------------------------------------------------------------------------
ifdef timing
MYFLAGS += -DTIMING
endif
ifdef mpi
LIBCOMM = lib_mpi.a
endif
ifdef omp
F90 += -openmp
MYFLAGS += -D_OPENMP
endif
ifdef shmem
LDFLAGS += -lsma
LIBCOMM = lib_shmem.a
endif
ifdef shmempi
LDFLAGS += -lsma
LIBCOMM = lib_shmempi.a
endif
ifdef debug
FFLAGS = -g -O0 $(FFLAGS_STD)
CFLAGS = -g -O0 $(CFLAGS_STD)
else
FFLAGS = -O2 $(FFLAGS_STD)
CFLAGS = -O2 $(CFLAGS_STD)
endif
ifeq ($(libd),1)
LIBD = libd.a
MYFLAGS += -DD3_BUFFER_VOL=1
endif
ifeq ($(libd),2)
LIBD = libd2.a
MYFLAGS += -DD3_BUFFER_VOL=1
endif
ifeq ($(libd),21)
LIBD = libd21.a
MYFLAGS += -DD3_BUFFER_VOL=1
endif
ifeq ($(libd),3)
LIBD = libd3.a
MYFLAGS += -DD3_BUFFER_VOL='$(d3_buffer_vol)'
endif
#===============================================================================
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2005, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! action.F90 - calculation of actions
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
REAL function sf(para, conf) ! returns S_f
use typedef_hmc
use module_p_interface
use module_vol
implicit none
type(hmc_para), intent(in) :: para
type(hmc_conf), intent(in) :: conf
P_SPINCOL_FIELD, save :: a, b
REAL, external :: dotprod, clover_action
integer :: iterations
external :: mtdagmt
if (para%kappa == ZERO) then
sf = ZERO
else
ALLOCATE_SC_FIELD(a)
ALLOCATE_SC_FIELD(b)
call flip_bc(conf%u)
call sc_copy(a, conf%phi) ! A = phi
call cg(mtdagmt, a, conf%phi, para, conf, iterations) ! A = inv(M~+ M~) Phi
call mtil(b, a, para, conf) ! B = M~ A
sf = dotprod(b, b, SIZE_SC_FIELD)
call flip_bc(conf%u)
endif
if (para%csw_kappa /= ZERO) sf = sf + clover_action(conf%b(1,1,ODD))
end
!-------------------------------------------------------------------------------
REAL function sg(u) ! returns S_g
use module_nn
use module_vol
implicit none
GAUGE_FIELD :: u
REAL :: plaq, global_sum, p
SU3 :: uuu
integer :: i, e, o, mu, nu, j1, j2
REAL, external :: Re_Tr_uu
TIMING_START(timing_bin_plaq)
plaq = 0
do mu = 1, DIM
do e = EVEN, ODD
o = EVEN + ODD - e
do nu = mu + 1, DIM
p = ZERO
!$omp parallel do reduction(+: p) private(j1, j2, uuu)
do i = 1, VOLH
! (j2,o) --<-- x nu
! | |
! v ^ ^
! | | |
! (i,e) -->-- (j1,o) x--> mu
j1 = nn(i, e, mu, FWD)
j2 = nn(i, e, nu, FWD)
uuu = 0
call uuu_fwd(uuu, u(1, 1, j1, o, nu), &
u(1, 1, j2, o, mu), &
u(1, 1, i, e, nu))
p = p + Re_Tr_uu(uuu, u(1, 1, i, e, mu))
enddo
!$omp end parallel do
plaq = plaq + p
enddo
enddo
enddo
plaq = global_sum(plaq)
sg = (6 * volume) - plaq / THREE
TIMING_STOP(timing_bin_plaq)
end
!-------------------------------------------------------------------------------
REAL function sp(p) ! returns action of momenta p
use module_vol
implicit none
REAL, external :: dotprod
GENERATOR_FIELD :: p
integer :: mu, eo
sp = ZERO
do mu = 1, DIM
do eo = EVEN, ODD
sp = sp + dotprod(p(1, 1, eo, mu), p(1, 1, eo, mu), NGEN * volh)
enddo
enddo
sp = sp * HALF
end
!===============================================================================
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2005, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! bqcd.F90 - main program and read/write of parameters
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
! JuBE
! use kernel_a as a subroutine in the qcd-bench, this was the main function
! in the original code
subroutine kernel_a()
use typedef_flags
use typedef_para
use module_input
use module_function_decl
implicit none
type(type_para) :: para
type(hmc_conf), dimension(MAX_TEMPER) :: conf
type(type_flags) :: flags
SECONDS :: time0, sekunden
integer :: kernel_number
kernel_number = 0
! JuBE
! call jube initial function
call jube_kernel_init(kernel_number)
! JuBE
! set the flags%input to the inputfile name: bqcd-input
flags%input = "kernel_A.input"
time0 = sekunden() ! start/initialize timer
TIMING_START(timing_bin_total)
call comm_init()
! JuBE
! there is no need for the following call, we ignore all cmd line arguments, non
! of them but the input file name (set above) is relevant for the benchmark
! call get_flags(flags)
call begin(UREC, "Job")
call input_read(flags%input)
call init_para(para, flags)
call init_counter(para, flags)
call init_ran(para, flags)
call init_cooling(input%measure_cooling_list)
call set_fmt_ensemble(para%n_temper)
call check_fmt(para%run, para%n_temper, para%maxtraj, para%L(4) - 1)
call init_common(para)
call init_modules()
call write_header(para)
call init_flip_bc()
call init_cg_para(para%cg_rest, para%cg_maxiter, para%cg_log)
call init_cg_stat()
call init_xbound()
call init_confs(para, conf)
call check_former(para%n_temper, conf)
! JuBE
! call jube kernel run function
call jube_kernel_run()
call mc(para, conf)
!!call xbound_test()
! JuBE
! call jube kernel finalize function
call jube_kernel_finalize()
call conf_write(.true., para, conf)
call write_counter(para%maxtraj)
call write_ran()
TIMING_STOP(timing_bin_total)
call write_footer(time0)
call end_A(UREC, "Job")
call comm_finalize()
! JuBE
! call jube kernel end function
call jube_kernel_end()
end subroutine kernel_a
!-------------------------------------------------------------------------------
subroutine init_para(para, flags)
! initialises module_para, module_switches and module_mre
use typedef_flags
use typedef_para
use module_bqcd
use module_input
use module_mre
use module_switches
implicit none
type(type_para) :: para
type(type_flags) :: flags
integer :: i
logical :: quenched, dynamical, clover, h_ext
quenched = .false.
dynamical = .false.
clover = .false.
h_ext = .false.
para%run = input%run
para%L = input%lattice
para%NPE = input%processes
para%bc_fermions = input%boundary_conditions_fermions
para%gamma_index = input%gamma_index
para%n_temper = input%ensembles
para%nstd = input%tempering_steps_without
para%nforce = input%hmc_accept_first
para%ntraj = input%mc_steps
para%maxtraj = input%mc_total_steps
para%nsave = input%mc_save_frequency
para%c_cg_rest = input%solver_rest
para%cg_maxiter = input%solver_maxiter
para%cg_log = input%solver_ignore_no_convergence
mre_n_vec = input%solver_mre_vectors
call check_bc_fermions(para%bc_fermions, para%gamma_index)
read(para%c_cg_rest, *) para%cg_rest
if (para%n_temper <= 0) call die("init_para(): n_temper <= 0")
if (para%n_temper > MAX_TEMPER) call die("init_para(): n_temper > max_temper")
do i = 1, para%n_temper
para%c_hmc(i)%beta = input%beta(i)
para%c_hmc(i)%kappa = input%kappa(i)
para%c_hmc(i)%csw = input%csw(i)
para%c_hmc(i)%h = input%h(i)
para%c_hmc(i)%traj_length = input%hmc_trajectory_length(i)
para%c_hmc(i)%ntau = input%hmc_steps(i)
para%c_hmc(i)%rho = input%hmc_rho(i)
para%c_hmc(i)%m_scale = input%hmc_m_scale(i)
para%info_file(i) = input%start_info_file(i)
read(para%c_hmc(i)%beta, *) para%hmc(i)%beta
read(para%c_hmc(i)%kappa, *) para%hmc(i)%kappa
read(para%c_hmc(i)%csw, *) para%hmc(i)%csw
read(para%c_hmc(i)%h, *) para%hmc(i)%h
read(para%c_hmc(i)%traj_length,*) para%hmc(i)%traj_length
read(para%c_hmc(i)%ntau, *) para%hmc(i)%ntau
read(para%c_hmc(i)%rho, *) para%hmc(i)%rho
read(para%c_hmc(i)%m_scale, *) para%hmc(i)%m_scale
if (para%hmc(i)%kappa == ZERO .and. para%hmc(i)%csw /= ZERO) then
para%hmc(i)%csw_kappa = para%hmc(i)%csw
para%c_hmc(i)%csw = "-1 (infinity)"
para%hmc(i)%csw = -1
else
para%hmc(i)%csw_kappa = para%hmc(i)%csw * para%hmc(i)%kappa
call check_csw(para%hmc(i)%beta, para%hmc(i)%csw)
endif
para%hmc(i)%tau = para%hmc(i)%traj_length / para%hmc(i)%ntau
write(para%c_hmc(i)%csw_kappa, "(f20.15)") para%hmc(i)%csw_kappa
write(para%c_hmc(i)%tau, "(f20.15)") para%hmc(i)%tau
if (para%hmc(i)%kappa == ZERO .and. para%hmc(i)%csw == ZERO) then
quenched = .true.
else
dynamical = .true.
endif
if (para%hmc(i)%csw /= ZERO) clover = .true.
if (para%hmc(i)%h /= ZERO) h_ext = .true.
para%hmc(i)%model = input%hmc_model
if (para%hmc(i)%model == "A" .and. para%hmc(i)%rho /= ZERO) then
call warn("init_para(): model == A but rho /= 0")
endif
if (para%hmc(i)%model /= "A" .and. para%hmc(i)%rho == ZERO) then
call warn("init_para(): model /= A but rho == 0")
endif
enddo
select case (input%start_configuration)
case ("hot"); para%start = START_HOT
case ("cold"); para%start = START_COLD
case ("file"); para%start = START_FILE
case default
call die("init_para(): start_configuration must be {hot|cold|file}")
end select
select case (input%start_random)
case ("random"); para%seed = -1
case ("default"); para%seed = 0
case default; read(input%start_random, *) para%seed
end select
select case (input%tempering_swap_sequence)
case ("up"); para%swap_seq = SWAP_UP
case ("down"); para%swap_seq = SWAP_DOWN
case ("random"); para%swap_seq = SWAP_RANDOM
case default
call die("init_para(): tempering_swap_sequence must be {up|down|random}")
end select
if (quenched .and. dynamical) call die("init_para(): quenched/dynamical mixed")
if (para%nforce < 0) call die("init_para(): nforce < 0")
if (flags%continuation_job) para%start = START_CONT
switches%quenched = quenched
switches%dynamical = dynamical
switches%clover = clover
switches%h_ext = h_ext
switches%hasenbusch = (input%hmc_model /= "A")
if (quenched) switches%hasenbusch = .false.
switches%tempering = .false.
switches%measure_polyakov_loop = .false.
switches%measure_traces = .false.
if (input%ensembles > 1) switches%tempering = .true.
if (input%measure_polyakov_loop /= 0) switches%measure_polyakov_loop = .true.
if (input%measure_traces /= 0) switches%measure_traces = .true.
if (input%hmc_test == 0) then
switches%hmc_test = .false.
else
switches%hmc_test = .true.
endif
end subroutine init_para
!-------------------------------------------------------------------------------
subroutine init_counter(para, flags)
use typedef_flags
use typedef_para
use module_counter
use module_function_decl
implicit none
type(type_para) :: para
type(type_flags) :: flags
FILENAME, external :: count_file, stop_file
if (f_exist(stop_file())) then
call die("init_counter(): found stop file " // stop_file())
endif
counter%run = para%run
counter%j_traj = 0
if (flags%continuation_job) then
open(UCOUNT, file = count_file(), action = "read", status = "old")
read(UCOUNT, *) counter%run
read(UCOUNT, *) counter%job
read(UCOUNT, *) counter%traj
close(UCOUNT)
if (counter%run /= para%run) call die("init_counter(): RUN inconsistent")
counter%job = counter%job + 1
else
counter%run = para%run
counter%job = 1
counter%traj = -para%nforce
endif
end subroutine init_counter
!-------------------------------------------------------------------------------
subroutine write_counter(maxtraj)
use module_counter
use module_function_decl
implicit none
integer :: maxtraj
FILENAME, external :: count_file, stop_file
if (my_pe() /= 0) return
open(UCOUNT, file = count_file(), action = "write")
write(UCOUNT, *) counter%run, " run"
write(UCOUNT, *) counter%job, " job"
write(UCOUNT, *) counter%traj, " traj"
close(UCOUNT)
if (counter%traj >= maxtraj) then
open(UCOUNT, file = stop_file(), status = "unknown")
close(UCOUNT)
endif
end subroutine write_counter
!-------------------------------------------------------------------------------
subroutine write_header(para)
use typedef_para
use module_bqcd
use module_counter
use module_function_decl
use module_input
use module_mre
use module_thread
implicit none
type(type_para) :: para
integer :: i
character(len = 50) :: fmt
character(len = 4), external :: format_ensemble
if (my_pe() == 0) then
fmt = "(1x,a," // format_ensemble() // ",2a)"
call begin(UREC, "Header")
if (input%comment /= "") then
write(UREC, 405) "Comment", trim(input%comment)
endif
write(UREC, 400) "Program", prog_name, prog_version
write(UREC, *) "Version_of_D ", version_of_d()
write(UREC, *) "Communication ", trim(comm_method())
write(UREC, *) "Run ", para%run
write(UREC, *) "Job ", counter%job
write(UREC, 405) "Host", rechner()
write(UREC, 400) "Date", datum(), uhrzeit()
write(UREC, 410) "L ", para%L
write(UREC, 410) "NPE ", para%NPE
write(UREC, 410) "bc_fermions", para%bc_fermions
write(UREC, 410) "gamma_index", para%gamma_index
write(UREC, *) "Threads ", n_thread
write(UREC, *) "Start ", para%start
if (para%start == START_FILE) then
do i = 1, para%n_temper
write(UREC, fmt) "StartConf_", i, " ", trim(para%info_file(i))
enddo
endif
write(UREC, *) "Seed ", para%seed
write(UREC, *) "Swap_seq", para%swap_seq
write(UREC, *) "N_force ", para%nforce
write(UREC, *) "N_traj ", para%ntraj
write(UREC, *) "N_save ", para%nsave
write(UREC, *) "N_temper", para%n_temper
do i = 1, para%n_temper
write(UREC, fmt) "beta_", i, " ", trim(para%c_hmc(i)%beta)
write(UREC, fmt) "kappa_", i, " ", trim(para%c_hmc(i)%kappa)
write(UREC, fmt) "csw_", i, " ", trim(para%c_hmc(i)%csw)
write(UREC, fmt) "csw_kappa_", i, " ", trim(para%c_hmc(i)%csw_kappa)
write(UREC, fmt) "h_", i, " ", trim(para%c_hmc(i)%h)
write(UREC, fmt) "tau_", i, " ", trim(para%c_hmc(i)%tau)
write(UREC, fmt) "N_tau_", i, " ", trim(para%c_hmc(i)%ntau)
write(UREC, fmt) "traj_length_", i, " ", trim(para%c_hmc(i)%traj_length)
write(UREC, fmt) "rho_", i, " ", trim(para%c_hmc(i)%rho)
write(UREC, fmt) "m_scale_", i, " ", trim(para%c_hmc(i)%m_scale)
enddo
write(UREC, *) "HMC_model ", para%hmc(1)%model
write(UREC, *) "REAL_kind ", RKIND
write(UREC, 405) "CG_rest ", trim(para%c_cg_rest)
write(UREC, *) "MRE_vectors ", mre_n_vec
call end_A(UREC, "Header")
400 format (3(1x,a))
405 format (2(1x,a))
410 format (1x,a,4i3)
endif
end subroutine write_header
!-------------------------------------------------------------------------------
subroutine write_footer(time0)
use module_function_decl
use module_thread
implicit none
SEED :: seed
SECONDS :: time0, sekunden
call ranget(seed)
call begin(UREC, "Footer")
if (my_pe() == 0) then
write(UREC, 400) "Date", datum(), uhrzeit()
write(UREC, *) "Seed", seed
write(UREC, 410) "CPU-Time", &
sekunden() - time0, "s on", num_pes() * n_thread, "CPUs"
endif
400 format (3(1x,a))
410 format (1x,a,1x,f8.1,1x,a,1x,i5,1x,a)
TIMING_WRITE(UREC)
call end_A(UREC, "Footer")
end subroutine write_footer
!-------------------------------------------------------------------------------
subroutine get_flags(flags)
use typedef_cksum
use typedef_flags
use module_bqcd
use module_function_decl
use module_input
implicit none
type(type_flags), intent(out) :: flags
integer :: iarg, length, stat, narg
integer, external :: ipxfargc
character(len = 2) :: opt
flags%continuation_job = .false.
flags%show_version = .false.
narg = ipxfargc()
iarg = 1
do while (iarg <= narg)
call pxfgetarg(iarg, opt, length, stat)
if (opt(1:1) == "-") then
if (length > 2) call usage()
select case (opt(2:2))
case ("c")
flags%continuation_job = .true.
iarg = iarg + 1
case ("I")
call input_dump(6)
call comm_finalize()
stop
case ("V")
flags%show_version = .true.
iarg = iarg + 1
case default
call usage
end select
else
exit
endif
enddo
if (flags%show_version) then
call version()
call comm_finalize()
stop
endif
call take_arg(iarg, flags%input, narg)
if (narg >= iarg) call usage
CONTAINS
subroutine usage()
implicit none
call die("Usage: " // prog_name // " [-c] [-I] [-V] input")
end subroutine usage
subroutine version()
implicit none
if (my_pe() == 0) then
write(6,*) "This is ", prog_name, " ", prog_version
write(6,*) " input format: ", input_version
write(6,*) " conf info format:", conf_info_version
write(6,*) " MAX_TEMPER: ", MAX_TEMPER
write(6,*) " real kind: ", RKIND
write(6,*) " version of D: ", version_of_d()
write(6,*) " D3: buffer vol: ", get_d3_buffer_vol()
write(6,*) " communication: ", trim(comm_method())
endif
end subroutine version
subroutine take_arg(iarg, arg, narg)
implicit none
integer, intent(inout) :: iarg
character(len = *), intent(out) :: arg
integer, intent(in) :: narg
integer :: length, stat
if (iarg > narg) call usage
call pxfgetarg(iarg, arg, length, stat)
if (length > len(arg)) then
call die("get_flags(): " // arg // ": argument too long")
endif
iarg = iarg + 1
end subroutine take_arg
end subroutine get_flags
!===============================================================================
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2005, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! cg.F90
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
module module_cg
type type_cg_para
real :: rest
integer :: maxiter
integer :: log
end type type_cg_para
type(type_cg_para), save :: cg_para
type type_cg_stat
integer :: niter
integer :: niter_max
integer :: niter_tot
integer :: ncall
end type type_cg_stat
type(type_cg_stat), save :: cg_stat
integer, save :: cg_iterations_total = 0 ! used in timing.F90
end
!-------------------------------------------------------------------------------
subroutine cg(matrix_mult, x, b, para, conf, iterations)
! solves "matrix_mult * x = b" and returns number of iterations
use module_cg
use module_function_decl
use module_p_interface
use module_vol
use typedef_hmc
implicit none
external :: matrix_mult
SPINCOL_OVERINDEXED, intent(out) :: x
SPINCOL_OVERINDEXED, intent(in) :: b
type(hmc_para), intent(in) :: para
type(hmc_conf), intent(in) :: conf
integer, intent(out) :: iterations
P_SPINCOL_OVERINDEXED, save :: r, p, aap
REAL :: ak, bk, rtr, rtrold, paap
integer :: i, niter
character(72) :: msg
TIMING_START(timing_bin_cg)
ALLOCATE_SC_OVERINDEXED(r)
ALLOCATE_SC_OVERINDEXED(p)
ALLOCATE_SC_OVERINDEXED(aap)
call matrix_mult(r, x, para, conf)
rtrold = ZERO
!$omp parallel do reduction(+: rtrold)
do i = 1, size_sc_field
r(i) = b(i) - r(i)
p(i) = r(i)
rtrold = rtrold + r(i)**2
enddo
rtrold = global_sum(rtrold)
do niter = 1, cg_para%maxiter
call matrix_mult(aap, p, para, conf)
paap = sc_dot(p, aap)
paap = global_sum(paap)
ak = rtrold / paap
rtr = ZERO
!$omp parallel do reduction(+: rtr)
do i = 1, size_sc_field
x(i) = x(i) + ak * p(i)
r(i) = r(i) - ak * aap(i)
rtr = rtr + r(i)**2
enddo
rtr = global_sum(rtr)
if (rtr <= cg_para%rest) goto 9999
bk = rtr / rtrold
rtrold = rtr
call sc_xpby(p, r, bk) ! p = r + bk * p
enddo
niter = niter - 1
if (cg_para%log /= 2) then
write(msg, *) "cg(): no convergence; rtr = ", rtr
call die(msg)
endif
9999 continue
cg_stat%ncall = cg_stat%ncall + 1
cg_stat%niter = niter
cg_stat%niter_max = max(cg_stat%niter_max, niter)
cg_stat%niter_tot = cg_stat%niter_tot + niter
cg_iterations_total = cg_iterations_total + niter
iterations = niter
TIMING_STOP(timing_bin_cg)
end
!-------------------------------------------------------------------------------
subroutine init_cg_para(rest, maxiter, log)
use module_cg
implicit none
real rest
integer maxiter, log
cg_para%rest = rest
cg_para%maxiter = maxiter
cg_para%log = log
end
!-------------------------------------------------------------------------------
subroutine init_cg_stat()
use module_cg
implicit none
cg_stat%ncall = 0
cg_stat%niter_max = 0
cg_stat%niter_tot = 0
end
!-------------------------------------------------------------------------------
subroutine get_cg_stat(ncall, niter_max, niter_tot)
use module_cg
implicit none
integer ncall, niter_max, niter_tot
ncall = cg_stat%ncall
niter_max = cg_stat%niter_max
niter_tot = cg_stat%niter_tot
end
!===============================================================================
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 2003, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! checks.F90
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
subroutine check_csw(beta, csw)
implicit none
REAL, intent(in) :: beta, csw
REAL :: g, c
if (beta == ZERO) return
if (csw == ZERO) return
g = SIX / beta
c = ONE - 0.454 * g - 0.175 * g**2 + 0.012 * g**3 + 0.045 * g**4
c = c / (ONE - 0.720 * g)
if (abs(c - csw) > 0.00005) then
call warn("check_csw(): c_sw differs more than 0.00005 from ALPHA value")
endif
end
!-------------------------------------------------------------------------------
subroutine check_bc_fermions(bc_fermions, gamma_index)
! warns if the number of anti-periodic fermionic b.c. is 1 and
! the anti-periodic direction is not the gamma_4 direction
implicit none
integer, dimension(DIM), intent(in) :: bc_fermions, gamma_index
integer :: i, i_anti, count
count = 0
do i = 1, DIM
if (bc_fermions(i) < 0) then
count = count + 1
i_anti = i
endif
enddo
if (count == 1 .and. gamma_index(i_anti) /= 4) then
call warn("check_bc_fermions(): anti-periodic b.c. not in gamma_4 direction")
endif
end
!===============================================================================
/*
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2001, Hinnerk Stueben, Zuse Institute Berlin
!
!-------------------------------------------------------------------------------
!
! Adopted from:
!
! ==================
! QCD SF/T3E PROGRAM
! ==================
!
! Calculate a modified cyclic redundancy check (CRC is specified by the
! POSIX.2 standard). Modification is necessary since I do not know how
! to handle a unsigned long in FORTRAN. Solution: CKSUM returns negative
! number if result is > LONG_MAX.
!
! CKSUM_GET() returns always positive numbers < LONG_MAX. (H.S.)
!
!
! Parts of the source come from cksum.c of the GNU text utilities (version
! 1.19) written by Q. Frank Xia.
!
! $Log: cksum.c,v $
! Revision 1.1 2007/11/14 13:10:15 mallalen
! *** empty log message ***
!
! Revision 1.1 1997/12/04 10:31:09 pleiter
! Initial writing attempt
!
!-----------------------------------------------------------------------------
*/
#ifdef NamesToLower_
# define CKSUM_INIT cksum_init_
# define CKSUM_ADD cksum_add_
# define CKSUM_GET cksum_get_
#endif
#ifdef NamesToLower
# define CKSUM_INIT cksum_init
# define CKSUM_ADD cksum_add
# define CKSUM_GET cksum_get
#endif
#ifdef LongLong
# define INT8 long long
#else
# define INT8 long
#endif
void CKSUM_INIT(void);
void CKSUM_ADD(void *, INT8 *);
void CKSUM_GET(INT8 *, INT8 *);
static unsigned INT8 the_crc = 0;
static INT8 the_bytes = 0;
static unsigned INT8 const crctab[256] =
{
0x0,
0x04C11DB7, 0x09823B6E, 0x0D4326D9, 0x130476DC, 0x17C56B6B,
0x1A864DB2, 0x1E475005, 0x2608EDB8, 0x22C9F00F, 0x2F8AD6D6,
0x2B4BCB61, 0x350C9B64, 0x31CD86D3, 0x3C8EA00A, 0x384FBDBD,
0x4C11DB70, 0x48D0C6C7, 0x4593E01E, 0x4152FDA9, 0x5F15ADAC,
0x5BD4B01B, 0x569796C2, 0x52568B75, 0x6A1936C8, 0x6ED82B7F,
0x639B0DA6, 0x675A1011, 0x791D4014, 0x7DDC5DA3, 0x709F7B7A,
0x745E66CD, 0x9823B6E0, 0x9CE2AB57, 0x91A18D8E, 0x95609039,
0x8B27C03C, 0x8FE6DD8B, 0x82A5FB52, 0x8664E6E5, 0xBE2B5B58,
0xBAEA46EF, 0xB7A96036, 0xB3687D81, 0xAD2F2D84, 0xA9EE3033,
0xA4AD16EA, 0xA06C0B5D, 0xD4326D90, 0xD0F37027, 0xDDB056FE,
0xD9714B49, 0xC7361B4C, 0xC3F706FB, 0xCEB42022, 0xCA753D95,
0xF23A8028, 0xF6FB9D9F, 0xFBB8BB46, 0xFF79A6F1, 0xE13EF6F4,
0xE5FFEB43, 0xE8BCCD9A, 0xEC7DD02D, 0x34867077, 0x30476DC0,
0x3D044B19, 0x39C556AE, 0x278206AB, 0x23431B1C, 0x2E003DC5,
0x2AC12072, 0x128E9DCF, 0x164F8078, 0x1B0CA6A1, 0x1FCDBB16,
0x018AEB13, 0x054BF6A4, 0x0808D07D, 0x0CC9CDCA, 0x7897AB07,
0x7C56B6B0, 0x71159069, 0x75D48DDE, 0x6B93DDDB, 0x6F52C06C,
0x6211E6B5, 0x66D0FB02, 0x5E9F46BF, 0x5A5E5B08, 0x571D7DD1,
0x53DC6066, 0x4D9B3063, 0x495A2DD4, 0x44190B0D, 0x40D816BA,
0xACA5C697, 0xA864DB20, 0xA527FDF9, 0xA1E6E04E, 0xBFA1B04B,
0xBB60ADFC, 0xB6238B25, 0xB2E29692, 0x8AAD2B2F, 0x8E6C3698,
0x832F1041, 0x87EE0DF6, 0x99A95DF3, 0x9D684044, 0x902B669D,
0x94EA7B2A, 0xE0B41DE7, 0xE4750050, 0xE9362689, 0xEDF73B3E,
0xF3B06B3B, 0xF771768C, 0xFA325055, 0xFEF34DE2, 0xC6BCF05F,
0xC27DEDE8, 0xCF3ECB31, 0xCBFFD686, 0xD5B88683, 0xD1799B34,
0xDC3ABDED, 0xD8FBA05A, 0x690CE0EE, 0x6DCDFD59, 0x608EDB80,
0x644FC637, 0x7A089632, 0x7EC98B85, 0x738AAD5C, 0x774BB0EB,
0x4F040D56, 0x4BC510E1, 0x46863638, 0x42472B8F, 0x5C007B8A,
0x58C1663D, 0x558240E4, 0x51435D53, 0x251D3B9E, 0x21DC2629,
0x2C9F00F0, 0x285E1D47, 0x36194D42, 0x32D850F5, 0x3F9B762C,
0x3B5A6B9B, 0x0315D626, 0x07D4CB91, 0x0A97ED48, 0x0E56F0FF,
0x1011A0FA, 0x14D0BD4D, 0x19939B94, 0x1D528623, 0xF12F560E,
0xF5EE4BB9, 0xF8AD6D60, 0xFC6C70D7, 0xE22B20D2, 0xE6EA3D65,
0xEBA91BBC, 0xEF68060B, 0xD727BBB6, 0xD3E6A601, 0xDEA580D8,
0xDA649D6F, 0xC423CD6A, 0xC0E2D0DD, 0xCDA1F604, 0xC960EBB3,
0xBD3E8D7E, 0xB9FF90C9, 0xB4BCB610, 0xB07DABA7, 0xAE3AFBA2,
0xAAFBE615, 0xA7B8C0CC, 0xA379DD7B, 0x9B3660C6, 0x9FF77D71,
0x92B45BA8, 0x9675461F, 0x8832161A, 0x8CF30BAD, 0x81B02D74,
0x857130C3, 0x5D8A9099, 0x594B8D2E, 0x5408ABF7, 0x50C9B640,
0x4E8EE645, 0x4A4FFBF2, 0x470CDD2B, 0x43CDC09C, 0x7B827D21,
0x7F436096, 0x7200464F, 0x76C15BF8, 0x68860BFD, 0x6C47164A,
0x61043093, 0x65C52D24, 0x119B4BE9, 0x155A565E, 0x18197087,
0x1CD86D30, 0x029F3D35, 0x065E2082, 0x0B1D065B, 0x0FDC1BEC,
0x3793A651, 0x3352BBE6, 0x3E119D3F, 0x3AD08088, 0x2497D08D,
0x2056CD3A, 0x2D15EBE3, 0x29D4F654, 0xC5A92679, 0xC1683BCE,
0xCC2B1D17, 0xC8EA00A0, 0xD6AD50A5, 0xD26C4D12, 0xDF2F6BCB,
0xDBEE767C, 0xE3A1CBC1, 0xE760D676, 0xEA23F0AF, 0xEEE2ED18,
0xF0A5BD1D, 0xF464A0AA, 0xF9278673, 0xFDE69BC4, 0x89B8FD09,
0x8D79E0BE, 0x803AC667, 0x84FBDBD0, 0x9ABC8BD5, 0x9E7D9662,
0x933EB0BB, 0x97FFAD0C, 0xAFB010B1, 0xAB710D06, 0xA6322BDF,
0xA2F33668, 0xBCB4666D, 0xB8757BDA, 0xB5365D03, 0xB1F740B4
};
void CKSUM_INIT(void)
{
the_crc = 0;
the_bytes = 0;
}
void CKSUM_ADD(void *memptr, INT8 *nbytes)
{
register unsigned INT8 crc;
register INT8 i;
register unsigned char *cp;
crc = the_crc;
cp = (unsigned char *) memptr;
for (i=0; i < *nbytes; i++)
crc = (crc << 8) ^ crctab[((crc >> 24) ^ *(cp++)) & 0xFF];
the_crc = crc;
the_bytes += *nbytes;
}
void CKSUM_GET(INT8 *total_crc, INT8 *total_bytes)
{
register unsigned INT8 crc;
register INT8 i;
crc = the_crc;
for (i = the_bytes; i > 0; i >>= 8)
crc = (crc << 8) ^ crctab[((crc >> 24) ^ i) & 0xFF];
crc = (~crc & 0xFFFFFFFF);
*total_crc = (INT8) crc;
*total_bytes = the_bytes;
}
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2001, Hinnerk Stueben, Zuse Institute Berlin
!
!-------------------------------------------------------------------------------
!
! cksum_dummy.F90 - dummy routines that can replace the real C routines
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
subroutine cksum_init()
return
end
!-------------------------------------------------------------------------------
subroutine cksum_add(i, j)
integer i(*)
CHECK_SUM j
return
end
!-------------------------------------------------------------------------------
subroutine cksum_get(sum, bytes)
CHECK_SUM sum, bytes
sum = 0
bytes = 0
end
!===============================================================================
#===============================================================================
#
# BQCD -- Berlin Quantum ChromoDynamics programme
#
# Author: Hinnerk Stueben <stueben@zib.de>
#
# Copyright (C) 1998-2006, Hinnerk Stueben, Zuse-Institut Berlin
#
#-------------------------------------------------------------------------------
#
# clover/Makefile
#
#===============================================================================
include ../Makefile.defs
fpp = $(FPP) -I.. $(FPPFLAGS)
MODULES_DIR = ../modules
.SUFFIXES:
.SUFFIXES: .a .o .F90
.F90.o:
$(fpp) $< > $*.f90
$(F90) -c $(FFLAGS) -I$(MODULES_DIR) $*.f90
OBJS = \
clover_action.o \
clover_allocate.o \
clover_bsa.o \
clover_d.o \
clover_f_mu_nu.o \
clover_init.o \
clover_inv.o \
clover_mult_a.o \
clover_mult_ao.o \
clover_mult_b.o \
clover_t_init.o \
clover_ts.o \
clover_uuu.o \
clover_uuuu.o
OBJS_CTEST = \
ctest.o \
clover_inv.o \
clover_mult_a.o \
clover_mult_ao.o \
clover_mult_b.o
$(LIBCLOVER):
libclover.a: $(OBJS)
$(AR) $(ARFLAGS) $@ $(OBJS)
$(RANLIB) $@
fast:
$(FAST_MAKE)
ctest: $(OBJS_CTEST)
f90 -o $@ $(OBJS_CTEST)
clean:
rm -f *.[Tiod] *.f90 *.mod core work.pc work.pcl
clobber: clean
rm -f libclover.a ctest
#ifdef CLOVER_AS_COMPLEX_ARRAY
# define A11 Re(a(1,J,i))
# define A22 Im(a(1,J,i))
# define A33 Re(a(11,J,i))
# define A44 Im(a(11,J,i))
# define A55 Re(a(17,J,i))
# define A66 Im(a(17,J,i))
# define A12 a(2,J,i)
# define A13 a(3,J,i)
# define A14 a(4,J,i)
# define A15 a(5,J,i)
# define A16 a(6,J,i)
# define A23 a(7,J,i)
# define A24 a(8,J,i)
# define A25 a(9,J,i)
# define A26 a(10,J,i)
# define A34 a(12,J,i)
# define A35 a(13,J,i)
# define A36 a(14,J,i)
# define A45 a(15,J,i)
# define A46 a(16,J,i)
# define A56 a(18,J,i)
# define B11 Re(b(16,J,i))
# define B22 Im(b(16,J,i))
# define B33 Re(b(17,J,i))
# define B44 Im(b(17,J,i))
# define B55 Re(b(18,J,i))
# define B66 Im(b(18,J,i))
# define B21 b(1,J,i)
# define B31 b(2,J,i)
# define B32 b(3,J,i)
# define B41 b(4,J,i)
# define B42 b(5,J,i)
# define B43 b(6,J,i)
# define B51 b(7,J,i)
# define B52 b(8,J,i)
# define B53 b(9,J,i)
# define B54 b(10,J,i)
# define B61 b(11,J,i)
# define B62 b(12,J,i)
# define B63 b(13,J,i)
# define B64 b(14,J,i)
# define B65 b(15,J,i)
#else
# define A11 a%i11
# define A22 a%i22
# define A33 a%i33
# define A44 a%i44
# define A55 a%i55
# define A66 a%i66
# define A12 a%i12
# define A13 a%i13
# define A14 a%i14
# define A15 a%i15
# define A16 a%i16
# define A23 a%i23
# define A24 a%i24
# define A25 a%i25
# define A26 a%i26
# define A34 a%i34
# define A35 a%i35
# define A36 a%i36
# define A45 a%i45
# define A46 a%i46
# define A56 a%i56
# define B11 b%i11
# define B22 b%i22
# define B33 b%i33
# define B44 b%i44
# define B55 b%i55
# define B66 b%i66
# define B21 b%i21
# define B31 b%i31
# define B32 b%i32
# define B41 b%i41
# define B42 b%i42
# define B43 b%i43
# define B51 b%i51
# define B52 b%i52
# define B53 b%i53
# define B54 b%i54
# define B61 b%i61
# define B62 b%i62
# define B63 b%i63
# define B64 b%i64
# define B65 b%i65
#endif
# define A21 conjg(A12)
# define A31 conjg(A13)
# define A41 conjg(A14)
# define A51 conjg(A15)
# define A61 conjg(A16)
# define A32 conjg(A23)
# define A42 conjg(A24)
# define A52 conjg(A25)
# define A62 conjg(A26)
# define A43 conjg(A34)
# define A53 conjg(A35)
# define A63 conjg(A36)
# define A54 conjg(A45)
# define A64 conjg(A46)
# define A65 conjg(A56)
# define B12 conjg(B21)
# define B13 conjg(B31)
# define B23 conjg(B32)
# define B14 conjg(B41)
# define B24 conjg(B42)
# define B34 conjg(B43)
# define B15 conjg(B51)
# define B25 conjg(B52)
# define B35 conjg(B53)
# define B45 conjg(B54)
# define B16 conjg(B61)
# define B26 conjg(B62)
# define B36 conjg(B63)
# define B46 conjg(B64)
# define B56 conjg(B65)
# define SC1 1, 1
# define SC2 1, 2
# define SC3 1, 3
# define SC4 2, 1
# define SC5 2, 2
# define SC6 2, 3
# define SC7 3, 1
# define SC8 3, 2
# define SC9 3, 3
# define SC10 4, 1
# define SC11 4, 2
# define SC12 4, 3
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2003, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! clover_action.F90 - calculates: -2 Tr(log(T_oo))
!
!-------------------------------------------------------------------------------
# include "defs.h"
# include "clover.h"
!-------------------------------------------------------------------------------
REAL function clover_action(b)
use typedef_clover
use module_vol
implicit none
type(type_clover_b) :: b(2, volh)
integer :: i
REAL :: s, global_sum
s = ZERO
!$omp parallel do reduction(+: s)
do i = 1, volh
s = s + log(det(b(1, i)) * det(b(2, i)))
enddo
clover_action = TWO * global_sum(s)
CONTAINS
REAL function det(b) ! returns (1 / det)
type(type_clover_b) :: b
det = B11 * B22 * B33 * B44 * B55 * B66
end function det
end
!===============================================================================
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2003, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! clover_allocate.F90 - allocation of clover arrays
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
subroutine allocate_clover_field_a(a)
use typedef_clover
use module_vol
implicit none
P_CLOVER_FIELD_A :: a
allocate(a(2, volh, EVEN:ODD))
end
!-------------------------------------------------------------------------------
subroutine allocate_clover_field_b(b)
use typedef_clover
use module_vol
implicit none
P_CLOVER_FIELD_B :: b
allocate(b(2, volh, EVEN:ODD))
end
!===============================================================================
!===============================================================================
!
! BQCD -- Berlin Quantum ChromoDynamics program
!
! Author: Hinnerk Stueben <stueben@zib.de>
!
! Copyright (C) 1998-2003, Hinnerk Stueben, Zuse-Institut Berlin
!
!-------------------------------------------------------------------------------
!
! clover_bsa.F90 - calculates "B sigma A"
!
!-------------------------------------------------------------------------------
# include "defs.h"
!-------------------------------------------------------------------------------
subroutine clover_bsa(mu, nu, w, b, a) ! w = transposed(conjg(b) sigma_mu_nu a)
use module_vol
implicit none
integer :: mu, nu
SU3_FIELD :: w
SPINCOL_FIELD :: b, a
if (mu == 1) then
if (nu == 2) then ; call clover_bsa_12(w, b, a)
elseif (nu == 3) then ; call clover_bsa_13(w, b, a)
elseif (nu == 4) then ; call clover_bsa_14(w, b, a) ; endif
elseif (mu == 2) then
if (nu == 1) then ; call clover_bsa_21(w, b, a)
elseif (nu == 3) then ; call clover_bsa_23(w, b, a)
elseif (nu == 4) then ; call clover_bsa_24(w, b, a) ; endif
elseif (mu == 3) then
if (nu == 1) then ; call clover_bsa_31(w, b, a)
elseif (nu == 2) then ; call clover_bsa_32(w, b, a)
elseif (nu == 4) then ; call clover_bsa_34(w, b, a) ; endif
elseif (mu == 4) then
if (nu == 1) then ; call clover_bsa_41(w, b, a)
elseif (nu == 2) then ; call clover_bsa_42(w, b, a)
elseif (nu == 3) then ; call clover_bsa_43(w, b, a) ; endif
endif
end
!-------------------------------------------------------------------------------
subroutine clover_bsa_12(w, b, a)
# include "clover_bsa_head.h90"
a1 = -a(1, ca, i)
a2 = a(2, ca, i)
a3 = -a(3, ca, i)
a4 = a(4, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_21(w, b, a)
# include "clover_bsa_head.h90"
a1 = a(1, ca, i)
a2 = -a(2, ca, i)
a3 = a(3, ca, i)
a4 = -a(4, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_13(w, b, a)
# include "clover_bsa_head.h90"
a1 = -i_times(a(2, ca, i))
a2 = i_times(a(1, ca, i))
a3 = -i_times(a(4, ca, i))
a4 = i_times(a(3, ca, i))
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_31(w, b, a)
# include "clover_bsa_head.h90"
a1 = i_times(a(2, ca, i))
a2 = -i_times(a(1, ca, i))
a3 = i_times(a(4, ca, i))
a4 = -i_times(a(3, ca, i))
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_14(w, b, a)
# include "clover_bsa_head.h90"
a1 = a(4, ca, i)
a2 = a(3, ca, i)
a3 = a(2, ca, i)
a4 = a(1, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_41(w, b, a)
# include "clover_bsa_head.h90"
a1 = -a(4, ca, i)
a2 = -a(3, ca, i)
a3 = -a(2, ca, i)
a4 = -a(1, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_23(w, b, a)
# include "clover_bsa_head.h90"
a1 = -a(2, ca, i)
a2 = -a(1, ca, i)
a3 = -a(4, ca, i)
a4 = -a(3, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_32(w, b, a)
# include "clover_bsa_head.h90"
a1 = a(2, ca, i)
a2 = a(1, ca, i)
a3 = a(4, ca, i)
a4 = a(3, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_24(w, b, a)
# include "clover_bsa_head.h90"
a1 = -i_times(a(4, ca, i))
a2 = i_times(a(3, ca, i))
a3 = -i_times(a(2, ca, i))
a4 = i_times(a(1, ca, i))
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_42(w, b, a)
# include "clover_bsa_head.h90"
a1 = i_times(a(4, ca, i))
a2 = -i_times(a(3, ca, i))
a3 = i_times(a(2, ca, i))
a4 = -i_times(a(1, ca, i))
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_34(w, b, a)
# include "clover_bsa_head.h90"
a1 = a(3, ca, i)
a2 = -a(4, ca, i)
a3 = a(1, ca, i)
a4 = -a(2, ca, i)
# include "clover_bsa_tail.h90"
!-------------------------------------------------------------------------------
subroutine clover_bsa_43(w, b, a)
# include "clover_bsa_head.h90"
a1 = -a(3, ca, i)
a2 = a(4, ca, i)
a3 = -a(1, ca, i)
a4 = a(2, ca, i)
# include "clover_bsa_tail.h90"
!===============================================================================
use module_vol
implicit none
SU3_FIELD :: w
SPINCOL_FIELD :: b, a
COMPLEX :: a1, a2, a3, a4
integer :: i, ca, cb
! statement function:
COMPLEX :: i_times, c
i_times(c) = cmplx(-aimag(c), real(c))
!$omp parallel do private(ca, cb, a1, a2, a3, a4)
do i = 1, volh
do ca = 1, NCOL