Skip to content

Commit

Permalink
Update ActuatorLine to master state. Remove bounded procedures of AL.
Browse files Browse the repository at this point in the history
Fix LES comilation error in gfortran.
  • Loading branch information
oscarmarino committed Sep 18, 2024
1 parent d9ffb94 commit 844b4a3
Show file tree
Hide file tree
Showing 5 changed files with 566 additions and 496 deletions.
13 changes: 9 additions & 4 deletions Solver/src/NavierStokesSolver/SpatialDiscretization.f90
Original file line number Diff line number Diff line change
Expand Up @@ -382,7 +382,7 @@ END SUBROUTINE ComputeTimeDerivativeIsolated
subroutine TimeDerivative_ComputeQDot( mesh , particles, t)
use WallFunctionConnectivity
use TripForceClass, only: randomTrip
use ActuatorLine, only: farm
use ActuatorLine, only: farm, ForcesFarm
use SpongeClass, only: sponge
implicit none
type(HexMesh) :: mesh
Expand Down Expand Up @@ -586,16 +586,20 @@ subroutine TimeDerivative_ComputeQDot( mesh , particles, t)
!
! Add physical source term
! ************************
!!! $acc parallel loop gang
!$omp do schedule(runtime) private(i,j,k)
do eID = 1, mesh % no_of_elements
associate ( e => mesh % elements(eID) )
! the source term is reset to 0 each time Qdot is calculated to enable the possibility to add source terms to
! different contributions and not accumulate each call
e % storage % S_NS = 0.0_RP
do k = 0, e % Nxyz(3) ; do j = 0, e % Nxyz(2) ; do i = 0, e % Nxyz(1)
!!!$acc loop vector collapse(3)
do k = 0, e % Nxyz(1) ; do j = 0, e % Nxyz(2) ; do i = 0, e % Nxyz(1)

call UserDefinedSourceTermNS(e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), t, e % storage % S_NS(:,i,j,k), thermodynamics, dimensionless, refValues)
call randomTrip % getTripSource( e % geom % x(:,i,j,k), e % storage % S_NS(:,i,j,k) )
call farm % ForcesFarm(e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), e % storage % S_NS(:,i,j,k), t)
! call farm % ForcesFarm(e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), e % storage % S_NS(:,i,j,k), t)
call ForcesFarm(farm, e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), e % storage % S_NS(:,i,j,k), t)
end do ; end do ; end do
end associate
end do
Expand Down Expand Up @@ -824,7 +828,8 @@ subroutine TimeDerivative_ComputeQDotIsolated( mesh , t )
do k = 0, e % Nxyz(3) ; do j = 0, e % Nxyz(2) ; do i = 0, e % Nxyz(1)
call UserDefinedSourceTermNS(e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), t, e % storage % S_NS(:,i,j,k), thermodynamics, dimensionless, refValues)
call randomTrip % getTripSource( e % geom % x(:,i,j,k), e % storage % S_NS(:,i,j,k) )
call farm % ForcesFarm(e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), e % storage % S_NS(:,i,j,k), t)
! call farm % ForcesFarm(e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), e % storage % S_NS(:,i,j,k), t)
call ForcesFarm(farm, e % geom % x(:,i,j,k), e % storage % Q(:,i,j,k), e % storage % S_NS(:,i,j,k), t)
end do ; end do ; end do
end associate
end do
Expand Down
6 changes: 3 additions & 3 deletions Solver/src/libs/mesh/HexMesh.f90
Original file line number Diff line number Diff line change
Expand Up @@ -891,7 +891,7 @@ subroutine HexMesh_ProlongSolutionToFaces(self, nEqn)
integer :: fIDs(6)
integer :: eID

!$omp do schedule(runtime)
!!$omp do schedule(runtime)
!!$acc parallel loop gang present(self) private(fIDs)
!do eID = 1, size(self % elements)
! fIDs = self % elements(eID) % faceIDs
Expand All @@ -904,7 +904,7 @@ subroutine HexMesh_ProlongSolutionToFaces(self, nEqn)
! self % faces(fIDs(6)) )
!end do
!!$acc end parallel loop
!$omp end do
!!$omp end do

end subroutine HexMesh_ProlongSolutionToFaces

Expand Down Expand Up @@ -4712,4 +4712,4 @@ subroutine HexMesh_Assign (to, from)


end subroutine HexMesh_Assign
END MODULE HexMeshClass
END MODULE HexMeshClass
102 changes: 54 additions & 48 deletions Solver/src/libs/physics/navierstokes/LESModels.f90
Original file line number Diff line number Diff line change
Expand Up @@ -34,35 +34,36 @@ module LESModels
logical :: requiresWallDistances = .FALSE.
integer :: WallModel
integer :: LESModel
real(kind=RP) :: C
contains
procedure :: Initialize => LESModel_Initialize
procedure, private :: ComputeWallEffect => LESModel_ComputeWallEffect
procedure :: ComputeViscosity => LESModel_ComputeViscosity
! procedure, private :: ComputeWallEffect => LESModel_ComputeWallEffect
! procedure :: ComputeViscosity => LESModel_ComputeViscosity
procedure :: Describe => LESModel_Describe
end type LESModel_t

type, extends(LESModel_t) :: Smagorinsky_t
real(kind=RP) :: CS
! real(kind=RP) :: CS
contains
procedure :: Initialize => Smagorinsky_Initialize
procedure :: Describe => Smagorinsky_Describe
procedure :: ComputeViscosity => Smagorinsky_ComputeViscosity
! procedure :: ComputeViscosity => Smagorinsky_ComputeViscosity
end type Smagorinsky_t

type, extends(LESModel_t) :: WALE_t
real(kind=RP) :: Cw
! real(kind=RP) :: Cw
contains
procedure :: Initialize => WALE_Initialize
procedure :: Describe => WALE_Describe
procedure :: ComputeViscosity => WALE_ComputeViscosity
! procedure :: ComputeViscosity => WALE_ComputeViscosity
end type WALE_t

type, extends(LESModel_t) :: Vreman_t
real(kind=RP) :: C
! real(kind=RP) :: C
contains
procedure :: Initialize => Vreman_Initialize
procedure :: Describe => Vreman_Describe
procedure :: ComputeViscosity => Vreman_ComputeViscosity
! procedure :: ComputeViscosity => Vreman_ComputeViscosity
end type Vreman_t


Expand Down Expand Up @@ -188,14 +189,15 @@ subroutine LESModel_Initialize(self, controlVariables)

end subroutine LESModel_Initialize

pure real(kind=RP) function LESModel_ComputeWallEffect (self,LS,dWall)
pure real(kind=RP) function LESModel_ComputeWallEffect (LS,dWall,WallModel)
!$acc routine seq
implicit none
class(LESModel_t), intent(in) :: self
! class(LESModel_t), intent(in) :: self
real(kind=RP) , intent(in) :: LS
real(kind=RP) , intent(in) :: dWall
integer , intent(in) :: WallModel

select case (self % WallModel)
select case (WallModel)
case (LINEAR_WALLMODEL)
LESModel_ComputeWallEffect = min(LS, dWall * K_VONKARMAN)
case (NO_WALLMODEL)
Expand All @@ -205,20 +207,20 @@ pure real(kind=RP) function LESModel_ComputeWallEffect (self,LS,dWall)

end function LESModel_ComputeWallEffect

pure subroutine LESModel_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
!$acc routine seq
implicit none
!-arguments---------------------------------------------
class(LESModel_t), intent(in) :: this
real(kind=RP), intent(in) :: delta
real(kind=RP), intent(in) :: dWall
real(kind=RP), intent(in) :: Q(NCONS)
real(kind=RP), intent(in) :: Q_x(NGRAD)
real(kind=RP), intent(in) :: Q_y(NGRAD)
real(kind=RP), intent(in) :: Q_z(NGRAD)
real(kind=RP), intent(out) :: mu

end subroutine LESModel_ComputeViscosity
!pure subroutine LESModel_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
! !$acc routine seq
! implicit none
! !-arguments---------------------------------------------
! class(LESModel_t), intent(in) :: this
! real(kind=RP), intent(in) :: delta
! real(kind=RP), intent(in) :: dWall
! real(kind=RP), intent(in) :: Q(NCONS)
! real(kind=RP), intent(in) :: Q_x(NGRAD)
! real(kind=RP), intent(in) :: Q_y(NGRAD)
! real(kind=RP), intent(in) :: Q_z(NGRAD)
! real(kind=RP), intent(out) :: mu

!end subroutine LESModel_ComputeViscosity

subroutine LESModel_Describe(self)
implicit none
Expand All @@ -242,13 +244,13 @@ subroutine LESModel_Selector(this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)

select case (this % LESModel)
case (0) !none
call LESModel_ComputeViscosity(this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
! call LESModel_ComputeViscosity(delta, dWall, Q, Q_x, Q_y, Q_z, mu)
case (1) !Smagorinsky
call Smagorinsky_ComputeViscosity(this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
call Smagorinsky_ComputeViscosity(delta, dWall, Q, Q_x, Q_y, Q_z, mu, this % C, this % WallModel)
case (2) !WALE
call WALE_ComputeViscosity(this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
call WALE_ComputeViscosity(delta, dWall, Q, Q_x, Q_y, Q_z, mu, this % C)
case (3) !Vreman
call Vreman_ComputeViscosity(this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
call Vreman_ComputeViscosity(delta, dWall, Q, Q_x, Q_y, Q_z, mu, this % C)
end select

end subroutine LESModel_Selector
Expand All @@ -272,33 +274,35 @@ subroutine Smagorinsky_Initialize(self, controlVariables)
self % active = .true.

if ( controlVariables % containsKey(LESIntensityKey) ) then
self % CS = controlVariables % doublePrecisionValueForKey(LESIntensityKey)
self % C = controlVariables % doublePrecisionValueForKey(LESIntensityKey)

else
self % CS = 0.2_RP
self % C = 0.2_RP

end if

!$acc enter data copyin(self)
!$acc enter data copyin(self % CS)
!$acc enter data copyin(self % C)
!$acc enter data copyin(self % WallModel)
!$acc enter data copyin(self % LESModel)

end subroutine Smagorinsky_Initialize
!
!///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
!
pure subroutine Smagorinsky_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
pure subroutine Smagorinsky_ComputeViscosity (delta, dWall, Q, Q_x, Q_y, Q_z, mu, C, WallModel)
!$acc routine seq
implicit none
!-arguments---------------------------------------------
class(Smagorinsky_t), intent(in) :: this
! class(Smagorinsky_t), intent(in) :: this
real(kind=RP), intent(in) :: delta
real(kind=RP), intent(in) :: dWall
real(kind=RP), intent(in) :: Q(NCONS)
real(kind=RP), intent(in) :: Q_x(NGRAD)
real(kind=RP), intent(in) :: Q_y(NGRAD)
real(kind=RP), intent(in) :: Q_z(NGRAD)
real(kind=RP), intent(in) :: C
integer, intent(in) :: WallModel
real(kind=RP), intent(out) :: mu
!-local-variables---------------------------------------
real(kind=RP) :: S(NDIM, NDIM)
Expand Down Expand Up @@ -330,8 +334,8 @@ pure subroutine Smagorinsky_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q
!
! Compute viscosity and thermal conductivity
! ------------------------------------------
LS = this % CS * delta
LS = LESModel_ComputeWallEffect(this,LS,dWall)
LS = C * delta
LS = LESModel_ComputeWallEffect(LS,dWall, WallModel)
mu = Q(IRHO) * POW2(LS) * normS

end subroutine Smagorinsky_ComputeViscosity
Expand All @@ -347,7 +351,7 @@ subroutine Smagorinsky_Describe(self)
write(STD_OUT,*)
call SubSection_Header("LES Model")
write(STD_OUT,'(30X,A,A30,A)') "->","LES model: ","Smagorinsky"
write(STD_OUT,'(30X,A,A30,F10.3)') "->","LES model intensity: ", self % CS
write(STD_OUT,'(30X,A,A30,F10.3)') "->","LES model intensity: ", self % C

select case (self % WallModel)
case(LINEAR_WALLMODEL)
Expand Down Expand Up @@ -378,32 +382,33 @@ subroutine WALE_Initialize(self, controlVariables)
self % active = .true.

if ( controlVariables % containsKey(LESIntensityKey) ) then
self % Cw = controlVariables % doublePrecisionValueForKey(LESIntensityKey)
self % C = controlVariables % doublePrecisionValueForKey(LESIntensityKey)

else
self % Cw = 0.325_RP
self % C = 0.325_RP

end if

!$acc enter data copyin(self)
!$acc enter data copyin(self % Cw)
!$acc enter data copyin(self % C)
!$acc enter data copyin(self % WallModel)
!$acc enter data copyin(self % LESModel)

end subroutine WALE_Initialize

pure subroutine WALE_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
pure subroutine WALE_ComputeViscosity (delta, dWall, Q, Q_x, Q_y, Q_z, mu, C)
!$acc routine seq
implicit none
!-arguments---------------------------------------------
class(WALE_t), intent(in) :: this
! class(WALE_t), intent(in) :: this
real(kind=RP), intent(in) :: delta
real(kind=RP), intent(in) :: dWall
real(kind=RP), intent(in) :: Q(NCONS)
real(kind=RP), intent(in) :: Q_x(NGRAD)
real(kind=RP), intent(in) :: Q_y(NGRAD)
real(kind=RP), intent(in) :: Q_z(NGRAD)
real(kind=RP), intent(out) :: mu
real(kind=RP), intent(in) :: C
!-local-variables---------------------------------------
real(kind=RP) :: S(NDIM, NDIM)
real(kind=RP) :: gradV2(NDIM, NDIM), gradV(NDIM,NDIM)
Expand Down Expand Up @@ -444,7 +449,7 @@ pure subroutine WALE_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
Sd(3,3) = Sd(3,3) - 1.0_RP / 3.0_RP * divV2

normSd = sum(Sd*Sd)
LS = this % Cw * delta
LS = C * delta

mu = Q(IRHO) * POW2(LS) * (normSd**(3.0_RP / 2.0_RP) / (normS**(5.0_RP / 2.0_RP)+normSd**(5.0_RP / 4.0_RP)))

Expand All @@ -461,7 +466,7 @@ subroutine WALE_Describe(self)
write(STD_OUT,*)
call SubSection_Header("LES Model")
write(STD_OUT,'(30X,A,A30,A)') "->","LES model: ","Wale"
write(STD_OUT,'(30X,A,A30,F10.3)') "->","LES model intensity: ", self % Cw
write(STD_OUT,'(30X,A,A30,F10.3)') "->","LES model intensity: ", self % C

select case (self % WallModel)
case(NO_WALLMODEL)
Expand Down Expand Up @@ -507,18 +512,19 @@ subroutine Vreman_Initialize(self, controlVariables)

end subroutine Vreman_Initialize

pure subroutine Vreman_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, mu)
pure subroutine Vreman_ComputeViscosity (delta, dWall, Q, Q_x, Q_y, Q_z, mu, C)
!$acc routine seq
implicit none
!-arguments---------------------------------------------
class(Vreman_t), intent(in) :: this
! class(Vreman_t), intent(in) :: this
real(kind=RP), intent(in) :: delta
real(kind=RP), intent(in) :: dWall
real(kind=RP), intent(in) :: Q(NCONS)
real(kind=RP), intent(in) :: Q_x(NGRAD)
real(kind=RP), intent(in) :: Q_y(NGRAD)
real(kind=RP), intent(in) :: Q_z(NGRAD)
real(kind=RP), intent(out) :: mu
real(kind=RP), intent(in) :: C
!-local-variables---------------------------------------
real(kind=RP) :: G__ij(NDIM, NDIM)
real(kind=RP) :: gradV(NDIM, NDIM)
Expand Down Expand Up @@ -556,7 +562,7 @@ pure subroutine Vreman_ComputeViscosity (this, delta, dWall, Q, Q_x, Q_y, Q_z, m
& - G__ij(1,3) * G__ij(1,3)

if(alpha>1.0e-10_RP) then
mu = Q(IRHO) * this % C * sqrt (abs(Bbeta)/alpha)
mu = Q(IRHO) * C * sqrt (abs(Bbeta)/alpha)
else
mu = 0.0_RP
end if
Expand All @@ -582,4 +588,4 @@ subroutine Vreman_Describe(self)
end select

end subroutine Vreman_Describe
end module LESModels
end module LESModels
Loading

0 comments on commit 844b4a3

Please sign in to comment.