diff --git a/sammy/src/clm/CrystalLatticeBroadening_M.f90 b/sammy/src/clm/CrystalLatticeBroadening_M.f90
new file mode 100644
index 0000000000000000000000000000000000000000..ed0d5757cefc97bbe63f7f319f86c21b991e3299
--- /dev/null
+++ b/sammy/src/clm/CrystalLatticeBroadening_M.f90
@@ -0,0 +1,49 @@
+module CrystalLatticeBroadening_M
+use FreeGasDopplerBroadening_M
+use DopplerBroadening_M
+use DerivativeHandler_M
+use GridData_M
+use, intrinsic :: ISO_C_BINDING
+
+
+implicit none
+
+type, extends(FreeGasDopplerBroadening) :: CrystalLatticeBroadening
+contains
+procedure, pass(this) :: initialize => CrystalLatticeBroadening_initialize
+procedure, pass(this) :: broaden => CrystalLatticeBroadening_broaden
+procedure, pass(this) :: destroy => CrystalLatticeBroadening_destroy
+end type CrystalLatticeBroadening
+
+contains
+subroutine CrystalLatticeBroadening_initialize(this, hand, list, work)
+ implicit none
+ class(CrystalLatticeBroadening) :: this
+ class(DerivativeHandler)::hand
+ class(GridDataList)::list
+ class(GridDataList)::work
+ call FreeGasDopplerBroadening_initialize(this, hand, list, work)
+end subroutine
+
+subroutine CrystalLatticeBroadening_destroy(this)
+ implicit none
+ class(CrystalLatticeBroadening) :: this
+ call FreeGasDopplerBroadening_destroy(this)
+end subroutine
+
+
+
+subroutine CrystalLatticeBroadening_broaden(this)
+ class(CrystalLatticeBroadening) :: this
+
+ integer::ndatb
+ type(DerivativeHandler)::data
+
+ ! reserve the data
+ ndatb = this%getNumEnergyBroadened()
+ call this%getData(data)
+ call data%nullify()
+ call data%reserve(ndatb*data%getNnnsig(), this%numPar+1)
+end subroutine
+
+end module CrystalLatticeBroadening_M
diff --git a/sammy/src/cro/mcro5.f90 b/sammy/src/cro/mcro5.f90
index 831388cbadfc114108dffa29dd377e4a5056fc5f..afa5c3abac0ae59c79f5d52f41d7b0a985943dc4 100644
--- a/sammy/src/cro/mcro5.f90
+++ b/sammy/src/cro/mcro5.f90
@@ -5,6 +5,9 @@ contains
use DopplerAndResolutionBroadener_M
use SumIsoAndConvertToTrans_M
use array_sizes_common_m, only : calcData, calcDataSelf
+ use normalize_and_background, only : Nnneta
+ use broad_common_m, only : dopplerOption
+ use exploc_common_m, only : I_Iflmsc
use ifwrit_m, only : Kcros
use fixedi_m, only : numUsedPar
use, intrinsic :: ISO_C_BINDING
@@ -14,10 +17,12 @@ contains
class(DopplerAndResolutionBroadener)::broadener
type(SumIsoAndConvertToTrans)::summer
-
call summer%initialize(broadener, calcDataSelf, Kcros.eq.8)
summer%Need_Isotopes=needIso
summer%Another_Process_Will_Happen=moreBroadening
+ if (.not.moreBroadening.and.Kcros.EQ.6.and.dopplerOption%bType.EQ.0) then
+ CALL Nnneta (I_Iflmsc, calcData)
+ end if
call summer%sumAndConvert(numUsedPar)
call summer%destroy
end subroutine sumAndConvertAfterBroadening
diff --git a/sammy/src/cro/mnrm1.f90 b/sammy/src/cro/mnrm1.f90
index ef2d17e2e1f686da4998222227eaaf69cab58905..3509409454926a720982e774eb8ef3d52a8bd09b 100644
--- a/sammy/src/cro/mnrm1.f90
+++ b/sammy/src/cro/mnrm1.f90
@@ -359,7 +359,7 @@ contains
!
! --------------------------------------------------------------
!
- SUBROUTINE Nnneta (Iflmsc, deriv, irow, Niso)
+ SUBROUTINE Nnneta (Iflmsc, deriv)
!
! *** Purpose -- Generate derivative with respect to NU
!
@@ -369,7 +369,7 @@ contains
use DerivativeHandler_M
IMPLICIT None
!
- integer::irow, Niso
+ integer::irow
integer::Iflmsc(*)
real(kind=8)::val
type(DerivativeHandler)::deriv
@@ -380,10 +380,12 @@ contains
IF (Kjetan.LE.0) RETURN
IF (Iflmsc(Kjetan).EQ.0) RETURN
N = Iflmsc(Kjetan)
- DO Iso=1,Niso
- val = deriv%getDataNs(irow, 1, 0, iso)
- val = val/Etanuu
- call deriv%addDataNs(irow, 1, N, iso, val)
+ Do irow = 1, deriv%getLength()
+ DO Iso=1,deriv%getUsedIsotopes()
+ val = deriv%getDataNs(irow, 1, 0, iso)
+ val = val/Etanuu
+ call deriv%addDataNs(irow, 1, N, iso, val)
+ end do
end do
RETURN
END
diff --git a/sammy/src/dbd/HighEnergyFreeGas_m.f90 b/sammy/src/dbd/HighEnergyFreeGas_m.f90
index 1b525f5510edead6d7f5e9aedfb3ec25790fe6b3..293fbcd84ea30bc89b769bf6655229c31931418a 100644
--- a/sammy/src/dbd/HighEnergyFreeGas_m.f90
+++ b/sammy/src/dbd/HighEnergyFreeGas_m.f90
@@ -35,6 +35,19 @@ end subroutine
subroutine HighEnergyFreeGas_broaden(this)
class(HighEnergyFreeGas) :: this
+
+ integer::ndatb
+ type(DerivativeHandler)::data
+
+ ! reserve the data
+ ndatb = this%getNumEnergyBroadened()
+ write(0,*)" I am called !!!!"
+ call this%getData(data)
+ call data%nullify()
+ call data%reserve(ndatb*data%getNnnsig(), this%numPar+1)
+
+ call this%setXoefWeights(2)
+ call this%nullifyWorkGrid(1)
end subroutine
end module HighEnergyFreeGas_m
diff --git a/sammy/src/dbd/mdbd0.f90 b/sammy/src/dbd/mdbd0.f90
index 6ebd13f78b887565399730ea6a0deba37bbf3a55..e85877da0293bcf125d461c9d72a5cc2bb0b4c3a 100644
--- a/sammy/src/dbd/mdbd0.f90
+++ b/sammy/src/dbd/mdbd0.f90
@@ -8,78 +8,34 @@ contains
! *** purpose -- perform High-Energy-Gaussian Approximation version of
! Doppler broadening
!
- use fixedi_m, only : K2reso, Kkkdex, Kkkrsl, Ndatd, Nudwhi, Numorr, Numorr, Numrpi
- use ifwrit_m, only : Kvers7, Ndatb
- use exploc_common_m
- use array_sizes_common_m
- use oopsch_common_m, only : Nowwww, Segmen
- use brdd_common_m, only : Weights
- use lbro_common_m, only : Debug, Yresol, Yssmsc
- use AllocateFunctions_m
- use rsl7_m
- use mxct27_m
+ use fixedi_m, only : numUsedPar
+ use ifwrit_m, only : Kdebug
use Dopplr_m
IMPLICIT None
!
- real(kind=8),allocatable,dimension(:)::A_Iwts
class(HighEnergyFreeGas)::broadener
- integer::Nm, Kdatb
+ integer::Kdatb
!
WRITE (6,99999)
99999 FORMAT (' *** SAMMY-DBD 15 Nov 07 ***')
- Segmen(1) = 'D'
- Segmen(2) = 'B'
- Segmen(3) = 'D'
- Nowwww = 0
!
CALL Initix
+
+ broadener%numPar = numUsedPar
+ broadener%debugOutput = .false.
+ if (Kdebug.ne.0) broadener%debugOutput = .true.
!
- Kdatb = Ndatd
- IF (Kdatb.EQ.0) Kdatb = Ndatb
-!
+ call broadener%broaden()
!
-! *** GUESSTIMATE SIZE OF ARRAY NEEDED FOR SAMDBD
- CALL Estdbd (Kdatb)
-!
- Nm = Kdatb
-! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -<
-! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - <
-! *** Xcoef generates coefficients to be used in broadening
- CALL Xcoef (Nm)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - >
!
- call allocate_real_data(A_Iwts, Kdatb)
! *** Dopplr PERFORMS DOPPLER BROADENING OPERATION (HEGA)
- CALL Dopplr(broadener, A_Idpiso , I_Iflmsc , Weights , A_Iwts)
-!
- deallocate(A_Iwts)
+ CALL Dopplr_broaden(broadener)
+
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ->
!
CALL Write_Commons_Many
RETURN
!
END
-!
-!
-! __________________________________________________________________
-!
- SUBROUTINE Estdbd (Kdatb)
-!
-! *** purpose -- estimate array size for this segment
-!
- IMPLICIT none
- integer::Kdatb
- logical Need_Isotopes
- integer::I, K, Idimen
- external Idimen
-!
-!x CALL Figure_Kws_1 (Kone)
-!x k = Kone + 2*Kdatb
- K = 2*Kdatb
- K = Idimen (K, 1, 'K, 1')
- I = Idimen (K, -1, 'K, -1')
- I = Idimen (0, 0, '0, 0')
-!
- RETURN
- END
end module samdbd_m
diff --git a/sammy/src/dbd/mdbd1.f90 b/sammy/src/dbd/mdbd1.f90
index 10a4d79afb9a26246f9501beb780b5444b979ad2..8a87ae529871ccb98a3a3147a239873325f7e818 100644
--- a/sammy/src/dbd/mdbd1.f90
+++ b/sammy/src/dbd/mdbd1.f90
@@ -1,65 +1,50 @@
module Dopplr_m
use HighEnergyFreeGas_m
-public Dopplr
+use DerivativeHandler_M
+use ifwrit_m, only : Kvtemp, Nolowb
+use fixedr_m, only : Emaxs, Emins
+use broad_common_m, only : Brdlim, Dopple
+use exploc_common_m, only : A_Idpiso, I_Iflmsc
+use, intrinsic :: ISO_C_BINDING
+implicit none
+public Dopplr_broaden
contains
+ SUBROUTINE Dopplr_broaden(calc)
+ class(HighEnergyFreeGas)::calc
+ call Dopplr(calc, A_Idpiso , I_Iflmsc)
+ end SUBROUTINE
!
!
! --------------------------------------------------------------
!
- SUBROUTINE Dopplr (calc, Dopwid, Iflmsc, Weight, Wts)
+ SUBROUTINE Dopplr (calc, Dopwid, Iflmsc)
!
! *** PURPOSE -- FORM DOPPLER-BROADENED CROSS SECTION and derivatives
!
- use fixedi_m, only : Ktruet,Ndats, numUsedPar
- use ifwrit_m, only : Kcros, Kdebug, Kfinit, Kvthck, Ndat, Nolowb, Nonu, Kvtemp
- use fixedr_m, only : Emaxs, Emins, Thick
- use broad_common_m, only : Brdlim, Dopple
- use brdd_common_m, only : Ipnts, Iup, Kc
- use lbro_common_m, only : Yresol, Yssmsc, Ytrans
- use xct2_m
- use mxct27_m
- use broad_common_m
- use Kountd_Dbd_m
- use EndfData_common_m, only : expData
- use AuxGridHelper_M, only : setAuxGridOffset,setAuxGridRowMax, getNumAuxGridPoints
- use SammyGridAccess_M
- use array_sizes_common_m, only : calcData, calcDataSelf
- use convert_to_transmission_m
- use normalize_and_background
- IMPLICIT none
- logical Another_Process_will_Happen
!
class(HighEnergyFreeGas)::calc
- type(SammyGridAccess)::grid, auxGrid
- real(kind=8):: Weight(*), Wts(*), Dopwid(*)
+ real(kind=8):: Dopwid(*)
integer::Iflmsc(*)
real(kind=8)::Two, val
real(kind=8)::Ddo, e1, eKdatb, Elow, EM, Eup
real(kind=8)::Vv, Wdop
- integer::Iffy, Iso, Isox, J, Jwhich, nauxStart, Ipar, N
+ integer::Iffy, Iso, Isox, J, Jwhich, nauxStart, Ipar, N, iposVel
integer::Kkkdat, Kkkkkk, Kkkmin, Now, numEl, nauxMax,insig
integer(C_SIZE_T)::isoC
+ type(DerivativeHandler) :: calcData
+ logical(C_BOOL)::oneExtra
+ integer::Ipnts, Kc, Ndats
! *** Brdlim is now set in segment DAT
DATA Two /2.0d0/
- call grid%initialize()
- call grid%setToExpGrid(expData)
- numEl = grid%getNumEnergies(expData)
- call auxGrid%initialize()
- call auxGrid%setToAuxGrid(expData)
- call setAuxGridOffset(1) ! reset the auxillary grid starting point
-
-!
-!
-! *** FORM DOPPLER-BROADENED CROSS SECTION, STORE IN Wg..
-! *** ALSO DOPPLER-BROADEN THE PARTIAL DERIVATIVES.
-!
- Another_Process_will_Happen = Yresol.OR.Yssmsc
- nauxMax = getNumAuxGridPoints()
+ oneExtra = .false.
+ numEl = calc%getNumEnergyBroadened()
+ nauxMax = calc%getNumEnergyUnbroadened()
!
+ call calc%getdata(calcData)
call calcData%nullify()
- call calcData%reserve(nauxMax*calcData%getNnnsig(),numUsedPar+1)
+ call calcData%reserve(nauxMax*calcData%getNnnsig(),calc%numPar+1)
Now = 0
!
!
@@ -74,8 +59,6 @@ contains
Isox = iso
!
Iffy = 0
- Kc = 1
- Iup = 0
Ddo = Dopple
!
! *** start of major loop over energy-points
@@ -83,7 +66,7 @@ contains
Kkkmin = 0
Kkkdat = 0
nauxStart = 0
- DO 70 J=1,nauxMax
+ DO 70 J=1,numEl
IF ((J/10000)*10000.EQ.J) WRITE (6,10000) J
10000 FORMAT (' *** on data point number', i10)
!q
@@ -95,15 +78,8 @@ contains
!q
!qqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqqq
!
- IF (Another_Process_will_Happen) THEN
- EM = auxGrid%getEnergy(J, expData)
- ELSE
- IF (J.GT.Ndat) THEN
- GO TO 80
- ELSE
- Em = grid%getEnergy(J, expData)
- END IF
- END IF
+ Em = calc%getEnergyBroadened(J) ! desired energy
+ iposVel = calc%getEmInUnbroadened(Em, iposVel) ! find position (iposVel) in un-broadened grid and get the velocity
!
IF (Em.LT.Emins) nauxStart = J
IF (Em.LT.Emins) GO TO 70
@@ -115,101 +91,90 @@ contains
! Wdop=DOPPLER WIDTH AT Energy EM
!
Eup = Em + Brdlim*Wdop
- eKdatb = auxGrid%getEnergy(nauxMax, expData)
+ eKdatb = calc%getEnergyUnbroadened(nauxMax)
IF (Eup.GT.eKdatb) Eup = eKdatb
Elow = Em - Brdlim*Wdop
- e1 = auxGrid%getEnergy(1, expData)
+ e1 = calc%getEnergyUnbroadened(1)
IF (Elow.LT.e1) Elow = e1
!
! ********* find how many and which points are in integral for j
- CALL Kountd_Dbd (Elow, Eup)
- IF (Iup.GT.nauxMax) Iup = nauxMax
+ call calc%calcIntegralSpan(Elow, Eup, oneExtra)
+ Kc = calc%getIntegralStart()
+ Ipnts = calc%getIntegralSpan()
!
Kkkkkk = Kkkkkk + 1
IF ((Nolowb.NE.0 .AND. calc%Elowbr.GT.Em) .OR. Ipnts.le.5) THEN
-! ********* IF (don't want low-E brdng, or) too few points, do not broaden
- CALL Which_Dbd (Em, Ipnts, Jwhich, Kc)
+! ********* IF (don't want low-E brdng, or) too few points, do not broaden
Now = Now + 1
- do ipar = 0, numUsedPar
+ do ipar = 0, calc%numPar
DO N=1, calcData%getNnnsig()
- val = calcData%getDataNsOld(Jwhich, N, Ipar, Isox)
+ val = calcData%getDataNsOld(iposVel, N, Ipar, Isox)
if ( val.eq.0.0d0) cycle
call calcData%addDataNs(Kkkkkk, N, Ipar, Isox, val)
end do
end do
- IF (numUsedPar.gt.0.and.Kvtemp.GT.0) THEN
+ IF (calc%numPar.gt.0.and.Kvtemp.GT.0) THEN
DO N=1, calcData%getNnnsig()
call calcData%addDataNs(Kkkkkk, N, Kvtemp, Isox, 0.0d0)
end do
END IF
ELSE
! ********* Doppler-broaden
- CALL Xdoppl (auxGrid, Weight, Wts,Em, Wdop, Iffy, Isox, Kkkkkk)
+ CALL Xdoppl (calc, Em, Wdop, Iffy, Isox, Kkkkkk)
END IF
!
!
-
- IF (.NOT. Another_Process_will_Happen) THEN
-! ********* eta?
- IF (Kcros.EQ.6) CALL Nnneta (Iflmsc, calcData, Kkkkkk, 1)
- END IF
!
70 CONTINUE
! ****** end of energy do-loop
!
80 CONTINUE
END DO
+
+
! *** End of do loop on nuclides
!
nauxStart = nauxStart + 1
- call setAuxGridRowMax(Kkkdat)
IF (Iffy.EQ.0) WRITE (21,99998)
99998 FORMAT (/,' ** NOTE -- NO DOPPLER BROADENING ACTUALLY OCCURED **')
IF (Now.NE.0) WRITE (21,99997) Now, Kkkdat*calcData%getUsedIsotopes()
- IF (Now.NE.0 .AND. Kdebug.NE.0) WRITE (06,99997) Now, Kkkdat*calcData%getUsedIsotopes()
+ IF (Now.NE.0 .AND. calc%debugOutput) WRITE (06,99997) Now, Kkkdat*calcData%getUsedIsotopes()
99997 FORMAT (' No Doppler broadening occured', I8, ' times of a possible', I8)
+ Ndats = 1
+ Ndats = calc%getEmInUnbroadened(Emins, Ndats)
IF (Kkkdat.LT.Ndats-2) WRITE (6,10100) Kkkdat, Ndats
10100 FORMAT (' Kkkdat, ndats =', 2I10)
!
- call grid%destroy()
- call auxGrid%destroy()
- call dopplerOption%highEnergyFreeGas%updateBroadenedOffset(nauxStart)
- call dopplerOption%highEnergyFreeGas%setLength(Kkkdat)
- call setAuxGridOffset(nauxStart)
+ call calc%updateBroadenedOffset(nauxStart)
+ call calc%setLength(Kkkdat)
RETURN
END
!
!
! --------------------------------------------------------------
!
- SUBROUTINE Xdoppl ( grid, Weight, Wts, EM, Wdop, Iffy, Isox, irow)
+ SUBROUTINE Xdoppl ( calc, EM, Wdop, Iffy, Isox, irow)
!
! *** Perform integration for Doppler-Broadening from point number
! *** Kc to point number Iup
!
- use fixedi_m, only : numUsedPar
- use ifwrit_m, only : Kdebug, Kvtemp
- use fixedr_m, only : Temp
- use brdd_common_m, only : Ipnts, Iup, Kc
- use EndfData_common_m, only : expData
- use SammyGridAccess_M
- use, intrinsic :: ISO_C_BINDING
- use array_sizes_common_m, only : calcData
- IMPLICIT None
- type(SammyGridAccess)::grid
+ class(HighEnergyFreeGas)::calc
+ type(DerivativeHandler) :: calcData
!
- real(kind=8)::Weight(*), Wts(*)
real(kind=8)::Em, Wdop, val
integer::Iffy, Isox, Ipos, irow
!
- real(kind=8)::Sigpls(100), Sigmns(100), Derdop(100)
- real(kind=8)::Zero, One
+ real(kind=8)::Sigpls(100), Sigmns(100), Derdop(100)
real(kind=8)::Del, Derdp1, Derdp2, ewhich, Sigt, W
integer::I, Ikc, Ipar, Iso, Jwhich, K, KS, N, Iparx
logical(C_BOOL)::accu
- DATA Zero /0.0d0/, One /1.0d0/
+ real(kind=8),parameter::Zero=0.0d0, One=1.0d0
+ integer::Kc, Ipnts
+ real(kind=8)::Wts
+
+ call calc%getData(calcData)
IF (calcData%getNnnsig().GT.100) then
STOP '[STOP in Xdoppl in dbd/mdbd1.f]'
end if
@@ -218,14 +183,16 @@ contains
Iffy = Iffy + 1
KS = 1
Iso = Isox
+ Kc = calc%getIntegralStart()
+ Ipnts = calc%getIntegralSpan()
!
Del = 0.02d0
Derdop = 0d0
- IF (numUsedPar.gt.0 .AND. Kvtemp.GT.0) THEN
+ IF (calc%numPar.gt.0 .AND. Kvtemp.GT.0) THEN
W = Wdop*(One+Del)
- CALL Fundop (W, EM, grid, Weight, Wts, Sigpls, Isox)
+ CALL Fundop (calc, W, EM, Sigpls, Isox)
W = Wdop*(One-Del)
- CALL Fundop (W, EM, grid, Weight, Wts, Sigmns, Isox)
+ CALL Fundop (calc, W, EM, Sigmns, Isox)
DO I=1,calcData%getNnnsig()
Derdop(I) = 0.5d0*(Sigpls(I)-Sigmns(I))/(Del*Wdop)
END DO
@@ -233,21 +200,22 @@ contains
!
! *** Simpson rule -- figure Wts = 0.5 * gaussian *
! *** [ Energy(I+1)-Energy(I-1) ]
- IF (Ipnts.LE.12) CALL Simpsn (EM, grid, Wts, Wdop)
+ IF (Ipnts.LE.12) CALL Simpsn (calc, EM, Wdop)
!
! *** 4-point quadrature rule Wts = Weight * gaussian
- IF (Ipnts.GT.12) CALL Dopbrd (EM, grid, Weight, Wts, Wdop)
+ IF (Ipnts.GT.12) CALL Dopbrd (calc, EM, Wdop)
!
accu = .true.
call calcData%setAccumulate(accu)
DO N=1,calcData%getNnnsig()
! update cross section (ipar=0) as well as derivatives
- Do Ipar = 0, numUsedPar
+ Do Ipar = 0, calc%numPar
val = 0.0
DO I=1,Ipnts
- IF (Wts(I).eq.Zero) cycle
+ Wts = calc%getWorkData(1, 1, I)
+ IF (Wts.eq.Zero) cycle
Ikc = I + Kc - 1
- val = val + calcData%getDataNsOld(Ikc,N, Ipar, Iso)*Wts(I)
+ val = val + calcData%getDataNsOld(Ikc,N, Ipar, Iso)*Wts
end do
if (val.eq.0) cycle
call calcData%addDataNs(irow, N, Ipar, Iso, val)
@@ -256,15 +224,15 @@ contains
accu = .false.
call calcData%setAccumulate(accu)
!
- IF (Kvtemp.GT.0.and.numUsedPar.GT.0) THEN
+ IF (Kvtemp.GT.0.and.calc%numPar.GT.0) THEN
DO N=1,calcData%getNnnsig()
- val = Derdop(N)*Wdop*0.5d0/Temp
+ val = Derdop(N)*Wdop*0.5d0/calc%getTemperature()
call calcData%addDataNs(irow, N, Kvtemp, Iso, val)
END DO
END IF
!
- IF (Kdebug.EQ.0 .OR. Kvtemp.LE.0) RETURN
- IF (numUsedPar.EQ.0) RETURN
+ IF (.not.calc%debugOutput.OR. Kvtemp.LE.0) RETURN
+ IF (calc%numPar.EQ.0) RETURN
!
! *** here for debug info only
Sigt = calcData%getDataNs(irow, 1, 0, Iso)
@@ -279,8 +247,7 @@ contains
RETURN
!
60 CONTINUE
- Jwhich = (Kc+Iup)/2
- ewhich = grid%getEnergy(Jwhich, expData)
+ ewhich = em
WRITE (6,99999) ewhich, Derdp1, Derdop(1), Derdp2
99999 FORMAT (' ** WARNING ON Dopplr DERIV--E,D1,D,D2=', F8.1, &
3(1PG11.3))
@@ -290,37 +257,39 @@ contains
!
! --------------------------------------------------------------
!
- SUBROUTINE Dopbrd (Em, grid, Weight, Wts, Widgau)
+ SUBROUTINE Dopbrd (calc, Em, Widgau)
!
! *** FORM THE GAUSSIAN DOPPLER WeightS AND NORMALIZE THEM
!
- use brdd_common_m, only : Ipnts, Kc
- use EndfData_common_m, only : expData
- use SammyGridAccess_M
- IMPLICIT None
-!
- real(kind=8)::Weight(*), Wts(*), Em, Widgau
- type(SammyGridAccess)::grid ! this is possible in fixed form g77 as in same compile unit
+ real(kind=8):: Em, Widgau
real(kind=8)::e1, S, Y, Z
integer::I
+ class(HighEnergyFreeGas)::calc
+ integer::Ipnts, Kc
+ real(kind=8)::Wts
+
!
Y = 85.0D0
+ Ipnts = calc%getIntegralSpan()
+ Kc = calc%getIntegralStart()
DO I=1,Ipnts
- e1 = grid%getEnergy(I+Kc-1, expData)
+ e1 = calc%getEnergyUnbroadened(I+Kc-1)
Z = (e1-Em)/Widgau
Z = Z*Z
IF (Z.GT.Y) Z = Y
Z = dEXP(-Z)
- Wts(I) = Weight(I+kc-1)*Z
+ Wts = calc%getWorkData(2,1, I+kc-1)*Z
+ call calc%setWorkData(1, 1, I, Wts)
END DO
!
S = 0.0D0
DO I=1,Ipnts
- S = S + Wts(I)
+ S = S + calc%getWorkData(1,1, I)
END DO
S = 1.0D0/S
DO I=1,Ipnts
- Wts(I) = Wts(I)*S
+ Wts = calc%getWorkData(1,1, I)*S
+ call calc%setWorkData(1, 1, I, Wts)
END DO
!
RETURN
@@ -329,59 +298,59 @@ contains
!
! --------------------------------------------------------------
!
- SUBROUTINE Simpsn (Em, grid, Wts, Widgau)
+ SUBROUTINE Simpsn (calc, Em, Widgau)
!
! *** Form the Gaussian Doppler Weights via Simpson's rule, and
! *** normalize them
!
- use brdd_common_m, only : Ipnts, Kc
- use EndfData_common_m, only : expData
- use AuxGridHelper_M, only : getNumAuxGridPoints
- use SammyGridAccess_M
- IMPLICIT None
!
- real(kind=8)::Wts(*)
+ real(kind=8)::Wts
real(kind=8)::Em, Widgau
- type(SammyGridAccess)::grid ! this is possible in fixed form g77 as in same compile unit
+ class(HighEnergyFreeGas)::calc
real(kind=8)::eval1, eval2, eval3, S, Y, Z
- integer::I, J, nauxMax
+ integer::I, J, nauxMax, Ipnts, Kc
!
+ Ipnts = calc%getIntegralSpan()
+ Kc = calc%getIntegralStart()
J = Ipnts
- nauxMax = getNumAuxGridPoints()
+ nauxMax = calc%getNumEnergyUnbroadened()
IF (Kc+Ipnts.GT.nauxMax) J = Ipnts - 1
Ipnts = J
!
Y = 85.0D0
I = 1
- eval1 = grid%getEnergy(I+Kc-1, expData)
+ eval1 = calc%getEnergyUnbroadened(I+Kc-1)
Z = (eval1-Em)/Widgau
Z = Z*Z
IF (Z.GT.Y) Z = Y
Z = dEXP(-Z)
- eval2 = grid%getEnergy(I+Kc, expData)
- Wts(I) = (eval2-eval1)*Z
+ eval2 = calc%getEnergyUnbroadened(I+Kc)
+ Wts = (eval2-eval1)*Z
+ call calc%setWorkData(1,1, I, Wts)
! *** note that we must treat i=1 case separately to avoid
! *** trying to use index at 0
DO I=2,Ipnts
- eval1 = grid%getEnergy(I+Kc-1, expData)
+ eval1 = calc%getEnergyUnbroadened(I+Kc-1)
Z = (eval1-Em)/Widgau
Z = Z*Z
IF (Z.GT.Y) Z = Y
Z = dEXP(-Z)
- eval2 = grid%getEnergy(I+Kc, expData)
- eval3 = grid%getEnergy(I+Kc - 2, expData)
- Wts(I) = (eval2-eval3)*Z
+ eval2 = calc%getEnergyUnbroadened(I+Kc)
+ eval3 = calc%getEnergyUnbroadened(I+Kc - 2)
+ Wts = (eval2-eval3)*Z
+ call calc%setWorkData(1,1, I, Wts)
! ****** note there should be a 1/2 factor, but we'll normalize anyway
END DO
!
S = 0.0D0
DO I=1,Ipnts
- S = S + Wts(I)
+ S = S + calc%getWorkData(1,1,I)
END DO
S = 1.0D0/S
DO I=1,Ipnts
- Wts(I) = Wts(I)*S
+ Wts = calc%getWorkData(1,1,I)*S
+ call calc%setWorkData(1,1, I, Wts)
END DO
!
RETURN
@@ -390,29 +359,30 @@ contains
!
! --------------------------------------------------------------
!
- SUBROUTINE Fundop (Wdop, Em, grid, Weight, Wts, Sum, Isox)
+ SUBROUTINE Fundop (calc, Wdop, Em, Sum, Isox)
!
- use brdd_common_m, only : Ipnts, Kc
- use EndfData_common_m, only : expData
- use SammyGridAccess_M
- use array_sizes_common_m, only : calcData
- IMPLICIT None
- type(SammyGridAccess)::grid ! this is possible in fixed form g77 as in same compile unit
- real(kind=8):: Wts(*), Sum(*), weight(*), Wdop, Em
+ type(DerivativeHandler) :: calcData
+ real(kind=8):: Sum(*), Wdop, Em
integer::isox
- real(kind=8)::Zero, val
- integer::I, Ikc, N
- DATA Zero /0.0d0/
-!
- IF (Ipnts.LE.12) CALL Simpsn (Em, grid, Wts, Wdop)
- IF (Ipnts.GT.12) CALL Dopbrd (Em, grid, Weight, Wts, Wdop)
- CALL Zero_Array (Sum, calcData%getNnnsig())
+ real(kind=8)::val
+ integer::I, Ikc, N, Ipnts, Kc
+ class(HighEnergyFreeGas)::calc
+ real(kind=8),parameter :: Zero = 0.0d0
+ real(kind=8)::Wts
+!
+ Ipnts = calc%getIntegralSpan()
+ Kc = calc%getIntegralStart()
+ call calc%getdata(calcData)
+ IF (Ipnts.LE.12) CALL Simpsn (calc, Em, Wdop)
+ IF (Ipnts.GT.12) CALL Dopbrd (calc, Em, Wdop)
+ Sum(1:calcData%getNnnsig()) = Zero
DO I=1,Ipnts
- IF (Wts(I).NE.Zero) THEN
+ Wts = calc%getWorkData(1,1,I)
+ IF (Wts.NE.Zero) THEN
Ikc = I + Kc - 1
DO N=1,calcData%getNnnsig()
val = calcData%getDataNsOld(Ikc, N, 0, Isox)
- Sum(N) = Sum(N) + Wts(I)*val
+ Sum(N) = Sum(N) + Wts*val
END DO
END IF
END DO
diff --git a/sammy/src/dbd/mdbd3.f90 b/sammy/src/dbd/mdbd3.f90
deleted file mode 100644
index 9ee969a84c52e43924e3a0dee523ade4e552dbff..0000000000000000000000000000000000000000
--- a/sammy/src/dbd/mdbd3.f90
+++ /dev/null
@@ -1,71 +0,0 @@
-module Kountd_Dbd_m
-public Kountd_Dbd, Which_Dbd
-contains
-!
-!
-! --------------------------------------------------------------
-!
- SUBROUTINE Kountd_Dbd (Elow, Eup)
-!
-! *** Purpose -- FIND Kc AND Iup, ALSO Ipnts
-!
- use brdd_common_m, only : Ipnts, Iup, Kc
- use EndfData_common_m, only : expData
- use SammyGridAccess_M
- IMPLICIT none
- type(SammyGridAccess)::grid
- real(kind=8)::Elow, Eup
- real(kind=8)::elowGrid, eupGrid
- integer::Ilow
-
- call grid%initialize()
- call grid%setToAuxGrid(expData)
-
- Ilow = Kc - 1
- 10 Ilow = Ilow + 1
- elowGrid = grid%getEnergy(Ilow, expData)
- IF (elowGrid.LE.Elow) GO TO 10
- 20 Iup = Iup + 1
- eupGrid = grid%getEnergy(Iup, expData)
- IF (eupGrid.LT.Eup) GO TO 20
- Ilow = Ilow - 2
- IF (Ilow.LT.0) Ilow = 0
- Ipnts = Iup - Ilow
- Kc = Ilow + 1
- Iup = Iup - 1
- call grid%destroy()
- RETURN
- END
-!
-!
-! --------------------------------------------------------------
-!
- SUBROUTINE Which_Dbd (EM, Ipnts, Jwhich, Kc)
- use EndfData_common_m, only : expData
- use SammyGridAccess_M
- IMPLICIT None
- type(SammyGridAccess)::grid
- integer::Ipnts, Jwhich, Kc
- real(kind=8)::EM
- real(kind=8)::A, Ddd, e1, e2
- integer::I
-
- call grid%initialize()
- call grid%setToAuxGrid(expData)
-!
- Jwhich = 0
- e1 = grid%getEnergy(Ipnts+Kc-1, expData)
- e2 = grid%getEnergy(kc, expData)
- Ddd = e1 - e2
- DO I=1,Ipnts
- e1 = grid%getEnergy(I+Kc-1, expData)
- A = dABS(e1-EM)
- IF (A.LE.Ddd) Jwhich = I
- IF (A.LE.Ddd) Ddd = A
- END DO
- IF (Jwhich.EQ.0) Jwhich = 1
- Jwhich = Kc + Jwhich - 1
- call grid%destroy()
- RETURN
- END
-end module Kountd_Dbd_m
diff --git a/sammy/src/dop/LealHwangBroadening_M.f90 b/sammy/src/dop/LealHwangBroadening_M.f90
new file mode 100644
index 0000000000000000000000000000000000000000..e9c4279309ce3dc5df2123251e89e47b0aacb6a7
--- /dev/null
+++ b/sammy/src/dop/LealHwangBroadening_M.f90
@@ -0,0 +1,49 @@
+module LealHwangBroadening_M
+use FreeGasDopplerBroadening_M
+use DopplerBroadening_M
+use DerivativeHandler_M
+use GridData_M
+use, intrinsic :: ISO_C_BINDING
+
+
+implicit none
+
+type, extends(FreeGasDopplerBroadening) :: LealHwangBroadening
+contains
+procedure, pass(this) :: initialize => LealHwangBroadening_initialize
+procedure, pass(this) :: broaden => LealHwangBroadening_broaden
+procedure, pass(this) :: destroy => LealHwangBroadening_destroy
+end type LealHwangBroadening
+
+contains
+subroutine LealHwangBroadening_initialize(this, hand, list, work)
+ implicit none
+ class(LealHwangBroadening) :: this
+ class(DerivativeHandler)::hand
+ class(GridDataList)::list
+ class(GridDataList)::work
+ call FreeGasDopplerBroadening_initialize(this, hand, list, work)
+end subroutine
+
+subroutine LealHwangBroadening_destroy(this)
+ implicit none
+ class(LealHwangBroadening) :: this
+ call FreeGasDopplerBroadening_destroy(this)
+end subroutine
+
+
+
+subroutine LealHwangBroadening_broaden(this)
+ class(LealHwangBroadening) :: this
+
+ integer::ndatb
+ type(DerivativeHandler)::data
+
+ ! reserve the data
+ ndatb = this%getNumEnergyBroadened()
+ call this%getData(data)
+ call data%nullify()
+ call data%reserve(ndatb*data%getNnnsig(), this%numPar+1)
+end subroutine
+
+end module LealHwangBroadening_M
diff --git a/sammy/src/fgm/mfgm1.f90 b/sammy/src/fgm/mfgm1.f90
index cea3f5bfae3d21cffedf73ee437acd66e30bdc94..b976d457ea820e88b29eaeb21a6f091914626855 100644
--- a/sammy/src/fgm/mfgm1.f90
+++ b/sammy/src/fgm/mfgm1.f90
@@ -1,4 +1,14 @@
module mfgm1_m
+use ifwrit_m, only : Ksitmp, Nolowb, Kvtemp
+use fixedr_m, only : Emax, Emin, Emins, Sitemp, Emaxs
+use EndfData_common_m, only : resparData
+use mfgm3_M
+use mfgm4_m
+use FreeGasDopplerBroadening_M
+IMPLICIT none
+
+public Dopfgm
+
contains
!
!
@@ -8,14 +18,6 @@ contains
!
! *** PURPOSE -- FORM DOPPLER-BROADENED CROSS SECTION AND DERIVATIVES
!
- use ifwrit_m, only : Ksitmp, Nolowb, Kvtemp
- use fixedr_m, only : Emax, Emin, Emins, Sitemp, Emaxs
- use EndfData_common_m, only : resparData
- use mfgm3_M
- use mfgm4_m
- use FreeGasDopplerBroadening_M
-
- IMPLICIT None
type(FreeGasDopplerBroadening)::broadener
type(DerivativeHandler)::calcData
diff --git a/sammy/src/fgm/mfgm4.f90 b/sammy/src/fgm/mfgm4.f90
index 8ceca181e1c13854cad46da4b44538db1a33ca33..58cac2ccca933101213828622891e7235b745444 100644
--- a/sammy/src/fgm/mfgm4.f90
+++ b/sammy/src/fgm/mfgm4.f90
@@ -3,6 +3,7 @@ module mfgm4_m
use DerivativeHandler_M
use FreeGasDopplerBroadening_M
use, intrinsic :: ISO_C_BINDING
+ IMPLICIT None
public Xdofgm
@@ -18,7 +19,6 @@ module mfgm4_m
! *** PERFORM INTEGRATION FOR DOPPLER-BROADENING FROM POINT NUMBER
! *** calc%Kc TO POINT NUMBER Iup
!
- IMPLICIT None
integer::Iffy, Ktempx, Isox, Ngtvx
real(kind=8)::Velcty(*), Wts(*), Em, Vv, Ddo, Tempx
@@ -180,7 +180,6 @@ module mfgm4_m
!
SUBROUTINE Funfgm (calc, Ddo, Vv, Velcty, Wts, Sum, Isox, derivs)
!
- IMPLICIT None
real(kind=8)::Velcty(*), Wts(*), Sum(*)
real(kind=8)::Ddo, Vv
integer::Isox, Nb
diff --git a/sammy/src/sammy/CMakeLists.txt b/sammy/src/sammy/CMakeLists.txt
index 948e59bf5af9fc97df2e278f38b3d61926b1c44f..109b3532b6d213685806655aa6a93298aff047fb 100644
--- a/sammy/src/sammy/CMakeLists.txt
+++ b/sammy/src/sammy/CMakeLists.txt
@@ -60,6 +60,7 @@ APPEND_SET(SAMMY_SOURCES
../clm/mclm7.f
../clm/mclm8.f
../clm/mclm8a.f90
+ ../clm/CrystalLatticeBroadening_M.f90
../clq/mclq0.f90
../clq/ArtificalCross_M.f90
@@ -94,7 +95,6 @@ APPEND_SET(SAMMY_SOURCES
../dbd/mdbd0.f90
../dbd/mdbd1.f90
../dbd/mdbd2.f90
- ../dbd/mdbd3.f90
../dbd/HighEnergyFreeGas_m.f90
../dex/mdex0.f
@@ -103,6 +103,7 @@ APPEND_SET(SAMMY_SOURCES
../dop/mdop0.f90
../dop/mdop1.f90
../dop/mdop2.f90
+ ../dop/LealHwangBroadening_M.f90
../end/mend0.f
../end/mend1.f