Problem 1366

Summary: Bug in G4NeutronHPContAngularPar::Sample()
Product: Geant4 Reporter: Liam Russell <lfr.russell>
Component: processes/hadronic/models/neutron_hpAssignee: tkoi
Status: RESOLVED FIXED    
Severity: normal    
Priority: P5    
Version: 9.5   
Hardware: All   
OS: All   

Description Liam Russell 2012-10-19 20:55:08 CEST
The code for G4NeutronHPContAngularPar::Sample() has a bug located in the for loop at line 550 in G4NeutronHPContAngularPar.cc (version 4.9.5p01).  

This code is only used when dealing with outgoing angles that are saved in the G4NDL format as a table of cos(theta) values, with corresponding probabilities, that are indexed by the incoming neutron energy.  Using the default G4NDL data (G4NDL 3.14 or 4.0), this is not a problem because the outgoing angle is described by legendre polynomials.  However, when using different neutron data, the following code causes a runtime error:


	for(i=0; i<nAngularParameters; i++)
	{
		theBuff1.SetX(i, theAngular[it-1].GetValue(i));
		theBuff1.SetY(i, theAngular[it-1].GetValue(i+1));
		theBuff2.SetX(i, theAngular[it].GetValue(i));
		theBuff2.SetY(i, theAngular[it].GetValue(i+1));
		i++;
	}


The first problem is that "i" increments by two for every iteration of the loop.  Therefore, the SetX() and SetY() functions try to set elements [0,2,4,...] of the buffer objects "theBuff1" and "theBuff2".  This causes the runtime error since elements are skipped in these buffers.  Additionally, the angular data, "theAngular", has the following format

	theAngular = {P_{energy},cos_1(theta),P_1,cos_2(theta),P_2,...}

where "P_{energy}" is the probability of the outgoing neutron energy, "cos_i(theta)" is the ith outgoing angle, and "P_i" is the probability of the ith outgoing angle.  Therefore, the GetValue() function should start at 1, not 0.

Fixing this bug is relatively simple, and the solution I used is shown below


	for(i=0, j=1; i<nAngularParameters; i++, j+=2)
	{
		theBuff1.SetX(i, theAngular[it-1].GetValue(j));
		theBuff1.SetY(i, theAngular[it-1].GetValue(j+1));
		theBuff2.SetX(i, theAngular[it].GetValue(j));
		theBuff2.SetY(i, theAngular[it].GetValue(j+1));
	}


In the fixed code, a separate index "j" is used to get the cos(theta) data from theAngular.
Comment 1 tkoi 2012-10-20 04:03:22 CEST
Thank you for reporting this problem.

Could you kindly tell us which data library (like JENDL3.3 JEF3.1) you got this problem
and if possilbe let us know isotope name(s) as well

Thank you very much!

Tatsumi
Comment 2 tkoi 2012-11-07 18:16:56 CET
Hi 

The codes are edited along with your suggestion. This will be available coming release of v9.6 (2012Dec)
However, even G4NDL4.2 the part is not really used, therefore if you have unexpected behavior, then please let us know.