Skip to content

Commit 5d977c3

Browse files
Hyeoksu LeeHyeoksu Lee
authored andcommitted
integrating mixture routines
1 parent 77ddd1e commit 5d977c3

File tree

5 files changed

+41
-244
lines changed

5 files changed

+41
-244
lines changed

src/common/m_checker_common.fpp

Lines changed: 2 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -312,30 +312,22 @@ contains
312312
integer :: i
313313

314314
bub_fac = 0
315-
if (bubbles_euler .and. (num_fluids == 1)) bub_fac = 1
315+
if (bubbles_euler) bub_fac = 1
316316

317-
do i = 1, num_fluids
317+
do i = 1, num_fluids + bub_fac
318318
call s_int_to_str(i, iStr)
319319
@:PROHIBIT(.not. f_is_default(fluid_pp(i)%gamma) .and. fluid_pp(i)%gamma <= 0._wp, &
320320
"fluid_pp("//trim(iStr)//")%gamma must be positive")
321321

322322
@:PROHIBIT(model_eqns == 1 .and. (.not. f_is_default(fluid_pp(i)%gamma)), &
323323
"model_eqns = 1 does not support fluid_pp("//trim(iStr)//")%gamma")
324324

325-
@:PROHIBIT((i <= num_fluids + bub_fac .and. fluid_pp(i)%gamma <= 0._wp) .or. &
326-
(i > num_fluids + bub_fac .and. (.not. f_is_default(fluid_pp(i)%gamma))), &
327-
"for fluid_pp("//trim(iStr)//")%gamma")
328-
329325
@:PROHIBIT(.not. f_is_default(fluid_pp(i)%pi_inf) .and. fluid_pp(i)%pi_inf < 0._wp, &
330326
"fluid_pp("//trim(iStr)//")%pi_inf must be non-negative")
331327

332328
@:PROHIBIT(model_eqns == 1 .and. (.not. f_is_default(fluid_pp(i)%pi_inf)), &
333329
"model_eqns = 1 does not support fluid_pp("//trim(iStr)//")%pi_inf")
334330

335-
@:PROHIBIT((i <= num_fluids + bub_fac .and. fluid_pp(i)%pi_inf < 0._wp) .or. &
336-
(i > num_fluids + bub_fac .and. (.not. f_is_default(fluid_pp(i)%pi_inf))), &
337-
"for fluid_pp("//trim(iStr)//")%pi_inf")
338-
339331
@:PROHIBIT(fluid_pp(i)%cv < 0._wp, &
340332
"fluid_pp("//trim(iStr)//")%cv must be positive")
341333
end do

src/common/m_variables_conversion.fpp

Lines changed: 38 additions & 224 deletions
Original file line numberDiff line numberDiff line change
@@ -33,8 +33,6 @@ module m_variables_conversion
3333
s_initialize_mv, &
3434
s_convert_to_mixture_variables, &
3535
s_convert_mixture_to_mixture_variables, &
36-
s_convert_species_to_mixture_variables_bubbles, &
37-
s_convert_species_to_mixture_variables_bubbles_acc, &
3836
s_convert_species_to_mixture_variables, &
3937
s_convert_species_to_mixture_variables_acc, &
4038
s_convert_conservative_to_primitive_variables, &
@@ -93,11 +91,7 @@ contains
9391
call s_convert_mixture_to_mixture_variables(q_vf, i, j, k, &
9492
rho, gamma, pi_inf, qv)
9593

96-
else if (bubbles_euler) then
97-
call s_convert_species_to_mixture_variables_bubbles(q_vf, i, j, k, &
98-
rho, gamma, pi_inf, qv, Re_K)
99-
else
100-
! Volume fraction model
94+
else ! Volume fraction model
10195
call s_convert_species_to_mixture_variables(q_vf, i, j, k, &
10296
rho, gamma, pi_inf, qv, Re_K, G_K, G)
10397
end if
@@ -230,129 +224,6 @@ contains
230224

231225
end subroutine s_convert_mixture_to_mixture_variables
232226

233-
!> This procedure is used alongside with the gamma/pi_inf
234-
!! model to transfer the density, the specific heat ratio
235-
!! function and liquid stiffness function from the vector
236-
!! of conservative or primitive variables to their scalar
237-
!! counterparts. Specifically designed for when subgrid bubbles_euler
238-
!! must be included.
239-
!! @param q_vf primitive variables
240-
!! @param j Cell index
241-
!! @param k Cell index
242-
!! @param l Cell index
243-
!! @param rho density
244-
!! @param gamma specific heat ratio
245-
!! @param pi_inf liquid stiffness
246-
!! @param qv fluid reference energy
247-
subroutine s_convert_species_to_mixture_variables_bubbles(q_vf, j, k, l, &
248-
rho, gamma, pi_inf, qv, Re_K)
249-
250-
type(scalar_field), dimension(sys_size), intent(in) :: q_vf
251-
252-
integer, intent(in) :: j, k, l
253-
254-
real(wp), intent(out), target :: rho
255-
real(wp), intent(out), target :: gamma
256-
real(wp), intent(out), target :: pi_inf
257-
real(wp), intent(out), target :: qv
258-
259-
real(wp), optional, dimension(2), intent(out) :: Re_K
260-
261-
integer :: i, q
262-
real(wp), dimension(num_fluids) :: alpha_rho_K, alpha_K
263-
264-
! Constraining the partial densities and the volume fractions within
265-
! their physical bounds to make sure that any mixture variables that
266-
! are derived from them result within the limits that are set by the
267-
! fluids physical parameters that make up the mixture
268-
do i = 1, num_fluids
269-
alpha_rho_K(i) = q_vf(i)%sf(j, k, l)
270-
alpha_K(i) = q_vf(advxb + i - 1)%sf(j, k, l)
271-
end do
272-
273-
if (mpp_lim) then
274-
275-
do i = 1, num_fluids
276-
alpha_rho_K(i) = max(0._wp, alpha_rho_K(i))
277-
alpha_K(i) = min(max(0._wp, alpha_K(i)), 1._wp)
278-
end do
279-
280-
alpha_K = alpha_K/max(sum(alpha_K), 1.e-16_wp)
281-
282-
end if
283-
284-
! Performing the transfer of the density, the specific heat ratio
285-
! function as well as the liquid stiffness function, respectively
286-
287-
if (model_eqns == 4) then
288-
rho = q_vf(1)%sf(j, k, l)
289-
gamma = fluid_pp(1)%gamma !qK_vf(gamma_idx)%sf(i,j,k)
290-
pi_inf = fluid_pp(1)%pi_inf !qK_vf(pi_inf_idx)%sf(i,j,k)
291-
qv = fluid_pp(1)%qv
292-
else if ((model_eqns == 2) .and. bubbles_euler) then
293-
rho = 0._wp; gamma = 0._wp; pi_inf = 0._wp; qv = 0._wp
294-
295-
if (mpp_lim .and. (num_fluids > 2)) then
296-
do i = 1, num_fluids
297-
rho = rho + q_vf(i)%sf(j, k, l)
298-
gamma = gamma + q_vf(i + E_idx)%sf(j, k, l)*fluid_pp(i)%gamma
299-
pi_inf = pi_inf + q_vf(i + E_idx)%sf(j, k, l)*fluid_pp(i)%pi_inf
300-
qv = qv + q_vf(i)%sf(j, k, l)*fluid_pp(i)%qv
301-
end do
302-
else if (num_fluids == 2) then
303-
rho = q_vf(1)%sf(j, k, l)
304-
gamma = fluid_pp(1)%gamma
305-
pi_inf = fluid_pp(1)%pi_inf
306-
qv = fluid_pp(1)%qv
307-
else if (num_fluids > 2) then
308-
!TODO: This may need fixing for hypo + bubbles_euler
309-
do i = 1, num_fluids - 1 !leave out bubble part of mixture
310-
rho = rho + q_vf(i)%sf(j, k, l)
311-
gamma = gamma + q_vf(i + E_idx)%sf(j, k, l)*fluid_pp(i)%gamma
312-
pi_inf = pi_inf + q_vf(i + E_idx)%sf(j, k, l)*fluid_pp(i)%pi_inf
313-
qv = qv + q_vf(i)%sf(j, k, l)*fluid_pp(i)%qv
314-
end do
315-
! rho = qK_vf(1)%sf(j,k,l)
316-
! gamma_K = fluid_pp(1)%gamma
317-
! pi_inf_K = fluid_pp(1)%pi_inf
318-
else
319-
rho = q_vf(1)%sf(j, k, l)
320-
gamma = fluid_pp(1)%gamma
321-
pi_inf = fluid_pp(1)%pi_inf
322-
qv = fluid_pp(1)%qv
323-
end if
324-
end if
325-
326-
#ifdef MFC_SIMULATION
327-
! Computing the shear and bulk Reynolds numbers from species analogs
328-
if (viscous) then
329-
if (num_fluids == 1) then ! need to consider case with num_fluids >= 2
330-
do i = 1, 2
331-
332-
Re_K(i) = dflt_real; if (Re_size(i) > 0) Re_K(i) = 0._wp
333-
334-
do q = 1, Re_size(i)
335-
Re_K(i) = (1 - alpha_K(Re_idx(i, q)))/fluid_pp(Re_idx(i, q))%Re(i) &
336-
+ Re_K(i)
337-
end do
338-
339-
Re_K(i) = 1._wp/max(Re_K(i), sgm_eps)
340-
341-
end do
342-
end if
343-
end if
344-
#endif
345-
346-
! Post process requires rho_sf/gamma_sf/pi_inf_sf/qv_sf to also be updated
347-
#ifdef MFC_POST_PROCESS
348-
rho_sf(j, k, l) = rho
349-
gamma_sf(j, k, l) = gamma
350-
pi_inf_sf(j, k, l) = pi_inf
351-
qv_sf(j, k, l) = qv
352-
#endif
353-
354-
end subroutine s_convert_species_to_mixture_variables_bubbles
355-
356227
!> This subroutine is designed for the volume fraction model
357228
!! and provided a set of either conservative or primitive
358229
!! variables, computes the density, the specific heat ratio
@@ -417,34 +288,43 @@ contains
417288
end do
418289

419290
alpha_K = alpha_K/max(sum(alpha_K), 1.e-16_wp)
420-
421291
end if
422292

423293
! Calculating the density, the specific heat ratio function, the
424294
! liquid stiffness function, and the energy reference function,
425295
! respectively, from the species analogs
426296
rho = 0._wp; gamma = 0._wp; pi_inf = 0._wp; qv = 0._wp
427297

428-
do i = 1, num_fluids
429-
rho = rho + alpha_rho_K(i)
430-
gamma = gamma + alpha_K(i)*gammas(i)
431-
pi_inf = pi_inf + alpha_K(i)*pi_infs(i)
432-
qv = qv + alpha_rho_K(i)*qvs(i)
433-
end do
298+
if (bubbles_euler) then
299+
rho = alpha_rho_K(1)
300+
gamma = gammas(1)
301+
pi_inf = pi_infs(1)
302+
qv = qvs(1)
303+
else
304+
do i = 1, num_fluids
305+
rho = rho + alpha_rho_K(i)
306+
gamma = gamma + alpha_K(i)*gammas(i)
307+
pi_inf = pi_inf + alpha_K(i)*pi_infs(i)
308+
qv = qv + alpha_rho_K(i)*qvs(i)
309+
end do
310+
end if
311+
434312
#ifdef MFC_SIMULATION
435313
! Computing the shear and bulk Reynolds numbers from species analogs
436-
do i = 1, 2
314+
if (viscous) then
315+
do i = 1, 2
437316

438-
Re_K(i) = dflt_real; if (Re_size(i) > 0) Re_K(i) = 0._wp
317+
Re_K(i) = dflt_real; if (Re_size(i) > 0) Re_K(i) = 0._wp
439318

440-
do j = 1, Re_size(i)
441-
Re_K(i) = alpha_K(Re_idx(i, j))/fluid_pp(Re_idx(i, j))%Re(i) &
442-
+ Re_K(i)
443-
end do
319+
do j = 1, Re_size(i)
320+
Re_K(i) = alpha_K(Re_idx(i, j))/fluid_pp(Re_idx(i, j))%Re(i) &
321+
+ Re_K(i)
322+
end do
444323

445-
Re_K(i) = 1._wp/max(Re_K(i), sgm_eps)
324+
Re_K(i) = 1._wp/max(Re_K(i), sgm_eps)
446325

447-
end do
326+
end do
327+
end if
448328
#endif
449329

450330
if (present(G_K)) then
@@ -504,15 +384,21 @@ contains
504384
end do
505385

506386
alpha_K = alpha_K/max(alpha_K_sum, sgm_eps)
507-
508387
end if
509388

510-
do i = 1, num_fluids
511-
rho_K = rho_K + alpha_rho_K(i)
512-
gamma_K = gamma_K + alpha_K(i)*gammas(i)
513-
pi_inf_K = pi_inf_K + alpha_K(i)*pi_infs(i)
514-
qv_K = qv_K + alpha_rho_K(i)*qvs(i)
515-
end do
389+
if (bubbles_euler) then
390+
rho_K = alpha_rho_K(1)
391+
gamma_K = gammas(1)
392+
pi_inf_K = pi_infs(1)
393+
qv_K = qvs(1)
394+
else
395+
do i = 1, num_fluids
396+
rho_K = rho_K + alpha_rho_K(i)
397+
gamma_K = gamma_K + alpha_K(i)*gammas(i)
398+
pi_inf_K = pi_inf_K + alpha_K(i)*pi_infs(i)
399+
qv_K = qv_K + alpha_rho_K(i)*qvs(i)
400+
end do
401+
end if
516402

517403
if (present(G_K)) then
518404
G_K = 0._wp
@@ -525,7 +411,6 @@ contains
525411
end if
526412

527413
if (viscous) then
528-
529414
do i = 1, 2
530415
Re_K(i) = dflt_real
531416

@@ -537,77 +422,12 @@ contains
537422
end do
538423

539424
Re_K(i) = 1._wp/max(Re_K(i), sgm_eps)
540-
541425
end do
542426
end if
543427
#endif
544428

545429
end subroutine s_convert_species_to_mixture_variables_acc
546430

547-
subroutine s_convert_species_to_mixture_variables_bubbles_acc(rho_K, &
548-
gamma_K, pi_inf_K, qv_K, &
549-
alpha_K, alpha_rho_K, Re_K)
550-
$:GPU_ROUTINE(function_name='s_convert_species_to_mixture_variables_bubbles_acc', &
551-
& parallelism='[seq]', cray_inline=True)
552-
553-
real(wp), intent(inout) :: rho_K, gamma_K, pi_inf_K, qv_K
554-
555-
real(wp), dimension(num_fluids), intent(in) :: alpha_K, alpha_rho_K !<
556-
!! Partial densities and volume fractions
557-
558-
real(wp), dimension(2), intent(out) :: Re_K
559-
560-
integer :: i, j !< Generic loop iterators
561-
562-
#ifdef MFC_SIMULATION
563-
rho_K = 0._wp
564-
gamma_K = 0._wp
565-
pi_inf_K = 0._wp
566-
qv_K = 0._wp
567-
568-
if (mpp_lim .and. (model_eqns == 2) .and. (num_fluids > 2)) then
569-
do i = 1, num_fluids
570-
rho_K = rho_K + alpha_rho_K(i)
571-
gamma_K = gamma_K + alpha_K(i)*gammas(i)
572-
pi_inf_K = pi_inf_K + alpha_K(i)*pi_infs(i)
573-
qv_K = qv_K + alpha_rho_K(i)*qvs(i)
574-
end do
575-
else if ((model_eqns == 2) .and. (num_fluids > 2)) then
576-
do i = 1, num_fluids - 1
577-
rho_K = rho_K + alpha_rho_K(i)
578-
gamma_K = gamma_K + alpha_K(i)*gammas(i)
579-
pi_inf_K = pi_inf_K + alpha_K(i)*pi_infs(i)
580-
qv_K = qv_K + alpha_rho_K(i)*qvs(i)
581-
end do
582-
else
583-
rho_K = alpha_rho_K(1)
584-
gamma_K = gammas(1)
585-
pi_inf_K = pi_infs(1)
586-
qv_K = qvs(1)
587-
end if
588-
589-
if (viscous) then
590-
if (num_fluids == 1) then ! need to consider case with num_fluids >= 2
591-
592-
do i = 1, 2
593-
Re_K(i) = dflt_real
594-
595-
if (Re_size(i) > 0) Re_K(i) = 0._wp
596-
597-
do j = 1, Re_size(i)
598-
Re_K(i) = (1._wp - alpha_K(Re_idx(i, j)))/Res_vc(i, j) &
599-
+ Re_K(i)
600-
end do
601-
602-
Re_K(i) = 1._wp/max(Re_K(i), sgm_eps)
603-
604-
end do
605-
end if
606-
end if
607-
#endif
608-
609-
end subroutine s_convert_species_to_mixture_variables_bubbles_acc
610-
611431
!> The computation of parameters, the allocation of memory,
612432
!! the association of pointers and/or the execution of any
613433
!! other procedures that are necessary to setup the module.
@@ -883,9 +703,6 @@ contains
883703
if (elasticity) then
884704
call s_convert_species_to_mixture_variables_acc(rho_K, gamma_K, pi_inf_K, qv_K, alpha_K, &
885705
alpha_rho_K, Re_K, G_K, Gs_vc)
886-
else if (bubbles_euler) then
887-
call s_convert_species_to_mixture_variables_bubbles_acc(rho_K, gamma_K, pi_inf_K, qv_K, &
888-
alpha_K, alpha_rho_K, Re_K)
889706
else
890707
call s_convert_species_to_mixture_variables_acc(rho_K, gamma_K, pi_inf_K, qv_K, &
891708
alpha_K, alpha_rho_K, Re_K)
@@ -1495,9 +1312,6 @@ contains
14951312
call s_convert_species_to_mixture_variables_acc(rho_K, gamma_K, pi_inf_K, qv_K, &
14961313
alpha_K, alpha_rho_K, Re_K, &
14971314
G_K, Gs_vc)
1498-
else if (bubbles_euler) then
1499-
call s_convert_species_to_mixture_variables_bubbles_acc(rho_K, gamma_K, &
1500-
pi_inf_K, qv_K, alpha_K, alpha_rho_K, Re_K)
15011315
else
15021316
call s_convert_species_to_mixture_variables_acc(rho_K, gamma_K, pi_inf_K, qv_K, &
15031317
alpha_K, alpha_rho_K, Re_K)

src/simulation/m_cbc.fpp

Lines changed: 1 addition & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -813,11 +813,7 @@ contains
813813
adv_local(i) = q_prim_rs${XYZ}$_vf(0, k, r, E_idx + i)
814814
end do
815815
816-
if (bubbles_euler) then
817-
call s_convert_species_to_mixture_variables_bubbles_acc(rho, gamma, pi_inf, qv, adv_local, alpha_rho, Re_cbc)
818-
else
819-
call s_convert_species_to_mixture_variables_acc(rho, gamma, pi_inf, qv, adv_local, alpha_rho, Re_cbc)
820-
end if
816+
call s_convert_species_to_mixture_variables_acc(rho, gamma, pi_inf, qv, adv_local, alpha_rho, Re_cbc)
821817
822818
$:GPU_LOOP(parallelism='[seq]')
823819
do i = 1, contxe

0 commit comments

Comments
 (0)