Skip to content

Commit 50d75cd

Browse files
committed
fix stress reflective BC
1 parent 314c05a commit 50d75cd

File tree

4 files changed

+203
-0
lines changed

4 files changed

+203
-0
lines changed

src/common/m_boundary_conditions.fpp

Lines changed: 42 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -354,6 +354,13 @@ contains
354354
q_prim_vf(i)%sf(j - 1, k, l)
355355
end do
356356

357+
if (elasticity) then
358+
do i = 1, shear_BC_flip_num
359+
q_prim_vf(shear_BC_flip_indices(1, i))%sf(-j, k, l) = &
360+
-q_prim_vf(shear_BC_flip_indices(1, i))%sf(j - 1, k, l)
361+
end do
362+
end if
363+
357364
if (hyperelasticity) then
358365
q_prim_vf(xibeg)%sf(-j, k, l) = &
359366
-q_prim_vf(xibeg)%sf(j - 1, k, l)
@@ -404,6 +411,13 @@ contains
404411
q_prim_vf(i)%sf(m - (j - 1), k, l)
405412
end do
406413

414+
if (elasticity) then
415+
do i = 1, shear_BC_flip_num
416+
q_prim_vf(shear_BC_flip_indices(1, i))%sf(m + j, k, l) = &
417+
-q_prim_vf(shear_BC_flip_indices(1, i))%sf(m - (j - 1), k, l)
418+
end do
419+
end if
420+
407421
if (hyperelasticity) then
408422
q_prim_vf(xibeg)%sf(m + j, k, l) = &
409423
-q_prim_vf(xibeg)%sf(m - (j - 1), k, l)
@@ -457,6 +471,13 @@ contains
457471
q_prim_vf(i)%sf(l, j - 1, k)
458472
end do
459473

474+
if (elasticity) then
475+
do i = 1, shear_BC_flip_num
476+
q_prim_vf(shear_BC_flip_indices(2, i))%sf(l, -j, k) = &
477+
-q_prim_vf(shear_BC_flip_indices(2, i))%sf(l, j - 1, k)
478+
end do
479+
end if
480+
460481
if (hyperelasticity) then
461482
q_prim_vf(xibeg + 1)%sf(l, -j, k) = &
462483
-q_prim_vf(xibeg + 1)%sf(l, j - 1, k)
@@ -504,6 +525,13 @@ contains
504525
q_prim_vf(i)%sf(l, n - (j - 1), k)
505526
end do
506527

528+
if (elasticity) then
529+
do i = 1, shear_BC_flip_num
530+
q_prim_vf(shear_BC_flip_indices(2, i))%sf(l, n + j, k) = &
531+
-q_prim_vf(shear_BC_flip_indices(2, i))%sf(l, n - (j - 1), k)
532+
end do
533+
end if
534+
507535
if (hyperelasticity) then
508536
q_prim_vf(xibeg + 1)%sf(l, n + j, k) = &
509537
-q_prim_vf(xibeg + 1)%sf(l, n - (j - 1), k)
@@ -556,6 +584,13 @@ contains
556584
q_prim_vf(i)%sf(k, l, j - 1)
557585
end do
558586

587+
if (elasticity) then
588+
do i = 1, shear_BC_flip_num
589+
q_prim_vf(shear_BC_flip_indices(3, i))%sf(k, l, -j) = &
590+
-q_prim_vf(shear_BC_flip_indices(3, i))%sf(k, l, j - 1)
591+
end do
592+
end if
593+
559594
if (hyperelasticity) then
560595
q_prim_vf(xiend)%sf(k, l, -j) = &
561596
-q_prim_vf(xiend)%sf(k, l, j - 1)
@@ -603,6 +638,13 @@ contains
603638
q_prim_vf(i)%sf(k, l, p - (j - 1))
604639
end do
605640

641+
if (elasticity) then
642+
do i = 1, shear_BC_flip_num
643+
q_prim_vf(shear_BC_flip_indices(3, i))%sf(k, l, p + j) = &
644+
-q_prim_vf(shear_BC_flip_indices(3, i))%sf(k, l, p - (j - 1))
645+
end do
646+
end if
647+
606648
if (hyperelasticity) then
607649
q_prim_vf(xiend)%sf(k, l, p + j) = &
608650
-q_prim_vf(xiend)%sf(k, l, p - (j - 1))

src/post_process/m_global_parameters.fpp

Lines changed: 51 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -146,6 +146,15 @@ module m_global_parameters
146146
type(int_bounds_info) :: bc_x, bc_y, bc_z
147147
!> @}
148148

149+
integer :: shear_num !! Number of shear stress components
150+
integer, dimension(3) :: shear_indices !<
151+
!! Indices of the stress components that represent shear stress
152+
integer :: shear_BC_flip_num !<
153+
!! Number of shear stress components to reflect for boundary conditions
154+
integer, dimension(3, 2) :: shear_BC_flip_indices !<
155+
!! Indices of shear stress components to reflect for boundary conditions.
156+
!! Size: (1:3, 1:shear_BC_flip_num) for (x/y/z, [indices])
157+
149158
logical :: parallel_io !< Format of the data files
150159
logical :: sim_data
151160
logical :: file_per_process !< output format
@@ -569,6 +578,27 @@ contains
569578
stress_idx%end = sys_size + (num_dims*(num_dims + 1))/2
570579
! number of distinct stresses is 1 in 1D, 3 in 2D, 6 in 3D
571580
sys_size = stress_idx%end
581+
582+
! shear stress index is 2 for 2D and 2,4,5 for 3D
583+
if (num_dims == 1) then
584+
shear_num = 0
585+
else if (num_dims == 2) then
586+
shear_num = 1
587+
shear_indices(1) = stress_idx%beg - 1 + 2
588+
shear_BC_flip_num = 1
589+
shear_BC_flip_indices(1:2, 1) = shear_indices(1)
590+
! Both x-dir and y-dir: flip tau_xy only
591+
else if (num_dims == 3) then
592+
shear_num = 3
593+
shear_indices(1:3) = stress_idx%beg - 1 + (/2, 4, 5/)
594+
shear_BC_flip_num = 2
595+
shear_BC_flip_indices(1, 1:2) = shear_indices((/1, 2/))
596+
shear_BC_flip_indices(2, 1:2) = shear_indices((/1, 3/))
597+
shear_BC_flip_indices(3, 1:2) = shear_indices((/2, 3/))
598+
! x-dir: flip tau_xy and tau_xz
599+
! y-dir: flip tau_xy and tau_yz
600+
! z-dir: flip tau_xz and tau_yz
601+
end if
572602
end if
573603

574604
if (hyperelasticity) then
@@ -610,6 +640,27 @@ contains
610640
stress_idx%end = sys_size + (num_dims*(num_dims + 1))/2
611641
! number of stresses is 1 in 1D, 3 in 2D, 6 in 3D
612642
sys_size = stress_idx%end
643+
644+
! shear stress index is 2 for 2D and 2,4,5 for 3D
645+
if (num_dims == 1) then
646+
shear_num = 0
647+
else if (num_dims == 2) then
648+
shear_num = 1
649+
shear_indices(1) = stress_idx%beg - 1 + 2
650+
shear_BC_flip_num = 1
651+
shear_BC_flip_indices(1:2, 1) = shear_indices(1)
652+
! Both x-dir and y-dir: flip tau_xy only
653+
else if (num_dims == 3) then
654+
shear_num = 3
655+
shear_indices(1:3) = stress_idx%beg - 1 + (/2, 4, 5/)
656+
shear_BC_flip_num = 2
657+
shear_BC_flip_indices(1, 1:2) = shear_indices((/1, 2/))
658+
shear_BC_flip_indices(2, 1:2) = shear_indices((/1, 3/))
659+
shear_BC_flip_indices(3, 1:2) = shear_indices((/2, 3/))
660+
! x-dir: flip tau_xy and tau_xz
661+
! y-dir: flip tau_xy and tau_yz
662+
! z-dir: flip tau_xz and tau_yz
663+
end if
613664
end if
614665

615666
if (hyperelasticity) then

src/pre_process/m_global_parameters.fpp

Lines changed: 51 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -120,6 +120,15 @@ module m_global_parameters
120120
type(int_bounds_info) :: bc_x, bc_y, bc_z !<
121121
!! Boundary conditions in the x-, y- and z-coordinate directions
122122

123+
integer :: shear_num !! Number of shear stress components
124+
integer, dimension(3) :: shear_indices !<
125+
!! Indices of the stress components that represent shear stress
126+
integer :: shear_BC_flip_num !<
127+
!! Number of shear stress components to reflect for boundary conditions
128+
integer, dimension(3, 2) :: shear_BC_flip_indices !<
129+
!! Indices of shear stress components to reflect for boundary conditions.
130+
!! Size: (1:3, 1:shear_BC_flip_num) for (x/y/z, [indices])
131+
123132
logical :: parallel_io !< Format of the data files
124133
logical :: file_per_process !< type of data output
125134
integer :: precision !< Precision of output files
@@ -659,6 +668,27 @@ contains
659668
stress_idx%end = sys_size + (num_dims*(num_dims + 1))/2
660669
! number of stresses is 1 in 1D, 3 in 2D, 6 in 3D
661670
sys_size = stress_idx%end
671+
672+
! shear stress index is 2 for 2D and 2,4,5 for 3D
673+
if (num_dims == 1) then
674+
shear_num = 0
675+
else if (num_dims == 2) then
676+
shear_num = 1
677+
shear_indices(1) = stress_idx%beg - 1 + 2
678+
shear_BC_flip_num = 1
679+
shear_BC_flip_indices(1:2, 1) = shear_indices(1)
680+
! Both x-dir and y-dir: flip tau_xy only
681+
else if (num_dims == 3) then
682+
shear_num = 3
683+
shear_indices(1:3) = stress_idx%beg - 1 + (/2, 4, 5/)
684+
shear_BC_flip_num = 2
685+
shear_BC_flip_indices(1, 1:2) = shear_indices((/1, 2/))
686+
shear_BC_flip_indices(2, 1:2) = shear_indices((/1, 3/))
687+
shear_BC_flip_indices(3, 1:2) = shear_indices((/2, 3/))
688+
! x-dir: flip tau_xy and tau_xz
689+
! y-dir: flip tau_xy and tau_yz
690+
! z-dir: flip tau_xz and tau_yz
691+
end if
662692
end if
663693

664694
if (hyperelasticity) then
@@ -699,6 +729,27 @@ contains
699729
stress_idx%end = sys_size + (num_dims*(num_dims + 1))/2
700730
! number of stresses is 1 in 1D, 3 in 2D, 6 in 3D
701731
sys_size = stress_idx%end
732+
733+
! shear stress index is 2 for 2D and 2,4,5 for 3D
734+
if (num_dims == 1) then
735+
shear_num = 0
736+
else if (num_dims == 2) then
737+
shear_num = 1
738+
shear_indices(1) = stress_idx%beg - 1 + 2
739+
shear_BC_flip_num = 1
740+
shear_BC_flip_indices(1:2, 1) = shear_indices(1)
741+
! Both x-dir and y-dir: flip tau_xy only
742+
else if (num_dims == 3) then
743+
shear_num = 3
744+
shear_indices(1:3) = stress_idx%beg - 1 + (/2, 4, 5/)
745+
shear_BC_flip_num = 2
746+
shear_BC_flip_indices(1, 1:2) = shear_indices((/1, 2/))
747+
shear_BC_flip_indices(2, 1:2) = shear_indices((/1, 3/))
748+
shear_BC_flip_indices(3, 1:2) = shear_indices((/2, 3/))
749+
! x-dir: flip tau_xy and tau_xz
750+
! y-dir: flip tau_xy and tau_yz
751+
! z-dir: flip tau_xz and tau_yz
752+
end if
702753
end if
703754

704755
if (hyperelasticity) then

src/simulation/m_global_parameters.fpp

Lines changed: 59 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -289,6 +289,17 @@ module m_global_parameters
289289

290290
!$acc declare create(sys_size, buff_size, E_idx, gamma_idx, pi_inf_idx, alf_idx, n_idx, stress_idx, b_size, tensor_size, xi_idx, species_idx)
291291

292+
integer :: shear_num !! Number of shear stress components
293+
integer, dimension(3) :: shear_indices !<
294+
!! Indices of the stress components that represent shear stress
295+
integer :: shear_BC_flip_num !<
296+
!! Number of shear stress components to reflect for boundary conditions
297+
integer, dimension(3, 2) :: shear_BC_flip_indices !<
298+
!! Indices of shear stress components to reflect for boundary conditions.
299+
!! Size: (1:3, 1:shear_BC_flip_num) for (x/y/z, [indices])
300+
301+
!$acc declare create(shear_num, shear_indices, shear_BC_flip_num, shear_BC_flip_indices)
302+
292303
! END: Simulation Algorithm Parameters
293304

294305
! Fluids Physical Parameters
@@ -895,6 +906,30 @@ contains
895906
stress_idx%end = sys_size + (num_dims*(num_dims + 1))/2
896907
! number of distinct stresses is 1 in 1D, 3 in 2D, 6 in 3D
897908
sys_size = stress_idx%end
909+
910+
! shear stress index is 2 for 2D and 2,4,5 for 3D
911+
print *, 'DEBUG DEBUG'
912+
if (num_dims == 1) then
913+
shear_num = 0
914+
else if (num_dims == 2) then
915+
print *, 'HERE'
916+
shear_num = 1
917+
shear_indices(1) = stress_idx%beg - 1 + 2
918+
shear_BC_flip_num = 1
919+
shear_BC_flip_indices(1:2, 1) = shear_indices(1)
920+
! Both x-dir and y-dir: flip tau_xy only
921+
else if (num_dims == 3) then
922+
shear_num = 3
923+
shear_indices(1:3) = stress_idx%beg - 1 + (/2, 4, 5/)
924+
shear_BC_flip_num = 2
925+
shear_BC_flip_indices(1, 1:2) = shear_indices((/1, 2/))
926+
shear_BC_flip_indices(2, 1:2) = shear_indices((/1, 3/))
927+
shear_BC_flip_indices(3, 1:2) = shear_indices((/2, 3/))
928+
! x-dir: flip tau_xy and tau_xz
929+
! y-dir: flip tau_xy and tau_yz
930+
! z-dir: flip tau_xz and tau_yz
931+
end if
932+
!$acc update device(shear_num, shear_indices, shear_BC_flip_num, shear_BC_flip_indices)
898933
end if
899934
900935
if (hyperelasticity) then
@@ -933,6 +968,30 @@ contains
933968
stress_idx%end = sys_size + (num_dims*(num_dims + 1))/2
934969
! number of stresses is 1 in 1D, 3 in 2D, 6 in 3D
935970
sys_size = stress_idx%end
971+
972+
! shear stress index is 2 for 2D and 2,4,5 for 3D
973+
print *, 'DEBUG DEBUG'
974+
if (num_dims == 1) then
975+
shear_num = 0
976+
else if (num_dims == 2) then
977+
print *, 'HERE'
978+
shear_num = 1
979+
shear_indices(1) = stress_idx%beg - 1 + 2
980+
shear_BC_flip_num = 1
981+
shear_BC_flip_indices(1:2, 1) = shear_indices(1)
982+
! Both x-dir and y-dir: flip tau_xy only
983+
else if (num_dims == 3) then
984+
shear_num = 3
985+
shear_indices(1:3) = stress_idx%beg - 1 + (/2, 4, 5/)
986+
shear_BC_flip_num = 2
987+
shear_BC_flip_indices(1, 1:2) = shear_indices((/1, 2/))
988+
shear_BC_flip_indices(2, 1:2) = shear_indices((/1, 3/))
989+
shear_BC_flip_indices(3, 1:2) = shear_indices((/2, 3/))
990+
! x-dir: flip tau_xy and tau_xz
991+
! y-dir: flip tau_xy and tau_yz
992+
! z-dir: flip tau_xz and tau_yz
993+
end if
994+
!$acc update device(shear_num, shear_indices, shear_BC_flip_num, shear_BC_flip_indices)
936995
end if
937996
938997
if (hyperelasticity) then

0 commit comments

Comments
 (0)