ADfer and SFCF. v1.1

input/FerFlow.in

@@ -17,6 +17,12 @@ eps      = 0.01
 ns       = 30
 mclength = -1 # Number of measurements. Negative to ignore
 
+[AD]
+nder    = 2
+
+[sfcf]
+meas    = true
+
 [Fermion1]
 kappa   = 0.1348
 theta   = 0.0

main.jl

@@ -15,6 +15,7 @@ reset_timer!()
 
 include("./src/io.jl")
 include("./src/meas.jl")
+include("./src/sfcf.jl")
 include("./src/mc.jl")
 
 @timeit "Input reading and space allocation" begin
@@ -37,6 +38,8 @@ while RUN_ON
 
     @timeit "One point function" Backflow_pt()
 
+    meas_sfcf()
+
     @timeit "Saving" save_data()
 
     check_run_status()

src/io.jl

@@ -20,6 +20,12 @@ function read_input()
         error("Output file already exists with this run name. Use flag -R to continue the run.")
     end
 
+    if ((params["Space"]["bc"] != 1) && (params["Space"]["bc"] != 2) && (params["sfcf"]["meas"]))
+        error("Sfcf only available with SF bc")
+    end
+
+
+
     return nothing
 end
 
@@ -66,18 +72,19 @@ function load_structs()
     CUDA.device!(parsed_args["G"])
     global lp = SpaceParm{4}(tuple(params["Space"]["size"]...), tuple(params["Space"]["blocks"]...),params["Space"]["bc"], (0,0,0,0,0,0))
     global gp = GaugeParm{Float64}(SU3{Float64},params["HMC"]["beta"],params["HMC"]["c0"],(params["HMC"]["cG"],0.0),(0.0,0.0),lp.iL);
-    global dws = DiracWorkspace(SU3fund,2,Float64,lp);
+    global dws = DiracWorkspace(SU3fund,params["AD"]["nder"]+1,Float64,lp);
     global ymws = YMworkspace(SU3,Float64,lp);
     global int = wfl_rk3(Float64, params["Frontflow"]["epsilon"], params["Backflow"]["tol"])
     global intsch = omf4(Float64,params["HMC"]["eps"], params["HMC"]["ns"]);
     global fernames = [replace(k, "Fermion" => "") for k in keys(params) if startswith(k, "Fermion")]
     # SU3 assumed in the constructor
-    global dpar = [DiracParam{Float64}(
-    Series{Float64,2}(( (1/(2*params["Fermion"*f]["kappa"])) - 4, 1.0 )),
-    Series{Float64,2}(( params["Fermion"*f]["csw"], 0.0 )),
+    global Nder = params["AD"]["nder"]+1
+    global dpar = [DiracParam{Float64,Nder}(
+    (1/(2*params["Fermion"*f]["kappa"])) - 4,
+    params["Fermion"*f]["csw"],
     ntuple(i -> exp(((i!=4)*im*params["Fermion"*f]["theta"]/lp.iL[i]) + ((i==4)*im*params["Fermion"*f]["theta_t"]/lp.iL[i])), 4),
-    Series{Float64,2}(( 0.0, 0.0 )),
-    Series{Float64,2}(( params["Fermion"*f]["ct"], 0.0 )) ) for f in fernames];
+    0.0,
+    params["Fermion"*f]["ct"], mder = (Nder > 1) ) for f in fernames];
 
     global flow_times = params["Backflow"]["Flow_times"]
     global MCid = 0
@@ -94,7 +101,7 @@ end
 
 function write_log()
 
-    println(log_file,"Running Ferflow.jl adfer v1.0 by ", params["Run"]["user"],". Name of the run: ",params["Run"]["name"])
+    println(log_file,"Running Ferflow.jl adfer v1.1 by ", params["Run"]["user"],". Name of the run: ",params["Run"]["name"])
     println(log_file,"")
 
     print(log_file,"Calling: ")
@@ -110,6 +117,7 @@ function write_log()
     println(log_file,"Parameters:")
     println(log_file,"Lattice size:        ", lp.iL)
     println(log_file,"Boundary conditions: ", params["Space"]["bc"])
+    println(log_file,"Number of derivatives: ", params["AD"]["nder"])
     for f in fernames
         println(log_file,"Fermion"*f*" parameters:")
         println(log_file,"kappa              = ", params["Fermion"*f]["kappa"])
@@ -137,14 +145,14 @@ end
 
 function save_data()
 
-    ihdr = [convert(Int32,1740655571)]
+    ihdr = [convert(Int32,1742298500)]
     fname = "./output/"*params["Run"]["name"]*".bdio"
 
     if isfile(fname)
         fb = BDIO_open(fname, "a")
         println(log_file,"\nAppending output to "*fname*"\n")
     else
-        fb = BDIO_open(fname, "w", "BDIO output from ferflow.jl adfer v1.0")
+        fb = BDIO_open(fname, "w", "BDIO output from ferflow.jl adfer v1.1")
         println(log_file,"Creating new BDIO output file "*fname)
 
         BDIO_start_record!(fb, BDIO_BIN_GENERIC, 14)
@@ -169,37 +177,33 @@ function save_data()
         BDIO_write!(fb, [convert(Int32,params["Backflow"]["tsource"])])
         length(flow_times) > 0 ? BDIO_write!(fb, flow_times) : nothing
 
+        BDIO_write!(fb, [convert(Int32,params["AD"]["nder"])])
+
+        if params["sfcf"]["meas"]
+            BDIO_write!(fb, [convert(Int32,1)])
+        else
+            BDIO_write!(fb, [convert(Int32,0)])
+        end
+
         BDIO_write_hash!(fb)
     end
 
 
     BDIO_start_record!(fb, BDIO_BIN_GENERIC, 8, true)
 
-    for f in 1:length(dpar)
-        for noi in 1:params["Frontflow"]["N_noise"]
-            BDIO_write!(fb, getindex.(getfield.(pp_corr_t0[:,noi,f], :c), 1))
-            BDIO_write!(fb, getindex.(getfield.(ap_corr_t0[:,noi,f], :c), 1))
-            BDIO_write!(fb, getindex.(getfield.(pphat_t0[:,noi,f], :c), 1))
-            BDIO_write!(fb, getindex.(getfield.(pptilde_t0[:,noi,f], :c), 1))
-            for fl in 1:params["Frontflow"]["nsteps"]+1
-                BDIO_write!(fb, getindex.(getfield.(pp_corr_t[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(ap_corr_t[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(pphat_t[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(pptilde_t[:,noi,fl,f], :c), 1))
-            end
-        end
-    end
-    for f in 1:length(dpar)
-        for noi in 1:params["Frontflow"]["N_noise"]
-            BDIO_write!(fb, getindex.(getfield.(pp_corr_t0[:,noi,f], :c), 2))
-            BDIO_write!(fb, getindex.(getfield.(ap_corr_t0[:,noi,f], :c), 2))
-            BDIO_write!(fb, getindex.(getfield.(pphat_t0[:,noi,f], :c), 2))
-            BDIO_write!(fb, getindex.(getfield.(pptilde_t0[:,noi,f], :c), 2))
-            for fl in 1:params["Frontflow"]["nsteps"]+1
-                BDIO_write!(fb, getindex.(getfield.(pp_corr_t[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(ap_corr_t[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(pphat_t[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(pptilde_t[:,noi,fl,f], :c), 2))
+    for der in 1:Nder
+        for f in 1:length(dpar)
+            for noi in 1:params["Frontflow"]["N_noise"]
+                BDIO_write!(fb, getindex.(getfield.(pp_corr_t0[:,noi,f], :c), der))
+                BDIO_write!(fb, getindex.(getfield.(ap_corr_t0[:,noi,f], :c), der))
+                BDIO_write!(fb, getindex.(getfield.(pphat_t0[:,noi,f], :c), der))
+                BDIO_write!(fb, getindex.(getfield.(pptilde_t0[:,noi,f], :c), der))
+                for fl in 1:params["Frontflow"]["nsteps"]+1
+                    BDIO_write!(fb, getindex.(getfield.(pp_corr_t[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(ap_corr_t[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(pphat_t[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(pptilde_t[:,noi,fl,f], :c), der))
+                end
             end
         end
     end
@@ -211,31 +215,32 @@ function save_data()
         BDIO_write!(fb, Qt[fl,:])
     end
 
-    for f in 1:length(dpar)
-        for noi in 1:params["Backflow"]["N_noise"]
-            for fl in 1:length(flow_times)
-                BDIO_write!(fb, getindex.(getfield.(pp_corr_tfl[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(ap_corr_tfl[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond_cfl[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond2[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond2_cfl[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(ChiDchi[:,noi,fl,f], :c), 1))
-                BDIO_write!(fb, getindex.(getfield.(ChiDchi_cfl[:,noi,fl,f], :c), 1))
+    for der in 1:Nder
+        for f in 1:length(dpar)
+            for noi in 1:params["Backflow"]["N_noise"]
+                for fl in 1:length(flow_times)
+                    BDIO_write!(fb, getindex.(getfield.(pp_corr_tfl[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(ap_corr_tfl[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(Quark_cond[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(Quark_cond_cfl[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(Quark_cond2[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(Quark_cond2_cfl[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(ChiDchi[:,noi,fl,f], :c), der))
+                    BDIO_write!(fb, getindex.(getfield.(ChiDchi_cfl[:,noi,fl,f], :c), der))
+                end
             end
         end
     end
-    for f in 1:length(dpar)
-        for noi in 1:params["Backflow"]["N_noise"]
-            for fl in 1:length(flow_times)
-                BDIO_write!(fb, getindex.(getfield.(pp_corr_tfl[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(ap_corr_tfl[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond_cfl[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond2[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(Quark_cond2_cfl[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(ChiDchi[:,noi,fl,f], :c), 2))
-                BDIO_write!(fb, getindex.(getfield.(ChiDchi_cfl[:,noi,fl,f], :c), 2))
+
+    if params["sfcf"]["meas"]
+        for der in 1:Nder
+            for f in 1:length(dpar)
+                BDIO_write!(fb, getindex.(getfield.(fp_sf[:,f], :c), der))
+                BDIO_write!(fb, getindex.(getfield.(fa_sf[:,f], :c), der))
+                BDIO_write!(fb, getindex.(getfield.(kv_sf[:,f], :c), der))
+                BDIO_write!(fb, getindex.(getfield.(kt_sf[:,f], :c), der))
+                BDIO_write!(fb, [getindex(getfield(f1_sf[f], :c), der),])
+                BDIO_write!(fb, [getindex(getfield(k1_sf[f], :c), der),])
             end
         end
     end