Skip to content
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.

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

Fix Benchmark & Add Monopole Checks & Add Monopole Tests #497

Merged
merged 7 commits into from
Jul 2, 2024
Merged
Show file tree
Hide file tree
Changes from 6 commits
Commits
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
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)%delay)) then
ChrisZYJ marked this conversation as resolved.
Show resolved Hide resolved
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
Loading