add support for performing 2-flavour and BSM inversions from external program via lib_wrapper

[kostrzewa]

@pittlerf @Marcogarofalo

I've adjusted where things are initialised, such that all BSM operators are now also usable when inversions are driven via tmLQCD_invert in wrapper/lib_wrapper.c. In the process I've also changed how the scalar field index is kept track of (at least in invert.c by externalizing this to a function in io/scalar_read.c).

I only touched contractions_BSM in a few places. Could you please check that it still works?

I was a bit surprised that you have zero-padding of the scalar field index in contractions_BSM, while in the past it was not zero-padded. Does this mean that we've agreed on zero-padding now?

[kostrzewa]

I'm not sure if I want to keep this interface like this. Perhaps it makes sense after all that the scalar index is computed externally and simply passed (just like for the gauge field). The reasoning for this was to keep the index calculation in the load_next_scalar_field function in io/scalar_read.c to avoid repetition, but I might just externalize this into a separate function which can then be called also by the external program.

[kostrzewa]

Yeah, I've switched away from this to a solution with two functions.

[kostrzewa]

This is a major problem unfortunately. I don't really want to call init_bispinor_field from every executable that we have (because it is a bona-fide dependency of the BSM operators). At the same time, this does not really fit into the way that tmLQCD uses the bispinor fields... Not quite sure yet what a good solution would be.

[kostrzewa]

I've now extended the support for the smeared gauge field in a more general fashion. When SmearedGaugeConfigInputFile is set, memory is allocated for a separate smeared gauge field. I've adjusted invert.c to make use of this, such that the BSM3 operator can actually be inverted by the invert program.

I still need to take care of interfacing the smearing of the scalar field.

[kostrzewa]

The two inversions performed in op_invert for the BSM operators need to be controlled via some kind of interface such that it's possible to get either D^{-1} or [D^\dagger]^{-1}. Of course, this was never intended in tmLQCD so it will take some creativity to get there.

It's also quite silly that op_invert contains the loop over source flavours. This was fine when the default was to write out the propagators to disk, but it's not a good idea at this stage...

[pittlerf]

op_invert

yeah, maybe a boolean inside the parameters, whether it should be daggered or not?

[pittlerf]

I actually have problems with the parts of
f( (i = read_gauge_field(conf_filename,g_smeared_gauge_field)) !=0)
in contractions_BSM

[kostrzewa]

I actually have problems with the parts of
f( (i = read_gauge_field(conf_filename,g_smeared_gauge_field)) !=0)
in contractions_BSM

can you explain?

[pittlerf]

I actually have problems with the parts of
f( (i = read_gauge_field(conf_filename,g_smeared_gauge_field)) !=0)
in contractions_BSM

can you explain?

yeah, sorry I got segmentation fault, when trying to read

Trying to read smeared gauge field from file conf_smeared.0000 in double precision.

Constructing LIME reader for file conf_smeared.0000 ...

WARNING: IO CHECKS HAVE BEEN DISABLED

[cyclamen:17843] *** Process received signal ***
[cyclamen:17843] Signal: Segmentation fault (11)
[cyclamen:17843] Signal code: Address not mapped (1)
[cyclamen:17843] Failing at address: (nil)
[cyclamen:17843] [ 0] /lib64/libpthread.so.0(+0xf5d0)[0x7f415c4215d0]
[cyclamen:17843] [ 1] ./contractions_BSM[0x6bc196]
[cyclamen:17843] [ 2] ./contractions_BSM[0x6b540b]
[cyclamen:17843] [ 3] ./contractions_BSM[0x4035ba]
[cyclamen:17843] [ 4] /lib64/libc.so.6(__libc_start_main+0xf5)[0x7f415c0673d5]
[cyclamen:17843] [ 5] ./contractions_BSM[0x403039]
[cyclamen:17843] *** End of error message ***
Segmentation fault

I have to look at it, maybe something in the initializiation changed

[kostrzewa]

You have to set SmearedGaugeConfigInputFile. When this has non-zero length, init_gauge_field is called with smeared = 1, which causes memory for the smeared gauge field to be allocated.

[pittlerf]

You have to set SmearedGaugeConfigInputFile. When this has non-zero length, init_gauge_field is called with smeared = 1, which causes memory for the smeared gauge field to be allocated.

thank you

[kostrzewa]

op_invert
yeah, maybe a boolean inside the parameters, whether it should be daggered or not?

sure, I just don't know yet how to integrate this elegantly with the other behaviours

[pittlerf]

yeah, basically in the other case we do not need the dagger at all

[kostrzewa]

I think there's a way to handle multiplying the source by the operator before the inversion in prepare_source, I'm just not 100% sure yet if it can be exposed nicely all the way to lib_wrapper. Maybe I'll have some time 2mw to implement it.

[kostrzewa]

@urbach I would like to modify the logic for dealing with two-flavour operators in op_invert. On the one hand, this will make it easier to deal with the needs of the BSM operators, but it will also translate to being able to drive ND inversions via lib_wrapper in a straightforward manner. We can then backport these changes into the quda and eventually master branches for doing ND inversions comfortably.

1. remove the loop over source flavours in op_invert
2. add source flavour loop in invert.c (if using invert)
3. add source flavour index argument to prepare_source

Would value your input if you see issues.

[kostrzewa]

In 29df069 I've implemented what I came up with to solve our problem in a (hopefully) minimally-invasive way. I've broken contractions_BSM in the process but will attempt to fix it today or 2mw. Basically I've moved the two loops over flavour and non-dagger, dagger out of op_invert to be done externally. To make this work, I've added a flav_idx argument to prepare_source and op_invert and a dagger argument to op_invert.

[kostrzewa]

@pittlerf @Marcogarofalo do I understand it correctly that the scalar field used for the propagators (the one stored in g_scalar_field) is already pre-smeared? Or in other words, it is smeared externally and then loaded via read_scalar_field[_parallel]?

[pittlerf]

Yes, hypercubically smeared.

[kostrzewa]

Thanks.

[kostrzewa]

@pittlerf I've got the crosscheck for any combination of gamma and tau structures at source and sink ready. What I don't have yet is the pseudo-Wilson density but I'll likely be able to add that in the coming weeks too. I'll send my correlators to the list in a few minutes.

[pittlerf]

Thank you very much

[kostrzewa]

@pittlerf sorry I just noticed that I had not actually included the axial vector ;) I'll give an update here soon.

[kostrzewa]

Current state of the cross-check

See also https://github.com/HISKP-LQCD/bsmsimul/commit/018a355e2fe1543f0e0ebc8df5fadfaac429467a and https://github.com/HISKP-LQCD/bsmsimul/commit/907e6e480da306ce8ec80a8777e40ebe30e9ad32

[kostrzewa]

@pittlerf what does "LOCALCURRENTTAU0P0" mean? More specifically, what is TAU0 ?

[Marcogarofalo]

I suspect that there was a typo that shifted the enumeration so
LOCALCURRENTTAU0P0 should be renamed LOCALCURRENTTAU1P1
LOCALCURRENTTAU1P1 should be renamed LOCALCURRENTTAU2P2
LOCALCURRENTTAU2P2 should be renamed LOCALCURRENTTAU3P3

[pittlerf]

I suspect that there was a typo that shifted the enumeration so LOCALCURRENTTAU0P0 should be renamed LOCALCURRENTTAU1P1 LOCALCURRENTTAU1P1 should be renamed LOCALCURRENTTAU2P2 LOCALCURRENTTAU2P2 should be renamed LOCALCURRENTTAU3P3

Hi, yes, this is a typo,thanks for noticing it, I am correcting now

[kostrzewa]

Thanks, I figured but wanted to make sure that I had not misunderstood (tau0 could be a shorthand for the identity to build the iscoscalar combination, for example).

This is the current state of the cross-check:

[pittlerf]

ah, ok, thanks there are problems with tau2 and tau3, but tau1 is correct

[kostrzewa]

ah, ok, thanks there are problems with tau2 and tau3, but tau1 is correct

it seems so but I wonder how that can happen

[kostrzewa]

ah, ok, thanks there are problems with tau2 and tau3, but tau1 is correct

it seems so but I wonder how that can happen

@pittlerf should we compare some intermediate fields or something like that?

[pittlerf]

ah, ok, thanks there are problems with tau2 and tau3, but tau1 is correct

it seems so but I wonder how that can happen

@pittlerf should we compare some intermediate fields or something like that?
Hi Batek, yes we could do that, I have tried to build nyom on cyclamen,
but I got an error regarding yamlcpp:
`CMake Error at /cyclamen/home/fpittler/build/yamlcpp/yaml-cpp-config.cmake:11 (include):
include could not find load file:

/cyclamen/home/fpittler/build/yamlcpp/yaml-cpp-targets.cmake

Call Stack (most recent call first):
CMakeLists.txt:107 (find_package)`

[kostrzewa]

I have tried to build nyom on cyclamen, but I got an error regarding yamlcpp:

yaml-cpp (https://github.com/jbeder/yaml-cpp) is a dependency and one either has to have it installed somewhere where CMake looks for cmake files by default (for example, in /usr, such that the cmake files are under /usr/share/cmake/yaml-cpp) or one sets -DCMAKE_SYSTEM_PREFIX_PATH=[...] to the directory where one has installed it locally. For example, I do:

-DCMAKE_SYSTEM_PREFIX_PATH="/nfshome/bartek/build/yaml-cpp/install_dir/share/cmake/yaml-cpp"

in my cmake call to configure nyom.

Hi Bartek, yes we could do that

So how should we proceed? I can relatively easily write out any intermediate quantity, but I'm not sure whether this makes so much sense or whether we should instead discuss in terms of formulae. I'm a bit surprised that we have compatible results for <P^2 P^2> (up to sign) but see a difference when the sink gamma matrix is changed to \gamma_0 \gamma_5 instead of just \gamma_5...