Skip to content

Commit

Permalink
Fix Benchmark & Add Monopole Checks & Add Monopole Tests (#497)
Browse files Browse the repository at this point in the history
  • Loading branch information
ChrisZYJ authored Jul 2, 2024
1 parent 0b1abfa commit 9630d35
Show file tree
Hide file tree
Showing 23 changed files with 1,876 additions and 32 deletions.
27 changes: 14 additions & 13 deletions benchmarks/viscous_weno5_sgb_mono/case.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@
procs = DICT["nodes"] * DICT["tasks_per_node"]
ncells = math.floor(ppg * procs * ARGS["gbpp"])
s = math.floor((ncells / 2.0) ** (1/3))
Nx, Ny, Nz = 2*s, s, s
Nx, Ny, Nz = s, s, 2*s

x0 = 10.E-04
y0 = 10.E-04
Expand Down Expand Up @@ -86,8 +86,8 @@
cfl = 0.01
Ldomain = 20.E-03
L = Ldomain/x0
dx = L/float(Nx)
dt = cfl*dx*c0/cact
dz = L/float(Nz)
dt = cfl*dz*c0/cact
Lpulse = 0.3*Ldomain
Tpulse = Lpulse/cact

Expand All @@ -98,12 +98,12 @@
# ==========================================================

# Computational Domain Parameters ==========================
'x_domain%beg' : -10.E-03/x0,
'x_domain%end' : 10.E-03/x0,
'x_domain%beg' : -5.E-03/x0,
'x_domain%end' : 5.E-03/x0,
'y_domain%beg' : -5.E-03/y0,
'y_domain%end' : 5.E-03/y0,
'z_domain%beg' : -5.E-03/z0,
'z_domain%end' : 5.E-03/z0,
'y_domain%end' : 5.E-03/y0,
'z_domain%beg' : -10.E-03/z0,
'z_domain%end' : 10.E-03/z0,
'stretch_x' : 'F',
'cyl_coord' : 'F',
'm' : Nx,
Expand Down Expand Up @@ -154,9 +154,9 @@
'patch_icpp(1)%x_centroid' : 0.,
'patch_icpp(1)%y_centroid' : 0.,
'patch_icpp(1)%z_centroid' : 0.,
'patch_icpp(1)%length_x' : 20.E-03/x0,
'patch_icpp(1)%length_x' : 10.E-03/x0,
'patch_icpp(1)%length_y' : 10.E-03/y0,
'patch_icpp(1)%length_z' : 10.E-03/z0,
'patch_icpp(1)%length_z' : 20.E-03/z0,
'patch_icpp(1)%vel(1)' : 0.0,
'patch_icpp(1)%vel(2)' : 0.0,
'patch_icpp(1)%vel(3)' : 0.0,
Expand All @@ -172,9 +172,9 @@
'patch_icpp(2)%x_centroid' : 0.,
'patch_icpp(2)%y_centroid' : 0.,
'patch_icpp(2)%z_centroid' : 0.,
'patch_icpp(2)%length_x' : 5.E-03/x0,
'patch_icpp(2)%length_x' : 10.E-03/x0,
'patch_icpp(2)%length_y' : 10.E-03/y0,
'patch_icpp(2)%length_z' : 10.E-03/z0,
'patch_icpp(2)%length_z' : 5.E-03/z0,
'patch_icpp(2)%alter_patch(1)' : 'T',
'patch_icpp(2)%vel(1)' : 0.0,
'patch_icpp(2)%vel(2)' : 0.0,
Expand Down Expand Up @@ -230,7 +230,8 @@
# Acoustic source ==========================================
'Monopole' : 'T',
'num_mono' : 1,
'Mono(1)%loc(1)' : -5.E-03/x0,
'Mono(1)%support' : 4,
'Mono(1)%loc(3)' : -5.E-03/x0,
'Mono(1)%npulse' : 1,
'Mono(1)%dir' : 1.,
'Mono(1)%pulse' : 1,
Expand Down
77 changes: 77 additions & 0 deletions src/simulation/m_checker.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,7 @@ contains
call s_check_inputs_riemann_solver
call s_check_inputs_time_stepping
call s_check_inputs_model_eqns
if (monopole) call s_check_inputs_monopole
if (hypoelasticity) call s_check_inputs_hypoelasticity
if (bubbles) call s_check_inputs_bubbles
if (adap_dt) call s_check_inputs_adapt_dt
Expand Down Expand Up @@ -159,6 +160,82 @@ contains
end if
end subroutine s_check_inputs_model_eqns

!> Checks constraints on monopole parameters
subroutine s_check_inputs_monopole
integer :: j
character(len=5) :: jStr

if (num_mono == dflt_int) then
call s_mpi_abort('num_mono must be specified for monopole. Exiting ...')
elseif (num_mono < 0) then
call s_mpi_abort('num_mono must be non-negative. Exiting ...')
end if

do j = 1, num_mono
call s_int_to_str(j, jStr)
if (mono(j)%support == dflt_int) then
call s_mpi_abort('mono('//trim(jStr)//')%support must be '// &
'specified. Exiting ...')
elseif (f_is_default(mono(j)%mag)) then
call s_mpi_abort('mono('//trim(jStr)//')%mag must be '// &
'specified. Exiting ...')
elseif (mono(j)%mag <= 0d0) then
call s_mpi_abort('mono('//trim(jStr)//')%mag must be '// &
'positive. Exiting ...')
end if

if (n == 0) then ! 1D
if (.not. any(mono(j)%support == (/0, 1/))) then ! undocumented support 0
call s_mpi_abort('Only Mono(i)support = 1 is allowed for '// &
'1D simulations. Exiting ...')
end if
if (mono(j)%support == 1 .and. f_is_default(mono(j)%loc(1))) then
call s_mpi_abort('mono_loc(1) must be specified for '// &
'Mono(i)support = 1. Exiting ...')
end if
elseif (p == 0) then ! 2D
if (.not. any(mono(j)%support == (/1, 2, 3, 4, 5/))) then
call s_mpi_abort('Only Mono(i)support = 1, 2, 3, 4, or 5 is '// &
'allowed for 2D simulations. Exiting ...')
end if
if (any(mono(j)%support == (/1, 2, 3, 5/)) .and. &
(f_is_default(mono(j)%loc(1)) .or. &
f_is_default(mono(j)%loc(2)))) then
call s_mpi_abort('mono('//trim(jStr)//')%loc(1:2) must be '// &
'specified for Mono(i)support = 1, 3, or 5. '// &
'Exiting ...')
elseif (mono(j)%support == 4 .and. f_is_default(mono(j)%loc(1))) then
call s_mpi_abort('mono('//trim(jStr)//')%loc(1) must be '// &
'specified for Mono(i)support = 4. Exiting ...')
end if
else ! 3D
if (.not. any(mono(j)%support == (/3, 4, 5, 6/))) then
call s_mpi_abort('Only Mono(i)support = 3, 4, 5, or 6 is '// &
'allowed for 3D simulations. Exiting ...')
elseif (mono(j)%support == 6 .and. (.not. cyl_coord)) then
call s_mpi_abort('Mono(i)support = 6 requires cyl_coord = true. '// &
'Exiting ...')
elseif (cyl_coord .and. mono(j)%support /= 6) then
call s_mpi_abort('cyl_coord = true requires Mono(i)support = 6. '// &
'Exiting ...')
end if
if (any(mono(j)%support == (/3, 5, 6/)) .and. &
(f_is_default(mono(j)%loc(1)) .or. &
f_is_default(mono(j)%loc(2)) .or. &
f_is_default(mono(j)%loc(3)))) then
call s_mpi_abort('mono('//trim(jStr)//')%loc(1:3) must be '// &
'specified for Mono(i)support = 3, 5, or 6. '// &
'Exiting ...')
elseif (mono(j)%support == 4 .and. &
(f_is_default(mono(j)%loc(3)))) then
call s_mpi_abort('mono('//trim(jStr)//')%loc(3) must be '// &
'specified for Mono(i)support = 4. Exiting ...')
end if
end if
end do

end subroutine s_check_inputs_monopole

!> Checks constraints on hypoelasticity parameters
subroutine s_check_inputs_hypoelasticity
if (riemann_solver /= 1) then
Expand Down
16 changes: 8 additions & 8 deletions src/simulation/m_global_parameters.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -595,7 +595,7 @@ contains

! Monopole source
monopole = .false.
num_mono = 1
num_mono = dflt_real

! Surface tension
sigma = dflt_real
Expand All @@ -619,16 +619,16 @@ contains
mono(j)%mag = dflt_real
mono(j)%length = dflt_real
mono(j)%delay = dflt_real
mono(j)%dir = 1.d0
mono(j)%npulse = 1.d0
mono(j)%pulse = 1
mono(j)%support = 1
mono(j)%dir = dflt_real
mono(j)%npulse = dflt_int
mono(j)%pulse = dflt_int
mono(j)%support = dflt_int
mono(j)%foc_length = dflt_real
mono(j)%aperture = dflt_real
mono(j)%support_width = dflt_real
! The author suggested the support width is typically on the order of
! the width of the characteristic cells. Here, we choose 2.5 cell width
! as the default value.
mono(j)%support_width = 2.5d0
! the width of the characteristic cells.
! The default value of support_width is 2.5 cell widths.
end do

fd_order = dflt_int
Expand Down
37 changes: 28 additions & 9 deletions src/simulation/m_monopole.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -69,19 +69,39 @@ contains
@:ALLOCATE_GLOBAL(mag(1:num_mono), support(1:num_mono), length(1:num_mono), npulse(1:num_mono), pulse(1:num_mono), dir(1:num_mono), delay(1:num_mono), loc_mono(1:3, 1:num_mono), foc_length(1:num_mono), aperture(1:num_mono), support_width(1:num_mono))

do i = 1, num_mono
do j = 1, 3
loc_mono(j, i) = mono(i)%loc(j)
end do
mag(i) = mono(i)%mag
support(i) = mono(i)%support
length(i) = mono(i)%length
npulse(i) = mono(i)%npulse
pulse(i) = mono(i)%pulse
dir(i) = mono(i)%dir
delay(i) = mono(i)%delay
foc_length(i) = mono(i)%foc_length
aperture(i) = mono(i)%aperture
support_width(i) = mono(i)%support_width
do j = 1, 3
loc_mono(j, i) = mono(i)%loc(j)
end do
if (mono(i)%npulse == dflt_int) then
npulse(i) = 1
else
npulse(i) = mono(i)%npulse
end if
if (mono(i)%pulse == dflt_int) then
pulse(i) = 1
else
pulse(i) = mono(i)%pulse
end if
if (f_is_default(mono(i)%dir)) then
dir(i) = 1d0
else
dir(i) = mono(i)%dir
end if
if (f_is_default(mono(i)%delay)) then
delay(i) = 0d0
else
delay(i) = mono(i)%delay
end if
if (f_is_default(mono(i)%support_width)) then
support_width(i) = 2.5d0
else
support_width(i) = mono(i)%support_width
end if
end do
!$acc update device(mag, support, length, npulse, pulse, dir, delay, foc_length, aperture, loc_mono, support_width)

Expand Down Expand Up @@ -332,7 +352,6 @@ contains
if (the_time > t0 .and. the_time < sigt) then
f_g = mag(nm)
end if
else
end if

end function f_g
Expand Down
Loading

0 comments on commit 9630d35

Please sign in to comment.