Skip to content

Lattice offset and disable offset #164

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
ChasingNeutrons wants to merge 7 commits into
base: main
Choose a base branch
from
Open

Lattice offset and disable offset #164

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
83df212
Fix to keff for MG scattering
ChasingNeutrons Jan 22, 2025
0791fe2
BEAVRS_HZP and geometry fixes
ChasingNeutrons Feb 6, 2025
8e51a1c
Merge branch 'CambridgeNuclear:main' into develop
ChasingNeutrons Feb 6, 2025
f20df29
Allow offset to be disabled in nested universes
ChasingNeutrons Jun 23, 2025
fe0f845
Add docs, small fix, and added example input
ChasingNeutrons Jun 24, 2025
3db8907
Added flux tallies to BEAVRS HZP
ChasingNeutrons Jun 24, 2025
387e6ce
Update Geometry/Universes/latUniverse_class.f90
ChasingNeutrons Aug 12, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
144 changes: 141 additions & 3 deletions DataStructures/dictionary_class.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,8 @@ module dictionary_class
! -> 1D array of shortInt - integer(shortInt) :: i(:)
! -> scalar character of length "charLen" defined below - character(charLen) :: c
! -> array of characters of length "charLen" defined below - character(charLen) :: c(:)
! -> boolean logical - logical(defBool) :: b
! -> array of boolean logicals - logical(defBool) :: b(:)
! -> another dictionary - type(dictionary) :: dict
!
! To add additional structure <type> to store it is necessary to :
Expand Down Expand Up @@ -100,7 +102,7 @@ module dictionary_class
! *** Note that dictionary is defined as pointer not allocatable
! *** This is because gfortran < 7.0 does not supports circular derived types with
! *** allocatable keyword. This line may change in the future
type(dictionary), pointer :: dict0_alloc => null()
type(dictionary), pointer :: dict0_alloc => null()

! dictContent type ID
integer(shortInt) :: type = empty
Expand Down Expand Up @@ -149,7 +151,9 @@ module dictionary_class
getCharArray_alloc_new,&
getCharArray_ptr_new,&
getDict_new,&
getBool_new
getBool_new, &
getBoolArray_alloc_new,&
getBoolArray_ptr_new

procedure,private :: getReal_new
procedure,private :: getRealArray_alloc_new
Expand All @@ -162,6 +166,8 @@ module dictionary_class
procedure,private :: getCharArray_ptr_new
procedure,private :: getDict_new
procedure,private :: getBool_new
procedure,private :: getBoolArray_alloc_new
procedure,private :: getBoolArray_ptr_new

generic :: getOrDefault => getOrDefault_real ,&
getOrDefault_realArray_alloc ,&
Expand All @@ -172,7 +178,9 @@ module dictionary_class
getOrDefault_char ,&
getOrDefault_charArray_alloc ,&
getOrDefault_charArray_ptr, &
getOrDefault_bool
getOrDefault_bool, &
getOrDefault_boolArray_alloc, &
getOrDefault_boolArray_ptr

procedure,private :: getOrDefault_real
procedure,private :: getOrDefault_realArray_alloc
Expand All @@ -184,6 +192,8 @@ module dictionary_class
procedure,private :: getOrDefault_charArray_alloc
procedure,private :: getOrDefault_charArray_ptr
procedure,private :: getOrDefault_bool
procedure,private :: getOrDefault_boolArray_alloc
procedure,private :: getOrDefault_boolArray_ptr

! Keys inquiry procedures
procedure :: keys
Expand Down Expand Up @@ -779,6 +789,83 @@ subroutine getBool_new(self,value,keyword)
call fatalError(Here,'Entry under keyword ' // keyword // ' is not an integer')
end select
end subroutine getBool_new

!!
!! Loads a boolean rank 1 from dictionary into provided variable
!! If keyword is associated with an integer which is not 1 or 0, it returns an error.
!! Variable needs to be allocatable. It will be deallocated before assignment
!!
subroutine getBoolArray_alloc_new(self,value,keyword)
class(dictionary), intent(in) :: self
logical(defBool),dimension(:),allocatable,intent(inout) :: value
character(*),intent(in) :: keyword
integer(shortInt) :: idx
integer(shortInt),dimension(:),allocatable :: tempInt
character(100),parameter :: Here='getBoolArray_alloc_new (dictionary_class.f90)'

idx = self % search(keyword, Here, fatal =.true.)

if(allocated(value)) deallocate(value)

select case (self % entries(idx) % getType())
case(arrInt)
tempInt = self % entries(idx) % int1_alloc
allocate(value(size(tempInt)))

if (any(tempInt < 0) .or. any(tempInt > 1)) then
call fatalError(Here,'Entry under keyword ' // keyword // ' has non-logical values (not 0 or 1).')
else
do idx = 1, size(tempInt)
value(idx) = tempInt(idx) == 1
end do
end if

case default
call fatalError(Here,'Entry under keyword ' // keyword // ' is not an int array.')

end select

end subroutine getBoolArray_alloc_new

!!
!! Loads a boolean rank 1 from dictionary into provided variable.
!! If keyword is associated with an integer which is not 1 or 0, it returns an error.
!! Variable needs to be pointer. It will be deallocated before assignment
!!
subroutine getBoolArray_ptr_new(self,value,keyword)
class(dictionary), intent(in) :: self
logical(defBool),dimension(:),pointer,intent(inout) :: value
character(*),intent(in) :: keyword
integer(shortInt) :: idx, N
integer(shortInt),dimension(:),allocatable :: tempInt
character(100),parameter :: Here='getBoolArray_ptr_new (dictionary_class.f90)'

idx = self % search(keyword, Here, fatal =.true.)

if(associated(value)) deallocate(value)

select case (self % entries(idx) % getType())
case(arrInt)

N = size(self % entries(idx) % int1_alloc)
allocate(value(N))
allocate(tempInt(N))
tempInt = self % entries(idx) % int1_alloc

if (any(tempInt < 0) .or. any(tempInt > 1)) then
call fatalError(Here,'Entry under keyword ' // keyword // ' has non-logical values (not 0 or 1).')
else
do idx = 1, N
value(idx) = tempInt(idx) == 1
Copy link
Member

@valeriaRaffuzzi valeriaRaffuzzi Aug 11, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should this be a pointer?

Copy link
Collaborator Author

@ChasingNeutrons ChasingNeutrons Aug 12, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think so because tempInt will be destroyed at the end of the subroutine, but I might be doing other bad things here. My Fortran isn't good enough for this... I think in other cases we equate value (pointer) with an allocatable array which I think is allowable? I might be doing something bad in this case. Wisdom from @Mikolaj-A-Kowalski would be appreciated...

end do
end if

case default
call fatalError(Here,'Entry under keyword ' // keyword // ' is not an int array.')

end select

end subroutine getBoolArray_ptr_new

!!
!! Loads a real rank 0 from a dictionary into provided variable
Expand Down Expand Up @@ -1044,6 +1131,57 @@ subroutine getOrDefault_bool(self,value,keyword,default)
end if
end subroutine getOrDefault_bool

!!
!! Loads a boolean rank 1 from a dictionary.
!! Variable needs to be allocatable. It will be deallocated before assignment
!!
subroutine getOrDefault_boolArray_alloc(self, value, keyword, default)
class(dictionary), intent(in) :: self
logical(defBool),dimension(:),allocatable,intent(inout) :: value
character(*),intent(in) :: keyword
logical(defBool),dimension(:),intent(in) :: default
integer(shortInt) :: idx
character(100),parameter :: Here='getOrDefault_boolArray_alloc (dictionary_class.f90)'

idx = self % search(keyword, Here, fatal =.false.)

if(allocated(value)) deallocate(value)

if (idx == targetNotFound) then
value = default
else
call self % get(value, keyword)
end if
end subroutine getOrDefault_boolArray_alloc


!!
!! Loads a boolean rank 1 from a dictionary.
!! Variable needs to be pointer. It will be deallocated before assignment
!!
!! For further details refer to doc of getOrDefault_real
!!
subroutine getOrDefault_boolArray_ptr(self, value, keyword, default)
class(dictionary), intent(in) :: self
logical(defBool),dimension(:),pointer,intent(inout) :: value
character(*),intent(in) :: keyword
logical(defBool),dimension(:),intent(in) :: default
integer(shortInt) :: idx
character(100),parameter :: Here='getOrDefault_boolArray_ptr (dictionary_class.f90)'

idx = self % search(keyword, Here, fatal =.false.)

if(associated(value)) deallocate(value)

if (idx == targetNotFound) then
allocate(value(size(default)))
value = default
else
call self % get(value, keyword)
end if

end subroutine getOrDefault_boolArray_ptr

!!
!! Returns an array of all keywords
!!
Expand Down
4 changes: 2 additions & 2 deletions Geometry/Surfaces/box_class.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -347,7 +347,7 @@ subroutine explicitBC(self, r, u)
r0 = r(ax) - self % origin(ax)

! Skip if particle is well inside the domain
if (abs(r0) <= self % halfwidth(ax) - self % surfTol()) cycle axis
if (abs(r0) <= self % halfwidth(ax) * (ONE - self % surfTol())) cycle axis

! Choose correct BC
if (r0 < ZERO) then
Expand Down Expand Up @@ -396,7 +396,7 @@ subroutine transformBC(self, r, u)

! Calculate halfwidth reduced by the surface_tolerance
! Necessary to capture particles at the boundary
a_bar = self % halfwidth - self % surfTol()
a_bar = self % halfwidth * (ONE - self % surfTol())

! Calculate distance (in # of transformations) in each direction
Ri = ceiling(abs(r - self % origin) / a_bar) / 2
Expand Down
4 changes: 2 additions & 2 deletions Geometry/Surfaces/squareCylinder_class.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -402,7 +402,7 @@ subroutine explicitBC(self, r, u)
r0 = r(ax) - self % origin(i)

! Skip if particle is well inside the domain
if (abs(r0) <= self % halfwidth(i) - self % surfTol()) cycle axis
if (abs(r0) <= self % halfwidth(i) * (ONE- self % surfTol())) cycle axis

! Choose correct BC
if (r0 < ZERO) then
Expand Down Expand Up @@ -452,7 +452,7 @@ subroutine transformBC(self, r, u)

! Calculate halfwidth reduced by the surface_tolerance
! Necessary to capture particles at the boundary
a_bar = self % halfwidth - self % surfTol()
a_bar = self % halfwidth * (ONE - self % surfTol())

! Calculate distance (in # of transformations) in each direction
Ri = ceiling(abs(r(self % plane) - self % origin) / a_bar) / 2
Expand Down
17 changes: 16 additions & 1 deletion Geometry/Universes/Tests/latUniverse_test.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -18,11 +18,15 @@ module latUniverse_test
character(*), parameter :: UNI1_DEF = &
"id 1; type latUniverse; origin (0.0 0.0 0.0); rotation (0.0 0.0 0.0); &
&pitch (1.0 2.0 3.0); shape (3 2 2); padMat void; &
&map ( 3 4 5 &
& map ( 3 4 5 &
& 7 4 8 &
& &
& 1 2 3 &
& 4 5 6); "
!& offsetMap ( 1 1 1 &
Copy link
Member

@valeriaRaffuzzi valeriaRaffuzzi Aug 11, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could this be confused with a boolArray?
This might happen also to a real input.

!& 0 1 1 &
!& 1 1 1 &
!& 1 1 1 ); &

character(*), parameter :: UNI2_DEF = &
"id 2; type latUniverse; pitch (1.0 2.0 0.0); shape (2 1 0); padMat u<1>; &
Expand Down Expand Up @@ -389,6 +393,17 @@ subroutine test_cellOffset()

ref = [0.0_defReal, 1.0_defReal, 1.5_defReal]
@assertEqual(ref, uni1 % cellOffset(pos), TOL)

!! Inside but at the position without an offset
!pos % r = [-1.5_defReal, 1.0_defReal, -0.5_defReal]
!pos % dir = [-ONE, ONE, -ONE]
!pos % dir = pos % dir / norm2(pos % dir)
!pos % uniIdx = 8
!pos % cellIdx = 0
!pos % localId = 4

!ref = [0.0_defReal, 0.0_defReal, 0.0_defReal]
!@assertEqual(ref, uni1 % cellOffset(pos), TOL)

! Outside
pos % r = [-7.0_defReal, 0.0_defReal, 0.5_defReal]
Expand Down
Loading