Skip to content

Commit 70f53b0

Browse files
albertjimenezramosalbertjimenezramos
authored
Merge pull request #297 from loganoz/296-acoustics-add-source-term
296 acoustics add source term
2 parents 566971c + 4a10d20 commit 70f53b0

5 files changed

Lines changed: 114 additions & 2 deletions

File tree

‎Solver/src/AcousticSolver/SpatialDiscretization.f90‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -209,6 +209,7 @@ END SUBROUTINE ComputeTimeDerivativeIsolated
209209
subroutine TimeDerivative_ComputeQDot( mesh , particles, t)
210210
! use ActuatorLine, only: farm
211211
use SpongeClass, only: sponge, addSourceSponge
212+
use APESourceClass, only: apeSource, addAPESource
212213
implicit none
213214
type(HexMesh) :: mesh
214215
type(Particles_t) :: particles
@@ -333,6 +334,9 @@ subroutine TimeDerivative_ComputeQDot( mesh , particles, t)
333334
!$omp end do
334335
! for the sponge, loops are in the internal subroutine as values are precalculated
335336
call addSourceSponge(sponge,mesh)
337+
338+
! add source term for the acoustics
339+
call addAPESource(apeSource, mesh)
336340
!
337341
! ***********************
338342
! Now add the source term

‎Solver/src/libs/mesh/StorageClass.f90‎

Lines changed: 11 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -95,7 +95,8 @@ module StorageClass
9595
#ifdef ACOUSTIC
9696
real(kind=RP), allocatable :: Qbase(:,:,:,:) ! Base flow State vector
9797
real(kind=RP), allocatable :: Lambbase(:,:,:,:) ! Lamb vector base
98-
real(kind=RP), allocatable :: Lamb_NS(:,:,:,:) ! Lamb vector from NS
98+
real(kind=RP), allocatable :: Lamb_NS(:,:,:,:) ! Lamb vector from NS
99+
real(kind=RP), allocatable :: grada2base(:,:,:,:) ! Gradient of the squared sound velocity
99100
#endif
100101
#ifdef CAHNHILLIARD
101102
real(kind=RP), dimension(:,:,:,:), allocatable :: c ! CHE concentration
@@ -199,6 +200,7 @@ module StorageClass
199200
real(kind=RP), dimension(:,:,:), allocatable :: Qbase ! Base flow State vector
200201
real(kind=RP), dimension(:,:,:), allocatable :: Lambbase ! Lamb vector base
201202
real(kind=RP), dimension(:,:,:), allocatable :: Lamb_NS ! Lamb vector
203+
real(kind=RP), dimension(:,:,:), allocatable :: grada2base ! Gradient of the squared sound velocity
202204
#endif
203205
#ifdef MULTIPHASE
204206
real(kind=RP), dimension(:,:), allocatable :: invMa2
@@ -819,6 +821,7 @@ elemental subroutine ElementStorage_Construct(self, Nx, Ny, Nz, computeGradients
819821
ALLOCATE( self % Qbase (NCONSB,0:Nx,0:Ny,0:Nz) )
820822
ALLOCATE( self % Lambbase (NDIM,0:Nx,0:Ny,0:Nz) )
821823
ALLOCATE( self % Lamb_NS (NDIM,0:Nx,0:Ny,0:Nz) )
824+
ALLOCATE( self % grada2base (NDIM,0:Nx,0:Ny,0:Nz) )
822825
#endif
823826
if (computeGradients) then
824827
ALLOCATE( self % U_xNS (NGRAD,0:Nx,0:Ny,0:Nz) )
@@ -901,6 +904,7 @@ elemental subroutine ElementStorage_Construct(self, Nx, Ny, Nz, computeGradients
901904
self % Qbase = 0.0_RP
902905
self % Lambbase = 0.0_RP
903906
self % Lamb_NS = 0.0_RP
907+
self % grada2base = 0.0_RP
904908
#endif
905909
if (computeGradients) then
906910
self % U_xNS = 0.0_RP
@@ -1012,6 +1016,7 @@ elemental subroutine ElementStorage_Assign(to, from)
10121016
to % Qbase = from % Qbase
10131017
to % Lambbase = from % Lambbase
10141018
to % Lamb_NS = from % Lamb_NS
1019+
to % grada2base = from % grada2base
10151020
#endif
10161021

10171022
#ifndef ACOUSTIC
@@ -1110,6 +1115,7 @@ elemental subroutine ElementStorage_Destruct(self)
11101115
safedeallocate(self % Qbase)
11111116
safedeallocate(self % Lambbase)
11121117
safedeallocate(self % Lamb_NS)
1118+
safedeallocate(self % grada2base)
11131119
#endif
11141120

11151121
if (self % computeGradients) then
@@ -1421,6 +1427,7 @@ pure subroutine FaceStorage_Construct(self, NDIM, Nf, Nel, computeGradients, ana
14211427
ALLOCATE( self % Qbase (NCONSB,0:Nf(1),0:Nf(2)) )
14221428
ALLOCATE( self % Lambbase (NDIM,0:Nf(1),0:Nf(2)) )
14231429
ALLOCATE( self % Lamb_NS (NDIM,0:Nf(1),0:Nf(2)) )
1430+
ALLOCATE( self % grada2base (NDIM,0:Nf(1),0:Nf(2)) )
14241431
#endif
14251432
! Biggest Interface flux memory size is u\vec{n}
14261433
! ----------------------------------------------
@@ -1486,6 +1493,7 @@ pure subroutine FaceStorage_Construct(self, NDIM, Nf, Nel, computeGradients, ana
14861493
self % Qbase = 0.0_RP
14871494
self % Lambbase = 0.0_RP
14881495
self % Lamb_NS = 0.0_RP
1496+
self % grada2base = 0.0_RP
14891497
#endif
14901498

14911499
self % rho = 0.0_RP
@@ -1586,6 +1594,7 @@ elemental subroutine FaceStorage_Destruct(self)
15861594
safedeallocate(self % Qbase )
15871595
safedeallocate(self % Lambbase )
15881596
safedeallocate(self % Lamb_NS )
1597+
safedeallocate(self % grada2base )
15891598
#endif
15901599

15911600
self % anJacobian = .FALSE.
@@ -1752,6 +1761,7 @@ elemental subroutine FaceStorage_Assign(to, from)
17521761
to % Qbase = from % Qbase
17531762
to % Lambbase = from % Lambbase
17541763
to % Lamb_NS = from % Lamb_NS
1764+
to % grada2base = from % grada2base
17551765
#endif
17561766
if (to % computeGradients) then
17571767
to % U_xNS = from % U_xNS

‎Solver/src/libs/sources/APESourceClass.f90‎

Lines changed: 93 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,93 @@
1+
!//////////////////////////////////////////////////////
2+
!
3+
!This class represents the source term for Acoustic Perturbation Equation
4+
5+
#include "Includes.h"
6+
Module APESourceClass !
7+
#if defined(ACOUSTIC)
8+
use SMConstants
9+
use HexMeshClass
10+
use SolutionFile
11+
use PhysicsStorage, only: NCONS
12+
use PhysicsStorage_CAA
13+
use StorageClass, only: ElementStorage_t
14+
Implicit None
15+
16+
private
17+
public apeSource, addAPESource, constructAPESource
18+
19+
! Acoustic Perturbation Equation number
20+
integer, parameter :: APEEQ_4 = 4
21+
22+
! Definition of apeSource class
23+
type apeSource_t
24+
integer :: apeeq
25+
contains
26+
end type apeSource_t
27+
28+
type(apeSource_t) :: apeSource
29+
30+
contains
31+
32+
!/////////////////////////////////////////////////////////////////////////
33+
! APE PROCEDURES --------------------------
34+
!/////////////////////////////////////////////////////////////////////////
35+
36+
subroutine constructAPESource(self, controlVariables)
37+
use FTValueDictionaryClass
38+
implicit none
39+
type(apeSource_t) :: self
40+
type(FTValueDictionary), intent(in) :: controlVariables
41+
42+
! Read APE equation from control file
43+
self % apeeq = controlVariables % getValueOrDefault("APE equation", APEEQ_4)
44+
if (self % apeeq .ne. APEEQ_4) then
45+
print *, "Only APE-4 equation is implemented."
46+
error stop
47+
end if
48+
49+
end subroutine constructAPESource
50+
!
51+
Subroutine addAPESource(self, mesh)
52+
type(apeSource_t) :: self
53+
type(HexMesh), intent(inout) :: mesh
54+
55+
select case (self % apeeq)
56+
case (APEEQ_4)
57+
call addAPE4Source(self, mesh)
58+
end select
59+
end subroutine addAPESource
60+
61+
!/////////////////////////////////////////////////////////////////////////
62+
! APE-4
63+
!/////////////////////////////////////////////////////////////////////////
64+
Subroutine addAPE4Source(self, mesh)
65+
type(apeSource_t) :: self
66+
type(HexMesh), intent(inout) :: mesh
67+
68+
! Local variables
69+
integer :: i, j, k, eID
70+
real(rp) :: qrho, qvel(NDIM), qp, aux(NDIM)
71+
72+
qrho = 0.0
73+
74+
!$omp do schedule(runtime) private(i,j,k,eID)
75+
do eID = 1, mesh % no_of_elements
76+
associate ( e => mesh % elements(eID) )
77+
do k = 0, e % Nxyz(3) ; do j = 0, e % Nxyz(2) ; do i = 0, e % Nxyz(1)
78+
qvel = e % storage % Lambbase(:,i,j,k) - e % storage % Lamb_NS(:,i,j,k)
79+
aux = e % storage % Qbase(IBRHO,i,j,k) * e % storage % Q(ICAAU:ICAAW,i,j,k) + &
80+
e % storage % Qbase(IBU:IBW,i,j,k) * e % storage % Q(ICAAP,i,j,k) / e % storage % Qbase(IBA2,i,j,k)
81+
qp = dot_product(aux, e % storage % grada2base(:,i,j,k))
82+
e % storage % S_NS(ICAARHO,i,j,k) = e % storage % S_NS(ICAARHO,i,j,k) + qrho
83+
e % storage % S_NS(ICAAU:ICAAW,i,j,k) = e % storage % S_NS(ICAAU:ICAAW,i,j,k) + qvel
84+
e % storage % S_NS(ICAAP,i,j,k) = e % storage % S_NS(ICAAP,i,j,k) + qp
85+
end do ; end do ; end do
86+
end associate
87+
end do
88+
!$omp end do
89+
90+
End Subroutine addAPE4Source
91+
!
92+
#endif
93+
End Module APESourceClass

‎Solver/src/libs/sources/Makefile‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -33,7 +33,8 @@ LIBOBJS_MU = ./build_mu/ControllerInterface.o \
3333

3434
LIBOBJS_CAA = ./build_caa/ControllerInterface.o \
3535
./build_caa/ActuatorLine.o \
36-
./build_caa/SpongeClass.o
36+
./build_caa/SpongeClass.o \
37+
./build_caa/APESourceClass.o
3738

3839
LIB = sources
3940

‎Solver/src/libs/timeintegrator/TimeIntegrator.f90‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -528,6 +528,9 @@ subroutine IntegrateInTime( self, sem, controlVariables, monitors, samplings, Co
528528
#if defined(NAVIERSTOKES) || defined(INCNS) || defined(MULTIPHASE)
529529
use ActuatorLine, only: farm, ConstructFarm, DestructFarm, UpdateFarm, WriteFarmForces, FarmUpdateControlVars, WriteControlVars, SaveControllerState
530530
#endif
531+
#if defined(ACOUSTIC)
532+
use APESourceClass, only: apeSource, constructAPESource
533+
#endif
531534

532535
use IBMClass
533536

@@ -617,6 +620,7 @@ subroutine IntegrateInTime( self, sem, controlVariables, monitors, samplings, Co
617620
#endif
618621
#if defined(ACOUSTIC)
619622
call sem % mesh % LoadLambVectorStatistics(controlVariables)
623+
call constructAPESource(apeSource, controlVariables)
620624
#endif
621625
!
622626
! ----------------------------------

0 commit comments

Comments
 (0)