Problems with JEFF-3.3 Cr50, Cr53, and Cr54

[HunterBelanger]

I recently processed all of JEFF-3.3, and upon trying a few benchmarks, I noticed that Cr50, Cr53, and Cr54 all had a similar problem in the generated ACE files. Each of these evaluations has an MT 5, for which the distribution is given in MT 6 as LAW 1 LANG 2 (Kalbach Mann). For each of these evaluations, I find that for at least one incident energy, the out going energy grid is not sorted. This also leads to the CDF for the outgoing energy to have negative values.

As an example, using Cr50, this problem occurs at an incident energy of 4MeV, for which I am seeing an out going energy grid of:

Eout = [0.001, 3.71661e-05, 4.4198e-05, 5.25604e-05, 6.25049e-05, 7.4331e-05, 8.83946e-05,
        0.000105119, 0.000125008, 0.00014866, 0.000176787, 0.000210236, 0.000250013,
        0.000297316, 0.000353569, 0.000420465, 0.000500018, 0.000594623, 0.000707127,
        0.000840917, 0.00100002, 0.00118923, 0.00141423, 0.0016818,
        0.002, 0.005, 0.01, 0.02, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 1.1, 1.2, 1.3, 1.4,
        1.5, 1.6, 1.7, 1.8, 1.9, 2, 2.2, 2.4, 2.6, 2.8, 3, 3.2, 3.4, 3.6, 3.8, 4]

From this grid, we see the that the exit energy grid isn't sorted due to what appears to be the very first energy point. Looking in the evaluation however, I found that the given outgoing energy grid is supposed to be

Eout = [0.001, 0.002, 0.005, 0.01, 0.02, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1,
        1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2, 2.2, 2.4, 2.6, 2.8, 3, 3.2, 3.4, 3.6, 3.8, 4]

It would appear as though the first point is the correct one, but that NJOY is adding a lot of extra points between first and second out going energy. This not only causes the out going energy grid to be unsorted, but also seems to be leading to the negative CDF values as well, as these added points also have weird PDFs.

The problems encountered with Cr53 and Cr54 appear to be nearly identical, happening with a Kalbach Mann distribution in MT 5, though the incident energy at which this occurs seems to vary.

[jchsublet]

I am afraid that this glitch will not only happen on your above examples, does not only belong to Cr JEFF-3.3 targets I mean. The outgoing energy grid is automatically extended below 0.001 in that case during NJOY2016 Acer processing following recommendations made and explained in:

Page 15 https://t2.lanl.gov/nis/publications/Bigten_Study.pdf
Page 30 https://t2.lanl.gov/nis/publications/200711a.pdf

The fact that its also does it between e1 and e2 is obviously wrong

It has already been raised or flagged in issue #131

When physics is concerned it does not make any sense to have an MT=5/MF=6 KM energy range going down to 1e-5, but when in ENDF-6 forms it apparently does ! Au royaume des aveugles les borgnes sont rois

[HunterBelanger]

Ah yes, I had not realized that this problem was related to the smoothing operation. Turning of smoothing indeed removes the exceptions being thrown by my library. As you mention in #131, it might not be wise to apply smoothing to MT5, but the addition of these points between E1 and E2 is certainly an independent bug, as you mention. I have never noticed this issue with other Kalbach Mann distributions, making me think that the problem is due to the histogram interpolation used in the mentioned Cr evaluations (LEP=1). Looking at acefc.f90 in the acelf6 subroutine, I think the issue occurs somewhere in this if statement:

NJOY2016/src/acefc.f90

Lines 7522 to 7561 in 14f1ac1

	! extend low histogram bins as sqrt(e) using log energy scale
	! only do this for outgoing neutrons with law=1, lang=2
	! only do this if there are no discrete data
	if (ismooth.gt.0.and.law.eq.1.and.lang.eq.2.and.&
	lep.eq.1.and.zap.eq.1.and.nd.eq.0) then
	fx=.8409
	ex=40
	cx=scr(7+ncyc)*scr(8)
	do while (n.gt.2)
	cxx=cx+scr(8+ncyc)*(scr(7+2*ncyc)-scr(7+ncyc))
	if (abs(cxx/scr(7+2*ncyc)**1.5-cx/scr(7+ncyc)**1.5)&
	.gt.cx/scr(7+ncyc)**1.5/50) exit
	scr(8)=(scr(8)*scr(7+ncyc)&
	+scr(8+ncyc)*(scr(7+2*ncyc)&
	-scr(7+ncyc)))/scr(7+2*ncyc)
	do ix=1,nx-2*ncyc
	scr(6+ix+ncyc)=scr(6+ix+2*ncyc)
	enddo
	cx=cxx
	nx=nx-ncyc
	n=n-1
	enddo
	write(nsyso,'('' extending histograms as sqrt(E) below'',&
	&1p,e10.2,'' MeV for E='',e10.2,'' MeV mt='',i3)')&
	scr(7+ncyc)/emev,ee,mt
	do while (scr(7+ncyc).gt.ex)
	do ix=nx,1,-1
	scr(6+ncyc+ix)=scr(6+ix)
	enddo
	scr(7+ncyc)=fx*scr(7+2*ncyc)
	scr(7+ncyc)=sigfig(scr(7+ncyc),6,0)
	val=scr(8)
	scr(8)=sqrt(fx)*val
	scr(8)=sigfig(scr(8),6,0)
	scr(8+ncyc)=(1-fx*sqrt(fx))*val/(1-fx)
	scr(8+ncyc)=sigfig(scr(8+ncyc),6,0)
	nx=nx+ncyc
	jscr=jscr+ncyc
	n=n+1
	enddo

This is the branch where histogram interpolation KM is smoothied, and it seems fairly different from the branch immediately bellow it, which handles linear interpolation. In the linear interpolation branch, there is mentioning of inserting the new points, and also re-normalizing the distribution. While I don't see comments for this in the histogram branch, that doesn't mean it isn't happening. J'avoue que je ne comprend pas trop bien ce qui se passe là dans les sources ; pour moi c'est un peu de la magie noire.