@ChrisWhite The terminal output tells you what the expectation value of the Polyakov loop, and the expectation of its absolute value are. Same is stored in the filename you input, and in #number.dat there are the values of each sweep for #number gridsize
@ACuriousMind Now I have very little understanding of what you're trying to do, but as everyone else is looking at the code, I had to too. The grid was very small: are you sure you are supposed to see a phase transition at these sizes?
I'm not fully convinced, especially since the reference papers he sent me report their data completely differently, talking about "renormalized values of the loop" and so on, but they also say their gridsize was a whooping 4, and I can't see any non-trivial behaviour in my simulation at all.
Also, I can't make heads or tails of their graphs...they label the axis by $\mathrm{tr}(W)$, which means it should lie between -2 and 2. Then the axis goes from 0 to 16...
@alarge Ah, yes, that's the 100 lines I said I could have less if I'd stored the links differently. I tested it with a fixed lattice configuration where I computed some of the staples by hand - at least in that test case, it did what it was supposed to
Also you might want to add "-O3 -march=native" to your compiler options (for speed, but apparently the algo is so fast this is not really necessary). Your code also does waaaay too much copying, but if performance is not really a concern, it's ok I suppose.
@alarge curious to look. missed where the code was posted...? maybe some of this can be chalked up to difficulty in replicating others results. aka "open science".
@alarge I'm sure experienced programmers would write this thing differently, but it's the first longer piece of code I've written in years, and never really in C before, so the efficient use of pointers is not to be expected ;)
@ACuriousMind In StapleInDir, are you not supposed to evaluate anything at "current lattice point"? I mean if STAPDIR is XDIR, you look up stuff at x+1 and x-1, but nothing at x.
@alarge The lookup for the thing at current lattice point is independent of the directions switches, and is done before the case...switch in Aux = GetLink(t,x,y,z,StapDIR); SF = MatHermConj(&Aux);
Also...you actually looked closer at that staple function? I really don't want to keep you from doing something...more fun
@ChrisWhite My dad learned programming on a punch-card machine at midnight (when the undergrads were allowed to run code). He says you can never have enough braces.