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

Revert "Fix hypoelastic instability" #776

Merged
merged 4 commits into from
Dec 27, 2024
Merged
Show file tree
Hide file tree
Changes from 2 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
173 changes: 114 additions & 59 deletions src/simulation/m_hypoelastic.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -8,14 +8,15 @@
module m_hypoelastic

use m_derived_types !< Definitions of the derived types

use m_global_parameters !< Definitions of the global parameters

use m_finite_differences
use m_helper
use m_mpi_proxy !< Message passing interface (MPI) module proxy
Copy link
Contributor

@ChrisZYJ ChrisZYJ Dec 26, 2024

Choose a reason for hiding this comment

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

Could we please remove the use m_mpi_proxy line? I think it wasn't present in the previous version.


implicit none

private; public :: s_initialize_hypoelastic_module, &
s_finalize_hypoelastic_module, &
s_compute_hypoelastic_rhs

real(wp), allocatable, dimension(:) :: Gs
Expand All @@ -29,11 +30,16 @@
real(wp), allocatable, dimension(:, :, :) :: rho_K_field, G_K_field
!$acc declare create(rho_K_field, G_K_field)

real(wp), allocatable, dimension(:, :) :: fd_coeff_x_h
real(wp), allocatable, dimension(:, :) :: fd_coeff_y_h
real(wp), allocatable, dimension(:, :) :: fd_coeff_z_h
!$acc declare create(fd_coeff_x_h,fd_coeff_y_h,fd_coeff_z_h)

contains

subroutine s_initialize_hypoelastic_module

integer :: i
integer :: i, k, r

@:ALLOCATE(Gs(1:num_fluids))
@:ALLOCATE(rho_K_field(0:m,0:n,0:p), G_K_field(0:m,0:n,0:p))
Expand All @@ -51,6 +57,29 @@
end do
!$acc update device(Gs)

@:ALLOCATE(fd_coeff_x_h(-fd_number:fd_number, 0:m))

Check warning on line 60 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L60

Added line #L60 was not covered by tests
if (n > 0) then
@:ALLOCATE(fd_coeff_y_h(-fd_number:fd_number, 0:n))

Check warning on line 62 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L62

Added line #L62 was not covered by tests
end if
if (p > 0) then
@:ALLOCATE(fd_coeff_z_h(-fd_number:fd_number, 0:p))

Check warning on line 65 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L65

Added line #L65 was not covered by tests
end if

! Computing centered finite difference coefficients
call s_compute_finite_difference_coefficients(m, x_cc, fd_coeff_x_h, buff_size, &
fd_number, fd_order)
!$acc update device(fd_coeff_x_h)
if (n > 0) then
call s_compute_finite_difference_coefficients(n, y_cc, fd_coeff_y_h, buff_size, &
fd_number, fd_order)
!$acc update device(fd_coeff_y_h)
end if
if (p > 0) then
call s_compute_finite_difference_coefficients(p, z_cc, fd_coeff_z_h, buff_size, &
fd_number, fd_order)
!$acc update device(fd_coeff_z_h)
end if

end subroutine s_initialize_hypoelastic_module

!> The purpose of this procedure is to compute the source terms
Expand All @@ -66,7 +95,7 @@

real(wp) :: rho_K, G_K

integer :: i, k, l, q !< Loop variables
integer :: i, k, l, q, r !< Loop variables
integer :: ndirs !< Number of coordinate directions

ndirs = 1; if (n > 0) ndirs = 2; if (p > 0) ndirs = 3
Expand All @@ -79,82 +108,91 @@
do q = 0, p
do l = 0, n
do k = 0, m
du_dx(k, l, q) = &
(q_prim_vf(momxb)%sf(k - 2, l, q) &
- 8._wp*q_prim_vf(momxb)%sf(k - 1, l, q) &
+ 8._wp*q_prim_vf(momxb)%sf(k + 1, l, q) &
- q_prim_vf(momxb)%sf(k + 2, l, q)) &
/(12._wp*dx(k))
du_dx(k, l, q) = 0._wp

Check warning on line 111 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L111

Added line #L111 was not covered by tests
end do
end do
end do
!$acc end parallel loop

!$acc parallel loop collapse(3) gang vector default(present)
do q = 0, p
do l = 0, n
do k = 0, m
!$acc loop seq
do r = -fd_number, fd_number
du_dx(k, l, q) = du_dx(k, l, q) &
+ q_prim_vf(momxb)%sf(k + r, l, q)*fd_coeff_x_h(r, k)

Check warning on line 124 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L123-L124

Added lines #L123 - L124 were not covered by tests
end do

end do
end do
end do
!$acc end parallel loop

if (ndirs > 1) then
!$acc parallel loop collapse(3) gang vector default(present)
do q = 0, p
do l = 0, n
do k = 0, m
du_dy(k, l, q) = &
(q_prim_vf(momxb)%sf(k, l - 2, q) &
- 8._wp*q_prim_vf(momxb)%sf(k, l - 1, q) &
+ 8._wp*q_prim_vf(momxb)%sf(k, l + 1, q) &
- q_prim_vf(momxb)%sf(k, l + 2, q)) &
/(12._wp*dy(l))
dv_dx(k, l, q) = &
(q_prim_vf(momxb + 1)%sf(k - 2, l, q) &
- 8._wp*q_prim_vf(momxb + 1)%sf(k - 1, l, q) &
+ 8._wp*q_prim_vf(momxb + 1)%sf(k + 1, l, q) &
- q_prim_vf(momxb + 1)%sf(k + 2, l, q)) &
/(12._wp*dx(k))
dv_dy(k, l, q) = &
(q_prim_vf(momxb + 1)%sf(k, l - 2, q) &
- 8._wp*q_prim_vf(momxb + 1)%sf(k, l - 1, q) &
+ 8._wp*q_prim_vf(momxb + 1)%sf(k, l + 1, q) &
- q_prim_vf(momxb + 1)%sf(k, l + 2, q)) &
/(12._wp*dy(l))
du_dy(k, l, q) = 0._wp; dv_dx(k, l, q) = 0._wp; dv_dy(k, l, q) = 0._wp

Check warning on line 137 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L137

Added line #L137 was not covered by tests
end do
end do
end do
!$acc end parallel loop

!$acc parallel loop collapse(3) gang vector default(present)
do q = 0, p
do l = 0, n
do k = 0, m
!$acc loop seq
do r = -fd_number, fd_number
du_dy(k, l, q) = du_dy(k, l, q) &
+ q_prim_vf(momxb)%sf(k, l + r, q)*fd_coeff_y_h(r, l)
dv_dx(k, l, q) = dv_dx(k, l, q) &
+ q_prim_vf(momxb + 1)%sf(k + r, l, q)*fd_coeff_x_h(r, k)
dv_dy(k, l, q) = dv_dy(k, l, q) &
+ q_prim_vf(momxb + 1)%sf(k, l + r, q)*fd_coeff_y_h(r, l)

Check warning on line 154 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L149-L154

Added lines #L149 - L154 were not covered by tests
end do
end do
end do
end do
!$acc end parallel loop

! 3D
if (ndirs == 3) then

!$acc parallel loop collapse(3) gang vector default(present)
do q = 0, p
do l = 0, n
do k = 0, m
du_dz(k, l, q) = &
(q_prim_vf(momxb)%sf(k, l, q - 2) &
- 8._wp*q_prim_vf(momxb)%sf(k, l, q - 1) &
+ 8._wp*q_prim_vf(momxb)%sf(k, l, q + 1) &
- q_prim_vf(momxb)%sf(k, l, q + 2)) &
/(12._wp*dz(q))
dv_dz(k, l, q) = &
(q_prim_vf(momxb + 1)%sf(k, l, q - 2) &
- 8._wp*q_prim_vf(momxb + 1)%sf(k, l, q - 1) &
+ 8._wp*q_prim_vf(momxb + 1)%sf(k, l, q + 1) &
- q_prim_vf(momxb + 1)%sf(k, l, q + 2)) &
/(12._wp*dz(q))
dw_dx(k, l, q) = &
(q_prim_vf(momxe)%sf(k - 2, l, q) &
- 8._wp*q_prim_vf(momxe)%sf(k - 1, l, q) &
+ 8._wp*q_prim_vf(momxe)%sf(k + 1, l, q) &
- q_prim_vf(momxe)%sf(k + 2, l, q)) &
/(12._wp*dx(k))
dw_dy(k, l, q) = &
(q_prim_vf(momxe)%sf(k, l - 2, q) &
- 8._wp*q_prim_vf(momxe)%sf(k, l - 1, q) &
+ 8._wp*q_prim_vf(momxe)%sf(k, l + 1, q) &
- q_prim_vf(momxe)%sf(k, l + 2, q)) &
/(12._wp*dy(l))
dw_dz(k, l, q) = &
(q_prim_vf(momxe)%sf(k, l, q - 2) &
- 8._wp*q_prim_vf(momxe)%sf(k, l, q - 1) &
+ 8._wp*q_prim_vf(momxe)%sf(k, l, q + 1) &
- q_prim_vf(momxe)%sf(k, l, q + 2)) &
/(12._wp*dz(q))
du_dz(k, l, q) = 0_wp; dv_dz(k, l, q) = 0_wp; dw_dx(k, l, q) = 0_wp;
dw_dy(k, l, q) = 0_wp; dw_dz(k, l, q) = 0_wp;

Check warning on line 169 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L168-L169

Added lines #L168 - L169 were not covered by tests
end do
end do
end do
!$acc end parallel loop

!$acc parallel loop collapse(3) gang vector default(present)
do q = 0, p
do l = 0, n
do k = 0, m
!$acc loop seq
do r = -fd_number, fd_number
du_dz(k, l, q) = du_dz(k, l, q) &
+ q_prim_vf(momxb)%sf(k, l, q + r)*fd_coeff_z_h(r, q)
dv_dz(k, l, q) = dv_dz(k, l, q) &
+ q_prim_vf(momxb + 1)%sf(k, l, q + r)*fd_coeff_z_h(r, q)
dw_dx(k, l, q) = dw_dx(k, l, q) &
+ q_prim_vf(momxe)%sf(k + r, l, q)*fd_coeff_x_h(r, k)
dw_dy(k, l, q) = dw_dy(k, l, q) &
+ q_prim_vf(momxe)%sf(k, l + r, q)*fd_coeff_y_h(r, l)
dw_dz(k, l, q) = dw_dz(k, l, q) &
+ q_prim_vf(momxe)%sf(k, l, q + r)*fd_coeff_z_h(r, q)

Check warning on line 190 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L181-L190

Added lines #L181 - L190 were not covered by tests
end do
end do
end do
end do
!$acc end parallel loop
end if
end if

Expand All @@ -171,7 +209,7 @@
G_K_field(k, l, q) = G_K

!TODO: take this out if not needed
if (G_K < 1000) then
if (G_K < verysmall) then
G_K_field(k, l, q) = 0
end if
end do
Expand Down Expand Up @@ -296,4 +334,21 @@

end subroutine s_compute_hypoelastic_rhs

subroutine s_finalize_hypoelastic_module()

@:DEALLOCATE(Gs)
@:DEALLOCATE(rho_K_field, G_K_field)
@:DEALLOCATE(du_dx)
@:DEALLOCATE(fd_coeff_x_h)

Check warning on line 342 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L339-L342

Added lines #L339 - L342 were not covered by tests
if (n > 0) then
@:DEALLOCATE(du_dy,dv_dx,dv_dy)
@:DEALLOCATE(fd_coeff_y_h)

Check warning on line 345 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L344-L345

Added lines #L344 - L345 were not covered by tests
if (p > 0) then
@:DEALLOCATE(du_dz, dv_dz, dw_dx, dw_dy, dw_dz)
@:DEALLOCATE(fd_coeff_z_h)

Check warning on line 348 in src/simulation/m_hypoelastic.fpp

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_hypoelastic.fpp#L347-L348

Added lines #L347 - L348 were not covered by tests
end if
end if

end subroutine s_finalize_hypoelastic_module

end module m_hypoelastic
1 change: 1 addition & 0 deletions src/simulation/m_start_up.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -1627,6 +1627,7 @@ contains
subroutine s_finalize_modules

call s_finalize_time_steppers_module()
if (hypoelasticity) call s_finalize_hypoelastic_module()
if (hyperelasticity) call s_finalize_hyperelastic_module()
call s_finalize_derived_variables_module()
call s_finalize_data_output_module()
Expand Down
Loading
Loading