Hi all.


CCTK_REAL masscenter_point type = GF Timelevels = 3 tags='Prolongation="none" tensortypealias="Scalar" '
{

  masscenterX_point
  masscenterY_point

}"mass center to be used in the global reductions"


CCTK_REAL global_masscenter type = SCALAR
{

 masscenterX
 masscenterY

}"mass center coordinates"


What I want to do is evaluate in every grid point the values for masscenterX_point and masscenterY_point (I do this in the function my_analysis_point_calc) and then integrate them (using CCTK_Reduce) over the entire domain and save the result in the scalar variables masscenterX and masscenterY (and I do it in the function my_analysis_reduce). So far everything seems to work well, if I plot the masscenter it has reasonable values. The problem happens when, in the same thorn, I want to use the values of masscenterX and masscenterY in another function ( my_analysis_modes_point_calc) to evaluate some other variable (let's call it mode_m1) point by point and then integrate it (in my_analysis_modes_reduce).


At the initial time I get a reasonable value for mode_m1, but the second time I print it ( after 12 timesteps) I get nan, the third time I print it (after 24 timesteps) I get a reasonable value and after that it seems to be fine for the rest of the simulation. When I print masscenterX and masscenterY they have reasonable values at every iteration.  This is the output:
"

#1:iteration 2:time 3:data
#data columns: 3:masscenterX 4:masscenterY
0 0 2.04233950307922e-07 1.43009709532618e-07
12 0.6 1.64952865679123e-07 1.15302233769628e-07
24 1.2 2.59441645250822e-08 1.8122269479532e-08

#1:iteration 2:time 3:data
#data columns: 3:mode_m1
0 0 5.09771568069084e-09
12 0.6 -nan
24 1.2 5.98971400034287e-10

"

I thought that the problematic line in my code could be this one:

            phi = atan((yy-*masscenterY)/(xx-*masscenterX));

            if(((xx-*masscenterX)<0.0)){
                phi = atan((yy-*masscenterY)/(xx-*masscenterX)) + 3.14159265358979323846;
            }

where I try to use the values of masscenterX and masscenterY inside the function "my_analysis_modes_point_calc". If I change that part not using masscenterX and masscenterY I don't get nans anymore in mode_m1.

After doing different tests, I actually found out that if I try to plot scalar output every 4 iterations I get also NaN in masscenter at time 0.2; this is the output in this case:
"
#1:iteration 2:time 3:data
#data columns: 3:masscenterX 4:masscenterY

0 0 2.04233950307922e-07 1.43009709532618e-07
4 0.2 -nan -nan
8 0.4 1.76472702563281e-07 1.23428120452789e-07
12 0.6 1.64952865679123e-07 1.15302233769628e-07
16 0.8 1.05057935424561e-07 7.33797284607643e-08
20 1 5.6677119219615e-08 3.93971526043537e-08
24 1.2 2.59441645250822e-08 1.8122269479532e-08
28 1.4 5.57004919317397e-08 3.90090013497981e-08
32 1.6 1.2454179781868e-08 8.57409404496548e-09
36 1.8 7.81397952778142e-08 5.4397631432516e-08

#1:iteration 2:time 3:data
#data columns: 3:mode_m1

0 0 5.09771568069084e-09
4 0.2 -nan
8 0.4 -nan
12 0.6 -nan
16 0.8 2.49105795013804e-09
20 1 -nan
24 1.2 5.98971400034287e-10
28 1.4 3.14686031011036e-10
32 1.6 2.57981958462494e-10
36 1.8 6.79137649228778e-10

"