Geometry error connected to Lattice

Geometry and Materials
1

Dear colleagues,
in my simulation I have an inhomogeneous hadronic calorimeter, made of different layers, each layer being divided in a scintillator part and an absorber part. I modeled quite accurately the first layer, then I make use of the LATTICE feature to replicate it.

After simulating ~ 10k primaries, the simulation crashed with the following error:

Geofar: Particle in region     9 (cell #   20)
  in position  -2.706968472E+01  2.306588177E+01  4.151000208E+02
  is now causing trouble, requesting a step of  4.821616320E+00 cm
  to direction  3.129836027E-01  9.497585294E-01 -1.808674501E-05
  end position -2.556059787E+01  2.764525300E+01  4.150999336E+02
  R2:  3.556406519E+01 R3:  4.166207268E+02 cm  error count:  0
  X*U   (2D):  1.343465050E+01 X*U   (3D):  1.342714269E+01 cm
  X*UOLD(2D):  1.736273063E+01 X*UOLD(3D): -2.990530501E+00 cm
  Kloop:  848497963, Irsave:      9, Irsav2:     12, error code: -33 Nfrom:  5000
  old direction  1.935300256E-01  9.798683452E-01 -4.903218527E-02, lagain, lstnew, lsense, lsnsct  F  F  F  F
  Particle index    7 total energy   1.440058354E-03 GeV  Nsurf    0
 Try again to establish the current region  moving the particle of a  1.218815943E-08 long step
 We succeeded in saving the particle:  current region is n.     9 (cell #   20)
 Geofar: Particle in region     9 (cell #   20)
...
...
 Abort called from FLKAG1 reason TOO MANY ERRORS IN GEOMETRY Run stopped!
 STOP TOO MANY ERRORS IN GEOMETRY

Indeed, the point ( -2.706968472E+01, 2.306588177E+01, 4.151000208E+02) where the error is reported is indeed at the interface between two lattice regions.
I tried to use FLAIR to check if any error was there, and after zooming it (a lot!) in that region, I saw something strange, as reported below.

image
image1532×932 156 KB

Zoom:

image
image1634×630 164 KB

As you see, there is a white region between Z=415.1 cm and Z=415.1+1E-8 cm, that I cannot understand.

I attach here the input file that I am currently using:

na64_2022B_h_50GeV.inp (324.0 KB)

I am currently using FLUKA 4.3.3 - I’ll soon update to the newest FLUKA version, but it seems to me, looking at the changelog, that this is not related to the version upgrade (maybe I am wrong?)

Thanks,
Bests,
Andrea

3

Dear colleagues,
I’d like to follow up on this. I spent some time to improve the geometry, so that the number of touching bodies / surfaces is drastically reduced. In particular, at this moment with FLAIR I do not get any warning regarding this for the hadronic calorimeter.

Still, I get the same error - note that the point where the error is reported is very close to the boundary between two lattice regions (Z=3.913000191E+02, the boundary is at Z=3.913E+02).

 Geofar: Particle in region     9 (cell #   13)
  in position  -4.059446966E+01  1.990311402E+00  3.913000191E+02
  is now causing trouble, requesting a step of  7.144723851E+00 cm
  to direction  5.441959219E-01  8.389581623E-01 -2.448448904E-05
  end position -3.670634007E+01  7.984435794E+00  3.912998442E+02
  R2:  4.064323198E+01 R3:  3.934051058E+02 cm  error count:  0
  X*U   (2D): -2.042155684E+01 X*U   (3D): -2.043113762E+01 cm
  X*UOLD(2D): -2.042155684E+01 X*UOLD(3D): -2.043113762E+01 cm
  Kloop:  187058584, Irsave:      9, Irsav2:      9, error code: -33 Nfrom:  1001
  old direction  5.441959219E-01  8.389581623E-01 -2.448448904E-05, lagain, lstnew, lsense, lsnsct  F  F  F  F
  Particle index    8 total energy   9.402900280E-01 GeV  Nsurf    0
 Try again to establish the current region  moving the particle of a  1.218815943E-08 long step
 We succeeded in saving the particle:  current region is n.     9 (cell #   13)
 Geofar: Particle in region     9 (cell #   13)
  in position  -4.059446965E+01  1.990311412E+00  3.913000191E+02
  is now causing trouble, requesting a step of  7.144723839E+00 cm
  to direction  5.441959219E-01  8.389581623E-01 -2.448448904E-05
  end position -3.670634007E+01  7.984435794E+00  3.912998442E+02
  R2:  4.064323198E+01 R3:  3.934051058E+02 cm  error count:  1
  X*U   (2D): -2.042155683E+01 X*U   (3D): -2.043113761E+01 cm
  X*UOLD(2D): -2.042155683E+01 X*UOLD(3D): -2.043113761E+01 cm
  Kloop:  187058593, Irsave:      9, Irsav2:      9, error code: -33 Nfrom:  1001
  old direction  5.441959219E-01  8.389581623E-01 -2.448448904E-05, lagain, lstnew, lsense, lsnsct  F  F  F  F
  Particle index    8 total energy   9.402900280E-01 GeV  Nsurf    0
 Try again to establish the current region  moving the particle of a  1.218815943E-08 long step
 We succeeded in saving the particle:  current region is n.     9 (cell #   13)
...
 Abort called from FLKAG1 reason TOO MANY ERRORS IN GEOMETRY Run stopped!
 STOP TOO MANY ERRORS IN GEOMETRY

I upload the new input file.
na64_2022B_h_50GeV.inp (310.3 KB)

Thanks,
Bests,
Andrea

4

Caro Andrea,

you are right.

I believe that the zoom above is just affected by a display artifact due to the excessively tiny scale.

Also, you geometry appears to be formally correct, already in its previous version.

The problem you meet is related to an intrinsic tracking limitation seemingly triggered by the double (region and lattice) boundary. I did not meet any geometry error after 20k primaries (for 50 GeV/c negative pions without custom source routine), but - in case the issue is not rare and systematically affects your simulation runs - I’d consider a couple of tentative measures (I did not test myself since I got no error anyway).

One is the use of the accuracy parameter in the GEOBEGIN card (WHAT(2)), which you could try tightening to 1E-05, from the default value of 1.22E-04 obtained in your case (as printed in the output file).
The other is the application of the inverse of the ROT-DEFI transformation used in your LATTICE card to the bodies defining the LATTICE region. What I mean is, instead of defining the H0PlZ002 plane at z = 391.3 cm, to define it at z = 384.5 cm (as in the HCL0BxIn RPP) and put it under a $Start_transform -tH0L002
directive (this is meant to increase accuracy through self-consistency).

5

Ciao Francesco,
indeed, the effect I see is quite “rare”. By this I mean that, running the code on CERN batch, and having 50 jobs each with 100k events, I see few of them failing (~5).

I tried with the first suggestion (set accuracy to 1E-05), but it looks this did not help - now ALL the job failed, again with the geometry error; in my out file I also see many, many lines such as:

IRLTGG=  17 SB=-0.28894757D+02-0.98811642D+00 0.38790066D+03 UB=-0.79718309D+00 0.60373746D+00-0.44476743D-03
 IR=   9 XB=-0.28894757D+02-0.98811641D+00 0.40830066D+03 WB=-0.79718309D+00 0.60373746D+00-0.44476743D-03 DIST= 0.14758661D+01

Unfortunately, the other suggestion is not easy for me to implement, since I have my own set of python scripts to produce the geometry, and it is not trivial to change them for this fix.

May I ask you if the new lines point to some fix/solution for the problem?

Thanks,
Andrea

6

Thanks for the quick feedback. I’d ask you to test then 1E-03 as accuracy parameter.
(A 10% fraction of aborting jobs is definitely disturbing, but in this case is not a sign of wrong geometry input).

7

Dear Francesco,
with the accuracy set to 1E-03, all the cycles work properly. May I ask you if I should then use this value for all simulations in which my geometry has an overall size similar to that in the example before, and similarly the “granularity” of the different elements is similar?

Is there any drawback in this change, from the default value of ~1E-04?

Thanks,
Bests,
Andrea

8

Good, go ahead with it then.
The need for this user’s adaptation is not just related to the input geometry size, rather to the lattice use and respective input specifications too.
This change is increasing the tolerance for (hardly avoidable) numerical inconsistencies in the boundary identification, and in general should be adopted with care.
In this case, it looks like a reasonable solution to the reported problem, considering that the accuracy value is in units of 1.E-6 cm (so resulting here in 0.01 nm).