[BUG] Nucleotide tracking: issue with particle_position observable

[SimonGrall]

Describe the bug I am tracking the ith nucleotide of a DNA strand with 2 approaches that, I thought, should give identical results.

1) I generate trajectory.dat and browse through the file with some code snippet to retrieve the XYZ position of the ith nucleotide at each time step.

2) I generate a custom observable in the input file tracking the ith nucleotide, adding the following to my input file (for the 40th nucleotide in a double strand DNA here) : data_output = { print_every = 100 name = TRACK.txt col_1 = { type = step } col_2 = { type = particle_position particle_id = 39 absolute = 1 } }

When looking at the results, I get reasonable measurements for the method 1, with a position roughly coherent with ~0.3 nm/nucleotides, as expected for DNA. But not for method 2, with about 2 times less nm/nucleotides and a different distribution in space at the end. I suspect I am messing up at some point with the indexing of particles, but could not get information on this. Oxdna was compiled in 2018 for method 1, and in 2024 for method 2, if that make any difference. Aside this, all parameters are identical.

I am sorry if this is a stupid request, by I am just out of things to check now. Thank you in advance, whoever ends up helping.

Simon.

[lorenzo-rovigatti]

I'm not sure what the problem could be without having a look at the input files, at the custom script you use and at your expected/actual output. Can you add some details?

[SimonGrall]

Hello Lorenzo, Here attached the input files, last conf files and ext force used, as well as the script retrieving the ith particle (40th in that case on a 20T double strand DNA). The file track_ushort.py has as input trajectory.dat, some index counter, the list of particle to track (here only 1 element) and some inputfile of our own storing some basic information.

Thank you for having a look at this. Simon. OXDNA2018.zip

OXDNA2024.zip

[lorenzo-rovigatti]

Thank you. I had a quick look and I noticed that the external forces in the two simulations are different. Would this be enough to explain the difference you observe?

[SimonGrall]

No, the only difference is the position of the harmonic trap, which is roughly the center of the simulation box's XY plane. That should not have an impact on it. I actually made small video snippets using cogli, and the DNA look alike in both, attached far away from the edges of the box.

[lorenzo-rovigatti]

OK. What do you want to see by tracking the position of a single nucleotide? I thought it had something to do with the average position of that specific nucleotide, which should be affected by the external forces, so that if those vary than you would expect to see some differences between the two simulations.

EDIT: the topology files are missing

[SimonGrall]

I updated the ZIP files to include the generated.top files.

[lorenzo-rovigatti]

Thank you, but there are still information missing. How do I analyse the results with track_ushort.py? What do I have to compare between the two simulations?

[SimonGrall]

Sorry for the oversights. Track ushort outputs a track of only the followed nucleotide. This should be compared with the output of the observable particle_position. track_ushort.py takes as inputs (obtained through sys.argv): trajectory.dat, i, [40], homemade_inputfile This will create a file track-Mol1-i.txt to be compared with the output of the observable, for the positions of the tracked nucleotide. (Here tracking the 40th nucleotide).

As for thet homemade_inputfile, it should contain the following:

nbbases=40 nbbasesOpposite=0 nbmolecules=1 nbmoleculesOpposite=0

[lorenzo-rovigatti]

The two simulations start from different initial configurations, use different external forces and don't set the seed of the random number generator. Is there any reason why you expect the two trajectories to be the same?

[SimonGrall]

They simulate the same DNA strand, with the same forces applied on it. I expect the average position over lets say 10^9 steps to be the same with respect to the position of the harmonic trap (i.e. where the DNA stand is attached). The only difference in forces is the 3D coordinate of the origin, i.e. the location of the DNA strand within the box. The boxes have different sizes, and should be roughly in the middle, far enough from the edges. As long as the simulation goes on correctly, I thought that the dynamic of a double stranded DNA attached at one end would be the same no matter the box coordinates, no?

[lorenzo-rovigatti]

I used your script on a trajectory generated with a simulation where I also used the particle_position observable and the two methods yield exactly the same positions. I can try to help you if you can manage to prepare an example where the bug can be reproduced with short-ish simulations, and if you tell me exactly what I have to do to reproduce it.