*cpl njoy,reconr,broadr,unresr,heatr,thermr,groupr,errorr *cpl covr,moder,ccccr,matxsr,acer,viewr,purr,gaspr,wimsr *cpl plotr,leapr *set sw *ident up1 */ reconr -- 04jul98 */ fix an error made in enhancing the equality tests in */ reconr when building njoy97. this error causes the code */ to miss discontinuities, such as those at the ends of */ resonance ranges, where it normally tries to put in points */ just below and just above the discontinuity. problem */ noticed by little (lanl). *d reconr.1456 if (abs(srnext-sr).lt.small*sr) go to 260 */ reconr -- 04jul98 */ change the wide constant to prevent infinite loops */ (this showed up for u-232). *d reconr.2808 data wide/3.01d0/ *d reconr.2814 data wide/3.01/ *d reconr.3018 data wide/3.01d0/ *d reconr.3024 data wide/3.01/ */ reconr -- 02jul98 */ fix up some mistakes in the use of the cmplx intrinsic */ function, fix some variables left as real when they should */ be double, and make sure that the complex variables go to */ double precision properly when requested. *d reconr.2582,2583 dimension qre(4),gc(4),sigc(4) *if sw real*8 k,kfac complex*16 one,eye,rfcnn,slsjnn,phi0,phi,r0,exphi *else real k,kfac complex one,eye,rfcnn,slsjnn,phi0,phi,r0,exphi *endif *d reconr.2588 data zero,uno/0.d0,1.d0/ *d reconr.2592 data zero,uno/0.,1./ *d reconr.2595,2596 *if sw eye=dcmplx(zero,uno) one=dcmplx(uno,zero) *else eye=cmplx(zero,uno) one=cmplx(uno,zero) *endif *d reconr.2668 *if sw phi0=dcmplx(phir,zero) *else phi0=cmplx(phir,zero) *endif *d reconr.2673 *if sw rfcnn=dcmplx(zero,zero) *else rfcnn=cmplx(zero,zero) *endif *d reconr.2729 *if sw r0=dcmplx(zero,zero) *else r0=cmplx(zero,zero) *endif *d reconr.2741 *if sw r0=dcmplx(rr0,ri0) *else r0=cmplx(rr0,ri0) *endif *d reconr.2758 *if sw phi=dcmplx(phr0,phri) *else phi=cmplx(phr0,phri) *endif *d reconr.2762 exphi=exp(-2*eye*phi) *d reconr.2765,2766 *if sw rpart=dble(slsjnn) *else rpart=real(slsjnn) *endif sabs=abs(slsjnn) sigp(1)=sigp(1)+2*pifac*gj*rpart */ reconr -- 02jul98 */ fix up the use of float. for double precision version, */ float will put out a real*4, which is not what we want */ when the user has asked for "*set sw". *d reconr.2040 spot=4*(2*ll+1)*pifac*sinsq *d reconr.2201 fl=ll *d reconr.2348 fl=ll *d reconr.3123 spot=4*pi*(2*ll+1)*(sin(ps)/k)**2 */ reconr -- 02jul98 */ fix some places where mixed integer and real variables */ could cause problems. *d reconr.1309 l=iscr+6+2*nint(a(iscr+4)) *d reconr.2686 lc(ii)=nint(a(inow+ii+4)) *d reconr.3100 jnow=nint(a(inow+4))+inow+6 *ident up2 */ broadr -- 4jul98 */ patches an error in the calculation of the thermal */ quantities involving fission, and then further */ improves the energy integrals so that the values */ computed in broadr match those that would be */ computed by groupr. we also print out the */ 1/v-equivalent version of k1, that is, the k1 */ integral times 2 over sqrt(pi). *i broadr.118 data pi/3.14159265d0/ *i broadr.120 data pi/3.14159265/ *i broadr.453 picon=2/sqrt(pi) *d broadr.473 2240 step=1.01 enow=tt(1) if (enow.gt.step*elast.and.elast.gt.0.) then enow=step*elast sf=slf+(tt(llf)-slf)*(enow-elast)/(tt(1)-elast) sc=slc+(tt(llc)-slc)*(enow-elast)/(tt(1)-elast) else sf=tt(llf) sc=tt(llc) endif *d broadr.477 ss=slf+(therm-elast)*(sf-slf)/(enow-elast) *d broadr.491 ss=slc+(therm-elast)*(sc-slc)/(enow-elast) *d broadr.499 ss=slf+(tev-elast)*(sf-slf)/(enow-elast) *d broadr.514 ss=slc+(tev-elast)*(sc-slc)/(enow-elast) *i broadr.533 fnow=enow*exp(-enow/tev) *d broadr.535 *d broadr.537 cint=cint+(fnow*sc+flast*slc)*(enow-elast)/2 *d broadr.539 *d broadr.542 *d broadr.544,551 fint=fint+(fnow*sf+flast*slf)*(enow-elast)/2 if (sf.gt.0..and.slf.gt.0.) 1 alint=alint+(fnow*sc/sf+flast*slc/slf) 2 *(enow-elast)/2 etint=etint+(fnow*fnu*sf/(sf+sc) 1 +flast*fnul*slf/(slf+slc))*(enow-elast)/2 v1int=v1int+(fnow*(fnu*sf-sf-sc) 1 +flast*(fnul*slf-slf-slc))*(enow-elast)/2 *d broadr.553 *i broadr.554 flast=fnow *d broadr.559,560 ss=slf+(eone-elast)*(sf-slf)/(enow-elast) ssf=(sf/enow+ss/elast)*(enow-eone)/2 *d broadr.563,564 ss=slc+(eone-elast)*(sc-slc)/(enow-elast) ssc=(sc/enow+ss/elast)*(enow-eone) *d broadr.568,571 if (llf.gt.0) ssf=ssf+(sf/enow+slf/elast)*(enow-elast)/2 if (llc.gt.0) ssc=ssc+(sc/enow+slc/elast)*(enow-elast)/2 *d broadr.576,577 if (llf.gt.0) slf=sf if (llc.gt.0) slc=sc if (enow.lt.tt(1)) go to 2240 *d broadr.586 2 picon*cint/ctev *d broadr.598 2 picon*fint/ftev *i broadr.604 write(nsyso,'( 1 '' equivalant k1:'',1p,e12.4)') 2 picon*v1int *ident up3 */ heatr -- 02jul98 */ fix some variables inadvertantly typed as real. */ they should become double if "*set sw" was used. *d heatr.3263 *if sw real*8 nc,nb *else real nc,nb *endif */ heatr -- 02jul98 */ fix some places where mixed integer and real variables */ could cause problems. *d heatr.1746,1747 nr=nint(c(lnow+4)) np=nint(c(lnow+5)) *d heatr.2749 l=7+2*nint(c(5)) */ heatr -- 02jul98 */ change the name of the 939 mev constant. it should */ an implicit floating type. *d heatr.654 data emc2/939.512d6/ *d heatr.659 data emc2/939.512e6/ *d heatr.732 tm=emc2*(awr+1) *d heatr.1245 data emc2/939.512d6/ *d heatr.1253 data emc2/939.512e6/ *d heatr.1265 tm=emc2*(awr+1) *d heatr.3444 data emc2/939.512d6/ *d heatr.3448 data emc2/939.512e6/ *d heatr.3453 ein=2*emc2*(awr+1) *d heatr.4108 data emc2/939.512d6/ *d heatr.4111 data emc2/939.512e6/ *d heatr.4120 tm=emc2*(awr+1) *d heatr.4546 data emc2/939.512d6/ *d heatr.4550 data emc2/939.512e6/ *d heatr.4555 ein=2*emc2*(awr+1) *d heatr.4605 data emc2/939.512d6/ *d heatr.4607 data emc2/939.512e6/ *d heatr.4611 tm=emc2*(awr+1) *ident up4 */ groupr -- 19jun98 */ fix a bad write format and a bad message. *d groupr.1917 1 '' gamma group structure......read in'')') *d groupr.3733 1 call error('stounr','storage exceeded.',' ') */ groupr -- 02jul98 */ fix problems with complex arithmetic when double */ precisioning is used. *d groupr.6846 *if sw complex*16 carg1,carg2,cs1,cs2 *else complex carg1,carg2,cs1,cs2 *endif *i groupr.6887 data zero,uno/0.d0,1.d0/ *i groupr.6925 data zero,uno/0.,1./ *d groupr.7024 *if sw cs1=dcmplx(c(8+np),c(9+np))/2 *else cs1=cmplx(c(8+np),c(9+np))/2 *endif *d groupr.7026 *if sw cs1=cs1+(2*it+1)*p(it+1)*dcmplx(c(8+np+2*it),c(9+np+2*it))/2 *else cs1=cs1+(2*it+1)*p(it+1)*cmplx(c(8+np+2*it),c(9+np+2*it))/2 *endif *d groupr.7028,7029 *if sw carg1=dcmplx(zero,uno)*eta*log((1-x)/2) sigi=(-2*eta/(1-x))*dble(cs1*exp(carg1)) *else carg1=cmplx(zero,uno)*eta*log((1-x)/2) sigi=(-2*eta/(1-x))*real(cs1*exp(carg1)) *endif *d groupr.7037 *if sw cs1=dcmplx(c(8+nt),c(9+nt))/2 *else cs1=cmplx(c(8+nt),c(9+nt))/2 *endif *d groupr.7041,7042 *if sw cs1=cs1+(2*it+1)*p(it+1)*dcmplx(c(8+nt+2*it),c(9+nt+2*it))/2 cs2=cs2+sgn*(2*it+1)*p(it+1)*dcmplx(c(8+nt+2*it),c(9+nt+2*it))/2 *else cs1=cs1+(2*it+1)*p(it+1)*cmplx(c(8+nt+2*it),c(9+nt+2*it))/2 cs2=cs2+sgn*(2*it+1)*p(it+1)*cmplx(c(8+nt+2*it),c(9+nt+2*it))/2 *endif *d groupr.7045,7048 *if sw carg1=dcmplx(zero,uno)*eta*log((1-x)/2) carg2=dcmplx(zero,uno)*eta*log((1+x)/2) sigi=(-2*eta/(1-x*x))*dble(cs1*(1+x)*exp(carg1) 1 +cs2*(1-x)*exp(carg2)) *else carg1=cmplx(zero,uno)*eta*log((1-x)/2) carg2=cmplx(zero,uno)*eta*log((1+x)/2) sigi=(-2*eta/(1-x*x))*real(cs1*(1+x)*exp(carg1) 1 +cs2*(1-x)*exp(carg2)) *endif */ groupr -- 02jul98 */ fix some variables inadvertantly typed as real. */ they should become double if desired. *d groupr.5413 *if sw real*8 nc,nb *else real nc,nb *endif */ groupr -- 02jul98 */ fix some places where mixed integer and real variables */ could cause problems. *d groupr.4551 npsx=nint(c(l+5)) *d groupr.4645 law=nint(c(l1+3)) *d groupr.6367 nbrag=nint(p(6)) *d groupr.6753,6754 il=iflg+nint(a(iloca+i-1)) im=iflg+nint(a(iloca+nipj-1)) *d groupr.7158,7160 lp=nint(a(l+2)) nr=nint(a(l+4)) np=nint(a(l+5)) *d groupr.7913 nktot=nint(c(ic+4)) */ groupr -- 02jul98 */ change the name of the mev constant. *d groupr.1803 data emev/1.d6/ *d groupr.1833 data emev/1.e6/ *d groupr.1854 egg(ig)=eg2(ig)*emev *d groupr.1861 egg(ig)=eg3(ig)*emev *d groupr.1868 egg(ig)=eg4(ig)*emev *d groupr.1875 egg(ig)=eg5(ig)*emev *d groupr.1882 egg(ig)=eg6(ig)*emev *d groupr.1898 egg(ig)=eg8(ig)*emev *d groupr.1903 egg(ig)=eg8(ig+1)*emev *ident up5 */ acer -- 04jul98 */ fix an uninitialized variable that always used to come */ through as zero on the type1 xsdir "stub" cards. *d acer.7902 lrec=0 */ acer -- 22may98 */ improve the printout, especially for secondary photons. */ there are no differences in calculated values caused */ by this change. requested by frankle, lanl/x-ci. *d acer.7061 1 '' intt ='',i2,'' nd ='',i4,'' np ='',i3// *d acer.7065 5 e2,intt,nd,nn,(xss(j+loci),xss(j+nn+loci), *d acer.7186 1 '' intt ='',i2,'' nd ='',i4,'' np ='',i3// *d acer.7192 7 e2,intt,nd,nn,(xss(j+loci),xss(j+nn+loci), *i acer.7586 mtrn=nint(xss(i-1+mtrp)) call mtname(mti,title(1)) nlaw=1 *d acer.7593 write(nsyso,'(// 1 '' energy distribution for secondary photons from '', 2 ''reaction '',a10,'' (mt='',i6,'') with '',i1,'' law(s)'')') 3 title(1),mtrn,nlaw write(nsyso,'(/ 1 '' law ='',i2,i5,''st of'',i2,'' laws''/)') 2 law,nlaw,nlaw *i acer.7631 if (iprint.eq.1.and.m.ne.0) write(nsyso,'(12x,''ne ='',i4)') ne *d acer.7647 1 '' intt ='',i2,'' nd ='',i3,'' np ='',i3// *d acer.7651 5 eg,intt,nd,n,(xss(j+loci),xss(j+n+loci),xss(j+2*n+loci),j=1,n) */ acer -- 22may98 */ fix reaction names on plots for inelastic cross sections. */ we're having trouble with quotes insides quotes, so we will */ leave the classical ace stars, e.g., (n,n*2). */ fix one bad y-axis label where the quotes were omitted. */ fix the line type for the resonance fission plots. *d acer.10844,10846 *d acer.10881,10883 *d acer.11019,11021 *d acer.11563,11565 *d acer.11882,11884 *d acer.12008,12010 *d acer.12103,12105 *d acer.12210,12212 *d acer.12266,12268 *d acer.12359,12361 *d acer.12728,12730 *d acer.11626 write(nout,'(a,''ross section (barns)'',a,''/'')') qu,qu *d acer.11366 *i acer.11367 write(nout,'(''0 0 1/'')') */ acer -- 19jun98 */ be sure to compute enough legendre polynomials for */ the required limits. *d acer.2452 call legndr(u,pl,nmax+1) */ acer -- 19jun98 */ fix a bad message format. *d acer.3282 1 '' ---message from summer---for distribution as per '', */ acer -- 02jul98 */ fix some variables inadvertantly typed as real. they */ should become double if desired. *d acer.8902 *if sw real*8 nc,nb *else real nc,nb *endif */ acer -- 02jul98 */ change the name of the mev constant. *d acer.4577 data emev/1.d6/ *d acer.4581 data emev/1.e6/ *d acer.4680 a(i+j)=sigfig(a(2*i+2*m+4+iscr)/emev,7,0) *d acer.4805 xss(lqr+ir)=sigfig(a(iscr+1)/emev,7,0) *d acer.4828 230 xss(ih+j)=sigfig(s/emev/xss(it+j),7,0) *d acer.4833 if (mth.eq.444) s=s/emev *d acer.4950 xss(ie+j)=sigfig(a(iscr+1)/emev,7,0) *d acer.5010,5011 xss(next+5)=sigfig(xss(esz+ja)/emev,7,0) xss(next+6)=sigfig(xss(esz+jb)/emev,7,0) *d acer.5053 xss(next-1+ie)=sigfig(a(jj)/emev,7,0) *d acer.5106 a(iegn-1+i)=sigfig(a(iscr+1)/emev,7,0) *d acer.5108 xss(next)=sigfig(a(iscr+5+j)/emev,7,0) *d acer.5126 xss(esz+i)=sigfig(xss(esz+i)/emev,9,0) *d acer.5157 data emev/1.d6/ *d acer.5159 data emev/1.e6/ *d acer.5201 xss(j+next)=sigfig(a(2*j+l)/emev,7,0) *d acer.5252 xss(next+j)=sigfig(e/emev,7,0) *d acer.5271,5272 xss(ki+nexd)=sigfig(a(iscr+6+2*ki)/emev,7,0) xss(ki+n+nexd)=sigfig(a(iscr+7+2*ki)*emev,7,0) *d acer.5310,5311 xss(j+next)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) xss(j+n+next)=sigfig(a(2*j+2*m+5+iscr)/emev,7,0) *d acer.5321 xss(j+next)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) *d acer.5346,5347 xss(j+next)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) xss(j+n+next)=sigfig(a(2*j+2*m+5+iscr)/emev,7,0) *d acer.5350 xss(next)=u/emeV *d acer.5364,5365 xss(j+next)=sigfig(a(2*j+4+iscr)/emev,7,0) xss(j+m+next)=sigfig(a(2*j+5+iscr)/emev,7,0) *d acer.5373,5374 xss(j+next)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) xss(j+n+next)=sigfig(a(2*j+2*m+5+iscr)/emev,7,0) *d acer.5386,5387 xss(j+next)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) xss(j+n+next)=sigfig(a(2*j+2*m+5+iscr)/emev,7,0) *d acer.5395,5396 xss(j+next)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) xss(j+n+next)=sigfig(a(2*j+2*m+5+iscr)*emev,7,0) *d acer.5399 xss(next)=sigfig(u/emev,7,0) *d acer.5428 xs(2)=emev *d acer.5459 xss(next+ie)=sigfig(tme/emev,7,0) *d acer.5466,5467 xss(ki+nexd)=sigfig(xxs(ki)/emev,7,0) xss(ki+js+nexd)=sigfig(yys(ki)*renorm*emev,7,0) *d acer.5505 data emev/1.d6/ *d acer.5507 data emev/1.e6/ *d acer.5631 xss(next+1+j)=sigfig(eyl/emev,7,0) *d acer.5683 xss(next+k)=sigfig(a(iscr+4+2*m+2*k)/emev,7,0) *d acer.5694 eyl=xss(igyl+1+j)*emev *d acer.5700 1 xss(lgyl+1),xss(lgyl+1+ngyl),eyl/emev,gyl,2) *d acer.5703 1 xss(lgyl),xss(lgyl+ngyl),eyl/emev,gyl,2) *d acer.5709 xss(next+j)=sigfig(eyl/emev,7,0) *d acer.5763,5764 450 xss(next+j)=sigfig(c2h/emev,7,0) ee=c2h/emev *d acer.5780 e=ee*emev *d acer.5785 3 6x,''patching...'')') mt,e/emev,ep/emev *d acer.5792 3 6x,''leaving it as is...'')') mt,e/emev,ep/emev *d acer.5794 xss(ki+nexd)=sigfig(a(iscr+6+ncyc*(ki-1))/emev,7,0) *d acer.5796,5797 if (ki.gt.nd) xss(ki+n+nexd) 1 =sigfig(a(iscr+7+ncyc*(ki-1))*emev,7,0) bzro=xss(nexd+ki+n)/emev *d acer.5898 ee=c2h/emev *d acer.5922,5923 xss(nexd+ki)=sigfig(a(iscr+8+2*(ki-1))/emev,7,0) xss(nexd+npep+ki)=sigfig(a(iscr+9+2*(ki-1))*emev,7,0) *d acer.6023 data emev/1.d6/ *d acer.6025 data emev/1.e6/ *d acer.6098 xss(nex+i)=sigfig(a(jscr+2*i)/emev,7,0) *d acer.6214 xss(nex+ie)=sigfig(c2h/emev,7,0) *d acer.6250,6251 xss(nex+5)=sigfig(ef/emev,7,0) xss(nex+6)=sigfig(el/emev,7,0) *d acer.6255 xss(nex+10)=sigfig(eg/emev,7,0) *d acer.6455 xss(nex+10)=sigfig(eg/emev,7,0) *d acer.6292 xss(j+nex)=sigfig(a(2*j+2*m+4+iscr)/emev,7,0) *d acer.6324 550 xss(nex-1+ie)=sigfig(c2h/emev,7,0) *d acer.6334,6335 xss(k+nexd)=sigfig(a(iscr+6+2*k)/emev,7,0) xss(k+n+nexd)=a(iscr+7+2*k)*emev *d acer.6375,6376 xss(nex+5)=sigfig(a(iscr+6+2*m)/emev,7,0) xss(nex+6)=sigfig(a(iscr+4+2*m+2*n)/emev,7,0) *d acer.6408 670 xss(nex-1+ie)=sigfig(c2h/emev,7,0) *d acer.6416 xss(ki+nexd)=sigfig(a(iscr+6+ncyc*(ki-1))/emev,7,0) *d acer.6418 if (ki.gt.nd) xss(ki+n+nexd)=a(iscr+7+ncyc*(ki-1))*emev *d acer.9353 data emev/1.d6/ *d acer.9356 data emev/1.e6/ *d acer.9411,9412 xss(i+itce)=a(2*(i-1)+iscr)/emev xss(i+nee+itce)=a(1+2*(i-1)+iscr)/emev/nmix *d acer.9425 xss(itce+i)=a(indx)/emev *d acer.9443 xss(itie+i)=a(iscr)/emev *d acer.9448 xss(indx+k)=a(iscr+1+k)/emev *d acer.9759 data emev/1.d6/ *d acer.9761 data emev/1.e6/ *d acer.9841 xss(l)=a(k)/emev *d acer.10095 data emev/1.d6/ *d acer.10105 data emev/1.e6/ *d acer.10175 xss(eszg+i)=xss(eszg+i)/emev */ acer -- 02jul98 */ fix some places where mixed integer and real variables */ could cause problems. *d acer.1204 ln=nint(b(4)) *d acer.1219 int=nint(a(iscr+5+2*i)) *d acer.2134 iso=nint(a(iscr+2)) *d acer.13187,13188 intmu=nint(a7(loci)) nmu=nint(a7(loci+1)) *d acer.13202,13203 210 intmu=nint(a7(loci)) nmu=nint(a7(loci+1)) */ acer -- 02jul98 */ fix up the use of float. for double precision version, */ float will put out a real*4, which is not what we want. */ these floats were left over from the old days. *d acer.2597 if (n.gt.1) afac=afac*(2*n-1)/n *d acer.2604,2605 aa=2*n+1 aa=aa/2 a(ixat+n-1)=cf*afac*aa/(n+jfi) *d acer.2610,2613 b=(n-ist)*(n-ifi) c=ivar*(2*n-ifi) d=n+jfi cf=cf*asign*b/c *d acer.2754 delmu1=(amumid+1)/mfmid *d acer.2760 delmu2=(1-amumid)/mfmid *d acer.2783 dprob=sum/mcoars *d acer.2794 ap=a(ipro+k-1)-dprob*k *d acer.2803 330 b=dprob*k *d acer.3010 dprob=sum/mcoars *ident up6 */ purr -- 26may98 */ add the capability of producing conditional probabilities */ for heating in the probability tables. these are best */ if partial heating cross sections are available for */ elastic (MT=302), fission (318), and capture (402). The */ coding will leave the heating unselfshielded if only the */ total heating (mt=301) is present. *i purr.20 c * a conditional probability for heating is added to the table. * c * if partial heating cross sections for elastic (302), fission * c * (318), and capture (402) are available from heatr, full * c * fluctuations will be provided for the total heating. * c * otherwise, the same value will be provided for each bin. * *i purr.59 external rdheat *i purr.184 c c ***read in the total and partial heating cross sections nw=ntemp*nunr*4 call reserv('heat',nw,iheat,a) call rdheat(a,a(iheat),a(ieunr),temp,ntemp,nunr,ihave,matd) if (ihave.eq.0) call mess('purr', 1 'no heating found on pendf', 2 'ur heating set to zero') if (ihave.eq.1) call mess('purr', 1 'no partial heating xsecs found on pendf', 2 'ur heating will not selfshield') *d purr.259 nc153=(1+6*nbin)*nunr *d purr.366 if (12+(1+6*nbin)*nunr.gt.maxscr) call error('purr', *d purr.381 a(n+10)=(1+6*nbin)*nunr *i purr.402 do 306 j=1,nbin n=n+1 a(n)=0 306 continue *i purr.415 k=iheat+4*nunr*(it-1)+4*(ie-1)-1 do 319 j=1,nbin l=n1+j+5*nbin if (ihave.eq.1) then if (lssf.eq.1) then a(l)=1 else a(l)=a(k+1)/sigu(1,1,1) endif else if (ihave.eq.2) then a(l)=a(l)+(a(k+1)-a(k+2)-a(k+3)-a(k+4)) h=a(k+2) if (lssf.eq.1) then h=h*a(n1+j+2*nbin) else if (sigu(2,1,1).ne.0) then h=h*a(n1+j+2*nbin)/sigu(2,1,1) endif a(l)=a(l)+h h=a(k+3) if (lssf.eq.1) then h=h*a(n1+j+3*nbin) else if (sigu(3,1,1).ne.0) then h=h*a(n1+j+3*nbin)/sigu(3,1,1) endif a(l)=a(l)+h h=a(k+4) if (lssf.eq.1) then h=h*a(n1+j+4*nbin) else if (sigu(4,1,1).ne.0) then h=h*a(n1+j+4*nbin)/sigu(4,1,1) endif a(l)=a(l)+h if (lssf.eq.1) then a(l)=a(l)/a(k+1) else a(l)=a(l)/sigu(1,1,1) endif endif 319 continue *d purr.419 nwds=6*nbin*nunr *i purr.897 c subroutine rdheat(a,heat,eunr,temp,ntemp,nunr,ihave,matd) c ****************************************************************** c read the total heating (mt=301) and partial heating cross c sections for elastic (302), fission (318), and capture (402) c from the pendf tape on the unresolved energy grid. c ****************************************************************** *if sw implicit real*8 (a-h,o-z) *endif common/mainio/nsysi,nsyso,nsyse,ntty common/units/nendf,nin,nout common/cont/c1h,c2h,l1h,l2h,n1h,n2h,math,mfh,mth,nsh,nsp,nsc dimension a(1) dimension heat(4,nunr,ntemp) dimension eunr(nunr),temp(ntemp) c c ***allocate storage call findex('scr',iscr,a) c c ***loop over the temperatures call repoz(nin) call tpidio(nin,0,0,a(iscr),nb,nw) itemp=0 110 itemp=itemp+1 call findf(matd,3,1,nin) ihave=0 120 call contio(nin,0,0,a(iscr),nb,nw) if (mth.eq.0) go to 150 ix=0 if (mth.eq.301) ix=1 if (mth.eq.302) ix=2 if (mth.eq.318) ix=3 if (mth.eq.402) ix=4 if (ix.gt.0) go to 130 call tosend(nin,0,0,a(iscr)) go to 120 130 e=0. call gety1(e,enext,idis,h,nin,a(iscr)) if (ix.eq.1) ihave=1 if (ix.gt.1) ihave=2 do 140 ie=1,nunr e=eunr(ie) call gety1(e,enext,idis,h,nin,a(iscr)) heat(ix,ie,itemp)=h 140 continue call tosend(nin,0,0,a(iscr)) go to 120 150 continue call tomend(nin,0,0,a(iscr)) if (itemp.lt.ntemp) go to 110 call repoz(nin) return end *ident up7 */ acer -- 26may98 */ adjust acer to accept the new probability tables that */ contain a conditional probability for heating. *d acer.5039 nurb=(nurb/nure-1)/6 *d acer.5052 jj=iurd+6+(ie-1)*(1+6*nurb) *d acer.5062 xss(ll+5*nurb+ib)=a(jj+5*nurb+ib)/emev *ident up8 */ acer -- 20jun98 */ fix an error introduced by up5. this problem damages the */ list of photon production yield mt values. *d up5.196 nex=nex+1 */ acer -- 6aug98 */ fix some problems with the ace files for photon angular */ distributions and with the corresponding output listing. */ noted by bob little (lanl). *d acer.6216 xss(nex+ne+ie)=lc-andp+2 *d acer.7552 if (k.gt.0) k=k+andp-1 *d acer.7559,7560 write(nsyso,'(6x,8(4x,a6,a4))') (ek,blank,ii=1,ic) write(nsyso,'(5x,1p,8e14.5)') (xss(ii+na),ii=iaa,ib) *d acer.7568 write(kk(m),'(1p,e14.5)') xss(ii) *d acer.7572 write(nsyso,'(1x,i4,8a14)') j,(kk(ii),ii=1,nkk) *ident up9 */ gaspr -- 14aug98 */ fix an error in the construction of the new file1 */ dictionary on the pendf tape after the insertion of the */ gas production reactions. this problem doesn't affect */ the numbers in file3. *d gaspr.628 nm=6*nx-j+k+5 do 400 i=1,nm *ident up10 */ reconr -- 15aug98 */ modify the reconstruction process in both lunion and resxs */ to use a tighter tolerance than err in the thermal range. */ this means that the .0253 ev cross section will be computed */ better, and the thermal parameters generated in broadr can */ be computed more accurately. requested by lubitz (kapl). *i reconr.1211 trange=.5d0 *i reconr.1218 trange=.5 *d reconr.1510 errn=err if (a(ix+i-2).lt.trange) errn=err/5 test=errn*ym+ssmall *i reconr.1608 data trange/.5d0/ *i reconr.1612 data trange/.5/ *d reconr.1736,1741 errn=err if (a(ix+i-2).lt.trange) errn=err/5 errm=errmax if (a(ix+i-2).lt.trange) errm=errmax/5 if (dm1.le.errn*sig(2).and. 1 dm2.le.errn*sig(3).and. 2 dm3.le.errn*sig(4)) go to 145 if (dm1.gt.errm*sig(2).or. 1 dm2.gt.errm*sig(3).or. 2 dm3.gt.errm*sig(4)) go to 175 */ reconr -- 20aug98 */ fix problem with reading in user's energy grid. this is */ a rarely used option. *d reconr.323 read(nsysi,*) (a(ienode+i-1),i=1,ngrid) *ident up11 */ broadr -- 15aug98 */ modify the reconstruction tests to use a tighter */ tolerance in the thermal range to be consistent with */ the changes to reconr in up10. change the upper limit */ for the thermal maxwellian integration to a more */ reasonable value. change the display of the thermal */ parameters so that it only appears if the temperature */ is close to 293.6 K (0.0253 ev). *d broadr.28,31 c * if the temperature is close to 293.6 K (.0253 ev), broadr * c * computes and displays thermal cross sections, maxwellian * c * integrals (one-group thermal cross sections), g-factors, * c * integral ratios (eta, alpha), the k1 integral and the * c * corresponding 1/v-equivalent, and resonance integrals. * *d broadr.177 8 '' starting material temperature ........ '',f10.1,''k''/ *d broadr.179,181 a '' max. energy .......................... '',1p,e10.3/ b '' errmax for thinning .................. '',1p,e10.3/ c '' errint for thinning .................. '',1p,e10.3)') *d broadr.185,186 1 '' final temperatures ................... '',1p,e10.3/ 2 (40x,e10.3))') (temp2(i),i=1,ntemp2) *d broadr.388 write(nsyse,'(1p,e11.3,'' deg'',54x,0p,f8.1,''s'')') 1 tempk,time *d broadr.446 c ***write out thermal quantities if t is close to 293.6 kelvin. *d broadr.449 if (abs(tev-therm).gt.therm/100) go to 2250 *d broadr.479 1 '' thermal fission xsec:'',1p,e12.4)') ss *d broadr.488 1 '' thermal fission nubar:'',1p,e12.4)') fnu *d broadr.493 1 '' thermal capture xsec:'',1p,e12.4)') ss *d broadr.501,502 *d broadr.510,511 *d broadr.516,517 *d broadr.532 break=1. *d up2.74 1 '' equivalent k1:'',1p,e12.4)') *i broadr.1002 data trange/.5d0/ *i broadr.1009 data trange/.5/ *i broadr.1087 errt=errthn if (es(is-1).lt.trange) errt=errt/5 errm=errmax if (es(is-1).lt.trange) errm=errm/5 *d broadr.1099,1100 if (dy.le.errt*abs(sn(i))+errmin) go to 145 if (dy.gt.errm*abs(sn(i))+errmin) go to 170 */ broadr -- 18aug98 */ modify the doppler-broadening method to extrapolate */ below the lowest energy point assuming a 1/v cross section */ rather than a constant. this is consistent with the newer */ versions of the sigma1 method from red cullen. otherwise, */ the 1/v cross sections at the lowest energies have a */ tendency to sag slightly from the 1/v shape. *d broadr.1245 c continue cross sections as 1/v to x=0. *d broadr.1248 aa=y *d broadr.1251 sbt(i)=sbt(i)-s(i,klow)*e(klow)*(oy**2*h(2)+oy*h(1)) *d broadr.1290 c extend cross section to x=0.0 as 1/v. *d broadr.1298 sbt(i)=sbt(i)-s(i,klow)*e(klow)*(oy**2*h(2)+oy*h(1)) *i broadr.1452 if (a.eq.alast) go to 120 *ident up12 */ unresr -- 18aug98 */ fix error made in rdunf3 when eliminating an arithmetic if. *d unresr.745,746 if (mtx(ix).eq.mth) go to 105 if (mtx(ix).gt.mth) go to 100 */ unresr -- 18aug98 */ fix the temperature format to allow for 293.6 degrees. *d unresr.107,108 1 '' temperatures ......................... '',1p,e10.3/ 2 (40x,e10.3))') (temp(i),i=1,ntemp) *d unresr.110,111 1 '' sigma zero values .................... '',1p,e10.3/ 2 (40x,e10.3))') (sigz(i),i=1,nsigz) *d unresr.120 1 '' mat = '',i4,3x,'' temp = '',1p,e10.3,37x,0p,f8.1,''s'')') */ unresr -- 24aug98 */ fix the wide parameter to match reconr. *d unresr.370 data wide/3.01d0/ *d unresr.376 data wide/3.01/ *ident up13 */ purr -- 18aug98 */ increase the storage space to make room for endf/b-vi pu-239. */ patch purr for elemental evaluations, and allow space for */ more spin sequences to handle the case of endf/b-vi mo-nat. *d purr.51 common/pustore/a(70000) *d purr.65 data namax/70000/, nidmax/24/, ipr/1/ *d purr.910,912 common/sigcon/e,t,cth(40),csz(40),cc2p(40),cs2p(40), 1 cgn(40),cgg(40),cgf(40),cgx(40),cgt(40),dbar(40), 2 spot,dbarin,sigi(4),ndfn(40),ndff(40),ndfx(40),nseq0 *d purr.1050 if (nseq0.gt.40) call error('unresx','too many sequences',' ') csz(nseq0)=abn*ab*gj *d purr.1070 temp=abn*pi*ab*gj*gnx/(2*dx) *d purr.1278,1280 common/sigcon/e,t,cth(40),csz(40),cc2p(40),cs2p(40), 1 cgn(40),cgg(40),cgf(40),cgx(40),cgt(40),dbar(40), 2 spot,dbarin,sigi(4),ndfn(40),ndff(40),ndfx(40),nseqz *d purr.1381,1383 common/sigcon/e,t,cth(40),csz(40),cc2p(40),cs2p(40), 1 cgn(40),cgg(40),cgf(40),cgx(40),cgt(40),dbar(40), 2 spot,dbarin,sigi(4),ndfn(40),ndff(40),ndfx(40),nseqz *d purr.1487,1489 common/sigcon/e,t,cth(40),csz(40),cc2p(40),cs2p(40), 1 cgn(40),cgg(40),cgf(40),cgx(40),cgt(40),dbar(40), 2 spot,dbarin,sigi(4),ndfn(40),ndff(40),ndfx(40),nseq0 */ purr -- 18aug98 */ change the temperature format to allow for 293.6. *d purr.114 1 '' temperatures ......................... '',1p,e10.3)') *d purr.116 if (ntemp.gt.1) write(nsyso,'(40x,1p,e10.3)') (temp(i),i=2,ntemp) *d purr.120 if (nsigz.gt.1) write(nsyso,'(40x,1p,e10.3)') (sigz(i),i=2,nsigz) *d purr.2020 write(nsyso,'(3x,1p,2e10.3,5e12.4)') *d purr.2071 write(nsyso,'(1x,a,1p,e10.3,10e11.3/(14x,10e11.3))') *d purr.2110 write(nsyso,'(3x,1p,2e10.3,5e12.4)') */ purr -- 24aug98 */ fix the wide parameter to match reconr. *d purr.483 data wide/3.01d0/ *d purr.489 data wide/3.01/ *ident up14 */ acer -- 24aug98 */ don't divide the heating by 1e6 ev when it is really */ a factor. this results in unresolved heating values too */ small by a factor of 1e6 for two endf/b-vi materials. */ found by little (los alamos). *d up7.10 if (lssf.eq.0) then xss(ll+5*nurb+ib)=a(jj+5*nurb+ib)/emev else xss(ll+5*nurb+ib)=a(jj+5*nurb+ib) endif */ acer -- 24aug98 */ fix some problems in the length of strings in common *d acer.11938 character hz*13,hd*10,hk*70,hm*10 *d acer.12141 character hz*13,hd*10,hk*70,hm*10 *d acer.12244 character hz*13,hd*10,hk*70,hm*10 *d acer.12334 character hz*13,hd*10,hk*70,hm*10 *d acer.12706 character hz*13,hd*10,hk*70,hm*10 */ acer -- 22sep98 */ fix the initialization of e for gety2 for subsections after */ the first. this was a byproduct of taking constants out of */ subroutine calls. it causes problems for materials with */ anisotropic photons, such as n and o from endf/b-vi. *i acer.3595 e=0 */ acer -- 22sep98 */ fix the logic for identifying which photon an angular */ distribution corresponds to. the old scheme only used */ energy, and if two different reactions produced photons */ with the same energy, the angular distribution would */ always be indexed to the first photon in the list with */ that energy. *d acer.6197 eps=1.e-6 mmm=nint(xss(mtrp+j-1)) mmm=mmm/10000 if (mmm.ne.mtd) go to 340 if (abs(eg-a(iphot+5*(j-1))).gt.eps*eg) go to 340 *d acer.6207 eps=1.e-6 mmm=nint(xss(mtrp+j-1)) mmm=mmm/10000 if (mmm.ne.mtd) go to 370 if (abs(eg-a(iphot+5*(j-1))).gt.eps*eg) go to 370 */ acer -- 22sep98 */ allow for more angles in the contour plots to handle the case */ of endf/b-vi oxygen. *d acer.12247 dimension loce(1200) *d acer.12271 if (ne.gt.1200) call error('aplof4', *ident up15 */ covr -- 23sep98 */ fix some problems with the covariance plots that arose when */ we changed to using quotes to delimit character strings. */ the reaction names containing "prime" do not print out */ right. we remove the prime, counting on the level subscript */ do tell the user that it is an inelastic reaction. also, */ change the blank labels to ".". this is a flag to viewr */ to use a blank label, but to go ahead with the numbering */ of the axis. these problems are visible in the file plot05 */ from the original njoy97 test set. noted by muir (iaea). *d covr.966 write(nplot,'(a,''.'',a,''/'')') qu,qu *d covr.968 write(nplot,'(a,''.'',a,''/'')') qu,qu *d covr.1022 write(nplot,'(a,''.'',a,''/'')') qu,qu *d covr.1024 write(nplot,'(a,''.'',a,''/'')') qu,qu *d covr.1557 data nmeb1/'(n,n'/, nmee1/'(delayed'/, nmef1/'(prompt'/ *d covr.1573 inamel=4 *ident up16 */ viewr -- 23sep98 */ when using quotes as character delimiters, it is hard to */ tell the difference between a blank label and a "default" */ label. we change the input so that a label of '.' results */ in a blank label on the plot with numbering on the axis */ ticks. a label of ' ' results in a "default" label, that */ is, no text label, and no axis numbering. these features */ are used for the covariance plots from covr (see up15). *i viewr.2438 if (label.eq.'.') return *ident up17 */ thermr -- 23sep98 */ modify output format to allow for 4-digit temperatures. *d thermr.2564 1 '' wrote thermal data for temp='',1p,e10.3,30x,0p,f8.1,''s'')') *d thermr.2567 1 '' wrote thermal data for temp='',1p,e10.3,30x,0p,f8.1,''s'')') *ident up18 */ acer -- 24sep98 */ move some of the legend blocks so that they don't interfere */ with the curves as much. this only works for curves with */ very predictable locations. *d acer.11134 xtag=5*xmin ytag=0.6*log10(ymin)+0.4*log10(ymax) ytag=10.**ytag write(nout,'(''4 0 2 1'',2e12.4,''/'')') xtag,ytag *d acer.11588,11610 270 xmin=0. xmax=xss(esz-1+nes) ymin=0. ymax=xss(esz+nes-1+nes) if (gpd.gt.0) then if (xss(gpd-1+nes).gt.ymax) ymax=xss(gpd-1+nes) endif ymax=2*ymax *d acer.11622 xtag=0.45*xmin+0.55*xmax ytag=0.1*ymin+0.9*ymax write(nout,'(''1 0 2 1'',2e12.4,''/'')') xtag,ytag *d acer.11987 xtag=0.3*xmin+0.7*xmax ytag=0.1*ymin+0.9*ymax write(nout,'(''1 0 2 1'',2e12.4,''/'')') xtag,ytag */ acer -- 24sep98 */ fix the expanded resonance plots to be more clever about */ the energy ranges plotted. check the number of points in */ each decade. if is is more than some cutoff, make a plot */ showing that decade. *d acer.11216,11230 maxii=120 e1=xss(esz) ii1=1 169 ii2=ii1 e2=10*e1 1169 ii2=ii2+1 if (ii2.ge.nes) go to 2169 if (xss(esz-1+ii2).ge.e2) go to 2169 go to 1169 2169 if (ii2-ii1.gt.maxii) go to 3169 if (ii2.ge.nes) go to 1174 e1=e2 ii1=ii2 go to 169 3169 e2=xss(esz-1+ii2) *i acer.11279 ii1=ii2+1 *i acer.11281 1174 continue *d acer.11284,11297 maxii=100 e1=xss(esz) ii1=1 168 ii2=ii1 e2=10*e1 1168 ii2=ii2+1 if (ii2.ge.nes) go to 2168 if (xss(esz-1+ii2).ge.e2) go to 2168 go to 1168 2168 if (ii2-ii1.gt.maxii) go to 3168 if (ii2.ge.nes) go to 183 e1=e2 ii1=ii2 go to 168 3168 e2=xss(esz-1+ii2) *i acer.11385 ii1=ii2+1 */ acer -- 24sep98 */ separate out the gas-production reactions onto their own page *i acer.12060 if (mt.lt.203.and.nint(xss(mtr+i)).ge.203) nlev=5 *i acer.12120 if (mt.lt.203.and.nint(xss(mtr+i)).ge.203) nlev=5 *d acer.12125 go to 210 */ acer -- 24sep98 */ fix another variable that is not defined except for type2 */ format. see up5 for the other one. *i up5.6 nern=0 */ acer -- 24sep98 */ fix a numerical problem that shows up for endf/b-vi h-2. *d acer.5976,5977 if (xx.gt.test) then xx=1 pn=0 else pn=sqrt(xx)*(1-xx)**(3*npsx/2-4) endif */ acer -- 24sep98 */ fix the temperature print format to allow for 293.6k. *d acer.252 2 '' temperature .......................... '',1p,e10.3/ *d acer.309 2 '' temperature .......................... '',1p,e10.3/ *ident up19 */ thermr -- 24sep98 */ fix a numerical problem, and allow for larger data sets */ (e.g., h in zrh from endf/b-vi rel. 3). *d thermr.94 common/tstore/a(80000) *d thermr.127 namax=80000 *d thermr.130 nwscr=40000 *d thermr.2057 data shade/.99999999d0/ *d thermr.2063 data shade/.99999999/ *ident up20 */ njoy -- 25sep98 */ fix some typing errors in njoy *d njoy.497 if (hdate(1:1).eq.' ') hdate(1:1)='0' *d njoy.510 if (htime(1:1).eq.' ') htime(1:1)='0' *ident up21 */ errorr -- 25sep98 */ permit user energies above 20 mev in errorr when iread=0 */ and correct a related minor code bug *d errorr.2986,2996 etop=a(iun+nunion) etop=sigfig(etop,ndig,0) if (a(iela+ng).ge.etop) go to 155 a(iela+ng)=etop write(strng,'(''since iread='',i2)') iread call mess('colaps', 1 'raising top of ngout structure to top of union structure',' ') 155 continue if (a(iela).gt.a(iun)) call error('colaps', 1 'ngout group structure does not span union grid.',' ') *d errorr.3147,3152 *i errorr.3164 if (ign.ne.1.and.ign.ne.19.and.egn(ngnp1).lt.a(ie+neni-1)) then call mess('uniong', 1 'raising top built-in multigroup boundary to top of file 33', 2 ' ') egn(ngnp1)=a(ie+neni-1) endif *d errorr.3181,3182 c c ***finished with endf energies in=i go to 170 c c ***treat endf energies below first group boundary *i errorr.3200 do 125 i=2,k if (a(iun+i-1).le.a(iun+i-2)) then write(strng,'(1p,e12.4,'' le '',1p,e12.4)') 1 a(iun+i-1),a(iun+i-2) call error('uniong','union energies out of order',strng) endif 125 continue *d errorr.3709,3717 call uniong(nendf,a) *ident up22 */ covr -- 25sep98 */ fix two additional places where a space needs to be a dot */ in order to get axis values labeled. See up15 for others. *d covr.1258 write(nplot,'(a,''.'',a,''/'')') qu,qu *d covr.1328 write(nplot,'(a,''.'',a,''/'')') qu,qu */ covr -- 25sep98 */ restore a deleted line and fix a typing error in covr *i covr.1750 ndig=nvf-6 *d covr.1838 test=1.e-10*abs(xctest) *ident up23 */ moder -- 21jun98 */ add new ltt=3 format. this format modification was introduced */ to help with adding high-energy data to endf evaluations. */ it allows for the use of legendre polynomials for elastic */ scattering at lower incident energies, with a change to */ tabulated distributions at higher energies. *d moder.595 if (lvt.eq.0.and.l1h.eq.1) go to 170 *d moder.618 go to 170 *i moder.632 170 if (ltt.lt.3) go to 200 call tab2io(nin,nout,nscr,a,nb,nw) ne=n2h do 190 ie=1,ne call tab1io(nin,nout,nscr,a,nb,nw) if (nb.eq.0) go to 190 180 call moreio(nin,nout,nscr,a,nb,nw) if (nb.gt.0) go to 180 190 continue *ident up24 */ heatr -- 21jun98 */ add new ltt=3 format to heatr. *d heatr.3517 save elo,ehi,nlo,nhi,flo,fhi,ltt,ltt3,lttn *i heatr.3532 ltt3=ltt if (ltt.eq.3) then ltt=1 lttn=1 endif *d heatr.3592,3593 if (nne.eq.ne.and.ltt3.eq.3.and.lttn.eq.1) then call tab2io(nin,0,0,c(1),nb,nwc) ne=nint(c(6)) nne=0 ir=1 ltt=2 lttn=2 else if (nne.eq.ne) then call error('hgtfle','desired energy above highest given.',' ') endif */ heatr -- 21jun98 */ make more room for lo=2 gammas and provide a protective */ error message. *i heatr.3827 lmax=100 *d heatr.3866,3868 call reserv('eg',lmax,ieg,a) call reserv('es',lmax,ies,a) call reserv('y',lmax,iy,a) *i heatr.3930 if (l.gt.lmax) call error('hconvr', 1 'too many lo=2 gammas',' ') *ident up25 */ groupr -- 21jun98 */ add new ltt=3 format to groupr. *i groupr.6142 save ltt3,lttn *i groupr.6165 ltt3=ltt if (ltt.eq.3) then ltt=1 lttn=1 endif *d groupr.6236,6237 if (nne.eq.ne.and.ltt3.eq.3.and.lttn.eq.1) then call tab2io(nin,0,0,c(1),nb,nwc) ne=nint(c(6)) nne=0 ir=1 ltt=2 lttn=2 else if (nne.eq.ne) then call error('getfle','desired energy above highest given.',' ') endif */ groupr -- 21jun98 */ make more room for lo=2 gammas and provide a protective */ error message. *i groupr.7243 lmax=100 *d groupr.7420,7422 call reserv('eg',lmax,ieg,a) call reserv('es',lmax,ies,a) call reserv('y',lmax,iy,a) *i groupr.7483 if (l.gt.lmax) call error('hconvr', 1 'too many lo=2 gammas',' ') *ident up26 */ njoy -- 25sep98 */ fix the typing errors in njoy again. *d up20.5 if (hdate(i:i).eq.' ') hdate(i:i)='0' *d up20.7 if (htime(i:i).eq.' ') htime(i:i)='0' *ident up27 */ acer -- 11oct98 */ add the word ptable to the xsdir cards when unresolved */ probability tables are present *d acer.7873 1 gpd,mtrp,lsigp,sigp,landp,andp,ldlwp,dlwp,yp,fis,end, 2 iurpt,jxsd(9) *d acer.8053,8055 if (iurpt.gt.0) then write(ndir, 1 '(a10,f12.6,'' filename route'',i2,i4,i7,2i6,1p,e10.3, 2 '' ptable'')') hz,aw0,itype,irec1,len2,lrec,nern,tz else write(ndir, 1 '(a10,f12.6,'' filename route'',i2,i4,i7,2i6,1p,e10.3)') 2 hz,aw0,itype,irec1,len2,lrec,nern,tz endif */ add more storage space to take care of the large probability */ tables in endf/b-vi pu-239. *d acer.184 common/astore/a(35000) *d acer.193 data namax/35000/, nidmax/27/ *d acer.402 common/astore/a(35000) *d acer.4570 common/astore/a(35000) *d acer.5151 common/astore/a(35000) *d acer.5499 common/astore/a(35000) *d acer.6019 common/astore/a(35000) *d acer.6625 common/astore/a(35000) *d acer.7346 common/astore/a(35000) *d acer.8999 common/astore/a(35000) *d acer.9347 common/astore/a(35000) *d acer.9747 common/astore/a(35000) *d acer.10072 common/astore/a(35000) *ident up28 */ purr -- 15oct98 */ fix some tests that are supposed to protect the sqrt */ function against small negative arguments. *d purr.1977 if (argt.lt.0.) argt=0 *d purr.1980 if (arge.lt.0.) arge=0 *d purr.1983 if (argf.lt.0.) argf=0 *d purr.1986 if (argc.lt.0.) argc=0 *ident up29 */ covr -- 5nov98 */ correct typos in 2 reaction names and add 7 new names. */ contributed by muir (iaea). *d covr.1542,1544 character*8 lnamel(2),hira(32),blank,nmea1,nmeb1,lnamer, 1 nmed1,nmee1,nmef1,nmeh1 dimension ira1(32),ira2(32) *d covr.1547 2 207,780,781,452,455,456,25,24,32,33,41,112,115,116/ *d covr.1549 1 5,4,7,6,3,3,2,5,3,1,6,4,2,2,2,6,4,3,5,5,5,5,5,5,6,6,3,3,4,5,3,3/ *d covr.1552,1554 2 'n]a<)','np)',')','cont.)',']g<)','p)','d)','t)', 3 'e3)',']a<)','2p)','e)',']a<0)',']a<1)',' ]n<)',' ]n<)', 4 ' ]n<)','3n]a<)','2n]a<)','nd)','nt)','2np)','p]a<)','pd)', 5 'pt)'/ *d covr.1556 1 nmed1/'(total'/ *d covr.1563 *d covr.1567 *d covr.1585,1587 */ covr -- 5nov98 */ correct error in color plot option (muir, iaea). *d covr.255 if (nplot.ne.0) write(nplot,'(''1 2 .22'',i3,''/'')') icolor */ covr -- 12nov98 */ correct diagnostic for empty plots and fix shading */ for color plots, b-w plots, and ndiv greater than 1. *d covr.1223,1230 call patlev(jpat,ilevel) if (jpat.ne.0) ishade=ishade+1 *d covr.1260,1318 c ***positive part tlow=0. do 420 ilevel=1,nlev call patlev(jpat,ilevel) c ***switch orientation for positive part of legend on b-w plots if(icolor.eq.0.and.jpat.ne.0) jpat=jpat+10 write(nplot,'(''0 0 0 0 0'',i3,''/'')') jpat write(nplot,'(''0'')') thi=a(ixlev-1+ilevel) write(nplot,'(''0. '',f6.3,''/'')') tlow write(nplot,'(''1. '',f6.3,''/'')') tlow write(nplot,'(''1. '',f6.3,''/'')') thi write(nplot,'(''0. '',f6.3,''/'')') thi write(nplot,'(''0. '',f6.3,''/'')') tlow write(nplot,'(''/'')') if (ilevel.lt.nlev) write(nplot,'(''2/'')') if (ilevel.lt.nlev) write(nplot,'(''/'')') tlow=thi 420 continue *d covr.1330,1388 c ***negative part tlow=0. do 440 iloop=1,nlev ilevel=-iloop call patlev(jpat,ilevel) write(nplot,'(''0 0 0 0 0'',i3,''/'')') jpat write(nplot,'(''0'')') thi=-a(ixlev-1+iloop) write(nplot,'(''0. '',f6.3,''/'')') tlow write(nplot,'(''1. '',f6.3,''/'')') tlow write(nplot,'(''1. '',f6.3,''/'')') thi write(nplot,'(''0. '',f6.3,''/'')') thi write(nplot,'(''0. '',f6.3,''/'')') tlow write(nplot,'(''/'')') if (iloop.lt.nlev) write(nplot,'(''2/'')') if (iloop.lt.nlev) write(nplot,'(''/'')') tlow=thi 440 continue *i covr.1432 subroutine patlev(jpat,ilevel) c ****************************************************************** c convert the correlation level to a color or b-w shading pattern c following the conventions of viewr c ****************************************************************** *if sw implicit real*8 (a-h,o-z) *endif common/cov0/icolor common/cov2/xcycle,nlev,ndiv,ixmin,ixmax, 1 izero,irelco,ismall,ishade c jpat=0 if (icolor.eq.0) then if (ilevel.gt.1) jpat=20-(nlev-ilevel)*2/ndiv if (ilevel.lt.-1) jpat=40-(nlev+ilevel)*2/ndiv else if (ilevel.gt.1) jpat=50-(nlev-ilevel)*2/ndiv if (ilevel.lt.-1) jpat=60-(nlev+ilevel)*2/ndiv endif return end *ident up30 */ njoy -- 5nov98 */ add a comment-card capability for njoy input decks. any */ lines starting with "-- " will be ignored (from muir, iaea). *i njoy.142 c c ***comment card c ***treat comment card as a null module else if (module.eq.'--') then go to 100 */ njoy -- 2nov98 */ increase size of field in storage-error diagnostic (muir, iaea) *d njoy.2517 1 '(''need '',i7,'' more words for id '',a4,''.'')') need,id *ident up31 */ reconr -- 5nov98 */ fix typographical error (noted by p.deleege, tu/delft) *d up1.73,77 *if sw phi=dcmplx(phro,phri) *else phi=cmplx(phro,phri) *endif *ident up32 */ unresr -- 5nov98 */ fix misplaced update from up12. Noticed by */ pelloni (psi) and deleege (tu/delft). *d up12.22 4 63000.,80000.,100000./ *d unresr.377 data wide/3.01/ *ident up33 */ purr -- 5nov98 */ fix misplaced update from up13. Noticed by */ pelloni (psi) and deleege (tu/delft). *d up13.50 4 63000.,80000.,100000./ *d purr.490 data wide/3.01/ *ident up34 */ reconr -- 5nov98 * allow for more unresolved-range grid energies to handle */ the case of fe-56 from eff-3.0 (fendl-2.0). it corrects */ a fatal error when the ur range exceeds 100 kev. this */ patch was provided by trkov (ijs/slovenia). *d reconr.491 dimension egridu(58) *d reconr.494 data ngridu/58/ *d reconr.501 5 63000.d0,80000.d0,100000.d0,125000.d0,150000.d0,200000.d0, 6 320000.d0,400000.d0,500000.d0,630000.d0,800000.d0,1.d6, 7 1.25d6,1.5d6,2.0d6,3.2d6,4.0d6,5.0d6,6.3d6,8.0d6/ *d reconr.512 5 63000.,80000.,100000.,125000.,150000.,200000., 6 320000.,400000.,500000.,630000.,800000.,1.e6, 7 1.25e6,1.5e6,2.0e6,3.2e6,4.0e6,5.0e6,6.3e6,8.0e6/ *ident up35 */ unresr -- 5nov98 */ corresponding changes to unresr (see up35) *d unresr.358 dimension egridu(58) *d unresr.362 data ngridu/58/ *d unresr.369 5 63000.d0,80000.d0,100000.d0,125000.d0,150000.d0,200000.d0, 6 320000.d0,400000.d0,500000.d0,630000.d0,800000.d0,1.d6, 7 1.25d6,1.5d6,2.0d6,3.2d6,4.0d6,5.0d6,6.3d6,8.0d6/ *d up32.6 5 63000.,80000.,100000.,125000.,150000.,200000., 6 320000.,400000.,500000.,630000.,800000.,1.e6, 7 1.25e6,1.5e6,2.0e6,3.2e6,4.0e6,5.0d6,6.3e6,8.0e6/ *ident up36 */ purr -- 5nov98 */ corresponding changes to purr (see up35) *d purr.472 dimension egridu(58) *d purr.475 data ngridu/58/ *d purr.482 5 63000.d0,80000.d0,100000.d0,125000.d0,150000.d0,200000.d0, 6 320000.d0,400000.d0,500000.d0,630000.d0,800000.d0,1.d6, 7 1.25d6,1.5d6,2.0d6,3.2d6,4.0d6,5.0d6,6.3d6,8.0d6/ *d up33.6 5 63000.,80000.,100000.,125000.,150000.,200000., 6 320000.,400000.,500000.,630000.,800000.,1.e6, 7 1.25e6,1.5e6,2.0e6,3.2e6,4.0e6,5.0d6,6.3e6,8.0e6/ *ident up37 */ groupr -- 5nov98 */ the addition of the logic to automatically determine the */ version of endf being used damaged the writing of the id */ record on the gendf file. in addition, the count for the */ title field was undefined from subroutine ruinb. this */ patch was modified from one provided by trkov (ijs/sovenia). *d groupr.292,298 *d groupr.316,317 math=1 mfh=0 mth=0 nsh=0 itend=0 nwds=17 text=' ' read(text,'(16a4,a2)') (a(iscr+i-1),i=1,nwds) if (ngout1.eq.0) call tpidio(0,ngout2,0,a(iscr),nb,nwds) *d groupr.749 common/gout1/ngout1,ngout2,nt */ groupr -- 5nov98 */ fix a problem with the input of a weighting function. */ it was introduced by the changeover from infree to read*. */ discovered by trkov (ijs/slovenia). *d groupr.2140,2142 *d groupr.2144,2148 nr=nint(a(iwght+4)) np=nint(a(iwght+5)) iw=6+2*nr+2*np *ident up38 */ acer -- 5nov98 */ previous versions of acer converted cm legendre polynomials */ in file6/law1 into kalbach distributions, but this does not */ always provide a good representation of the evaluator's */ intent. this patch converts these data into a laboratory */ file6/law7 form instead. it is needed for fe-56 from */ eff-3.1. the patch was contributed by trkov (ijs/slovenia). */ mcnp-4 will provide direct sampling from the original cm */ distribution, and this patch will no longer be effective. *d acer.2063 c set flag to convert file-6 Law-7 data representation to Law-7 if (lf.eq.1) new6=1 *d acer.2093 c ***convert mf6,law1 lab or cm distributions to law7 *d acer.3299,3300 c convert a section of file 6 using legendre or tabulated c angular distributions into law 7 format with 33 angles. *d acer.3317 ncos=33 *i acer.3321 awr=c2h lct=l2h c force laboratory coordinate system a(4)=1 *i acer.3352 ein=a(l1+1) *d acer.3380 *d acer.3398 csn=acos elb=ep drv=1 if (lct.eq.1) go to 296 if (ep.le.0) go to 296 aw1=awr+1 clb=csn c minimum lab cosine (= zero particle energy in C.M.) cmn=-1 qq=1-aw1*aw1*ep/ein if (qq.gt.0) cmn=sqrt(qq) if (clb.lt.cmn) then clb=cmn c zero distribution when far outside valid cosine range if(imu.lt.nmu.and.amu(imu+1).le.cmn) drv=0 endif c outgoing particle energy in the laboratory system qq=ep-ein*(1-clb*clb)/(aw1*aw1) if (qq.lt.0) qq=0 elb=clb*sqrt(ein)/aw1+sqrt(qq) elb=elb*elb c calculate corresponding cosine in the cm system csn=clb*sqrt(elb/ep)-sqrt(ein/ep)/aw1 c check the limits qq=-1 if (csn.lt.qq) csn=qq qq=1 if (csn.gt.qq) csn=qq 296 if (lang.gt.10) go to 301 c calculate the probability from Legendre polynomials call legndr(csn,p,nl) *i acer.3403 c calculate the probability from pointwise representation *d acer.3412,3413 305 a(j)=elb a(j+1)=fmu*drv *ident up39 */ acer -- 5nov98 */ in up14, the length of the hz string in common/mis1 was */ set to 13 in preparation for the new mcnpx format, but */ it also needs to be set in other routines. this was */ being done in the mcnp4/mcnpx jumbo acer patch (to be */ released later), but these changes are being moved here */ to keep the compilers happy while waiting for the */ jumbo patch. recommended by trkov (ijs/slovenia). *d acer.187 character hz*13,hd*10,hk*70,hm*10 *d acer.4553 character hz*13,hd*10,hk*70,hm*10 *d acer.6609 character hz*13,hd*10,hk*70,hm*10 *d acer.7869 character hz*13,hd*10,hk*70,hm*10 *d acer.9332 character hz*13,hd*10,hk*70,hm*10 *d acer.9467 character hz*13,hd*10,hk*70,hm*10 *d acer.9637 character hz*13,hd*10,hk*70,hm*10 *d acer.9742 character hz*13,hd*10,hk*70,hm*10 *d acer.9891 character hz*13,hd*10,hk*70,hm*10 *d acer.9971 character hz*13,hd*10,hk*70,hm*10 *d acer.10066 character hz*13,hd*10,hk*70,hm*10 *d acer.10384 character hz*13,hd*10,hk*70,hm*10 *d acer.10476 character hz*13,hd*10,hk*70,hm*10 *d acer.10558 character hz*13,hd*10,hk*70,hm*10 *d acer.10761 character hz*13,hd*10,hk*70,hm*10 *d acer.11081 character hz*13,hd*10,hk*70,hm*10 *ident up40 */ wimsr -- 5nov98 */ the following wimsr updates were provided by a.trkov, */ (ijs/slovenia). define nfid from rdfid, and make */ sure that some flags are consistent. *i wimsr.168 nfid=rdfid+0.001 *i wimsr.181 if (ifprod.gt.0) then ifprod=1 if (ires.gt.0) ifprod=2 endif *d wimsr.1861,1862 jfis=ifis if (ifprod.eq.1) jfis=-1 if (ifprod.eq.2) jfis=-2 */ wimsr -- 5nov98 */ extend absorption and multiple neutron emission reactions. *d wimsr.929 if (mth.eq.16.or.mth.eq.24) go to 236 if (mth.eq.17.or.mth.eq.25) go to 237 *i wimsr.982 c c ***n3n, n3n-alpha cross section 237 loc=in2n+jg-1 loca=l+lz+nl*nz a(loc)=a(loc)+2.*a(loca) go to 360 */ wimsr -- 5nov98 */ improve input instructions for burnup data. */ add delimiter records when burnup data are present. *d wimsr.79 c * the following cards 5 and 6 are for iburn gt 0 only * *i wimsr.81 c * for burnable materials ntis=2 * c * for fissile materials ntis>2 when fission product * c * yields are given. *d wimsr.84 c * card 6a * c * identa ident of capture product isotope * c * yield yield of product identa from capture * c * * c * card 6b * c * identa ident of decay product isotope (zero if stable) * c * lambda decay constant (s-1) * c * * c * card 6c (repeated ntis-2 times, if necessary) *i wimsr.1875 write(nout,'(2i15)') 999999999,3 *i wimsr.1879 write(nout,'(2i15)') 999999999,4 *d wimsr.1875 if (jcc.ge.8) ifisp(4)=jfis *ident up41 */ acer -- 05nov98 */ fix a problem with the dimensions of two arrays that was */ discovered by broeders (fzk). *d acer.5153 dimension xs(20),ys(20),xxs(200),yys(200) *d acer.5155 data ismax/20/, jsmax/200/ *d acer.5448 if (js.gt.jsmax) js=jsmax *ident up42 */ plotr -- 17nov98 */ fix some problems noted in plotr by m.mattes (u.stuttgart). */ the values of xtag and ytag are not being read in. factx and */ and facty are being applied to the x and y values in plotr, */ so they need to be reset to unity before going on to viewr. */ this is necessary to prevent the factors from being applied */ to the coordinates twice. *d plotr.452 read(nsysi,*) itype,jtype,igrid,ileg,xtag,ytag *d plotr.1218,1219 factx1=1 facty1=1 write(nplt,'(i4,i8,7f7.2,''/ 3d plot'')') iplot,iwcol, 1 factx1,facty1,xll,yll,ww,wh,wr *d plotr.1357,1358 factx=1 facty=1 write(nplt,'(i4,i8,7f7.2,''/ 3d plot'')') iplot,iwcol, 1 factx1,facty1,xll,yll,ww,wh,wr *d plotr.1693,1694 factx1=1 facty1=1 write(nplt,'(i4,i8,7f7.2,''/ 3d plot'')') iplot,iwcol, 1 factx1,facty1,xll,yll,ww,wh,wr *d plotr.1762,1763 610 factx1=1 facty1=1 write(nplt,'(i4,i8,7f7.2,''/ 2d plot'')') iplot,iwcol, 1 factx1,facty1,xll,yll,ww,wh,wr *ident up43 */ unresr -- 20nov98 */ the changes to the unresolved energy grid in up34 were */ not quite correct for the sw machines. trkov. *d unresr.365,368 1 63.d0,80.d0,100.d0,125.d0,150.d0,200.d0,250.d0,320.d0,400.d0, 2 500.d0,630.d0,800.d0,1000.d0,1250.d0,1500.d0,2000.d0,2500.d0, 3 3200.d0,4000.d0,5000.d0,6300.d0,8000.d0,10000.d0,12500.d0, 4 15000.d0,20000.d0,25000.d0,32000.d0,40000.d0,50000.d0, *ident up44 */ purr -- 20nov98 */ the changes to the unresolved energy grid in up35 were */ not quite correct for the sw machines. trkov. *d purr.478,481 1 63.d0,80.d0,100.d0,125.d0,150.d0,200.d0,250.d0,320.d0,400.d0, 2 500.d0,630.d0,800.d0,1000.d0,1250.d0,1500.d0,2000.d0,2500.d0, 3 3200.d0,4000.d0,5000.d0,6300.d0,8000.d0,10000.d0,12500.d0, 4 15000.d0,20000.d0,25000.d0,32000.d0,40000.d0,50000.d0, *ident up45 */ covr -- 20nov98 */ one change was left out of muir's covr patches in up29. */ it is needed to get correct labels for plots of the */ high threshold reactions. *d covr.1594 do 300 j=1,32 *ident up46 */ njoy -- 4dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision, and that any conversion */ of real numbers to integers is explicitly flagged for */ truncation or for nearest integer. *d njoy.1062 if (nin.eq.0.and.a(5).le.npage) nw=6+nint(a(5)) *d njoy.1213 zero=0.d0 *d njoy.1238 f=0 *d njoy.1246 n=int(log10(abs(xx))) *d njoy.1281 130 n=int(log10(abs(x))) *d njoy.1287 f=f/ten *d njoy.1594 z=0 *d njoy.1627 z=0 *d njoy.1661 z=0 *d njoy.1696 z=0 *d njoy.1831 x=j x=x/nl *d njoy.1834 inow=nint(xnow) *i njoy.1940 zero=0 *d njoy.1963 500 if (y1.eq.zero) go to 100 *d njoy.1968 600 if (y1.eq.zero) go to 100 *d njoy.1971 b=b/(1/t-1/sqrt(x2-thr6)) *d njoy.1993 f=0 *i njoy.2062 *if sw data shade/1.00001d0/ data xbig/1.d10/ data zero/0.d0/ *else data shade/1.00001e0/ data xbig/1.e10/ data zero/0.e0/ *endif *d njoy.2114,2116 150 if (x.lt.shade*a(jp)) go to 160 y=0 xnext=xbig *d njoy.2119,2120 xnext=xbig if (y.gt.zero) xnext=shade*a(jp) *d njoy.2124 170 y=0 *i njoy.2139 *if sw data break/0.1d0/ data zero/0.d0/ data one/1.d0/ *else data break/0.1e0/ data zero/0.e0/ data one/1.e0/ *endif *d njoy.2142 gral=0 *d njoy.2156 gral=(x2-x1)*(a+b*(x2+x1)/2) *d njoy.2160 300 if ((xl.le.zero).or.(xh.le.zero)) go to 200 *d njoy.2163,2165 if (abs(z).gt.break) go to 310 gral=(x2-x1)*(yl+b*log(x1/xl))+b*x1*z*z*(1+ 1 z*(-one/3+z*(one/6-z/10)))/2 *d njoy.2171 400 if ((yl.lt.zero).or.(yh.lt.zero)) go to 200 *d njoy.2175,2176 if (abs(z).gt.break) go to 410 gral=exp(a+b*x1)*(x2-x1)*(1+z*(one/2+z/6)) *d njoy.2178 410 gral=exp(a+b*x1)*(exp(z)-1)/b *d njoy.2182,2183 500 if ((xl.le.zero).or.(xh.le.zero))go to 400 if ((yl.le.zero).or.(yh.le.zero)) go to 300 *d njoy.2185,2187 z=(b+1)*log(x2/x1) if (abs(z).gt.break) go to 510 gral=yl*x1*((x1/xl)**b)*log(x2/x1)*(1+z*(one/2+z/6)) *d njoy.2189 510 gral=yl*x1*((x1/xl)**b)*(((x2/x1)**(b+1))-1)/(b+1) *i njoy.2206 *if sw data down/.99999d0/ data xbig/1.d10/ data zero/0.d0/ *else data down/.99999e0/ data xbig/1.e10/ data zero/0.e0/ *endif *d njoy.2210 if (x.gt.zero) go to 110 *d njoy.2227,2228 if (ylast.ne.zero) xlast=down*xlast if (ylast.ne.zero) return *d njoy.2230 if (a(loc+3).ne.zero) return *d njoy.2283 200 y1=0 *d njoy.2289,2290 xlast=xbig xnext=xbig *i njoy.2309 *if sw data down/.99999d0/ data zero/0.d0/ data xbig/1.d10/ *else data down/.99999e0/ data zero/0.e0/ data xbig/1.e10/ *endif *d njoy.2313 if (x.gt.zero) go to 110 *d njoy.2328 if (ylast.ne.zero) xlast=down*xlast *d njoy.2368 200 y2=0 *d njoy.2374,2375,2290 xlast=xbig xnext=xbig *d njoy.2408 a(i)=0 *i njoy.2436 zero=0 *d njoy.2455 if (a(3*i+4).ne.zero) nidtot=nidtot+1 *d njoy.2486 if (a(3*jd+4).eq.zero) go to 140 *i njoy.2655 data zero/0.d0/ *i njoy.2667 data zero/0.e0/ *d njoy.2669,2670 xx=0 if (x.eq.zero) go to 110 *d njoy.2673 if (aa.lt.zero) ipwr=ipwr-1 *d njoy.2672 ipwr=int(aa) *ident up47 */ reconr -- 4dec98 */ make sure this step factor is not a multiple of 2. */ without this patch, the if statement can be sensitive */ to differences in the 10th or 11th place, thus giving */ different results with different compilers. *i reconr.1606 data estp/4.1d0/ *i reconr.1610 data estp/4.1e0/ *d reconr.1753 145 est=estp*(a(ix+1-1)-res(1)) if (in.gt.3.and.dx.gt.est) go to 175 */ reconr -- 4dec98 */ at low energies in nonfissionable materials, the r-function */ in csrmat is obtained from a degenerate r matrix that is very */ close to unity. this leads to difficult numerics for low */ incident energies. this patch rearranges the terms to allow */ unit terms to cancel out, thus giving better results for the */ total and capture cross sections for materials in the */ fe,ni,cr range. the old approach could show differences */ between different machines starting in the 6th significant */ figure of the capture cross section. these differences are */ only significant when comparing results from different systems. *i reconr.2288 data small/3.d-4/ *i reconr.2294 data small/3.e-4/ *d reconr.2368 *i reconr.2418 r(1,1)=uno+r(1,1) r(2,2)=uno+r(2,2) r(3,3)=uno+r(3,3) *i reconr.2432 u11r=p1*(two*ri(1,1)-uno)+two*p2*si(1,1) u11i=p2*(uno-two*ri(1,1))+two*p1*si(1,1) termt=two*gj*(uno-u11r) termn=gj*((uno-u11r)**2+u11i**2) *d reconr.2435,2445 160 dd=r(1,1) rr=uno+dd ss=s(1,1) amag=rr**2+ss**2 rri=rr/amag ssi=-ss/amag uur=p1*(two*rri-uno)+two*p2*ssi uui=p2*(uno-two*rri)+two*p1*ssi if (abs(dd).lt.small.and.abs(phid).lt.small) then xx=2*dd xx=xx+2*(dd*dd+ss*ss+phid*phid+p2*ss) xx=xx-2*phid*phid*(dd*dd+ss*ss) xx=xx/amag termt=two*gj*xx termn=gj*(xx**2+uui**2) else termt=two*gj*(uno-uur) termn=gj*((uno-uur)**2+uui**2) endif termf=0 170 if (jj.gt.jjl.and.jj.lt.numj) then */ reconr -- 4dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. *i reconr.102 *if sw data elow,ehigh,elarge/1.d-5,20.d6,9.d9/ *else data elow,ehigh,elarge/1.e-5,20.e6,9.e9/ *endif *d reconr.193,194 zain=1 awin=1 *d reconr.198 zain=int(nsub/10) *d reconr.221,225 eresl=elow eresr=ehigh eresh=ehigh eresu=elarge eresm=elarge *i reconr.288 zero=0 *d reconr.293,295 tempr=0 errmax=-1 errint=-1 *d reconr.297 if (errmax.le.zero) errmax=20*err *d reconr.299 if (errint.le.zero) errint=err/10000 *d reconr.488 common/runits/nin,nout,nscrr,nscrs *d reconr.3754 common/runits/nin,nout,nscrr,nscrs *i reconr.506 data elow/1.d-5/ data elarge/9.d9/ data small/1.d-6/ data wide/3.d0/ *i reconr.517 data elow/1.e-5/ data elarge/9.e9/ data small/1.e-6/ data wide/3.e0/ *d reconr.530,536 eresl=elarge eresh=0 eresr=0 eresu=elarge eresm=elarge zero=0 *d reconr.569 if (abs(el-elow).gt.small) go to 107 *d reconr.666 if (enut.gt.ener+ener/100) go to 148 *d reconr.713 if (ascatl.eq.zero) ascatl=ascat *d reconr.759 a(jj+2)=0 *d reconr.776 if (enut.gt.ener+ener/100) go to 229 *d reconr.846 if (enut.gt.ener+ener/100) go to 339 *d reconr.871 hw=0 *d reconr.929 el=0 *d reconr.1009 se=0 *i reconr.1057 *if sw data big/1.d10/ *else data big/1.e10/ *endif *i reconr.1064 zero=0 *d reconr.1080 a(isunr+12)=big *d reconr.1083 e=0 *d reconr.1092 a(l+k)=0 *d reconr.1126 e=0 *d reconr.1131 if (e.lt.zero) go to 160 *i reconr.1205 elow=1.d-5 elim=.99d6 up=1.001d0 dn=0.999d0 *i reconr.1212 elow=1.e-5 elim=.99e6 up=1.001e0 dn=0.999e0 *i reconr.1233 zero=0 *d reconr.1277 if (mfh.eq.13) a(iscr+4)=1 *d reconr.1306 if (thr6.lt.0) thr6=0 *d reconr.1310 if (a(l+1).eq.zero) go to 183 *d reconr.1314 a(l+1)=0 *d reconr.1372 220 test=elow *d reconr.1376 if (sr.eq.zero) go to 225 *d reconr.1406 245 et=0 *d reconr.1475,1476 if (a(ix+i-1).gt.elim) go to 320 *d reconr.1479 if (a(ix+i-2).gt.up*xm.and.a(ix+i-1).lt.dn*xm) go to 315 *d reconr.1482 if (a(ix+i-2).gt.up*xm.and.a(ix+i-1).lt.dn*xm) go to 315 *d reconr.1485 if (a(ix+i-2).gt.up*xm.and.a(ix+i-1).lt.dn*xm) go to 315 *d reconr.1488 if (a(ix+i-2).gt.up*xm.and.a(ix+i-1).lt.dn*xm) go to 315 *d reconr.1497,1498 if (a(ix+i-1).gt.elim) go to 325 *d reconr.1552 egl=0 *d reconr.1556 if (eg.lt.zero) ngneg=ngneg+1 *d reconr.1631,1636 crej=0 cmax=0 cint=0 eint=0 erej=0 emax=0 zero=0 *d reconr.1649 e=0 *d reconr.1700,1702 cb=(c1+c2)/2 if (dc.lt.cb/2.and.i.ge.3) go to 120 crej=crej+dc*dx/(2*xm) *d reconr.1706 crej=crej+abs(cc*dx3/(6*dx2))/xm *d reconr.1710,1712 cb=(c1+c2)/2 if (dc.lt.cb/2.and.i.ge.3) go to 130 erej=erej+dc*dx/(2*xm) *d reconr.1716 erej=erej+abs(cc*dx3/(6*dx2))/xm *d reconr.1727 fr1=1-fr2 *d reconr.1730,1731 dm2=0 if (sig(3).eq.zero) go to 140 *d reconr.1748,1749 cmax=cmax+dm3*dx/(2*xm) emax=emax+dm1*dx/(2*xm) *d reconr.1770 cint=cint+(c1+c2)*dx/(2*xm) *d reconr.1773 eint=eint+(c1+c2)*dx/(2*xm) *d reconr.1778,1782 if (cint.eq.zero.or.eint.eq.zero) go to 160 cmax=100*cmax/cint emax=100*emax/eint crej=100*crej/cint erej=100*erej/eint *d reconr.1798,1803 cint=0 eint=0 cmax=0 emax=0 crej=0 erej=0 *d reconr.1879 sig(i)=0 *i reconr.1880 zero=0 *d reconr.1917 if (sigp(j).lt.zero) sigp(j)=0 *d reconr.1951 sigp(i)=0 *i reconr.1992,1993 zero=0 if (tempr.gt.zero) call findex('tr',itr,aa) if (tempr.gt.zero) call findex('ti',iti,aa) *d reconr.1997 sigp(i)=0 *d reconr.2030 pec=0 *d reconr.2041 if (tempr.eq.zero) go to 35 *d reconr.2157 zero=0 if (tempr.gt.zero) *d reconr.2163 sigp(i)=0 *d reconr.2192 pec=0 *d reconr.2200 sum=0 *d reconr.2204 nj=nint(ajmax-ajmin+1) *d reconr.2214 sigj(i,ii)=0 *d reconr.2222 j=nint(a(inow+1)-ajmin+1) *d reconr.2298 zero=0 if (tempr.gt.zero) *d reconr.2304 sigp(i)=0 *d reconr.2324,2326 gfa=0 gfb=0 gf=0 *d reconr.2337,2338 if (apl.ne.zero) rhoc=k*apl if (apl.ne.zero.and.naps.eq.1) rho=k*apl *d reconr.2351 numj=nint(ajmax-ajmin+1) *d reconr.2365,2366 s(j,i)=0 r(j,i)=0 *d reconr.2376,2377 if (abs(aj-ajc).gt.quar) go to 140 *d reconr.2387,2392 a2=0 if (apl.ne.zero) rhoc=k*apl if (apl.ne.zero.and.naps.eq.1) rho=k*apl a3=0 if (gfb.ne.zero) a3=sqrt(abs(gfb)) if (gfb.lt.zero) a3=-a3 *d reconr.2402 if (gfa.eq.zero.and.gfb.eq.zero) go to 140 *d reconr.2413 gf=1 *d reconr.2417 if (gf.eq.zero) go to 160 *d reconr.2558 c(i,j)=0 *d reconr.2599 if (tempr.gt.zero) call error('cshybr', *d reconr.2604 sigp(i)=0 *d reconr.2624 sigc(n)=0 *d reconr.2685 gc(ii)=0 *d reconr.2819 zero=0 if (e.gt.zero) go to 125 *d reconr.2830 elast=0 *d reconr.2832 elst(isect)=0 *d reconr.2858 sigp(i)=0 *d reconr.2891,2892 gfx=0 gxx=0 *d reconr.2915 mu=int(amun) *i reconr.3019 data small/1.d-8/ *i reconr.3025 data small/1.e-8/ *d reconr.3029 zero=0 if (e.gt.zero) go to 125 *d reconr.3042 elast=0 *d reconr.3044 elst(isect)=0 *d reconr.3068 sigp(i)=0 *d reconr.3115,3117 if (gxx.lt.small) gxx=0 if (gfx.lt.small) gfx=0 *d reconr.3213,3214 zero=0 if (tempr.gt.zero) call findex('tr',itr,aa) if (tempr.gt.zero) call findex('ti',iti,aa) *d reconr.3218 sigp(i)=0 *d reconr.3229 if (tempr.gt.zero) delta=1/sqrt(4*bk*tempr*e/awri) *d reconr.3240,3242 bakt=0 bakf=0 bakc=0 *d reconr.3277 if (tempr.gt.zero) go to 155 *i reconr.3319 *if sw data test/1.d-6/ *else data test/1.e-6/ *endif *d reconr.3327 *d reconr.3328 if ((phi/rho).lt.test) phi=0 *d reconr.3331 *d reconr.3332 if ((phi/rho).lt.test) phi=0 *d reconr.3338 *d reconr.3339 if ((phi/rho).lt.test) phi=0 *d reconr.3411 zero=0 s=0 *d reconr.3412,3417 if (galpha.le.zero) go to 1000 if (gamma.le.zero) go to 1000 if (gbeta.lt.zero) go to 1000 if (gbeta.gt.zero) go to 350 if (df.lt.zero) go to 1000 if (df.gt.zero) go to 250 *d reconr.3442,3443 350 if (df.lt.zero) go to 1000 if (df.gt.zero) go to 450 *i reconr.3519 data finity/.99d10/ data small/1.d-8/ *i reconr.3521 data finity/.99e10/ data small/1.e-8/ *d reconr.3535,3536 zero=0 *d reconr.3540 tot(i+1)=0 *d reconr.3560,3562 inn=in if (ig.eq.ngp.and.en.eq.eg) inn=-in call loada(inn,tot,ntot,iold,a(ibufo),nbuf) *d reconr.3591 e=0 *d reconr.3594 if (sn.ne.zero) thresh=sigfig(thresh,7,-1) *d reconr.3620 sn=0 *d reconr.3623 if (thresh.gt.1..and.abs(thresh-eg).lt.test) sn=0 *d reconr.3626 if (eg.ge.eresr.and.eg.lt.eresh.and.itype.gt.0) sn=0 *d reconr.3648,3650 inn=in if (ig.eq.ngo) inn=-in call loada(inn,tot,2,ngrid,a(ibufg),nbufg) *d reconr.3695 a(incs+nxc-1)=a(incs+nxc-1)+2+int((n2h+2)/3) *i reconr.3770 *if sw data emin/1.d-5/ *else data emin/1.e-5/ *endif *d reconr.3802,3803 a(iscr)=0 a(iscr+1)=0 *d reconr.3812,3816 a(iscr+1)=0 a(iscr+2)=0 a(iscr+3)=0 a(iscr+4)=int(10*zain) a(iscr+5)=0 *d reconr.3854,3855 a(j+1)=0 a(j+2)=0 *d reconr.3859 a(j+5)=3+int((np+2)/3) *d reconr.3900 a(iscr+1)=0 *d reconr.3902 a(iscr)=emin *d reconr.3937 a(iscr+2)=0 *d reconr.3951 a(iscr+1)=0 *d reconr.4055,4056 zero=0 aki=1 if (ax.lt.zero) aki=-1 *d reconr.4066,4067 10 ii=int(x*10) jj=int(y*10) *i reconr.4124 *if sw data tenth/.1d0/ *else data tenth/.1e0/ *endif *d reconr.4126,4129 x0=-tenth y0=-tenth dx=tenth dy=tenth *d reconr.4179,4180 zero=0 rew=0 aimw=0 *d reconr.4183,4185 if (abrez+aimz.ne.zero) go to 20 rew=1 aimw=0 *d reconr.4192,4195 if (abrez+brk1*aimz-brk2.gt.zero) go to 350 if (abrez+brk3*aimz-brk4.gt.zero) go to 350 if (r2+brk5*ai2-brk6.lt.zero) go to 340 if (r2+brk7*ai2-brk8.ge.zero) go to 340 if (aimz-brk9.ge.zero) go to 350 *d reconr.4200 if (aim1.ge.zero) go to 370 *d reconr.4209,4213 h=0 b=0 a=0 tempm=0 temel=0 *d reconr.4239,4246 c=0 b=0 ajsig=0 d=0 jsig=0 g=0 h=0 el=0 *d reconr.4287 if ((abs(tempm)+abs(temel)-up).lt.zero) go to 500 *d reconr.4297 500 if ((abs(tempm)+abs(temel)-dn).gt.zero) go to 520 *ident up48 */ broadr -- 11dec98 */ make sure this step factor is not a multiple of 2. */ without this patch, the if statement can be sensitive */ to differences in the 10th or 11th place, thus giving */ different results with different compilers. *i broadr.997 data estp/4.1d0/ *i broadr.1004 data estp/4.1e0/ *d broadr.1105 est=estp*(es(is)-tt(1)) if (j.gt.3.and.dx.gt.est) go to 170 */ broadr -- 11dec98 */ tighten tolerance for printing thermal quantities. */ fix problem in capture resonance integrals. *d up11.31 if (abs(tev-therm).gt.therm/1000) go to 2250 *d up2.60 ssc=(sc/enow+ss/elast)*(enow-eone)/2 */ broadr -- 11dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. *i broadr.119 fact=.99999d0 onemev=1.d6 tbreak=0.5d0 big=9.d9 step=1.01d0 *i broadr.122 fact=.99999e0 onemev=1.e6 tbreak=0.5e0 big=9.e9 step=1.01e0 *d broadr.130 *i broadr.131 zero=0 *d broadr.158 temp1=0 *d broadr.162.164 thnmx=onemev errmax=0 errint=0 *d broadr.166,167 if (errmax.eq.zero) errmax=20*errthn if (errint.eq.zero) errint=errthn/10000 *d broadr.244,245 emin=1 if (thnmax.lt.zero) emin=-thnmax *d broadr.276 test=1+temp1/1000 *d broadr.290 enow=0 *d broadr.303 if (en.gt.zero) sig=sig-sun *d broadr.306,307 if (enext.ge.big) kt=-kt *d broadr.311 if (enext.lt.big) go to 150 *d broadr.320 enow=0 *d broadr.340 if (en.gt.zero) sig=sig-sun *d broadr.352 if (thnmax.lt.zero) thnmax=-thnmax *d broadr.384 break=9999 *d broadr.401 if (mti(k).eq.0.and.tt(1+k).gt.zero) mti(k)=j *d broadr.428 if (mt.eq.1.or.mt.eq.3.or.mt.eq.19) a(iscr+i+3)=3+int((n2out+2)/3) *d broadr.435 a(iscr+i+3)=3+int((n2out-mt4i+3)/3) *d broadr.439 if (mt.eq.mtr(j)) a(iscr+i+3)=3+int((n2out-jjj+3)/3) *d broadr.447 eone=tbreak *d broadr.460,469 fint=0 cint=0 alint=0 etint=0 v1int=0 ssf=0 slf=0 ssc=0 slc=0 elast=0 flast=0 fnul=0 *d up2.17,19 2240 enow=tt(1) if (enow.gt.step*elast.and.elast.gt.zero) then *d broadr.521 fnu=0 *d up11.42 break=1 *d up2.44 if (sf.gt.zero.and.slf.gt.zero) *d broadr.584 if (ctev.ne.zero) write(nsyso,'( *d broadr.596 if (ftev.ne.zero) write(nsyso,'( *d broadr.651 275 ee=0 *d broadr.681 if (en.gt.zero) a(k)=a(k)+sun *d broadr.684 a(k)=0 *d broadr.694 if (en.gt.zero) a(k)=a(k)+sun *i broadr.763 zero=0 *d broadr.811,812 if (i.eq.1.or.e.lt.zero) eumax=abs(e) if (e.lt.zero) a(l)=-e *i broadr.834 *if sw data up/1.00001d0/ *else data up/1.00001e0/ *endif *d broadr.837 sun=0 *d broadr.855 enext=eumax*up *i broadr.883 zero=0 *d broadr.911 120 if (tempef.eq.zero) go to 130 *i broadr.1012 zero=0 *d broadr.1028,1029 if (abs(dn).lt.abs(s(i,k))/1000) dn=0 dl(i)=1 *d broadr.1031 if (dn.lt.zero) dl(i)=-1 *d broadr.1054 if (abs(dn).lt.abs(s(i,k))/1000) dn=0 *d broadr.1055,1056 if (dn.ge.zero) then dn=1 *d broadr.1058 dn=-1 *i broadr.997 one=1.d0 *i broadr.1004 one=1.e0 *d broadr.1092 test=one-one/100 *d broadr.1097 si=f*ss(i,is-1)+(1-f)*ss(i,is) *d broadr.1101 if (dy*dx/2.gt.errint*em) go to 170 *d broadr.1122 if (abs(dn).lt.abs(s(i,k))/1000) dn=0 *d broadr.1124 dl(i)=1 *d broadr.1126 dl(i)=-1 *d broadr.1203 aa=0 *d broadr.1208 sbt(i)=0 *d broadr.1246,1247 x=0 xx=0 *d broadr.1295 xx=0 *d broadr.1323 if (sbt(i).lt.sigmin) sbt(i)=0 *d broadr.1353 *if sw data stpmax/1.24d0/ *else data stpmax/1.24e0/ *endif *i broadr.1354 zero=0 *d broadr.1372,1373 if (denom.eq.zero) go to 150 denom=1/denom *i broadr.1508 data aerr,rerr/1.d-30,1.d-8/ *i broadr.1512 data aerr,rerr/1.e-30,1.e-8/ *d broadr.1519 sign=1 *d broadr.1537,1538 cm(1)=1 xk=1 *d broadr.1562 beta=0 *d broadr.1570 160 xk=1 *d broadr.1578 test=aerr+rerr*abs(s) *ident up49 */ unresr -- 11dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. */ the common with the unit numbers was changed to be the */ same in every subroutine to avoid confusion. *d unresr.187,190 bkgz(1)=0 bkgz(2)=0 bkgz(3)=0 bkgz(4)=0 *d unresr.224 if (nx.gt.0.and.new.gt.0) a(iscr+5)=a(iscr+5)+1 *d unresr.232 if (nx.gt.0.and.new.gt.0) a(iscr+5)=a(iscr+5)+1 *d unresr.351 common/uunits/nendf,nin,nout,nscr *d unresr.353 common/unen/nunr,matd,nsigz,nmtx,mtx(5),npnts(5) *i unresr.379 *if sw data onemev/1.d6/ data small/1.d-12/ data step/1.01d0/ *else data onemev/1.e6/ data small/1.e-12/ data step/1.01e0/ *endif *d unresr.388,389 zero=0 elr=0 ehr=onemev *d unresr.391 a(ieunr)=onemev *d unresr.400 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.407 115 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.417 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.420 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.422 call moreio(nendf,0,0,a(iscr),nb,nw) *d unresr.426 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.428 call moreio(nendf,0,0,a(iscr),nb,nw) *d unresr.430 126 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.433 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.435 call moreio(nendf,0,0,a(iscr),nb,nw) *d unresr.469 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.476 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.505 190 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.527 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.537 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.558 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.567 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.575 call listio(nendf,0,0,a(iscr),nb,nw) *d unresr.578 call moreio(nendf,0,0,a(loc),nb,nw) *d unresr.486 a(inow+5)=0 *d unresr.596 if (a(jnow+k).le.zero) a(jnow+k)=small *d unresr.629,630 if (enext.ge.onemev) go to 380 *d unresr.635 if (enut.gt.step*et) go to 363 *d unresr.648 en=0 *d unresr.653 if (et.lt.zero) iovlp=1 *d unresr.706 common/uunits/nendf,nin,nout,nscr *i unresr.710 *if sw data up,dn/1.00001d0,0.99999d0/ *else data up,dn/1.00001e0,0.99999e0/ *endif *d unresr.717 100 call contio(nendf,0,0,a(iscr),nb,nw) *d unresr.724 a(ie+ibase)=0 *d unresr.727 120 call tosend(nendf,0,0,a(iscr)) *d unresr.730,731 e=0 call gety1(e,enext,idis,bkg,nendf,a(iscr)) *d unresr.736,738 if (ie.eq.1) e=up*e if (ie.eq.nunr) e=dn*e call gety1(e,enext,idis,a(ibase+ie),nendf,a(iscr)) *d unresr.740 call tosend(nendf,0,0,a(iscr)) *i unresr.796 data small/1.d-8/ *i unresr.819 data small/1.e-8/ *d unresr.829,830 if (temp(1).lt.one) temp(1)=one *d unresr.835,836 spot=0 sint=0 *d unresr.900,901 amux=0 gxx=0 *d unresr.903 amuf=0 *d unresr.909 gfx=0 *d unresr.950 tk(ks,is0,itp)=0 *d unresr.952 t(kx,is0,itp)=0 *d unresr.957 *d unresr.959 if (gfx.le.small) mu=0 *d unresr.962 if (gxx.lt.small) lu=0 *d unresr.1022,1026 yy(1)=0 yy(2)=0 yy(3)=0 yj=0 yk=0 *d unresr.1028,1029 xk=0 xj=0 *i unresr.1118 *if sw data small/1.d-8/ *else data small/1.e-8/ *endif *d unresr.1131,1133 if (gxx.lt.small) gxx=0 if (gfx.lt.small) gfx=0 *i unresr.1156 zero=0 *d unresr.1158,1159 aki=1 if (ax.lt.zero) aki=-1 *d unresr.1169,1170 10 ii=int(x*10) jj=int(y*10) *d unresr.1196 a2=2*x*y *i unresr.1241 zero=0 *d unresr.1244 if (ep.gt.zero) go to 20 *d unresr.1249,1252 y=0 z=0 yk=0 zk=0 *d unresr.1254 rki=0 *d unresr.1268 x2=1/(1+b1/rew) *i unresr.1297 *if sw data tenth/0.1d0/ *else data tenth/0.1e0/ *endif *d unresr.1299,1302 x0=-tenth y0=-tenth dx=tenth dy=tenth *d unresr.1352,1353 zero=0 rew=0 aimw=0 *d unresr.1356,1358 if (abrez+aimz.ne.zero) go to 20 rew=1 aimw=0 *d unresr.1363,1367 if (abrez+brk1*aimz-brk2.gt.zero) go to 350 if (abrez+brk3*aimz-brk4.gt.zero) go to 350 if (r2+brk5*ai2-brk6.lt.zero) go to 340 if (r2+brk7*ai2-brk8.ge.zero) go to 340 if (aimz-brk9.ge.zero) go to 350 *d unresr.1373 if (aim1.ge.zero) go to 370 *d unresr.1381,1385 h=0 b=0 a=0 tempm=0 temel=0 *d unresr.1411,1418 c=0 b=0 ajsig=0 d=0 jsig=0 g=0 h=0 el=0 *d unresr.1459 if ((abs(tempm)+abs(temel)-up).lt.zero) go to 500 *d unresr.1469 500 if ((abs(tempm)+abs(temel)-dn).gt.zero) go to 520 *ident up50 */ heatr -- 11 dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. *d heatr.735,739 qs=0 q0=0 icon=0 ebar=0 dame=0 *d heatr.787 q0=0 *d heatr.832 yld=1 *d heatr.919 elst=0 *d heatr.921,922 if (icon.eq.0) yld=0 yld0=0 *d heatr.930 if (new6.eq.1) c(npkk)=0 *d heatr.985 ebal6=0 *d heatr.993,996 if (izap.eq.1) h=0 if (izap.eq.0.and.local.eq.0) h=0 if (izap.le.1) dame=0 ebal6=0 *i heatr.2528 zero=0 one=1 *d heatr.2596 if (y(i-1,1).ne.zero) f1=y(i-1,2)/y(i-1,1) *d heatr.2651 if (yy.ne.zero) b=b/yy *d heatr.2708 f=one/2 *d heatr.2823 zero=0 if (ep.gt.zero) go to 200 *d heatr.3035 zero=0 if (ep.gt.zero) go to 120 *d heatr.3044 if (cnow(lnow).le.zero) lnow=lnow+ncnow *d heatr.3090 if (t.lt.zero) t=0 *d heatr.3122 if (t.lt.zero) t=0 *d heatr.3180 zero=0 if (i.gt.zero) go to 120 *d heatr.3208 if (t.lt.zero) t=0 *d heatr.3234 if (t.lt.zero) t=0 *d heatr.3470 awc=awr+1 h=h+df(xr,z,awc,z,awr)*y *d heatr.3526 zero=0 if (e.gt.zero) go to 200 *d heatr.3612 320 fle(i)=0 *d heatr.3620 410 fle(1)=1 *d heatr.3622 fle(il)=0 *d heatr.3630 440 fle(1)=1 *d heatr.3633 fle(il)=0 *d heatr.3754 220 fl(il)=0 *d heatr.3762 fl(il)=0 *d heatr.3802 if (abs(fl(il)).gt.zero) nlz=il *d heatr.3871 a(ie)=0 *d heatr.3875,3877 a(ir+m)=0 if (i.eq.j) a(ir+m)=1 a(iaa+m)=0 *d heatr.3885 g=1 *d heatr.3905,3906 a(iaa+(k-1)*imax+j-1)=0 a(ir+(k-1)*imax+j-1)=0 *d heatr.3909 g=1 *d heatr.3921 ysum=0 *d heatr.3937 a(ie)=0 *d heatr.3973 a(iscr+1)=0 *d heatr.3989 a(iscr+2)=0 *d heatr.4123 e=0 *d heatr.4176 cerr=0 *d heatr.4179 e=0 *ident up51 */ thermr -- 15dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. */ in addition, any conversion of real numbers to integers is */ explicitly flagged for truncation or for nearest integer. *i thermr.114 data up/1.00001d0/ data zero/0.d0/ *i thermr.122 data up/1.00001e0/ data zero/0.e0/ *d thermr.163,164 eftemp(i)=0 eftmp2(i)=0 *d thermr.185,186 sz2=0 az2=0 *d thermr.192 if (az2.ne.zero) sb2=sz2*((az2+1)/az2)**2 *d thermr.198 if (sz2.eq.zero) go to 110 *d thermr.238 if (sz2.eq.zero) go to 115 *d thermr.328 e=0 *d thermr.336 if (e.gt.emax*(1+small)) ex(2)=0 *d thermr.342 if (abs(e-emax).lt.small*emax) enext=up*emax *d thermr.379 if (sz2.gt.zero) teff2=eftmp2(itemp) *d thermr.446 130 if (abs(tnow-temp).lt.temp/1000+5) return *i thermr.554 zero=0 *d thermr.557 if (eftemp(i).ne.zero) go to 110 *d thermr.559 jmat=nint(tabl(1,j)) *d thermr.561 test=5 *d thermr.565 if (eftemp(i).eq.zero) eftemp(i)=temp(i) *i thermr.601 zero=0 *d thermr.604 if (eftmp2(i).ne.zero) go to 110 *d thermr.606 jmat=nint(tabl(1,j)) *d thermr.608 test=5 *d thermr.612 if (eftmp2(i).eq.zero) eftmp2(i)=temp(i) *d thermr.627 common/tunits/nendf,nin,nout,nscr,nscr2 *i thermr.638 data tolmin/1.d-6/ data eps/3.d-5/ *i thermr.641 data tolmin/1.e-6/ data eps/3.e-5/ *d thermr.654 ex(nex+i)=0 *d thermr.659 *d thermr.668 e=0 *d thermr.712 if (a(istk+nx*(i-2))-a(istk+nx*(i-1)).lt.eps*xm) go to 160 *d thermr.740 if (a(istk+nx*(i-1)).gt.emax*(1+small)) ej(nex)=0 *d thermr.743 if (a(istk+nx*(i-1)).gt.emax*(1+small)) ej(nex+il)=0 *i thermr.844 data eps/0.05d0/ data zero/0.d0/ *i thermr.863 data eps/0.05e0/ data zero/0.e0/ *d thermr.868 if (e.gt.zero) go to 210 *d thermr.873 *d thermr.911 i1m=int(a*sqrt(phi)) *d thermr.916 i2m=int(half*(l1+sqrt(3*(a*a*phi-l1*l1)))) *d thermr.922 if (x.gt.zero) i3m=int(c*sqrt(x)) *d thermr.926 w1=2 *d thermr.928 w2=2 *d thermr.945 if (tsq.lt.b(ifl+2*i-2).or.tsq.ge.(1+eps)*b(ifl+2*i-2)) *d thermr.1003 200 k=int(b(ifl+5)) *d thermr.1005 blast=0 *d thermr.1022 s(il)=0 *d thermr.1029 if (e.gt.emax) f=0 *i thermr.1118 data up,dn/1.1d0,.9d0/ *i thermr.1125 data up,dn/1.1e0,.9e0/ *d thermr.1158 if (abs(temp-tt1).gt.temp/10) call error('iel', *d thermr.1164 if (temp.lt.dn*tt1.or.temp.gt.up*ttn) call error('iel', *d thermr.1213 a(iscr)=0 *d thermr.1237 if (iex.eq.ne) ej(nj)=0 *d thermr.1336 sum=0 *i thermr.1401 data zero/0.d0/ data tenth/0.1d0/ data up,dn/1.1d0,0.9d0/ data uumin/0.00001d0/ *i thermr.1418 data zero/0.e0/ data tenth/0.1e0/ data up,dn/1.1e0,0.9e0/ data uumin/0.00001e0/ *d thermr.1473,1474 sb2=0 az2=0 *d thermr.1476,1477 if (a(iscr+12).eq.zero) az2=a(iscr+14) if (az2.ne.zero) sb2=a(iscr+13)*((az2+1)/az2)**2 *d thermr.1480 if (a(iscr+18).eq.zero) *d thermr.1509 if (abs(t-temp).le.(t/50)) go to 130 *d thermr.1528 if (t.eq.t1) l=l+2*nint(a(loc+4))+1 *d thermr.1532 if (a(l).gt.zero) itrunc=1 *d thermr.1565 test=5 *d thermr.1587 if (abs(tempt-tt1).gt.tempt/10) call error('calcem', *d thermr.1596 if (tempt.lt.dn*tt1.or.tempt.gt.up*ttn) call error('calcem', *d thermr.1601 190 if (sb2.eq.zero) go to 300 *d thermr.1607 if (abs(tempt-tt1).gt.tempt/10) call error('calcem', *d thermr.1616 if (tempt.lt.dn*tt1.or.tempt.gt.up*ttn) call error('calcem', *d thermr.1630,1632 a(ibeta)=0 a(ibeta+1)=tenth a(ibeta+2)=25 *d thermr.1666 cliq=0 *d thermr.1678,1680 a(ixs+ie-1)=0 ubar(ie)=0 ep=0 *d thermr.1704 315 if (ep.gt.zero) go to 316 *d thermr.1726,1727 uu=0 uum=0 *d thermr.1737 test=2*tol*abs(uu)+uumin *d thermr.1745 uu=0 *d thermr.1768 if (ylast.ne.zero) jnz=j *d thermr.1769 ulast=0 *d thermr.1779 y(il,i)=0 *d thermr.1794 uu=0 *d thermr.1807 if (y(1,1).ne.zero) jnz=j *d thermr.1822 a(iscr)=0 *d thermr.1840 ex(ie)=0 *d thermr.1854 if (ie.eq.ne) xs=0 *d thermr.1864 if (ie.eq.ne) xs=0 *i thermr.1894 zero=0 *d thermr.1918 if (lasym.eq.1.and.bb.lt.zero) bbb=-b *d thermr.1950 sig=0 *d thermr.1953 if (sig.lt.sigmin) sig=0 *d thermr.1959 s=0 *d thermr.1963 if (sb2.eq.zero) go to 190 *d thermr.1966 s2=0 *d thermr.1969 190 if (sig.lt.sigmin) sig=0 *d thermr.1974 s=0 *d thermr.1978 if (sig.lt.sigmin) sig=0 *d thermr.1983 sig=0 *i thermr.2052 data zero/0.d0/ *i thermr.2058 data zero/0.e0/ *d thermr.2074 if (ep.ne.zero) seep=1/sqrt(e*ep) *d thermr.2078 sum=0 *d thermr.2083,2085 if (ep.eq.zero) x(2)=0 if (ep.ne.zero) x(2)=half*(e+ep-(s1bb-1)*az*tev)*seep if (abs(x(2)).gt.1-eps) x(2)=0 *d thermr.2117 s(il)=0 *d thermr.2127,2128 sum=0 gral=0 *d thermr.2130 s(il)=0 *d thermr.2139,2141 if (ep.eq.zero) x(2)=0 if (ep.ne.zero) x(2)=(e+ep-(s1bb-1)*az*tev)*seep*half if (abs(x(2)).gt.1-eps) x(2)=0 *d thermr.2189 if (disc.ge.zero) go to 178 *d thermr.2193,2194 178 if (f.gt.zero) xn=xl-(yl*rf)+sqrt(disc) if (f.lt.zero) xn=xl-(yl*rf)-sqrt(disc) *d thermr.2220,2221 sum=0 gral=0 *i thermr.2253 data etop/20.d6/ data zero/0.d0/ *i thermr.2255 data etop/20.e6/ data zero/0.e0/ *d thermr.2261 ex(i)=0 *d thermr.2318 a(i-1+idicn)=0 *d thermr.2426 a(iscr+1)=0 *d thermr.2440,2441 ex(1)=etop ex(ix)=0 *d thermr.2515 if (a(ixs+ie-1).eq.zero) 1 call error('tpend','cross section=0.',' ') *ident up52 */ groupr -- 14dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. */ in addition, any conversion of real numbers to integers is */ explicitly flagged for truncation or for nearest integer. *i groupr.271 zero=0 *d groupr.366 eps=tempin/10000 *d groupr.368 if (diff.gt.1+eps) go to 220 *d groupr.421 a(is+1)=0 *d groupr.445 if (ngg.eq.0) a(nw)=0 *d groupr.544 ee=0 *d groupr.548 test=-1 *d groupr.557 if (econst.eq.zero) go to 510 *d groupr.580 a(it)=0 *d groupr.614 if (econst.gt.zero) write(nsyso,'(/ *d groupr.647 a(iscr+1)=0 *i groupr.760 *if sw data sigzmx/1.d10/ data zero/0.d0/ *else data sigzmx/1.e10/ data zero/0.e0/ *endif *d groupr.798 title(1)=0 *d groupr.805 if (temp(1).eq.zero) write(nsyso,'( *d groupr.807 if (temp(1).ne.zero) write(nsyso,'( *d groupr.814 sigz(1)=sigzmx *d groupr.1390,1391 data u187a,u187b,u187c,e187d,e187e/-15.5d0,-14.125d0,-5.875d0, 1 2.6058d4,6.868d0/ *d groupr.1499,1500 data u187a,u187b,u187c,e187d,e187e/-15.5e0,-14.125e0,-5.875e0, 1 2.6058e4,6.868e0/ *d groupr.1601 eg(60)=e187e *d groupr.1839 egg(1)=0 *i groupr.2046 data small/1.d-10/ data zero/0.d0/ data onep5/1.5d0/ *i groupr.2106 data small/1.e-10/ data zero/0.e0/ data onep5/1.5e0/ *d groupr.2112 sigpot=0 *d groupr.2115,2119 alpha2=0 sam=0 beta=0 alpha3=0 gamma=0 *d groupr.2127,2128 if (alpha3.ne.zero.and.alpha2.eq.zero) alpha2=small if (alpha2.gt.zero) write(nsyso,'( *d groupr.2163 ac=1/(exp(-ec/tc)*ec**onep5) *i groupr.2261 data exmin/-89.d0/ data zero/0.d0/ *i groupr.2276 data exmin/-89.e0/ data zero/0.e0/ *d groupr.2308 if (wtf.eq.zero) go to 120 *d groupr.2351 cc=1 *d groupr.2363,2364 if (pow.gt.exmin) wtf=wtf+wt6h*exp(pow) *i groupr.2457 *if sw data up,dn/1.001d0,0.999d0/ data tenth/0.1d0/ data sigzmx/1.d10/ data etop/1.d10/ data zero/0.d0/ *else data up,dn/1.001e0,0.999d0/ data tenth/0.1e0/ data sigzmx/1.e10/ data etop/1.e10/ data zero/0.e0/ *endif *d groupr.2471 e=0 *d groupr.2473 if (enext.gt.up*egn(1)) *d groupr.2493,2494 if (alpha2.ne.zero) nalph=2 if (alpha3.ne.zero) nalph=3 *d groupr.2499,2500 if (elo.lt.tenth) elo=tenth *d groupr.2505 e=0 *d groupr.2517 e=0 *d groupr.2572 f4=1-f2*log(ejp/ej) *d groupr.2637 e=0 *d groupr.2644 e=0 *d groupr.2666 sigzj=sigzmx *d groupr.2675 e=0 *d groupr.2709,2710 if (en.gt.dn*etop) enext=etop *d groupr.2719 xtot=0 *d groupr.2721 if (en.lt.zero) go to 340 *d groupr.2726 fac=0 *d groupr.2737,2738 if (enext.gt.dn*etop) nett=-nett *d groupr.2781 zero=0 if (e.gt.zero) go to 200 *d groupr.2783,2788 el=0 wl=0 xl=0 en=0 wn=0 xn=0 *d groupr.2792 if (flag.lt.zero) go to 110 *d groupr.2801 if (flag.ge.zero) go to 130 *i groupr.2860 *if sw data ebig/1.d8/ *else data ebig/1.e8/ *endif *d groupr.2897 econst=0 *d groupr.2901,2903 if (mtd.eq.18.or.mtd.eq.19) econst=ebig if (mtd.eq.452.or.mtd.eq.455.or.mtd.eq.456) econst=ebig if (mtd.eq.102) econst=ebig *i groupr.2960 data zero/0.d0/ data smin/1.d-9/ data ttest/0.1d0/ *i groupr.2976 data zero/0.e0/ data smin/1.e-9/ data ttest/0.1e0/ *d groupr.2990 eqw=0 *d groupr.3027 114 if (sig(1,1).ne.zero) go to 125 *d groupr.3030,3031 ans(il,iz,1)=1 ans(il,iz,2)=0 *d groupr.3103,3106 if (elo.lt.ttest) go to 220 if (abs(sig(1,1)).gt.smin) go to 220 *d groupr.3109 ans(il,iz,2)=0 *i groupr.3152 data big/1.d10/ data zero/0.d0/ *i groupr.3155 data big/1.e10/ data zero/0.e0/ *d groupr.3228,3230 result(il)=0 if (ans(il,1,ig2).eq.zero) go to 230 if (ans(il,1,1).eq.zero) ans(il,1,1)=big *d groupr.3232,3233 if (result(il).ne.zero) then jjj=int(log10(result(1))) *d groupr.3262 if (ans(1,1,1).eq.zero) ans(1,1,1)=big *d groupr.3264,3265 if (abs(ans(1,1,j)).lt.smin/1000) ans(1,1,j)=0 if (ans(1,1,j).eq.zero) go to 260 *d groupr.3292,3293 result(iz)=0 if (ans(il,iz,ig2).eq.zero) go to 330 *d groupr.3296 if (ans(il,iz,locd).eq.zero) ans(il,iz,locd)=big *d groupr.3300 jjj=int(log10(result(1))) *d groupr.3304 iii=nint(fact*result(iz))+nint(ten**(ndig-11)) *d groupr.3416 zero=0 if (e.gt.zero) go to 105 *d groupr.3488 zero=0 if (e.gt.zero) go to 200 *d groupr.3584 zero=0 if (e.gt.zero) go to 300 *d groupr.3629 thresh=0 *d groupr.3651 if (en.lt.zero) return *d groupr.3670 sig(iz,il)=1 *d groupr.3750 if (abs(sigz(i)-a(iunr+5+j)).lt.sigz(i)/100) go to 170 *i groupr.3781 zero=0 *d groupr.3820 if (a(locl).lt.zero) iovl=1 *d groupr.3827 if (lssf.gt.0.and.sinf.ne.zero) sig(is)=sig(is)*s/sinf *d groupr.3835 if (a(loc).lt.0..and.a(loc1).lt.zero) iovl=1 *d groupr.3846 if (lssf.gt.0.and.sinf.ne.zero) sig(is)=sig(is)*s/sinf *i groupr.3860 *if sw data rlo,rhi/0.99d0,10.1d0/ *else data rlo,rhi/0.99e0,10.1e0/ *endif *d groupr.3862 sig=0 *d groupr.3866,3868 if (sigz.lt.sigs(ns).and.sm.gt.rhi*sigs(ns)) go to 110 if (sigz.lt.sigs(2).and.sm.gt.rhi*sigs(2)) go to 110 if (sigz.ge.sigs(2).and.sm.lt.rlo*sigs(2)) go to 110 *d groupr.3871 if (abs(sm-st).lt.sm/1000) go to 120 *d groupr.3873 if (st.gt.rhi*sigs(ns)) go to 120 *d groupr.3877 if (st.gt.rhi*sigs(2)) go to 120 *d groupr.3880 140 if (st.lt.rlo*sigs(2)) go to 120 *i groupr.3914 zero=0 *d groupr.3968 yld=1 *d groupr.3975 if (e.eq.zero) return *d groupr.3988,3989 if (e.eq.zero) ifirst=0 if (e.gt.zero) ifirst=ifirst+1 *d groupr.4012,4013 if (e.eq.zero) call reserv('sed',nwds,ised,a) if (e.ne.zero) call findex('sed',ised,a) *d groupr.4018,4019 if (e.gt.zero) go to 205 if (econst.eq.zero) return *d groupr.4026 if (jconst.le.2) econst=0 *d groupr.4065 if (ff(ik,jg).ne.zero) ighi=ig *d groupr.4080 if (e.gt.zero) go to 325 *d groupr.4088 if (econ.ne.zero) go to 305 *d groupr.4099 315 if (econst.eq.zero) go to 322 *d groupr.4106 if (jconst.le.2) econst=0 *d groupr.4124 if (a(ieyl-1+nyl).ne.zero) go to 330 *d groupr.4131 if (e.gt.zero) go to 330 *d groupr.4141 ff(il,ig)=0 *d groupr.4144,4145 if (a(ik-1+igyl).eq.zero) go to 350 if (a(ik-1+ieyl).eq.zero) go to 340 *d groupr.4149 if (a(ieyl+ik-1).eq.zero) tempo=tempo*a(ised+ig-1) *d groupr.4157 yld=0 *d groupr.4186 yld=0 *d groupr.4206 yld=1 *d groupr.4208 if (e.eq.zero) return *d groupr.4213,4214 if (abs(ff(il,ig)).le.small) go to 720 *i groupr.4397 data up,dn/1.00001d0,0.99999d0/ data eps/0.02d0/ data zero/0.d0/ *i groupr.4401 data up,dn/1.00001e0,0.99999e0/ data eps/0.02e0/ data zero/0.e0/ *d groupr.4405 if (ed.gt.zero) go to 200 *d groupr.4418 if (mfd.eq.18) zad=0 *d groupr.4552,4553 c(l)=0 c(l+1)=0 *d groupr.4670 199 enext=-1 *d groupr.4677 ans(l,j)=0 *d groupr.4681 yld=0 *d groupr.4696 if (pe.eq.zero) go to 450 *d groupr.4730 300 ep=0 *d groupr.4763,4764 if (el.lt.zero) el=ed*awr/(awr+1)-el if (eh.lt.zero) eh=ed*awr/(awr+1)-eh *d groupr.4769 if (j.eq.1) e1=0 *d groupr.4796,4798 700 if (ed.le.zero) return if (yld.eq.zero) return test=0 *d groupr.4804 if (q.gt.zero) go to 710 *d groupr.4807 ans(l,i)=0 *d groupr.4810,4811 if (test.eq.zero) return *d groupr.4813 1 .and.(ed.lt.up*elo.or.ed.gt.dn*ehi)) *d groupr.4851 c(jnow)=0 *d groupr.4856,4857 sum=0 x(2)=0 *d groupr.4966 zero=0 if (ep.gt.zero) go to 200 *d groupr.4973 elmax=0 *d groupr.4975 epmax=0 *d groupr.5001 term(l)=0 *d groupr.5007 term(l)=0 *d groupr.5048 if (un.lt.umin+tiny/10) un=umin *d groupr.5141 if (abs(term(l)).lt.test) term(l)=0 *i groupr.5175 data tomev/1.d-6/ data half/0.5d0/ data bach1,bach2/0.0236d0,0.0039d0/ data bach3,bach4/92.d0,90.d0/ data zero/0.d0/ *i groupr.5181 data tomev/1.e-6/ data half/0.5e0/ data bach1,bach2/0.0236e0,0.0039e0/ data bach3,bach4/92.e0,90.e0/ data zero/0.e0/ *d groupr.5186 if (ep.gt.zero) go to 120 *d groupr.5194 if (cnow(lnow).le.zero) lnow=lnow+ncnow *d groupr.5200,5201 efirst=0 f6ddx=0 *d groupr.5223 t=0 *d groupr.5229 t=0 *d groupr.5240 if (t.lt.zero) t=0 *d groupr.5245,5246 s=0 r=0 *d groupr.5257 t=half *d groupr.5259,5262 aa=bach1+bach2*l*(l+1) ss=l*(l+1) bb=bach3-bach4/sqrt(ss) sa=1/(1+exp(aa*(bb-ep*tomev))) tt=(2*l+1)*sa*p(l+1)/2 *d groupr.5272 if (t.lt.zero) t=0 *d groupr.5282 tii=0 *d groupr.5289 tjj=0 *i groupr.5326 *if sw data tomev/1.d-6/ data half/0.5d0/ data bach1,bach2/0.0236d0,0.0039d0/ data bach3,bach4/92.d0,90.d0/ data zero/0.d0/ *else data tomev/1.e-6/ data half/0.5e0/ data bach1,bach2/0.0236e0,0.0039e0/ data bach3,bach4/92.e0,90.e0/ data zero/0.e0/ *endif *d groupr.5330 if (i.gt.0) go to 120 *d groupr.5342 f6dis=0 *d groupr.5352 t=0 *d groupr.5358 if (t.lt.zero) t=0 *d groupr.5369 t=half *d groupr.5371,5374 aa=bach1+bach2*l*(l+1) ss=l*(l+1) bb=bach3-bach4/sqrt(ss) sa=1/(1+exp(aa*(bb-ep*tomev))) tt=(2*l+1)*sa*p(l+1)/2 *d groupr.5384 if (t.lt.zero) t=0 *d groupr.5391 t=0 *d groupr.5452 za=int(iza/1000) *d groupr.5454 zc=za+int(iza1/1000) *d groupr.5456 zb=zc-int(iza2/1000) *d groupr.5487,5488 fa=1 if (iza1.eq.2004) fa=0 *i groupr.5528 zero=0 *d groupr.5532 c(jnow)=0 *d groupr.5541 epnext=0 *d groupr.5566 term(l)=0 *d groupr.5582 sum=0 *d groupr.5587 if (term(1).ne.zero) fact=sum/term(1) *d groupr.5599 test=tol*abs(tmid)+small/100 if (abs(tmid-c(lnow-nl-1+i)).gt.test) go to 550 *i groupr.5646 data emin/1.d-5/ *i groupr.5655 data emin/1.e-5/ *d groupr.5659 zero=0 if (ep.gt.zero) go to 200 *d groupr.5675 if (clo(llo).le.zero) llo=llo+nclo *d groupr.5682 if (chi(lhi).le.zero) lhi=lhi+nchi *d groupr.5696 term(l)=0 *d groupr.5699 150 epnext=emin *d groupr.5735 term1(l)=0 *d groupr.5741 255 term2(l)=0 *d groupr.5751 265 term1(l)=0 *d groupr.5757 266 term2(l)=0 *d groupr.5774 term1(l)=0 *d groupr.5778 t0=0 *d groupr.5799 term2(l)=0 *d groupr.5803 t0=0 *i groupr.5915 zero=0 *d groupr.5946 ast=0 *d groupr.5961 ff(il,ng)=0 *d groupr.5980 nqp=int(npo*2*b) *d groupr.5987 flt(il)=0 *d groupr.6009 prob=0 *d groupr.6016 i2s=nint(2*spi) *d groupr.6019,6020 zt=int(izat/1000) zi=int(izap/1000) *d groupr.6023 sigc=0 *d groupr.6050 464 ff(1,1)=1 *d groupr.6069 600 if (e.gt.zero) go to 610 *d groupr.6080 if (thresh.gt.zero) go to 630 *d groupr.6072,6074 ecl=0 ech=0 ecn=0 *d groupr.6077 ecl=emax *d groupr.6090 if (disc.lt.zero) go to 620 *d groupr.6096 ech=emax *d groupr.6099 if (thresh.gt.zero) go to 650 *d groupr.6102 test=test/10 *d groupr.6111 if (disc.lt.zero) go to 640 *d groupr.6156 zero=0 if (e.gt.zero) go to 200 *d groupr.6257 320 fle(i)=0 *d groupr.6324 zero=0 if (e.gt.zero) go to 200 *d groupr.6358 ehi=0 *d groupr.6373 clast=0 *d groupr.6382 if ((cnow-clast).le.(eps*clast+eps)) go to 160 *d groupr.6439 aed(il,ig)=0 *d groupr.6445 fl(il)=0 *d groupr.6447 eg=0 *d groupr.6451,6453 egp=0 egp1=0 egp2=0 *d groupr.6479 if (ek2.gt.up*eg2) go to 385 *d groupr.6515 aed(il,i)=aed(il,i)+(fi(il)+fl(il))*(ei-eg)/2 *d groupr.6531 aed(il,ig)=0 *d groupr.6546 aed(il,i)=0 *d groupr.6613 630 if (aed(1,ig).eq.zero) go to 650 *d groupr.6618 650 aed(1,ig)=1 *d groupr.6641 fl(il)=0 *i groupr.6684 data up/1.1d0/ *i groupr.6687 data up/1.1e0/ *d groupr.6690 zero=0 if (ed.gt.zero) go to 200 *d groupr.6779 205 gfl(1,j)=1 *d groupr.6782 gfl(i,j)=0 *d groupr.6803 if (im.eq.np.and.ed.gt.up*ehi) call error('getgfl', *d groupr.6949 220 fl(il)=0 *d groupr.6961 240 fl(il)=0 *d groupr.6968 fl(il)=0 *d groupr.7004 356 y=0 *d groupr.7013,7014 zt=int(izat/1000) zi=int(izap/1000) *d groupr.7100 zero=0 if (ed.gt.zero) go to 200 *d groupr.7162 e=0 *d groupr.7168 if (econst.eq.zero) go to 190 *d groupr.7191 if (a(l).eq.zero) go to 230 *i groupr.7239 *if sw data etop/20.d6/ data zero/0.d0/ *else data etop/20.e6/ data zero/0.e0/ *endif *d groupr.7287,7288 312 za2=-1 if (mth.eq.2.or.(mth.ge.50.and.mth.lt.91)) za2=1 *d groupr.7425 a(ie)=0 *d groupr.7427,7428 a(i+ir)=0 a(i+iaa)=0 *d groupr.7431 a(imax*(i-1)+i-1+ir)=1 *d groupr.7458,7459 a(iaa+(k-1)*imax+j-1)=0 a(ir+(k-1)*imax+j-1)=0 *d groupr.7474 ysum=0 *d groupr.7482 if (y.eq.zero) go to 245 *d groupr.7490 a(ie)=0 *d groupr.7525,7526 a(iscr)=0 a(iscr+1)=0 *d groupr.7534 a(iscr+10)=etop *d groupr.7542 a(iscr+2)=0 *d groupr.7550 a(iscr+10)=etop *d groupr.7672 e1=0 *d groupr.7688 a(iscr+5+2*im)=sig(1,1)*(1-w)*a(iscr+5+2*im) *i groupr.7897 data ebig/1.d8/ data tenth/0.1d0/ data zero/0.d0/ *i groupr.7903 data ebig/1.e8/ data tenth/0.1e0/ data zero/0.e0/ *d groupr.7908 if (ed.gt.zero) go to 200 *d groupr.8002 if (ig.eq.1) e1=0 *d groupr.8004 if (ig.eq.ng) e2=ebig *d groupr.8043 econ=0 *d groupr.8055 sed(ik,ig)=0 *d groupr.8071 if (pe.eq.zero) go to 290 *d groupr.8082 if (ig.eq.1) e1=0 *d groupr.8084 if (ig.eq.ng) e2=ebig *d groupr.8121 sigup=0 *d groupr.8129 sed(ikt,ig)=0 *d groupr.8131 280 test=tenth *d groupr.8155 sed(1,ig)=0 *i groupr.8191 data etest/1.1d-5/ *i groupr.8199 data etest/1.1e-5/ *d groupr.8206,8207 if (ep1.lt.etest) new=1 *d groupr.8212 g=0 *d groupr.8219 zero=0 if (epl.gt.de*(1+small).or.de.le.zero) return *d groupr.8273 if (top.lt.small) top=0 *d groupr.8286 if (top.lt.small) top=0 *d groupr.8325 200 if (top.lt.small) top=0 *d groupr.8468,8469 cm(1)=1 xk=1 *i groupr.8416 data epsrel,epsabs/1.d-9,1.d-30/ data zero/0.d0/ *i groupr.8423 data epsrel,epsabs/1.e-9,1.e-30/ data zero/0.e0/ *d groupr.8432 if (aa.lt.zero) fa(1)=-fa(1) *d groupr.8434 if (bb.lt.zero) fb(1)=-fb(1) *d groupr.8450 200 sgn=1 *d groupr.8493 beta=0 *d groupr.8509 test=epsabs+epsrel*abs(s) *d groupr.8554 zero=0 if (ep.gt.zero) go to 200 *d groupr.8563 if (eimax.le.zero) eimax=100 *d groupr.8572 f6psp=0 *d groupr.8578 200 s=0 *ident up53 */ leapr -- 2dec98 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. */ in addition, any conversion of real numbers to integers is */ explicitly flagged for truncation or for nearest integer. */ an error in mixed-mode conversion was fixed that affects */ the statistical factor for cold hydrogen or deuterium. *d leapr.229,231 b7=0 aws=0 sps=0 *d leapr.318 call trans(itemp,a(issm),nalpha,nbeta,ntempr,a(ibeta), *d leapr.373 emax=5 *i leapr.420 data tiny/1.d-30/ *i leapr.423 data tiny/1.e-30/ *d leapr.431 sc=1 *d leapr.447,448 if (add.lt.tiny) add=0 *d leapr.473,474 if (add.lt.tiny) add=0 *d leapr.478 if (add.gt.ssm(k,j,itemp)/1000.and.j.lt.maxt(k)) maxt(k)=j *i leapr.511 bel=0 ff1l=0 ff2l=0 sum0=0 sum1=0 *d leapr.542.543 sum0=0 sum1=0 *d leapr.594 tau=1 tau=tau/2 *d leapr.658.662 v=1 an=1-2*mod(n,2) be=0 fs=0 w=1 *d leapr.685 i=int(be/delta) *d leapr.692 120 terpt=0 *i leapr.707 *if sw data tiny/1.d-30/ *else data tiny/1.e-30/ *endif *d leapr.709 ckk=0 *d leapr.711 tnext(k)=0 *d leapr.715 f1=0 *d leapr.719 f2=0 *d leapr.731,732 if (tnext(k).lt.tiny) tnext(k)=0 *i leapr.767 data tiny/1.d-30/ *i leapr.770 data tiny/1.e-30/ *d leapr.777 sc=1 *d leapr.820 s=0 *d leapr.836,837 if (s.lt.tiny) s=0 *i leapr.851 sum0=0 sum1=0 ff1l=0 ff2l=0 bel=0 *d leapr.867,868 sum0=0 sum1=0 *d leapr.929 be=0 *d leapr.935 c5=0 *d leapr.957 be=0 *d leapr.979,980 check0=0 check1=0 *i leapr.1006 *if sw data slim/-225.d0/ *else data slim/-225.e0/ *endif *d leapr.1008 i=int(be/delta) *d leapr.1014 st=slim *d leapr.1019 stp=slim *d leapr.1024,1026 terps=0 if (stt.gt.slim) terps=exp(stt) *d leapr.1028 120 terps=0 *i leapr.1042 *if sw data shade/1.00001d0/ data slim/-225.d0/ *else data shade/1.00001e0/ data slim/-225.e0/ *endif *d leapr.1053 if (j.eq.nbeta.and.b.lt.shade*beta(j)) go to 130 *d leapr.1065 st=slim *d leapr.1070 stm=slim *d leapr.1077,1079 sb(i)=0 if (arg.gt.slim) sb(i)=exp(arg) *d leapr.1081 135 sb(i)=0 *d leapr.1127 test=1 *i leapr.1170 data small/1.d-5/ data vsmall/1.d-10/ data tiny/1.d-20/ *i leapr.1172 data small/1.e-5/ data vsmall/1.e-10/ data tiny/1.e-20/ *d leapr.1180 sc=1 *d leapr.1184 dwt=0 *d leapr.1189 tsave=0 *d leapr.1228 sexpb(j)=0 *d leapr.1232,1233 ben(1)=0 wtn(1)=1 *d leapr.1248,1249 if (besn.gt.0..and.wtsn.lt.small) go to 160 *d leapr.1263,1264 if (wtsn.lt.small) go to 170 *d leapr.1279,1280 if (wtsn.lt.small) go to 185 *d leapr.1317,1318 if (wts(i).lt.small) go to 207 *d leapr.1326,1327 sumn=0 sumr=0 *d leapr.1344,1345 if (add.lt.tiny) go to 225 *d leapr.1354,1355 if (dwf.lt.vsmall) go to 265 *d leapr.1359 db=1000 *d leapr.1370,1371 if (add.lt.tiny) go to 260 *i leapr.1395 sum0=0 sum1=0 ff1l=0 ff2l=0 bel=0 *d leapr.1411,1412 sum0=0 sum1=0 *d leapr.1414 if (twt.eq.0.) sum0=sum0/(1-exp(-al*dwpix(itemp))) *i leapr.1457 data big/1.d10/ data tiny/1.d-30/ *i leapr.1469 data big/1.e10/ data tiny/1.e-30/ *d leapr.1495,1496 bn(imax)=0 bn(imax-1)=1 *d leapr.1500 if (bn(i).lt.big) go to 180 *d leapr.1502 bn(j)=bn(j)/big *d leapr.1516 bplus(i)=0 *d leapr.1518,1520 if (bplus(i).lt.tiny) bplus(i)=0 bminus(i)=0 *d leapr.1523,1524 if (bminus(i).lt.tiny) bminus(i)=0 *d leapr.1530 bplus(i)=0 *d leapr.1533,1535 if (bplus(i).lt.tiny) bplus(i)=0 bminus(i)=0 *d leapr.1538,1539 if (bminus(i).lt.tiny) bminus(i)=0 *i leapr.1553 *if sw data small/1.e-9/ *else data small/1.e-9/ *endif *d leapr.1560,1562 if (beta(1).gt.small) go to 110 bex(nbeta)=0 *i leapr.1589 *if sw data small/1.e-9/ *else data small/1.e-9/ *endif *d leapr.1596,1597 if (beta(1).gt.small) go to 110 *i leapr.1621 data slim/-225.d0/ *i leapr.1623 data slim/-225.e0/ *d leapr.1629 sv=0 *d leapr.1654 ss1=slim *d leapr.1659 ss3=slim *d leapr.1664,1666 sv=0 if (ex.gt.slim) sv=exp(ex) *d leapr.1688 *i leapr.1700 data cct/0.01381d0/ data deh/0.0147d0/ data ded/0.0074d0/ data amassh/3.465d0/ data amassd/6.69d0/ data disth/1.5d0/ data distd/0.71d0/ data sampch/0.356d0/ data sampcd/0.668d0/ data sampih/2.526d0/ data sampid/0.403d0/ data small/1.d-6/ *i leapr.1707 data cct/0.01381e0/ data deh/0.0147e0/ data ded/0.0074e0/ data amassh/3.465e0/ data amassd/6.69e0/ data disth/1.5e0/ data distd/0.71e0/ data sampch/0.356e0/ data sampcd/0.668e0/ data sampih/2.526e0/ data sampid/0.403e0/ data small/1.e-6/ *d leapr.1721,1722 t=cct*abs(temp) sc=1 *d leapr.1725,1726 de=deh if (law.eq.3) de=ded *d leapr.1728,1735 amassm=amassh if (law.gt.3) amassm=amassd bp=hbar/2*sqrt(2/disth/ev/pmass) if (law.gt.3) bp=hbar/2*sqrt(2/distd/ev/dmass) sampc=sampch if (law.gt.3) sampc=sampcd sampi=sampih if (law.gt.3) sampi=sampid *d leapr.1745 waven=sqrt(amassm*t*al)/hbar *d leapr.1750 if (nokap.eq.1) sk=1 *d leapr.1757,1764 120 if (law.eq.2) swe=sampi**2/3 if (law.eq.2) swo=sk*sampc**2+2*sampi**2/3 if (law.eq.3) swe=sk*sampc**2 if (law.eq.3) swo=sampi**2 if (law.eq.4) swe=sk*sampc**2+5*sampi**2/8 if (law.eq.4) swo=3*sampi**2/8 if (law.eq.5) swe=3*sampi**2/4 if (law.eq.5) swo=sk*sampc**2+sampi**2/4 *d leapr.1789,1790 sn=0 total=0 *d leapr.1804 snlg=0 *d leapr.1807,1810 betap=(-j*(j+1)+jp*(jp+1))*x/2 tmp=(2*jp+1)*pj*swe*4*sumh(j,jp,y) if (jj.eq.1.and.tmp.ge.small) then *d leapr.1826 snlk=0 *d leapr.1829 betap=(-j*(j+1)+jp*(jp+1))*x/2 *d leapr.1831,1832 if (jj.eq.1.and.tmp.ge.small) then *d leapr.1867 up=0 *i leapr.1874 sum0=0 bel=0 ff1l=0 ff2l=0 *d leapr.1888 sum0=0 *d leapr.1905 *if sw data half/0.5d0/ *else data half/0.5e0/ *endif c yy=half*j*(j+1) *d leapr.1907 b=0 *d leapr.1911 yy=half*k*(k+1) *d leapr.1934 120 sum1=0 *d leapr.1971,1972 s=0 fact=1 *d leapr.1974 zi=i *d leapr.1979,1980 s=0 fact=1 *d leapr.1982 zi=i *d leapr.1987,1988 s=0 fact=1 *d leapr.1990 zi=i *d leapr.1995,1996 s=0 fact=1 *d leapr.1998 zi=i *d leapr.2003,2004 s=0 b1=1 *d leapr.2006 zi=i *d leapr.2010,2011 s=0 b2=1 *d leapr.2013 zi=i *d leapr.2017,2018 s=0 b3=1 *d leapr.2020 zi=i *d leapr.2024,2025 s=0 b4=1 *d leapr.2027 zi=i *d leapr.2031,2033 rat=2*nn+1 rat=rat/(jj+ll+nn+1) iwign=(jj+ll-nn)/2 wign=(-1)**iwign wign=wign*sqrt(rat)*b1/a1*a2/b2*a3/b3*a4/b4 *d leapr.2035 190 wign=0 *i leapr.2048 *if sw data huge/1.d25/ data small/2.d-38/ data break1/3.d4/ data break2/7.d-4/ data break3/0.2d0/ data one/1.d0/ data ten/1.d1/ data hund/1.d2/ *else data huge/1.e25/ data small/2.e-38/ data break1/3.e4/ data break2/7.e-4/ data break3/0.2e0/ data one/1.e0/ data ten/1.e1/ data hund/1.e2/ *endif *d leapr.2051,2052 if (x.gt.break1) go to 240 *d leapr.2054,2056 if (x.gt.break2) go to 120 w=1 *d leapr.2060,2061 t1=3 t2=1 *d leapr.2064 t1=t1+2 *d leapr.2069,2070 120 if (x.ge.break3) go to 130 *d leapr.2072 w=1-y*(1-y/20)/6 *d leapr.2079,2080 160 if (x.lt.hund) go to 170 l=int(x/50+18) *d leapr.2082,2083 170 if (x.lt.ten) go to 180 l=int(x/10+10) *d leapr.2085,2086 180 if (x.le.one) go to 190 l=int(x/2+5) *d leapr.2093,2095 z=1/x t3=0 t2=small *d leapr.2098 t1=(2*k+3)*z*t2-t3 *d leapr.2101,2105 210 if (abs(t1).lt.huge) go to 220 t1=t1/huge t2=t2/huge sj=sj/huge *d leapr.2118 260 bessel=0 *d leapr.2133 i=int(be/delta) *d leapr.2140 120 terpk=1 *i leapr.2184 data toler/1.d-6/ data eps/.05d0/ *i leapr.2194 data toler/1.e-6/ data eps/.05/ *d leapr.2204 twopis=(2*pi)**2 *d leapr.2207,2209 tsqx=econ/20 *d leapr.2233 wint=0 *d leapr.2241 i1m=int(a*sqrt(phi)) *d leapr.2246 i2m=int((l1+sqrt(3*(a*a*phi-l1*l1)))/2) *d leapr.2252 if (x.gt.0.) i3m=int(c*sqrt(x)) *d leapr.2260 if (l1.eq.0.and.l2.eq.0) w2=1 if (l1.eq.0.and.l2.eq.0) w2=w2/2 *d leapr.2275 if (tsq.lt.b(ifl+2*i-2).or.tsq.ge.(1+eps)*b(ifl+2*i-2)) *d leapr.2296 if (tsq.lt.b(ifl+2*i-2).or.tsq.ge.(1+eps)*b(ifl+2*i-2)) *d leapr.2328 bel=-1 *d leapr.2333,2334 if (be-bel.ge.toler) go to 504 *i leapr.2360 data c1,c2,c3/7.54d0,4.24d0,11.31d0/ *i leapr.2362 data c1,c2,c3/7.54,4.24,11.31/ *d leapr.2379,2380 200 formf=(1+cos(2*pi*(2*l1+4*l2+3*l3)/6))*(c1+c2+c3*cos(3*pi*l3/4)) *i leapr.2409 *if sw data small/1.d-9/ data tiny/-999.d0/ data smin/2.d-38/ data tol/0.9d-7/ data up/1.01d0/ data therm/.0253d0/ *else data small/1.e-9/ data tiny/-999.e0/ data smin/2.e-38/ data tol/0.9e-7/ data up/1.01e0/ data therm/.0253e0/ *endif *d leapr.2425,2426 *d leapr.2433,2434 sb=spr*((1+awr)/awr)**2 if (aws.ne.0.) sbs=sps*((1+aws)/aws)**2 *d leapr.2493,2497 scr(1)=0 scr(2)=0 scr(3)=0 scr(4)=0 scr(5)=0 *d leapr.2500,2503 scr(1)=1 scr(2)=0 scr(3)=0 scr(4)=0 *d leapr.2507,2508 scr(1)=0 scr(2)=0 *d leapr.2532,2533 160 scr(1)=0 scr(2)=0 *d leapr.2540,2541 scr(ii+1)=0 scr(ii+2)=0 *d leapr.2551,2552 165 scr(ii+1)=0 scr(ii+2)=0 *d leapr.2578 scr(2)=0 *d leapr.2611,2612 sum=0 suml=0 *d leapr.2618 if (sum-suml.le.tol*sum) go to 190 *d leapr.2627 scr(2)=0 *d leapr.2637 sum=0 *d leapr.2639 if (nss.gt.0.and.b7.eq.0.) w=(dbw(i)+dbw1(i))/2 *d leapr.2658 scr(2)=0 *d leapr.2666 sum=0 *d leapr.2668 if (nss.gt.0.and.b7.eq.0.) w=(dbw(i)+dbw1(i))/2 *d leapr.2697,2698 scr(1)=0 scr(2)=0 *d leapr.2708 scr(11)=0 *d leapr.2714,2715 scr(16)=0 scr(17)=0 *d leapr.2718,2719 scr(1)=0 scr(2)=0 *d leapr.2732 sc=1 *d leapr.2751,2752 if (scr(8+2*j).ge.small) then *d leapr.2758 scr(8+2*j)=tiny *d leapr.2768,2769 if (scr(8+2*j).ge.small) then *d leapr.2775 scr(8+2*j)=tiny *d leapr.2784,2785 if (scr(8+2*j).ge.small) then *d leapr.2791 scr(8+2*j)=tiny *d leapr.2801,2802 if (scr(8+2*j).ge.small) then *d leapr.2808 scr(8+2*j)=tiny *d leapr.2818,2819 if (scr(8+2*j).ge.small) then *d leapr.2825 scr(8+2*j)=tiny *d leapr.2834,2835 if (scr(8+2*j).ge.small) then *d leapr.2841 scr(8+2*j)=tiny *d leapr.2849 if (ilog.eq.0.and.scr(8+2*j).lt.smin) scr(8+2*j)=0 *d leapr.2866,2867 if (scr(6+j).ge.small) then *d leapr.2873 scr(6+j)=0 *d leapr.2883,2884 if (scr(6+j).ge.small) then *d leapr.2890 scr(6+j)=tiny *d leapr.2899,2900 if (scr(6+j).ge.small) then *d leapr.2906 scr(6+j)=tiny *d leapr.2916,2917 if (scr(6+j).ge.small) then *d leapr.2923 scr(6+j)=tiny *d leapr.2933,2934 if (scr(6+j).ge.small) then *d leapr.2940 scr(6+j)=tiny *d leapr.2949,2950 if (scr(6+j).ge.small) then *d leapr.2956 scr(6+j)=tiny *d leapr.2964 if (ilog.eq.0.and.scr(6+j).lt.smin) scr(6+j)=0 *d leapr.2970,2971 scr(1)=0 scr(2)=0 *d leapr.2984 290 scr(2*i+7)=sigfig(up*scr(2*i+5),7,0) *d leapr.2988,2989 300 scr(1)=0 scr(2)=0 *d leapr.3002 305 scr(2*i+7)=sigfig(up*scr(2*i+5),7,0) *ident up54 */ gaspr -- 2dec98 */ increase the maximum energy grid size. *d gaspr.28 dimension egas(50000),sgas(5,50000) *d gaspr.33 *if sw data emax/1.d10/ *else data emax/1.e10/ *endif maxg=50000 */ fix various issues raised by ftnchek. these changes make */ sure that there are no 4-byte constants in statements that */ might lead to a loss in precision. we have even eliminated */ occurences of "0." for completeness, even though this value */ probably doesn't lead to problems with known compilers. */ in addition, any conversion of real numbers to integers is */ explicitly flagged for truncation or for nearest integer. */ also, add more fortran block structures, including do blocks. *d gaspr.74 zain=1 *d gaspr.78 zain=int(nsub/10) *d gaspr.81,84 do while (nb.ne.0) call moreio(nendf,0,0,a(1),nb,nw) enddo if (iverf.ne.4) nx=n2h *d gaspr.113,117 do while (nb.ne.0) call moreio(nendf,0,0,a(ll),nb,nw) ll=ll+nw enddo if (izap.eq.1001) then *d gaspr.119,121 do i=1,nw six(lsix+i-1)=a(i) enddo *d gaspr.125.127 do i=1,nw six(lsix+i-1)=a(i) enddo *d gaspr.131,133 do i=1,nw six(lsix+i-1)=a(i) enddo *d gaspr.137,139 do i=1,nw six(lsix+i-1)=a(i) enddo *d gaspr.143,145 do i=1,nw six(lsix+i-1)=a(i) enddo *d gaspr.157,162 do ie=1,ne call listio(nendf,0,0,a,nb,nw) do while (nb.ne.0) call moreio(nendf,0,0,a,nb,nw) enddo enddo *d gaspr.166.175 do ie=1,ne call tab2io(nendf,0,0,a,nb,nw) nmu=n2h do imu=1,nmu call tab1io(nendf,0,0,a,nb,nw) do while (nb.ne.0) call moreio(nendf,0,0,a,nb,nw) enddo enddo enddo *d gaspr.205 thrg=emax *d gaspr.210 en=0 *d gaspr.302 en=0 *d gaspr.312 if (enext.lt.emax) go to 260 *d gaspr.321 sgas(j,i)=0 *d gaspr.338 en=0 *d gaspr.342,512 y203=0 y204=0 y205=0 y206=0 y207=0 if (mth.eq.5.and.mf6mt5.eq.1) then if (l203.gt.0) y203=111 if (l204.gt.0) y204=111 if (l205.gt.0) y205=111 if (l206.gt.0) y206=111 if (l207.gt.0) y207=111 else if (mth.eq.11) then izr=izr-1004 y204=1 else if (mth.eq.16) then izr=izr-2 else if (mth.eq.17) then izr=izr-3 else if (mth.eq.22) then izr=izr-2005 y207=1 else if (mth.eq.23) then izr=izr-6013 y207=3 else if (mth.eq.24) then izr=izr-2006 y207=1 else if (mth.eq.25) then izr=izr-2007 y207=1 else if (mth.eq.28) then izr=izr-1002 y203=1 else if (mth.eq.29) then izr=izr-4009 y207=2 else if (mth.eq.30) then izr=izr-4010 y207=2 else if (mth.eq.32) then izr=izr-1003 y204=1 else if (mth.eq.33) then izr=izr-1004 y205=1 else if (mth.eq.34) then izr=izr-2004 y206=1 else if (mth.eq.35) then izr=izr-5011 y204=1 y207=2 else if (mth.eq.36) then izr=izr-5012 y205=1 y207=2 else if (mth.eq.37) then izr=izr-4 else if (mth.eq.41) then izr=izr-1003 y203=1 else if (mth.eq.42) then izr=izr-1004 y203=1 else if (mth.eq.44) then izr=izr-2003 y203=2 else if (mth.eq.45) then izr=izr-3006 y203=1 y207=1 else if (mth.ge.51.and.mth.le.91) then izr=izr-1 if (lr.eq.22) then izr=izr-2004 y207=1 else if (lr.eq.23) then izr=izr-6012 y207=3 else if (lr.eq.24) then izr=izr-2005 y207=1 else if (lr.eq.25) then izr=izr-2006 y207=1 else if (lr.eq.28) then izr=izr-1001 y203=1 else if (lr.eq.29) then izr=izr-4008 y207=2 else if (lr.eq.30) then izr=izr-4009 y207=2 else if (lr.eq.32) then izr=izr-1002 y204=1 else if (lr.eq.33) then izr=izr-1003 y205=1 else if (lr.eq.34) then izr=izr-2003 y206=1 else if (lr.eq.35) then izr=izr-5010 y204=1 y207=2 else if (lr.eq.36) then izr=izr-5011 y205=1 y207=2 else if (lr.eq.39) then izr=izr else if (lr.eq.40) then izr=izr endif else if (mth.eq.103) then izr=izr-1001 y203=1 else if (mth.eq.104) then izr=izr-1002 y204=1 else if (mth.eq.105) then izr=izr-1003 y205=1 else if (mth.eq.106) then izr=izr-2003 y206=1 else if (mth.eq.107) then izr=izr-2004 y207=1 else if (mth.eq.108) then izr=izr-4008 y207=2 else if (mth.eq.109) then izr=izr-6012 y207=3 else if (mth.eq.111) then izr=izr-2002 y203=2 else if (mth.eq.112) then izr=izr-3005 y203=1 y207=1 else if (mth.eq.113) then izr=izr-5011 y205=1 y207=2 else if (mth.eq.114) then izr=izr-5010 y204=1 y207=2 else if (mth.eq.115) then izr=izr-2003 y203=1 y204=1 else if (mth.eq.116) then izr=izr-2004 y203=1 y205=1 else if (mth.eq.117) then izr=izr-3007 y204=1 y207=1 endif if (izr.eq.1001) y203=y203+1 if (izr.eq.1002) y204=y204+1 if (izr.eq.1003) y205=y205+1 if (izr.eq.2003) y206=y206+1 if (izr.eq.2004) y207=y207+1 if (izr.eq.4008) y207=y207+2 *d gaspr.596,609 if (iverf.ge.5) then call contio(nscr1,0,0,a(k),nb,nw) k=k+nw endif if (iverf.ge.6) then call contio(nscr1,0,0,a(k),nb,nw) k=k+nw endif call hdatio(nscr1,0,0,a(k),nb,nw) k=k+nw if (iverf.ne.4) nx=n2h do while (nb.ne.0) call moreio(nscr1,0,0,a(k),nb,nw) k=k+nw enddo nw=nx *d gaspr.631,671 401 if (n203.ne.0) then a(j+1)=0 a(j+2)=0 a(j+3)=3 a(j+4)=203 a(j+5)=int((n203+2)/3) a(j+6)=1 j=j+6 endif if (n204.ne.0) then a(j+1)=0 a(j+2)=0 a(j+3)=3 a(j+4)=204 a(j+5)=int((n204+2)/3) a(j+6)=1 j=j+6 endif if (n205.ne.0) then a(j+1)=0 a(j+2)=0 a(j+3)=3 a(j+4)=205 a(j+5)=int((n205+2)/3) a(j+6)=1 j=j+6 endif if (n206.eq.0) then a(j+1)=0 a(j+2)=0 a(j+3)=3 a(j+4)=206 a(j+5)=int((n206+2)/3) a(j+6)=1 j=j+6 endif if (n207.ne.0) then a(j+1)=0 a(j+2)=0 a(j+3)=3 a(j+4)=207 a(j+5)=int((n207+2)/3) a(j+6)=1 j=j+6 endif nw=6 *d gaspr.676,688 if (iverf.ge.5) then call contio(0,noutp,0,a(k),nb,nw) k=k+nw endif if (iverf.ge.6) then call contio(0,noutp,0,a(k),nb,nw) k=k+nw endif if (iverf.ne.4) a(k+5)=nx+nsec-nold call hdatio(0,noutp,0,a(k),nb,nw) k=k+nw do while (nb.ne.0) call moreio(0,noutp,0,a(k),nb,nw) k=k+nw enddo nw=nx+nsec-nold *d gaspr.721,724 a(3)=0 a(4)=0 a(5)=0 a(6)=0 *d gaspr.728,731 a(1)=0 a(2)=0 a(3)=0 a(4)=0 *ident up55 */ reconr -- 13jan99 */ fix a typographical error in up47 *d up47.12 145 est=estp*(a(ix+i-1)-res(1)) */ reconr -- 13jan99 */ make trange incommensurable with normal grid energies */ to help in getting the same grid on different machines *d up10.9 trange=.4999d0 *d up10.11 trange=.4999 *d up10.17 data trange/.4999d0/ *d up10.19 data trange/.4999/ */ reconr -- 14jan99 */ fix errors made in up47 while eliminating real*4 */ constants from the code. from mattes (ike-stuttgart). *d up47.294,295 if (gfa.ne.zero) a2=sqrt(abs(gfa)) if (gfa.lt.zero) a2=-a2 */ reconr -- 15jan99 */ fix an incorrect comment line, and replace one */ last real*4 constant in an if (trkov, ijs slovenia). *d up47.408 one=1 if (thresh.gt.one.and.abs(thresh-eg).lt.test) sn=0 *d reconr.1494 c for energies less than elim, add points by dividing *d reconr.4191 *ident up56 */ broadr -- 13jan99 */ make trange incommensurable with normal grid energies */ to help in getting the same grid on different machines *d up11.46 data trange/.4999d0/ *d up11.48 data trange/.4999/ *ident up57 */ acer -- 13jan99 */ fix an error in up38. this update introduced a conversion */ to law=7 format for file 6 cm or lab tabulated or legendre */ data (needed for jef iron). unfortunately, the conversion */ is running for endf kalbach data also! *d up38.12,13 lang=l1h c set flag to convert file-6 law=1 data representation c to law-7 angle-energy format except for kalbach data. if (lf.eq.1.and.lang.ne.2) new6=1 *ident up58 */ heatr -- 14jan99 */ add initialization for zero. */ found by both deleege (dlft) and muir (iaea). *i heatr.3737 zero=0 *ident up59 */ groupr -- 14jan99 */ fix an error introduced by up52. the 64-bit versions of the */ constants bgam2 and tgam2 were clobbered by some of the */ other revised constants. this does not effect "sw" versions. */ noted by deleege (delft). *i groupr.1498 data bgam2,tgam2/27.63e0,-0.53063e0/ *ident up60 */ wimsr -- 15jan99 */ add additional space for fission products and fix the */ consistency between the container array size and dimensions. */ requested by trkov (ijs, slovenia). *d wimsr.101 common/wstore/a(400000) *d wimsr.110 common/burnup/yield(100),ifisp(100) *i wimsr.118 nymax=100 *i wimsr.202 if (ntis.gt.nymax) call error('wimsr', 1 'too many time-dependent isotopes',' ') *d wimsr.1832 common/wstore/a(400000) *d wimsr.1842 common/burnup/yield(100),ifisp(100) *ident up61 */ njoy -- 15jan99 */ make a change recommended by trkov (ijs slovenia) that */ prevents an infinite loop for some pc compilers. *d up46.59 if (y.gt.zero) xnext=shade*shade*a(jp) *ident up62 */ matxsr -- 15jan99 */ remove some obsolete coding. this shouldn't affect */ answers, but it stops a few annoying compiler messages. *d matxsr.1681 *d matxsr.1684 if (abs(a(nxsd)).lt.eps) go to 510 *d matxsr.1739,1757 *d matxsr.1851 *d matxsr.1859 call shufl(nscrr,irefr,ia(ijgll),ia(ijgll+noutg), *d matxsr.2007 subroutine shufl(nscr,iref,jg1lo,jg1hi,nout3,irec3,irite) *ident up63 */ njoy -- 27jan99 */ fix a typing error when adding double constants. */ reported by mattes, muir, panini. *d up46.148 *i njoy.2665 data zero/0.d0/ *ident up64 */ broadr -- 27jan99 */ fix a typing error when adding double constants. */ only affects 64-bit installations. from mattes. *d broadr.1200 data sigmin/1.e-15/ *ident up65 */ unresr -- 27jan99 */ fix a typing error when adding double constants. */ only affects 64-bit installations. from mattes. *d up35.15 7 1.25e6,1.5e6,2.0e6,3.2e6,4.0e6,5.0e6,6.3e6,8.0e6/ *ident up66 */ heatr -- 27jan99 */ fix a typing error when adding double constants. */ only affects 64-bit installations. from mattes. *d heatr.1598 data step/1.5e0/ *ident up67 */ groupr -- 27jan99 */ fix a typing error when adding double constants. */ only affects 64-bit installations. from mattes. *d up52.99 data up,dn/1.001e0,0.999e0/ *ident up68 */ purr -- 27jan99 */ fix a typing error when adding double constants. */ only affects 64-bit installations. from mattes. *d up36.15 7 1.25e6,1.5e6,2.0e6,3.2e6,4.0e6,5.0e6,6.3e6,8.0e6/ *ident up69 */ groupr -- 28jan99 */ add a missing save variable. this variable is an */ index to a scratch area, and the problem can be */ intermittent, depending on the machine. it only */ occurs for discrete levels (mt=51,52,...) when */ given in file 6. noticed by frankle (lanl). *d up25.5 save ltt3,lttn,ifls *ident up70 */ dtfr -- 28jan99 */ fix a problem with the graphics file for the */ neut-phot matrix from dtfr. it shows up in viewr. *d dtfr.1397 write(nplot,'(1x,a,''sec/e'',a,''/'')') qu,qu *ident up71 */ matxsr -- 28jan99 */ set uninitialized variables in mtxdat. *i matxsr.948 irec2=0 irec3=0 *ident up72 */ thermr -- 28jan99 */ fix a typing error in one of the data statements *d thermr.585 1 1095.d0,296.0d0,685.54d0,1095.d0,350.d0,712.02d0, *ident up73 */ reconr -- 1feb99 */ watch out for histogram gamma yields and add extra */ grid points above and below the mf12 energies. *d reconr.1264 if (mfl.eq.3.or.mfl.eq.10.or.mfl.eq.12.or.mfl.eq.13.or.mfl.eq.23) *i reconr.1269 if (mfh.eq.12) go to 140 *i reconr.1276 if (mfh.eq.12.and.nint(a(iscr+2)).ne.1) go to 150 if (mfh.eq.12) a(iscr+4)=1 *ident up74 */ heatr -- 1feb99 */ add two missing save variables. *i up24.5 save lct,awr *ident up75 */ leapr -- 5feb99 */ resolve some problems with constants for the cold hydrogen */ calculation. there is a data item, cct=.01381, which is */ actually bk*ev=.01380658, except for some factors of ten */ that were left out of the original constants because of */ cancellation. we take cct out of the equations. in */ addition, different values of the rotational level */ spacing were use to compute rotational energies */ (.0147 for h, .0074 for d) and the average bond length */ (1.5 and .71, with pesky factors of ten treated differently. */ We will change the average bond length calculation to */ agree with the level spacing values. *d up53.244 *d up53.249,250 *d leapr.1701 data ev/160.217733d0/ *d leapr.1708 data ev/160.217733/ *d up53.257 *d up53.262,263 *d up53.270 *d up53.278,279 bp=hbar/2*sqrt(2/deh/ev/pmass) if (law.gt.3) bp=hbar/2*sqrt(2/ded/ev/dmass) *d up53.285 waven=sqrt(amassm*tev*ev*al)/hbar *ident up76 */ heatr -- 5feb99 */ improve the treatment of the damage threshold. the atomic */ displacement energy "ed" appears in two places in a damage */ calculation: (1) in the damage function, displacement damage */ is zero for primary recoils with energy below "ed", and */ (2) in converting damage energy to dpa, you divide by twice */ the displacement energy "ed". the latter function is */ performed outside of njoy, but the former uses an internal */ default value of 25 ev. this update provides a table of */ default values for "ed" that varies from element to */ element, and it allows the user to enter a special value */ of the displacement energy, if desired. requested by */ ornl to improve the agreement between njoy and the current */ iron damage standard in the 1 keV to 1 MeV range. *d heatr.33.35 c * damage energy is computed using the lindhard electronic * c * screening damage function with a displacement threshold * c * from a table of default values for important elements * c * or a value provided by the user. * *i heatr.52 c * ed displacement energy for damage * c * (default from built-in table) * *i heatr.88 common/ed/break *i heatr.96 zero=0 *d heatr.126 break=0 read(nsysi,*) matd,npk,nqa,ntemp,local,iprint,break *i heatr.198 if (break.eq.zero) then write(nsyso,'( 1 '' damage displacement energy ........... default'')') else write(nsyso,'( 1 '' damage displacement energy ........... '',f10.1, 2 '' ev'')') break endif *i heatr.253 c c ***check on default damage displacement energy if (break.eq.zero) then iz=nint(za/1000) if (iz.eq.4) then break=31 else if (iz.eq.6) then break=31 else if (iz.eq.12) then break=25 else if (iz.eq.13) then break=27 else if (iz.eq.14) then break=25 else if (iz.eq.20) then break=40 else if (iz.ge.22.and.iz.le.29) then break=40 else if (iz.eq.40) then break=40 else if (iz.eq.41) then break=40 else if (iz.eq.42) then break=60 else if (iz.eq.47) then break=60 else if (iz.eq.73) then break=90 else if (iz.eq.74) then break=90 else if (iz.eq.79) then break=30 else if (iz.eq.82) then break=25 else break=25 endif write(nsyso,'(/'' default damage energy ='',f5.1, 1 '' ev'')') break endif *i heatr.1533 common/ed/break *d heatr.1544 *d heatr.1554 *ident up77 */ acer -- 13feb99 */ fix a problem with reading in file6,law7 as used for */ the endf 9be evaluation. *d acer.2076,2078 *d acer.2080 1116 if (nb.eq.0) go to 1118 *d acer.2083 *ident up78 */ wimsr -- 19feb99 */ fix a problem with reopening scratch files. the usage of scratch */ files varies on different machines, and this problem didn't */ show up for most systems tested. however, on dec alpha machines, */ test problem 11 is going to terminate with an end-of-file */ message from xseco. diagnosed by roberto orsi (enea, bologna). *d wimsr.1295,1298 *i wimsr.1317 c c ***set up scratch units for xsec output call openz(-nscr2,1) call openz(-nscr3,1) *ident up79 */ thermr -- 19feb99 */ the value of the neutron mass used for calculating bragg edges */ in sigcoh is currently 1.00866489 amu, but the endf standard is */ 1.0086652 amu. this change affects the bragg energies and the */ corresponding cross sections by only one digit in the last place */ for some of the values, and it has no practical consequences. */ making it now will improve consistency between njoy97 and */ future versions that use more consistent constants. similarly, */ we add more digits to pi. in some cases, it is necessary to hold */ 10 or 11 digits to get perfect agreement in the seventh place */ because of roundoff, e.g., x.xxxxxx4999 and x.xxxxxx5001 will */ differ by 1 in the seventh place when printed out. this happens! */ these changes only affect thermal coherent scatterers like */ graphite, be, and beo. *d thermr.837 data pi/3.14159265358979d0/ *d thermr.839 data amne/1.674929113d-24/ *d thermr.856 data pi/3.14159265358979e0/ *d thermr.858 data amne/1.674929113e-24/ *ident up80 */ heatr -- 24feb99 */ update the physical constants in heatr to agree with the */ standard constants to be used with the next version of njoy. */ this helps when making comparisons. emc2 uses the endf */ value for the neutron amu (1.0086652) and codata for the */ amu, c, and the ev conversion. *d up3.23 data emc2/939.56591826d6/ *d up3.25 data emc2/939.56591826e6/ *d up3.29 data emc2/939.56591826d6/ *d up3.31 data emc2/939.56591826e6/ *d up3.35 data emc2/939.56591826d6/ *d up3.37 data emc2/939.56591826e6/ *d up3.41 data emc2/939.56591826d6/ *d up3.43 data emc2/939.56591826e6/ *d up3.47 data emc2/939.56591826d6/ *d up3.49 data emc2/939.56591826e6/ *d up3.53 data emc2/939.56591826d6/ *d up3.55 data emc2/939.56591826e6/ *ident up81 */ reconr -- 24feb99 */ update the physical constants in reconr to agree with the */ standard constants to be used with the next version of njoy. */ this helps when making comparisons. cwaven uses the endf */ value for the neutron amu (1.0086652) and codata for the */ amu, hbar, and the ev conversion. *d reconr.503 data cwaven/2.1968074704d-3/ *d reconr.514 data cwaven/2.1968074704e-3/ *d reconr.1944 data forpi/12.5663706144d0/ *d reconr.1946 data forpi/12.5663706144e0/ *d reconr.1977 data pi/3.14159265358979d0/ *d reconr.1983 data pi/3.14159265358979e0/ *d reconr.1980 data cwaven/2.1968074704d-3/ *d reconr.1986 data cwaven/2.1968074704e-3/ *d reconr.2143 data pi/3.14159265358979d0/ *d reconr.2145 data cwaven/2.1968074704d-3/ *d reconr.2149 data pi/3.14159265358979e0/ *d reconr.2151 data cwaven/2.1968074704e-3/ *d reconr.2284 data pi/3.14159265358979d0/ *d reconr.2286 data cwaven/2.1968074704d-3/ *d reconr.2290 data pi/3.14159265358979e0/ *d reconr.2292 data cwaven/2.1968074704e-3/ *d reconr.2586 data pi/3.14159265358979d0/ *d reconr.2587 data croc/2.1968074704d-3/ *d reconr.2590 data pi/3.14159265358979e0/ *d reconr.2591 data croc/2.1968074704e-3/ *d reconr.2805 data pi/3.14159265358979d0/ *d reconr.2807 data cwaven/2.1968074704d-3/ *d reconr.2811 data pi/3.14159265358979e0/ *d reconr.2813 data cwaven/2.1968074704e-3/ *d reconr.3015 data pi/3.14159265358979d0/ *d reconr.3017 data cwaven/2.1968074704d-3/ *d reconr.3021 data pi/3.14159265358979e0/ *d reconr.3023 data cwaven/2.1968074704e-3/ *d reconr.3202 data rpi/1.772453850906d0/ *d reconr.3203 data cwaven/2.1968074704d-3/ *d reconr.3207 data rpi/1.772453850906e0/ *d reconr.3208 data cwaven/2.1968074704e-3/ *ident up82 */ unresr -- 24feb99 */ update the physical constants in unresr to agree with the */ standard constants to be used with the next version of njoy. */ this helps when making comparisons. cwaven uses the endf */ value for the neutron amu (1.0086652) and codata for the */ amu, hbar, and the ev conversion. *d unresr.791 data pi/3.14159265358979d0/ *d unresr.794 data cwaven/2.1968074704d-3/ *d unresr.814 data pi/3.14159265358979e0/ *d unresr.817 data cwaven/2.1968074704e-3/ *d unresr.1147 data rpi/1.772453850906d0/ *d unresr.1152 data rpi/1.772453850906e0/ *d unresr.1233 data pi2/1.57079632675d0/,sqpi/1.772453850906d0/ data eps/.0002d0/ *d unresr.1239 data pi2/1.57079632675e0/,sqpi/1.772453850906e0/ data eps/.0002e0/ *d unresr.1335 data rpi/1.772453850906d0/ *d unresr.1343 data rpi/1.772453850906e0/ *ident up83 */ thermr -- 24feb99 */ fix an indexing problem when writing out the thermal scattering */ data. it results in the highest e' point having a bad */ scattering probability value (should be zero). *d thermr.1800 a(1+jscr)=sigfig(y(1,i),7,0) *d thermr.1802 a(1+jscr)=sigfig(y(1,i),6,0) */ thermr -- 24feb99 */ improve the accuracy of a constant for better comparisons */ to newer versions of njoy. c is just the sqrt(4*pi). *d thermr.1885 data c/3.5449077018d0/, sigmin/1.d-10/ *d thermr.1890 data c/3.5449077018e0/, sigmin/1.e-10/ *ident up84 */ groupr -- 25feb99 */ fix an error in a data statement. this affects the */ epri-cell-lwr flux around .09 ev on short-word machines. *d groupr.1988 3 .921d0,.059d0,.918d0,.07d0,.892d0,.09d0,.799d0,.112d0,.686d0, */ groupr -- 25feb99 */ improve the accuracy of a constant for better comparisons */ to newer versions of njoy. rp4 is just the sqrt(pi/4). *d groupr.8185 data xmin/1.d-3/, rp4/.88622692544d0/ *d groupr.8193 data xmin/1.e-3/, rp4/.88622692544e0/ */ groupr -- 25feb99 */ the double and single versions of some of the quadrature sets */ are not quite consistent with each other. This can effect */ the high order digits for intermachine comparisons. *d groupr.5860 data qp4/-.8611363116d0,-.3399810436d0,.3399810436d0, 1 .8611363116d0/ *d groupr.5862,5863 data qp8/-.9602898565d0,-.7966664774d0,-.5255324099d0, 1 -.1834346425d0,.1834346425d0,.5255324099d0,.7966664774d0, 2 .9602898565d0/ *d groupr.5866,5868 data qp12/-.9815606342d0,-.9041172564d0,-.7699026742d0, 1 -.5873179543d0,-.3678314990d0,-.1252334085d0,.1252334085d0, 2 .3678314990d0,.5873179543d0,.7699026742d0,.9041172564d0, 3 .9815606342d0/ */ groupr -- 25feb99 */ set a sometimes unitialized variable. the bad values for */ this variable can show up on gendf files. *i groupr.467 izam=0 */ groupr -- 25feb99 */ change the approach used to reduce the precision of some small */ quantities. in some cases, small numbers (such as higher */ legendre moments) are computed as differences between larger */ numbers, and they don't have as much numerical precision as */ the larger numbers. njoy97 has some logic designed to */ reduce the number of significant figures in such numbers in */ order to make machines using different arithmetic produce */ consistent results. this is currently done at the cross-section */ level in displa. we are now moving it to the feed-function */ level, where the truncation can be done using methods more */ appropriate to each particular quantity. truncated cross */ sections are currently fairly obvious because they end in */ several zeros. this feature is lost, but comparisons between */ different machines should still be consistent. *d groupr.3151 *d groupr.3154 *i up52.193 ans(il,1,1)=sigfig(ans(il,1,1),7,0) *d up52.195,196 *i up52.198 ans(1,1,1)=sigfig(ans(1,1,1),7,0) *d groupr.3234,3237 result(il)=sigfig(result(il),7,0) *i groupr.3263 ans(1,1,j)=sigfig(ans(1,1,j),7,0) *i up52.203 if (ig2.eq.2) then ans(il,iz,1)=sigfig(ans(il,iz,1),7,0) if (irat.eq.1) ans(il,iz,3)=sigfig(ans(il,iz,3),7,0) endif *d groupr.3298,3307 result(iz)=sigfig(result(iz),7,0) *i groupr.3345 do 404 j=1,ng2 ans(1,1,j)=sigfig(ans(1,1,j),7,0) 404 continue igzero=1 *d groupr.3353 *i groupr.3357 do 429 j=1,ng2 do 428 il=1,nl ans(il,1,j)=sigfig(ans(il,1,j),7,0) 428 continue 429 continue *d groupr.3367 460 do 464 i=1,ng2 ans(1,1,i)=sigfig(ans(1,1,i),7,0) 464 continue igzero=1 if (iprint.ne.1) return do 465 i=1,ng2,6 *d groupr.3375 *d groupr.3379 470 igzero=1 do 474 i=1,ng2 ans(1,1,i)=sigfig(ans(1,1,i),7,0) 474 continue if (iprint.ne.1) return do 475 i=1,ng2,6 *d groupr.3386 *i groupr.5916 ten=10 *i groupr.6052 c reduce significant figures for small values do 463 il=1,nl ndig=7 fact=ten**ndig iii=nint(fact*ff(il,ng)+ten**(ndig-11)) ff(il,ng)=iii/fact 463 continue *d groupr.6945,6947 fl(il)=c(il+5) if (fl(il).ne.0.) nlz=il *d groupr.6957,6959 fl(il)=c(il+6) if (fl(il).ne.0.) nlz=il *d up52.695 *d up52.697 *d up52.699 200 continue */ groupr -- 25feb99 */ use sigfig to adjust the weight-function enext points to */ get the same grid on different machines. *d groupr.2316 go to 800 *d groupr.2304 go to 800 *d groupr.2317,2318 120 go to 800 *d groupr.2323 go to 800 *d groupr.2328 go to 800 *d groupr.2336 go to 800 *d groupr.2341 go to 800 *d groupr.2344 go to 800 *d groupr.2356 go to 800 *d groupr.2360 go to 800 *d groupr.2370 go to 800 *d groupr.2382 go to 800 *d groupr.2386 go to 800 *d groupr.2390 go to 800 *d groupr.2394 go to 800 *d groupr.2399 go to 800 *d groupr.2401 c c ***return enext on an even grid 800 enext=sigfig(enext,6,0) return */ groupr -- 26feb99 */ eliminate excess gendf send cards for reactions where the */ threshold is higher than the upper limit of the group structure. *d groupr.678 if (ig1.ne.0) call asend(ngout2,0) */ groupr -- 28feb99 */ add three more save variables. this can cause erratic */ behavior on some machines. *i groupr.5635 save elo,ehi,intt */ groupr -- 28feb99 */ truncate the quadrature energy to match the existing energy */ grid points. otherwise, additional panels sometimes occur. *d groupr.3062 170 eq=sigfig(eq,9,0) t1=(eq-elow)/(ehi-elow) *ident up85 */ njoy -- 26feb99 */ on the cray, sigfig(x,n,0) produces a number that is very */ slightly less than what is desired. for example, with x=.009, */ sigfig produces .00899999999999962. this can cause problems */ with searching for interpolation ranges. we change sigfig */ to make the results slightly larger than the target value. *i njoy.2664 data bias/1.0000000000001d0/ *i njoy.2667 data bias/1.0000000000001e0/ *i njoy.2681 xx=xx*bias *ident up86 */ reconr -- 26feb99 */ fix the tests in emerge. they should be fraction tests, not */ tests on 1e-10 barns difference! *d reconr.3621 if (thresh-eg.gt.test*thresh.and.itype.eq.0) go to 370 *d up55.28 if (thresh.gt.one.and.abs(thresh-eg).lt.test*thresh) sn=0 *d reconr.3627,3628 345 if (er-eg.gt.test*eg) go to 360 if (abs(eg-er).lt.test*eg) go to 355 *d reconr.3655 if (thresh-eg.gt.test*thresh.and.itype.eq.0) go to 440 */ reconr -- 26feb99 */ for slbw, the small elastic cross sections in the interference */ minima are computed as differences and lose precision. this */ change truncates the two parts to 8 significant figures before */ calculating the difference. the small differences then have */ fewer significant figures, and intermachine comparisions look */ better. there is no effect on practical applications. *d reconr.2119 50 sigp(2)=sigfig(sigp(2),8,0) spot=sigfig(spot,8,0) sigp(2)=sigp(2)+spot *ident up87 */ acer -- 28feb99 */ fix the threshold logic so that the damage threshold is */ at the minimum energy, even with an error in the heatr grid. *i acer.1839 if (mt.eq.444) thresh=1.e-5 */ acer -- 28feb99 */ change the tests used for thinning angular distributions */ to allow the thinning to be turned off with a small */ enough tolerance value. this makes it easier to do */ comparisions with newer versions of acer that don't thin. *d acer.2470 if (abs(sumi-sum).lt.(toler*abs(sum)+.000001)) go to 105 *d acer.2500 if (abs(p-s).lt.(toler*abs(p)+.000001))go to 240 *ident up88 */ heatr -- 28feb99 */ make sure that sixbar can find the highest energy point */ for the new representation of the energy with the */ new bias added to sigfig. *d heatr.2360 if (nne.eq.ne.and.e.lt.ehi*(1+small)) go to 400 */ heatr -- 28feb99 */ rationalize the numbers before printing to help make them */ come out the same on different machines. *i heatr.1007 ylpd=sigfig(yld,9,0) yp=sigfig(y,9,0) hp=sigfig(h,9,0) damep=sigfig(dame,9,0) ebal6p=sigfig(ebal6,9,0) *d heatr.1014,1015 if (idame.eq.0) then write(nsyso,'(1x,1p,7e14.4)') e,ebarp,yldp,yp,hp else write(nsyso,'(1x,1p,7e14.4)') e,ebarp,yldp,yp,hp,damep endif *d heatr.1017,1019 198 if (idame.eq.0) then write(nsyso,'(1x,1p,7e14.4)') e,q0,ebarp,yldp,yp,hp else write(nsyso,'(1x,1p,7e14.4)') e,q0,ebarp,yldp,yp,hp,damep endif *d heatr.1026 write(nsyso,'(53x,''ebal'',1p,e14.4)') ebal6p *d heatr.1028 202 write(nsyso,'(67x,''ebal'',1p,e14.4)') ebal6p *d heatr.4275 163 ebarp=sigfig(ebar,9,0) xp=sigfig(x,9,0) yp=sigfig(y,9,0) hp=sigfig(h,9,0) egamp=sigfig(egam,9,0) edamp=sigfig(edam,9,0) damep=sigfig(dame,9,0) cerr=sigfig(cerr,4,0) test=1.e-9 *d heatr.4278,4281 if (idame.eq.0) write(nsyso,'(1x,1p,8e14.4)') 1 e,ebarp,xp,yp,hp if (idame.gt.0) write(nsyso,'(1x,1p,8e14.4)') 1 e,ebar,egam,edamp,xp,yp,hp,damep *d heatr.4308,4310 ebarp=sigfig(ebar,9,0) xp=sigfig(x,9,0) yp=sigfig(y,9,0) hxp=sigfig(hx,9,0) hp=sigfig(h,9,0) write(nsyso,'(1x,1p,8e14.4)') e,ebarp,xp,yp,-hxp,hp if (ik.eq.nk.and.nk.gt.1) then subtot=-c(npkk) subtot=sigfig(subtot,9,0) write(nsyso,'(48x,''subtot '',1p,e14.4)') subtot endif *d heatr.4330,4332 ebarp=sigfig(ebar,9,0) yp=sigfig(y,9,0) hxp=sigfig(hx,9,0) hp=sigfig(h,9,0) write(nsyso,'(1x,1p,8e14.4)') e,ebarp,yp,-hxp,hp if (ik.eq.nk.and.nk.gt.1) then subtot=-c(npkk) subtot=sigfig(subtot,9,0) write(nsyso,'(48x,''subtot '',1p,e14.4)') subtot endif *d heatr.4701 do 111 j=1,npkk c(j)=sigfig(c(j),7,0) 111 continue e=c(1) *d heatr.4735 a(ibase+i)=c(inpk) *ident up89 */ broadr -- 29mar99 */ use sigfig on the input maximum energy for broadening and */ thinning. the old coding results in the break energy */ being included in the broadening for some machines, and */ not being included for others. the effects of this change */ will usually be in the seventh digit at the break energy. *i broadr.168 thnmx=sigfig(thnmx,7,0) *ident up90 */ acer -- 29mar99 */ the scratch tape containing the converted endf cross sections */ is not properly terminated. this results in a read error */ in some cases. *d acer.818 call totend(nine,noute,0,a(iscr)) */ acer -- 29mar99 */ reduce the precision of the some calculated values */ in order to get the same results using binary or ascii. *d acer.2825,2826 a(2*l-1+ibase)=sigfig(a(l-1+ixxmu),7,0) a(2*l+ibase)=sigfig(a(l-1+ipmu),7,0) *d acer.3183,3184 a(2*i-1+ibase)=sigfig(tbmu(i),7,0) a(2*i+ibase)=sigfig(pmu(i),7,0) *i acer.3554 y=sigfig(y,7,0) *ident up91 */ groupr -- 29mar99 */ use more digits for the lethargy constants for the gam-ii */ structure (suggested by muir, iaea). *d groupr.1389 data bgam2,tgam2/27.631021d0,-0.53062825d0/ *d up59.8 data bgam2,tgam2/27.631021e0,-0.53062825e0/ */ groupr -- 29mar99 */ protect against division by zero. this can occur for */ the recoil nucleus. noted by vladimirov (kurchatov). *d groupr.6093 if (af**2.eq.awr2) then ecl=emax else ecl=thresh/(1-af**2/awr2) endif *d groupr.6114 if (af**2.eq.awr2) then ech=emax else ech=thresh/(1-af**2/awr2) endif */ groupr -- 29mar99 */ fix groupr to handle the lct=3 option. this had been */ done before for heatr and acer, but groupr was omitted. */ noted by vladimirov (kurchatov). remember that lct=1 */ means lab, lct=2 means cm, and lct means cm for particles */ from alpha down, and lab for heavier particles. *b groupr.4403 alight=5 *d groupr.4538 147 if (lct.eq.1) go to 155 if (lct.eq.3.and.awp.gt.alight) go to 155 *i groupr.4704 if (law.eq.1.and.lct.eq.3.and.awp.lt.alight) langn=1 *i groupr.4751 if (lct.eq.3.and.awp.lt.alight) go to 450 *ident up92 */ viewr -- 29mar99 */ reduce the size of the legend block to make it easier */ to fit on plots without obscuring the curves. also, */ fix the axis labeling to allow non-integer values for */ the lower and upper limits without mislabeling (e.g., */ 1.5 to 4.5 by steps of 1.0). from muir (iaea). *d viewr.800 hleg=2*hlab/3 *d viewr.2106 *d viewr.2122 *i viewr.2132 ifracs=0 if (abs(vv-iv).gt..01) ifracs=1 *d viewr.2134 write(num,'(f5.1)') vv *ident up93 */ heatr -- 29mar99 */ use 6 significant figures for heating values larger */ than 1e10 to make the binary and ascii pendf files */ equivalent. helps with diffs in acer. *i heatr.4645 data big/1.d10/ *i heatr.4647 data big/1.e10/ *d up88.73 if (c(j).gt.big) then c(j)=sigfig(c(j),6,0) else c(j)=sigfig(c(j),7,0) endif *ident up94 */ njoy -- 30mar99 */ when a tab1 starts with a zero followed by a discontinuity, */ such as the kind of things found at resonance range boundaries */ or when photon production starts artificially as in MF13, shade */ xnext down slightly so that the threshold is not missed. */ for an example, check the threshold of mt51 in out08. */ found by frankle and little (lanl). *i njoy.2229 if (a(loc+3).ne.zero.and.ip.gt.1) xnext=down*xnext *ident up95 */ unresr -- 30mar99 */ reorder some arithmetic to give better values for */ the infinitely-dilute transport cross section. *d unresr.1050 term1=fact*(yk-yj)/(1-yk) *ident up96 */ njoy -- 14apr99 */ the bias factor must be real*8 *d njoy.2664 real*8 x,xx,aa,ten,bias,sigfig *ident up97 */ heatr -- 14apr99 */ fix a typo from up88. it only affects the listing. *d up88.12 ebarp=sigfig(ebar,9,0) yldp=sigfig(yld,9,0) *ident up98 */ groupr -- 14apr99 */ this attempt to improve the numerics causes some machines */ to go into a infinite loop. restore the old way. *d up84.154 800 continue */ groupr -- 15apr99 */ add the xmas group structure. requested by the nea. */ this structure contains about 80 thermal groups and is */ used in Europe for thermal reactor neutronics calculations. *i groupr.129 c * 18 xmas 172-group structure *i groupr.1250 c 18 xmas 172-group structure *i groupr.1264 dimension eg18(173) *i groupr.1386 data eg18/ &1.96403d+7,1.73325d+7,1.49182d+7,1.38403d+7,1.16183d+7,1.00000d+7, &8.18731d+6,6.70320d+6,6.06531d+6,5.48812d+6,4.49329d+6,3.67879d+6, &3.01194d+6,2.46597d+6,2.23130d+6,2.01897d+6,1.65299d+6,1.35335d+6, &1.22456d+6,1.10803d+6,1.00259d+6,9.07180d+5,8.20850d+5,6.08101d+5, &5.50232d+5,4.97871d+5,4.50492d+5,4.07622d+5,3.01974d+5,2.73237d+5, &2.47235d+5,1.83156d+5,1.22773d+5,1.11090d+5,8.22975d+4,6.73795d+4, &5.51656d+4,4.08677d+4,3.69786d+4,2.92830d+4,2.73944d+4,2.47875d+4, &1.66156d+4,1.50344d+4,1.11378d+4,9.11882d+3,7.46586d+3,5.53084d+3, &5.00451d+3,3.52662d+3,3.35463d+3,2.24867d+3,2.03468d+3,1.50733d+3, &1.43382d+3,1.23410d+3,1.01039d+3,9.14242d+2,7.48518d+2,6.77287d+2, &4.53999d+2,3.71703d+2,3.04325d+2,2.03995d+2,1.48625d+2,1.36742d+2, &9.16609d+1,7.56736d+1,6.79041d+1,5.55951d+1,5.15780d+1,4.82516d+1, &4.55174d+1,4.01690d+1,3.72665d+1,3.37201d+1,3.05113d+1,2.76077d+1, &2.49805d+1,2.26033d+1,1.94548d+1,1.59283d+1,1.37096d+1,1.12245d+1, &9.90555d+0,9.18981d+0,8.31529d+0,7.52398d+0,6.16012d+0,5.34643d+0, &5.04348d+0,4.12925d+0,4.00000d+0,3.38075d+0,3.30000d+0,2.76792d+0, &2.72000d+0,2.60000d+0,2.55000d+0,2.36000d+0,2.13000d+0,2.10000d+0, &2.02000d+0,1.93000d+0,1.84000d+0,1.75500d+0,1.67000d+0,1.59000d+0, &1.50000d+0,1.47500d+0,1.44498d+0,1.37000d+0,1.33750d+0,1.30000d+0, &1.23500d+0,1.17000d+0,1.15000d+0,1.12535d+0,1.11000d+0,1.09700d+0, &1.07100d+0,1.04500d+0,1.03500d+0,1.02000d+0,9.96000d-1,9.86000d-1, &9.72000d-1,9.50000d-1,9.30000d-1,9.10000d-1,8.60000d-1,8.50000d-1, &7.90000d-1,7.80000d-1,7.05000d-1,6.25000d-1,5.40000d-1,5.00000d-1, &4.85000d-1,4.33000d-1,4.00000d-1,3.91000d-1,3.50000d-1,3.20000d-1, &3.14500d-1,3.00000d-1,2.80000d-1,2.48000d-1,2.20000d-1,1.89000d-1, &1.80000d-1,1.60000d-1,1.40000d-1,1.34000d-1,1.15000d-1,1.00001d-1, &9.50000d-2,8.00000d-2,7.70000d-2,6.70000d-2,5.80000d-2,5.00000d-2, &4.20000d-2,3.50000d-2,3.00000d-2,2.50000d-2,2.00000d-2,1.50000d-2, &1.00000d-2,6.90000d-3,5.00000d-3,3.00000d-3,1.00001d-5/ *i groupr.1496 data eg18/ &1.96403e+7,1.73325e+7,1.49182e+7,1.38403e+7,1.16183e+7,1.00000e+7, &8.18731e+6,6.70320e+6,6.06531e+6,5.48812e+6,4.49329e+6,3.67879e+6, &3.01194e+6,2.46597e+6,2.23130e+6,2.01897e+6,1.65299e+6,1.35335e+6, &1.22456e+6,1.10803e+6,1.00259e+6,9.07180e+5,8.20850e+5,6.08101e+5, &5.50232e+5,4.97871e+5,4.50492e+5,4.07622e+5,3.01974e+5,2.73237e+5, &2.47235e+5,1.83156e+5,1.22773e+5,1.11090e+5,8.22975e+4,6.73795e+4, &5.51656e+4,4.08677e+4,3.69786e+4,2.92830e+4,2.73944e+4,2.47875e+4, &1.66156e+4,1.50344e+4,1.11378e+4,9.11882e+3,7.46586e+3,5.53084e+3, &5.00451e+3,3.52662e+3,3.35463e+3,2.24867e+3,2.03468e+3,1.50733e+3, &1.43382e+3,1.23410e+3,1.01039e+3,9.14242e+2,7.48518e+2,6.77287e+2, &4.53999e+2,3.71703e+2,3.04325e+2,2.03995e+2,1.48625e+2,1.36742e+2, &9.16609e+1,7.56736e+1,6.79041e+1,5.55951e+1,5.15780e+1,4.82516e+1, &4.55174e+1,4.01690e+1,3.72665e+1,3.37201e+1,3.05113e+1,2.76077e+1, &2.49805e+1,2.26033e+1,1.94548e+1,1.59283e+1,1.37096e+1,1.12245e+1, &9.90555e+0,9.18981e+0,8.31529e+0,7.52398e+0,6.16012e+0,5.34643e+0, &5.04348e+0,4.12925e+0,4.00000e+0,3.38075e+0,3.30000e+0,2.76792e+0, &2.72000e+0,2.60000e+0,2.55000e+0,2.36000e+0,2.13000e+0,2.10000e+0, &2.02000e+0,1.93000e+0,1.84000e+0,1.75500e+0,1.67000e+0,1.59000e+0, &1.50000e+0,1.47500e+0,1.44498e+0,1.37000e+0,1.33750e+0,1.30000e+0, &1.23500e+0,1.17000e+0,1.15000e+0,1.12535e+0,1.11000e+0,1.09700e+0, &1.07100e+0,1.04500e+0,1.03500e+0,1.02000e+0,9.96000e-1,9.86000e-1, &9.72000e-1,9.50000e-1,9.30000e-1,9.10000e-1,8.60000e-1,8.50000e-1, &7.90000e-1,7.80000e-1,7.05000e-1,6.25000e-1,5.40000e-1,5.00000e-1, &4.85000e-1,4.33000e-1,4.00000e-1,3.91000e-1,3.50000e-1,3.20000e-1, &3.14500e-1,3.00000e-1,2.80000e-1,2.48000e-1,2.20000e-1,1.89000e-1, &1.80000e-1,1.60000e-1,1.40000e-1,1.34000e-1,1.15000e-1,1.00001e-1, &9.50000e-2,8.00000e-2,7.70000e-2,6.70000e-2,5.80000e-2,5.00000e-2, &4.20000e-2,3.50000e-2,3.00000e-2,2.50000e-2,2.00000e-2,1.50000e-2, &1.00000e-2,6.90000e-3,5.00000e-3,3.00000e-3,1.00001e-5/ *d groupr.1509 1 1200,1000,1300,1300,1400), ign *i groupr.1687 c c ***xmas 172-group structure 1400 ng=172 do 1410 ig=1,173 eg(ig)=eg18(173-ig) 1410 continue go to 1900 *i groupr.1730 if (ign.eq.18) write(nsyso,'(/ 1 '' neutron group structure......xmas 172-group'')') *ident up99 */ reconr -- 12may99 */ don't extrapolate for the unresolved cross section. this */ can occur for elemental evaluations when the upper limit */ of the resolved ranges is greater than the upper limit of */ the unresolved ranges. an example is jeff zr-nat. *d reconr.1174 if (e.ge.e1.and.e.le.e2) then call terp1(e1,a(l1+i),e2,a(l2+i),e,sigu(i),intunr) else sigu(i)=0 endif *ident up100 */ matxsr -- 13may99 */ set uninitialized variable. reported by mattes. */ causes a storage exceeded message for binary mode. *i matxsr.889 nmax=2000 *ident up101 */ ccccr -- 13may99 */ fix up the pointers for delayed-neutron processing. */ patch provided by broeders and kiefhaber, fzk. *i ccccr.2897 lzero=next *d ccccr.2904,2909 lzeroh=(lzero-1)/mult+1 read(hname,'(a6)') ha(lzeroh) ha(lzeroh+1)=huse(1) ha(lzeroh+2)=huse(2) ia(3*mult+lzero)=ivers call rite(idlay,irec,a(lzero),nwds,0) *d ccccr.2914,2918 ia(lzero)=ngroup ia(lzero+1)=nisod ia(lzero+2)=nfam ia(lzero+3)=nfam call rite(idlay,irec,a(lzero),nwds,0) *d ccccr.2923 call rite(idlay,irec,a(l1),nwds,0) *d ccccr.2926 l3=l1+nwds *d ccccr.2929 l4=l1+nisod*mult+nfam+ngroup*nfam+ngroup *d ccccr.2993 next=lzero+nwds+4 *d ccccr.3227 c ***delayed neutron precursor yield data *ident up102 */ groupr -- 29aug97 */ enhance activation processing and provide an option */ to automatically process all the mats on the endf input tape. *d groupr.42,43 c * if ngout=0, matb<0 is an option to automatically * c * process all the mats on the endf input tape. * c * otherwise, matb<0 is a flag to add mts to and/or * c * replace individual mts on gout1. * *d groupr.189,192 c * 1zzzaaam nuclide production for zzzaaam * c * subsection from file 3 * c * 2zzzaaam nuclide production for zzzaaam * c * subsection from file 6 * c * 3zzzaaam nuclide production for zzzaaam * c * subsection from file 9 * c * 4zzzaaam nuclide production for zzzaaam * c * subsection from file 10 * *d groupr.206 c * 10/ do all isotope productions using mf8 * *d groupr.238,239 common/rlist/mf4(50),mf6(50),mf12(50),mf13(50),mf18(50), 1 mf4r(6,50),mf6p(6,50),mf10f(200),mf10s(200),mf10i(200) *d groupr.341 if (iaddmt.lt.0) then call contio(nendf,0,0,a(iscr),nb,nw) if (math.lt.0) go to 650 matb=math call skiprz(nendf,-1) else call findf(matb,1,451,nendf) endif *d groupr.483,484 if (mfd.gt.36.and.mfd.le.10000000) go to 381 *d groupr.470,472 call nextr(iauto,matd,mfd,mtdp,a(iscr)) izam=0 if (mfd.gt.10000000) izam=mod(mfd,10000000) *d groupr.496,498 call nextr(iauto,matd,mfd,mtdp) izam=0 if (mfd.gt.10000000) izam=mod(mfd,10000000) *d groupr.509,521 izam=mod(mfd,10000000) *d groupr.523 if (iaddmt.le.0) go to 490 *i groupr.701 call tomend(nendf,0,0,a(iscr)) *i groupr.703 if (iaddmt.lt.0) go to 160 *i groupr.706 if (matb.lt.0.and.ngout1.eq.0) iaddmt=-1 *i groupr.773 if (matb.lt.0.and.ngout1.eq.0) iaddmt=-1 *d groupr.837,838 common/rlist/mf4(50),mf6(50),mf12(50),mf13(50),mf18(50), 1 mf4r(6,50),mf6p(6,50),mf10f(200),mf10s(200),mf10i(200) *d groupr.877 c nuclide production (first time) *d groupr.881,883 mfd=mf10i(ir) if (mf10f(ir).eq.3) mfd=mfd+10000000 if (mf10f(ir).eq.6) mfd=mfd+20000000 if (mf10f(ir).eq.9) mfd=mfd+30000000 if (mf10f(ir).eq.10) mfd=mfd+40000000 mtd=mf10s(ir) *d groupr.960,976 c nuclide production (after first 10/ entry) 430 if (mfd.le.10000000) go to 280 if (mf10f(ir).eq.0) go to 280 mfd=mf10i(ir) if (mf10f(ir).eq.3) mfd=mfd+10000000 if (mf10f(ir).eq.6) mfd=mfd+20000000 if (mf10f(ir).eq.9) mfd=mfd+30000000 if (mf10f(ir).eq.10) mfd=mfd+40000000 mtd=mf10s(ir) *d groupr.1181,1182 common/rlist/mf4(50),mf6(50),mf12(50),mf13(50),mf18(50), 1 mf4r(6,50),mf6p(6,50),mf10f(200),mf10s(200),mf10i(200) *d groupr.1214,1215 common/rlist/mf4(50),mf6(50),mf12(50),mf13(50),mf18(50), 1 mf4r(6,50),mf6p(6,50),mf10f(200),mf10s(200),mf10i(200) *d groupr.3492 if (mft.gt.40000000) mft=10 if (mft.gt.30000000) mft=9 if (mft.gt.20000000) mft=6 *i groupr.3494 nk=nint(a(iyld+4)) ik=0 *d groupr.3498 if (mf.gt.10000000) lfs=mod(mf,10) iza=0 if (mf.gt.10000000) iza=mod(mf/10,1000000) *d groupr.3513,3525 130 ik=ik+1 loc=iyld call tab1io(itape,0,0,a(loc),nb,nw) 140 loc=loc+nw if (nb.eq.0) go to 145 call moreio(itape,0,0,a(loc),nb,nw) go to 140 145 if (mft.eq.9.or.mft.eq.10) izn=nint(a(iyld+2)) if (mft.eq.6) izn=nint(a(iyld)) lfn=nint(a(iyld+3)) if (iza.eq.izn.and.lfs.eq.lfn) go to 180 if (izn.eq.0.and.lfs.eq.lfn) go to 180 if (mft.eq.9.or.mft.eq.10) go to 170 law=nint(a(iyld+3)) if (law.eq.0.or.law.eq.2.or.law.eq.3) go to 170 loc=iyld call tab2io(itape,0,0,a(loc),nb,nw) ne=nint(a(iyld+5)) do 160 ie=1,ne call listio(itape,0,0,a(loc),nb,nw) 150 if (nb.eq.0) go to 160 call moreio(itape,0,0,a(loc),nb,nw) go to 150 160 continue 170 if (ik.lt.nk) go to 130 call error('getyld','requested nuclide not found',' ') *d groupr.3589,3591 if (mfd.gt.10000000) mf=3 if (mfd.gt.20000000) mf=3 if (mfd.gt.30000000) mf=3 if (mfd.gt.40000000) mf=10 if (mfd.gt.10000000) lfs=mod(mfd,10) *d groupr.3922,3923 if (mfd.gt.10000000) go to 100 *d groupr.3969 if (mfd.eq.12.or.mfd.gt.20000000) *d groupr.7228,7229 common/rlist/mf4(50),mf6(50),mf12(50),mf13(50),mf18(50), 1 mf4r(6,50),mf6p(6,50),mf10f(200),mf10s(200),mf10i(200) *ident up103 */ acer -- 13may99 */ the complex procedure used to break the two subsections of */ the f-19 (n,2n) reaction up into two separate reactions */ (mt6 and mt46) is not really needed to make proper mcnp runs. */ we are removing it for simplicity. *d acer.407 external openz,mess,reserv,releas,closz *d acer.429,447 *d acer.492,947 *d acer.2089,2091 *d acer.2097 call fix6(nin,nsix,ik) *d acer.3297 subroutine fix6(nin,nout,nk) *i acer.3301 c convert just the first nk subsections of this mf6 section. *d acer.3322 a(5)=nk *i acer.3324 if (nk.gt.1) & write(nsyso,'(i6,'' neutron subsections found'')') nk *d acer.3328 if (ik.gt.nk) go to 340 */ acer -- 13may99 */ add an additional consistency check for the gamma production */ sum. make sure it is equal to the gamma production (gpd). *i acer.11016 write(nsyso,'(/'' check photon production sum'')') do 1950 ie=1,nes e=xss(esz+ie-1) gsum=0 do 1945 i=1,ntrp loct=nint(xss(i-1+lsigp)+sigp-1) loc1=nint(xss(i-1+ldlwp)+dlwp-1) mftype=nint(xss(loct)) if (mftype.eq.13) go to 1941 loct=loct+1 mtmult=nint(xss(loct)) do ii=1,ntr if (mtmult.eq.nint(xss(mtr+ii-1))) then k=nint(xss(lsig+ii-1)+sig-1) iaa=nint(xss(k)) naa=nint(xss(k+1)) endif enddo if (ie.lt.iaa.or.ie.ge.iaa+naa) go to 1945 loct=loct+1 locc=loct m=nint(xss(loct)) loct=loct+1+2*m n=nint(xss(loct)) law=nint(xss(loc1+1)) if (law.ne.2) go to 1941 if (m.ne.0.or.n.gt.2) go to 1941 if (xss(loct+3).ne.xss(loct+4)) go to 1941 gsum=gsum+xss(loct+3)*xss(2+k+ie-iaa) go to 1945 1941 if (mftype.eq.13) go to 1942 call terpc(e,y,xss(locc)) gsum=gsum+y*xss(2+k+ie-iaa) go to 1945 1942 loct=loct+1 call terpc(e,ss,xss(loct)) gsum=gsum+ss 1945 continue gg=xss(gpd+ie-1) if (abs(gg-gsum).gt.gg/10000) then nerr=nerr+1 write(6,'(1p,3e16.6)') e,gg,gsum endif 1950 continue c *i acer.11071 c subroutine terpc(x,y,a) c ***************************************************************** c interpolate in an ace table c ***************************************************************** *if sw implicit real*8 (a-h,o-z) *endif dimension a(*) i=1 m=nint(a(i)) i=i+1 if (m.ne.0) i=i+2*m mm=2 if (m.gt.0) then jj=1 mm=nint(a(1+jj+m)) endif n=nint(a(i)) do j=1,n-1 if (m.gt.0.and.j.ge.nint(a(1+jj))) then jj=jj+1 mm=nint(a(1+jj+m)) endif if (a(i+j).le.x.and.a(i+j+1).ge.x) then call terp1(a(i+j),a(i+n+j),a(i+j+1),a(i+j+1+n),x,y,mm) endif enddo return end */ acer -- 13may99 */ provide additional scratch space for endf f-19. *d acer.1999 data nwmaxn/20000/ *ident up104 */ groupr -- 28jan00 */ fix a typographical error in up102 */ this causes groupr to fail for automatic input *d up102.40 call nextr(iauto,matd,mfd,mtdp,a(iscr)) *ident up105 */ heatr -- 30jan00 */ add missing external statement *i heatr.3823 external error *ident up106 */ groupr -- 30jan00 */ add missing external statement *i groupr.7239 external error */ groupr -- 30jan00 */ fix problem with finding the correct interpolation range */ in getunr. causes glitch at resolved-unresolved boundary */ on some machines. see the 200 eV region in test problem 2. *i groupr.3778 data del/1.d-10/ *i groupr.3780 data del/1.e-10/ *d groupr.3794 if (e.lt.e1-del*e1) return *d groupr.3796 if (e.gt.enunr+del*enunr) return *d groupr.3815 elo=abs(a(loc)) elo=elo-del*elo ehi=abs(a(loc1)) if (e.gt.elo.and.e.lt.ehi) go to 120 *ident up107 */ wimsr -- 30jan00 */ add missing external statement *i wimsr.112 external error *ident up108 */ groupr -- 28mar00 */ fix a problem introduced with the activation patch. */ add some statements to fix problems when running */ without static mode under g77 on linux. fix a problem */ identified by pelloni (psi) with tabulated delayed */ neutron data. these problems show up when doing */ delayed-neutron spectra with various evaluations. *d up102.88 nk=0 if (mf.gt.99) nk=nint(a(iyld+4)) *d up102.102 145 if (nk.eq.0) go to 180 if (mft.eq.9.or.mft.eq.10) izn=nint(a(iyld+2)) *i groupr.8183 save loct *i groupr.8250 np=nint(a(loct+6)) *i groupr.8123 if (mtd.eq.455) go to 280 */ groupr -- 28 mar 00 */ fix the indexing for the 172-group xmas structure *d up98.86 eg(ig)=eg18(174-ig) *ident up109 */ acer -- 03apr00 */ increase the available storage to handle the very large */ mf6/mt16 tabulation in JEFF-3 Be-9. *d up27.21 common/astore/a(80000) *d up27.23 data namax/80000/, nidmax/27/ *d up27.25 common/astore/a(80000) *d up103.108 data nwmaxn/65000/ *d up27.27 common/astore/a(80000) *d up27.29 common/astore/a(80000) *d up27.31 common/astore/a(80000) *d up27.33 common/astore/a(80000) *d up27.35 common/astore/a(80000) *d up27.37 common/astore/a(80000) *d up27.39 common/astore/a(80000) *d up27.41 common/astore/a(80000) *d up27.43 common/astore/a(80000) *d up27.45 common/astore/a(80000) */ acer -- 03apr00 */ adjust the fix6 routine to handle cm distributions that do not */ span the entire range of lab angles. *d acer.3395,3396 *d up38.34 *d up38.52 if (ep.gt.0.) csn=clb*sqrt(elb/ep)-sqrt(ein/ep)/aw1 *d up38.64 305 continue if (j.gt.l2+8.and.elb.le.a(j-2)) go to 306 a(j)=elb *d acer.3415,3426 306 if (j.ge.999) call error('fix6','storage in a exceeded',' ') if (iep.lt.nep) go to 290 nnep=(j-(l2+8))/2 if (nnep.eq.1) then a(l2+10)=a(l2+8)*2 a(l2+11)=0 nnep=2 endif a(l2+5)=nnep a(l2+6)=nnep j2=l2 call tab1io(0,nout,ndebug,a(j2),nb,nw) do while (nb.ne.0) j2=j2+nw call moreio(0,nout,ndebug,a(j2),nb,nw) enddo *ident up110 */ leapr -- 04apr00 */ fix an error made in up53 that destroys the calculations for */ cold parahydrogen, orthodeuterium, and paradeuterium. The */ evidence is too large cross sections at high energies (such */ as 0.5 eV). It was just a typing problem, but it didn't get */ caught until we tried the calculations with NJOY99. *d up53.274 if (law.gt.3) de=ded *ident up111 */ acer -- 30may00 */ more work on handling file 6 cm distributions by conversion */ to the lab angle-energy representation. add the jacobian */ for the transformation, and fix the check that makes sure */ that the outgoing energies are increasing. this is not a */ perfect way to do this conversion because it doesn't */ reconstruct the energy grid to provide a faithful */ representation in the lab. this patch is needed for several */ jef-2.2 evaluations. *d up38.44 c if(imu.lt.nmu.and.amu(imu+1).le.cmn) drv=0 *i up38.57 c add jacobian for the cm-to-lab transformation if (ep.ne.zero) drv=drv*sqrt(elb/ep) *d up109.42 if (j.gt.l2+8) go to 307 *i acer.3414 307 if (elb.le.a(j-2)) go to 306 a(j)=elb a(j+1)=fmu*drv j=j+2 */ fix the unit for photon production consistency messages *d up103.70 write(nsyso,'(1p,3e16.6)') e,gg,gsum *ident up112 */ viewr -- 30may00 */ increase the allowed size of 3d plots *d viewr.310 dimension a(1200),z(15),aa(75000) *d viewr.319 data maxaa/75000/ *d viewr.611 if (l+5000.ge.maxaa) then *d viewr.1168 dimension x(2000),y(2000),z(2000) *d viewr.1170 kmax=1999 *ident up113 */ acer -- 27jun00 */ add missing line for the sequential (n,2n) reactions for be-9. */ this line was accidentally removed in njoy 94.19 (jan96). the */ error was continued through njoy 97 and 99. discovered by */ bob little (lanl). *i acer.4832 if (mth.ge.46.and.mth.le.49) s=s/2 */ acer -- 26jun00 */ patches provides by bokyn seo (kaeri) to add a missing */ external statement and to fix an uninitialized variable in */ our recent update up111. *i acer.3312 external error */ fix an uninitialized variable *i up111.14 zero=0 */ acer -- 26jun00 */ fix an error in the reaction naming. it affects mt=44 (n,n2p) and */ mt=45 (n,npa). this problem was introduced in may of 1995. */ examples of cases that use these reactions are later releases of */ endf/b-vi al-27 and o-16. *d acer.7703,7704 a '(n,x) ', '(n,2np) ', '(n,3np) ', '(n,x) ', b '(n,n2p) ', '(n,npa) ', '(n,2/2*1) ', '(n,2/2*2) ', *ident up114 */ leapr -- 27jun00 */ fix two incorrect constants in leapr. one affects cases with */ diffusive effects, and it has been incorrect since njoy97.0 */ (oct 97). the other affects cold hydrogen calculations, and it */ has been incorrect since njoy97.53 (dec98). *d leapr.1097 data c0/.125d0/ *d leapr.1112 data c0/.125/ *d up53.247 data amassh/3.3465d0/ *d up53.260 data amassh/3.3465e0/ *ident up115 */ reconr -- 03sep00 */ add changes provided by Nancy Larson at Oak Ridge for the */ new convention added to the Reich-Moore format that uses the */ sign of AJ to indicate which of the two channel spins is */ to be used. *i reconr.2361 c take care of new format for indicating channel spin do 175 kchanl=1,2 inow=inowb kpstv=0 kngtv=0 *d reconr.2374 aj=abs(a(inow+1)) *i up47.291 if (a(inow+1).lt.zero)kngtv=kngtv+1 if (a(inow+1).gt.zero)kpstv=kpstv+1 if (kchanl.eq.1.and.a(inow+1).lt.zero) go to 140 if (kchanl.eq.2.and.a(inow+1).gt.zero)go to 140 *i reconr.2416 c take care of extra channel spin as define by sign of aj c kkkkkk = 0 => do not add anything in here c kkkkkk = 1 => add resonance contribution but not extra hard-sphere c kkkkkk = 2 => add resonance plus hard-sphere phase shift contrib kkkkkk = 0 if (kchanl.eq.1) then if (kpstv.gt.0) then if (kngtv.eq.0) then if (jj.gt.jjl.and.jj.lt.numj) then kkkkkk = 2 else kkkkkk = 1 endif else if (kngtv.gt.0) then kkkkkk = 1 endif else if (kpstv.eq.0) then if (kngtv.eq.0) then if (jj.gt.jjl.and.jj.lt.numj) then kkkkkk = 2 else kkkkkk = 1 endif else if (kngtv.gt.0) then kkkkkk = 0 endif endif else if (kchanl.eq.2) then if (kpstv.gt.0) then if (kngtv.eq.0) then else if (kngtv.gt.0) then kkkkkk = 1 endif else if (kpstv.eq.0) then if (kngtv.eq.0) then else if (kngtv.gt.0) then if (jj.gt.jjl.and.jj.lt.numj) then kkkkkk = 2 else kkkkkk = 1 endif endif endif endif if (kkkkkk.eq.0) go to 175 *d up47.60 170 if (kkkkkk.eq.2) then *i reconr.2453 175 continue *ident vers */ update the version name and date */ to reflect the date of the latest modifications *d njoy.8,9 c * version 97.115 * c * 03 sep 00 * *d njoy.279 data vers/'97.115 '/