505 type(
scalar_field),
dimension(sys_size),
intent(in) :: q_prim_vf
506 integer,
intent(in) :: t_step
509# 178 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
510 real(wp),
dimension(num_fluids) :: alpha
511 real(wp),
dimension(num_vels) :: vel
512# 181 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
520 real(wp),
dimension(2) :: re
522 real(wp) :: icfl_max_loc, icfl_max_glb
523 real(wp) :: vcfl_max_loc, vcfl_max_glb
524 real(wp) :: ccfl_max_loc, ccfl_max_glb
525 real(wp) :: rc_min_loc, rc_min_glb
526 real(wp) :: icfl, vcfl, ccfl, rc
532 rc_min_loc = huge(1.0_wp)
535# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
537# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
538#if defined(MFC_OpenACC)
539# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
541# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
543# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
544#elif defined(MFC_OpenMP)
545# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
547# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
549# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
551# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
553# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
555# 202 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
557# 205 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
561 call s_compute_enthalpy(q_prim_vf, pres, rho, gamma, pi_inf, re, h, alpha, vel, vel_sum, qv,
j,
k,
l)
567 do fl = 1, num_fluids
574 re(1) = 1._wp/max(re(1),
sgm_eps)
577 call s_compute_stability_from_dt(vel, c, rho, re,
j,
k,
l, icfl, vcfl, rc, ccfl)
579 icfl_max_loc = max(icfl_max_loc, icfl)
580 vcfl_max_loc = max(vcfl_max_loc, merge(vcfl, 0.0_wp, viscous))
581 ccfl_max_loc = max(ccfl_max_loc, merge(ccfl, 0.0_wp, surface_tension))
582 rc_min_loc = min(rc_min_loc, merge(rc, huge(1.0_wp), viscous))
587# 233 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
588#if defined(MFC_OpenACC)
589# 233 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
591# 233 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
592#elif defined(MFC_OpenMP)
593# 233 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
595# 233 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
597# 233 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
602 call s_mpi_reduce_stability_criteria_extrema(icfl_max_loc, vcfl_max_loc, rc_min_loc,
n_el_bubs_loc, icfl_max_glb, &
603 & vcfl_max_glb, rc_min_glb,
n_el_bubs_glb, ccfl_max_loc, ccfl_max_glb)
605 icfl_max_glb = icfl_max_loc
606 if (viscous) vcfl_max_glb = vcfl_max_loc
607 if (viscous) rc_min_glb = rc_min_loc
608 if (surface_tension) ccfl_max_glb = ccfl_max_loc
614 if (surface_tension)
then
624 write (3,
'(13X,I9,13X,F10.6,13X,F10.6,13X,F10.6)', advance=
"no") t_step, dt,
mytime, icfl_max_glb
626 if (surface_tension)
then
627 write (3,
'(13X,F10.6)', advance=
"no") ccfl_max_glb
631 write (3,
'(13X,F10.6,13X,ES16.6)', advance=
"no") vcfl_max_glb, rc_min_glb
634 if (bubbles_lagrange)
then
641 call s_mpi_abort(
'ICFL is NaN. Exiting.')
642 else if (icfl_max_glb > 1._wp)
then
643 print *,
'icfl', icfl_max_glb
644 call s_mpi_abort(
'ICFL is greater than 1.0. Exiting.')
649 call s_mpi_abort(
'VCFL is NaN. Exiting.')
650 else if (vcfl_max_glb > 1._wp)
then
651 print *,
'vcfl', vcfl_max_glb
652 call s_mpi_abort(
'VCFL is greater than 1.0. Exiting.')
656 if (bubbles_lagrange)
then
658 call s_mpi_abort(
'No Lagrangian bubbles remain in the domain. Exiting.')
672 type(
scalar_field),
dimension(sys_size),
intent(inout) :: q_prim_vf
673 integer,
intent(in) :: t_step
675 type(
integer_field),
dimension(1:num_dims,-1:1),
intent(in) :: bc_type
676 character(LEN=path_len + 2*name_len) :: t_step_dir
677 character(LEN=path_len + 3*name_len) :: file_path
678 logical :: file_exist
679 character(LEN=15) :: fmt
680 integer :: i,
j,
k,
l, r
681 real(wp) :: gamma, lit_gamma, pi_inf, qv
683 write (t_step_dir,
'(A,I0,A,I0)') trim(case_dir) //
'/p_all'
684 write (t_step_dir,
'(a,i0,a,i0)') trim(case_dir) //
'/p_all/p',
proc_rank,
'/', t_step
686 file_path = trim(t_step_dir) //
'/.'
691 file_path = trim(t_step_dir) //
'/x_cb.dat'
693 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
694 write (2)
x_cb(-1:m);
close (2)
697 file_path = trim(t_step_dir) //
'/y_cb.dat'
699 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
700 write (2)
y_cb(-1:n);
close (2)
703 file_path = trim(t_step_dir) //
'/z_cb.dat'
705 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
706 write (2)
z_cb(-1:p);
close (2)
711 write (file_path,
'(A,I0,A)') trim(t_step_dir) //
'/q_cons_vf', i,
'.dat'
713 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
715 write (2)
q_cons_vf(i)%sf(0:m,0:n,0:p);
close (2)
720 if (bubbles_lagrange)
then
721 write (file_path,
'(A,I0,A)') trim(t_step_dir) //
'/q_cons_vf', sys_size + 1,
'.dat'
723 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
725 write (2) beta%sf(0:m,0:n,0:p);
close (2)
728 if (qbmm .and. .not. polytropic)
then
731 write (file_path,
'(A,I0,A)') trim(t_step_dir) //
'/pb', sys_size + (i - 1)*
nnode + r,
'.dat'
733 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
735 write (2)
pb_ts(1)%sf(0:m,0:n,0:p,r, i);
close (2)
741 write (file_path,
'(A,I0,A)') trim(t_step_dir) //
'/mv', sys_size + (i - 1)*
nnode + r,
'.dat'
743 open (2, file=trim(file_path), form=
'unformatted', status=
'new')
745 write (2)
mv_ts(1)%sf(0:m,0:n,0:p,r, i);
close (2)
756 lit_gamma = gs_min(1)
766 write (t_step_dir,
'(A,I0,A,I0)') trim(case_dir) //
'/D'
767 file_path = trim(t_step_dir) //
'/.'
769 inquire (file=trim(file_path), exist=file_exist)
773 if ((prim_vars_wrt .or. (n == 0 .and. p == 0)) .and. (.not. igr))
then
777# 411 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
778#if defined(MFC_OpenACC)
779# 411 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
781# 411 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
782#elif defined(MFC_OpenMP)
783# 411 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
785# 411 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
791 q_prim_vf(eqn_idx%bub%beg)%sf = 1._wp
795 if (n == 0 .and. p == 0)
then
798 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/prim.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
800 open (2, file=trim(file_path))
803 if (((i >= eqn_idx%adv%beg) .and. (i <= eqn_idx%adv%end)))
then
806 write (2, fmt)
x_cb(
j), q_prim_vf(i)%sf(
j, 0, 0)
814 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/cons.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
816 open (2, file=trim(file_path))
823 if (qbmm .and. .not. polytropic)
then
826 write (file_path,
'(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/pres.', i,
'.', r,
'.',
proc_rank, &
827 &
'.', t_step,
'.dat'
829 open (2, file=trim(file_path))
838 write (file_path,
'(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/mv.', i,
'.', r,
'.',
proc_rank, &
839 &
'.', t_step,
'.dat'
841 open (2, file=trim(file_path))
857 if ((n > 0) .and. (p == 0))
then
859 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/cons.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
860 open (2, file=trim(file_path))
870 if (
present(beta))
then
871 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/beta.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
872 open (2, file=trim(file_path))
882 if (qbmm .and. .not. polytropic)
then
885 write (file_path,
'(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/pres.', i,
'.', r,
'.',
proc_rank, &
886 &
'.', t_step,
'.dat'
888 open (2, file=trim(file_path))
899 write (file_path,
'(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/mv.', i,
'.', r,
'.',
proc_rank, &
900 &
'.', t_step,
'.dat'
902 open (2, file=trim(file_path))
913 if (prim_vars_wrt .and. (.not. igr))
then
915 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/prim.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
917 open (2, file=trim(file_path))
921 if (((i >= eqn_idx%cont%beg) .and. (i <= eqn_idx%cont%end)) .or. ((i >= eqn_idx%adv%beg) &
922 & .and. (i <= eqn_idx%adv%end)))
then
925 write (2, fmt)
x_cb(
j),
y_cb(
k), q_prim_vf(i)%sf(
j,
k, 0)
943 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/cons.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
944 open (2, file=trim(file_path))
957 if (
present(beta))
then
958 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/beta.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
959 open (2, file=trim(file_path))
972 if (qbmm .and. .not. polytropic)
then
975 write (file_path,
'(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/pres.', i,
'.', r,
'.',
proc_rank, &
976 &
'.', t_step,
'.dat'
978 open (2, file=trim(file_path))
991 write (file_path,
'(A,I0,A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/mv.', i,
'.', r,
'.',
proc_rank, &
992 &
'.', t_step,
'.dat'
994 open (2, file=trim(file_path))
1007 if (prim_vars_wrt .and. (.not. igr))
then
1009 write (file_path,
'(A,I0,A,I2.2,A,I6.6,A)') trim(t_step_dir) //
'/prim.', i,
'.',
proc_rank,
'.', t_step,
'.dat'
1011 open (2, file=trim(file_path))
1016 if (((i >= eqn_idx%cont%beg) .and. (i <= eqn_idx%cont%end)) .or. ((i >= eqn_idx%adv%beg) &
1017 & .and. (i <= eqn_idx%adv%end)) .or. ((i >= eqn_idx%species%beg) &
1018 & .and. (i <= eqn_idx%species%end)))
then
1039 integer,
intent(in) :: t_step
1041 type(
integer_field),
dimension(1:num_dims,-1:1),
intent(in) :: bc_type
1045 integer :: ifile, ierr, data_size
1046 integer,
dimension(MPI_STATUS_SIZE) :: status
1047 integer(kind=MPI_OFFSET_kind) :: disp
1048 integer(kind=MPI_OFFSET_kind) :: m_mok, n_mok, p_mok
1049 integer(kind=MPI_OFFSET_kind) :: wp_mok, var_mok, str_mok
1050 integer(kind=MPI_OFFSET_kind) :: nvars_mok
1051 integer(kind=MPI_OFFSET_kind) :: mok
1052 character(LEN=path_len + 2*name_len) :: file_loc
1053 logical :: file_exist, dir_check
1054 character(len=10) :: t_step_string
1058 integer :: m_ds, n_ds, p_ds
1059 integer :: m_glb_ds, n_glb_ds, p_glb_ds
1060 integer :: m_glb_save, n_glb_save, p_glb_save
1062 if (down_sample)
then
1066 if (
present(beta))
then
1067 alt_sys = sys_size + 1
1072 if (file_per_process)
then
1075 if (down_sample)
then
1076 call s_initialize_mpi_data_ds(m_ds, n_ds, p_ds)
1087 file_loc = trim(case_dir) //
'/restart_data/lustre_' // trim(t_step_string)
1089 if (dir_check .neqv. .true.)
then
1094 call s_mpi_barrier()
1099 write (file_loc,
'(I0,A,i7.7,A)') t_step,
'_',
proc_rank,
'.dat'
1100 file_loc = trim(case_dir) //
'/restart_data/lustre_' // trim(t_step_string) // trim(mpiiofs) // trim(file_loc)
1101 inquire (file=trim(file_loc), exist=file_exist)
1102 if (file_exist .and.
proc_rank == 0)
then
1103 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1105 call mpi_file_open(mpi_comm_self, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1107 if (down_sample)
then
1108 data_size = (m_ds + 3)*(n_ds + 3)*(p_ds + 3)
1109 m_glb_save = m_glb_ds + 1
1110 n_glb_save = n_glb_ds + 1
1111 p_glb_save = p_glb_ds + 1
1113 data_size = (m + 1)*(n + 1)*(p + 1)
1114 m_glb_save =
m_glb + 1
1115 n_glb_save =
n_glb + 1
1116 p_glb_save =
p_glb + 1
1119 m_mok = int(m_glb_save + 1, mpi_offset_kind)
1120 n_mok = int(n_glb_save + 1, mpi_offset_kind)
1121 p_mok = int(p_glb_save + 1, mpi_offset_kind)
1122 wp_mok = int(storage_size(0._stp)/8, mpi_offset_kind)
1123 mok = int(1._wp, mpi_offset_kind)
1124 str_mok = int(
name_len, mpi_offset_kind)
1125 nvars_mok = int(sys_size, mpi_offset_kind)
1127 if (bubbles_euler)
then
1129 var_mok = int(i, mpi_offset_kind)
1131 call mpi_file_write_all(ifile,
mpi_io_data%var(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1133 if (qbmm .and. .not. polytropic)
then
1134 do i = sys_size + 1, sys_size + 2*
nb*
nnode
1135 var_mok = int(i, mpi_offset_kind)
1137 call mpi_file_write_all(ifile,
mpi_io_data%var(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1141 if (down_sample)
then
1143 var_mok = int(i, mpi_offset_kind)
1145 call mpi_file_write_all(ifile,
q_cons_temp_ds(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1149 var_mok = int(i, mpi_offset_kind)
1151 call mpi_file_write_all(ifile,
mpi_io_data%var(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1156 call mpi_file_close(ifile, ierr)
1165 else if (
present(beta))
then
1171 write (file_loc,
'(I0,A)') t_step,
'.dat'
1172 file_loc = trim(case_dir) //
'/restart_data' // trim(mpiiofs) // trim(file_loc)
1173 inquire (file=trim(file_loc), exist=file_exist)
1174 if (file_exist .and.
proc_rank == 0)
then
1175 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1177 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1179 data_size = (m + 1)*(n + 1)*(p + 1)
1181 m_mok = int(
m_glb + 1, mpi_offset_kind)
1182 n_mok = int(
n_glb + 1, mpi_offset_kind)
1183 p_mok = int(
p_glb + 1, mpi_offset_kind)
1184 wp_mok = int(storage_size(0._stp)/8, mpi_offset_kind)
1185 mok = int(1._wp, mpi_offset_kind)
1186 str_mok = int(
name_len, mpi_offset_kind)
1187 nvars_mok = int(alt_sys, mpi_offset_kind)
1189 if (bubbles_euler)
then
1191 var_mok = int(i, mpi_offset_kind)
1193 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1195 call mpi_file_set_view(ifile, disp, mpi_p,
mpi_io_data%view(i),
'native', mpi_info_int, ierr)
1196 call mpi_file_write_all(ifile,
mpi_io_data%var(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1198 if (qbmm .and. .not. polytropic)
then
1199 do i = sys_size + 1, sys_size + 2*
nb*
nnode
1200 var_mok = int(i, mpi_offset_kind)
1202 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1204 call mpi_file_set_view(ifile, disp, mpi_p,
mpi_io_data%view(i),
'native', mpi_info_int, ierr)
1205 call mpi_file_write_all(ifile,
mpi_io_data%var(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1210 var_mok = int(i, mpi_offset_kind)
1212 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1214 call mpi_file_set_view(ifile, disp, mpi_p,
mpi_io_data%view(i),
'native', mpi_info_int, ierr)
1215 call mpi_file_write_all(ifile,
mpi_io_data%var(i)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1219 if (
present(beta))
then
1220 var_mok = int(sys_size + 1, mpi_offset_kind)
1222 disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
1224 call mpi_file_set_view(ifile, disp, mpi_p,
mpi_io_data%view(sys_size + 1),
'native', mpi_info_int, ierr)
1225 call mpi_file_write_all(ifile,
mpi_io_data%var(sys_size + 1)%sf, data_size*mpi_io_type, mpi_io_p, status, ierr)
1228 call mpi_file_close(ifile, ierr)
1341 integer,
intent(in) :: t_step
1344 character(LEN=path_len + 2*name_len) :: file_loc
1345 integer(kind=MPI_OFFSET_KIND) :: disp
1346 integer(kind=MPI_OFFSET_KIND) :: WP_MOK
1347 integer :: ifile, ierr
1348 integer,
dimension(MPI_STATUS_SIZE) :: status
1349 logical :: file_exist, dir_check
1350 integer :: i, ib_idx
1351 integer,
parameter :: NFIELDS_PER_IB = 20
1352 real(wp) :: ib_buf(NFIELDS_PER_IB)
1353 integer :: file_unit
1354 character(len=10) :: t_step_string
1357 integer :: ib_start, ib_end, nibs_per_rank, remainder
1359 wp_mok = int(storage_size(0._wp)/8, mpi_offset_kind)
1361 if (file_per_process)
then
1365 file_loc = trim(case_dir) //
'/restart_data/lustre_' // trim(t_step_string)
1368 call s_mpi_barrier()
1371 write (file_loc,
'(A,I0,A,i7.7,A)')
'ib_state_', t_step,
'_',
proc_rank,
'.dat'
1372 file_loc = trim(case_dir) //
'/restart_data/lustre_' // trim(t_step_string) //
'/' // trim(file_loc)
1374 inquire (file=trim(file_loc), exist=file_exist)
1375 if (file_exist)
then
1376 open (newunit=file_unit, file=trim(file_loc), form=
'unformatted', access=
'stream', status=
'replace')
1378 open (newunit=file_unit, file=trim(file_loc), form=
'unformatted', access=
'stream', status=
'new')
1385 ib_buf(2:4) = patch_ib(ib_idx)%force(1:3)
1386 ib_buf(5:7) = patch_ib(ib_idx)%torque(1:3)
1387 ib_buf(8:10) = patch_ib(ib_idx)%vel(1:3)
1388 ib_buf(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1389 ib_buf(14:16) = patch_ib(ib_idx)%angles(1:3)
1390 ib_buf(17) = patch_ib(ib_idx)%x_centroid
1391 ib_buf(18) = patch_ib(ib_idx)%y_centroid
1392 ib_buf(19) = patch_ib(ib_idx)%z_centroid
1393 ib_buf(20) = patch_ib(ib_idx)%radius
1395 write (file_unit) patch_ib(ib_idx)%gbl_patch_id
1396 write (file_unit) ib_buf
1404 call s_mpi_barrier()
1406 write (file_loc,
'(A,I0,A)')
'/restart_data/ib_state_', t_step,
'.dat'
1407 file_loc = trim(case_dir) // trim(file_loc)
1409 inquire (file=trim(file_loc), exist=file_exist)
1410 if (file_exist .and.
proc_rank == 0)
then
1411 call mpi_file_delete(file_loc, mpi_info_int, ierr)
1413 call s_mpi_barrier()
1415 call mpi_file_open(mpi_comm_world, file_loc, ior(mpi_mode_wronly, mpi_mode_create), mpi_info_int, ifile, ierr)
1420 ib_buf(2:4) = patch_ib(ib_idx)%force(1:3)
1421 ib_buf(5:7) = patch_ib(ib_idx)%torque(1:3)
1422 ib_buf(8:10) = patch_ib(ib_idx)%vel(1:3)
1423 ib_buf(11:13) = patch_ib(ib_idx)%angular_vel(1:3)
1424 ib_buf(14:16) = patch_ib(ib_idx)%angles(1:3)
1425 ib_buf(17) = patch_ib(ib_idx)%x_centroid
1426 ib_buf(18) = patch_ib(ib_idx)%y_centroid
1427 ib_buf(19) = patch_ib(ib_idx)%z_centroid
1428 ib_buf(20) = patch_ib(ib_idx)%radius
1431 disp = int(patch_ib(ib_idx)%gbl_patch_id - 1, mpi_offset_kind)*int(nfields_per_ib, mpi_offset_kind)*wp_mok
1433 call mpi_file_write_at(ifile, disp, ib_buf, nfields_per_ib, mpi_p, status, ierr)
1436 call mpi_file_close(ifile, ierr)
1539 integer,
intent(in) :: t_step
1541 real(wp),
dimension(0:m,0:n,0:p),
intent(in) :: accel_mag
1542 real(wp),
dimension(-1:m) :: distx
1543 real(wp),
dimension(-1:n) :: disty
1544 real(wp),
dimension(-1:p) :: distz
1548 real(wp) :: lit_gamma, nbub
1550 real(wp),
dimension(num_vels) :: vel
1556 real(wp),
dimension(num_fluids) :: alpha
1561 real(wp) :: m00, m10, m01, m20, m02
1562 real(wp) :: varr, varv
1563 real(wp),
dimension(Nb) :: nr, r, nrdot, rdot
1566 real(wp) :: int_pres
1567 real(wp) :: max_pres
1568 real(wp),
dimension(2) :: re
1569 real(wp),
dimension(6) :: tau_e
1571 real(wp) :: dyn_p, t
1572 real(wp) :: damage_state
1573 integer :: i,
j,
k,
l, s, d
1574 real(wp) :: nondim_time
1577 real(wp) :: rad, thickness
1579 real(wp) :: rhoyks(1:num_species)
1583 if (time_stepper == 23)
then
1587 nondim_time = real(t_step + t_step_old, wp)*dt
1589 nondim_time = real(t_step, wp)*dt
1593 do i = 1, num_probes
1604 nr = 0._wp; r = 0._wp
1605 nrdot = 0._wp; rdot = 0._wp
1612 varr = 0._wp; varv = 0._wp
1614 do s = 1, (num_dims*(num_dims + 1))/2
1617 damage_state = 0._wp
1620 if ((probe(i)%x >=
x_cb(-1)) .and. (probe(i)%x <=
x_cb(m)))
then
1622 distx(s) =
x_cb(s) - probe(i)%x
1623 if (distx(s) < 0._wp) distx(s) = 1000._wp
1625 j = minloc(distx, 1)
1631 do d = 1, num_species
1632 rhoyks(d) =
q_cons_vf(eqn_idx%species%beg + d - 1)%sf(
j - 2,
k,
l)
1637 if (hypoelasticity)
then
1638 call s_convert_to_mixture_variables(
q_cons_vf,
j - 2,
k,
l, rho, gamma, pi_inf, qv, re, g_local, &
1644 vel(s) =
q_cons_vf(eqn_idx%cont%end + s)%sf(
j - 2,
k,
l)/rho
1646 do s = 1, num_fluids
1647 alpha(s) =
q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(
j - 2,
k,
l)
1650 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
1652 if (hypoelasticity)
then
1653 if (cont_damage)
then
1654 damage_state =
q_cons_vf(eqn_idx%damage)%sf(
j - 2,
k,
l)
1655 g_local = g_local*max((1._wp - damage_state), 0._wp)
1659 & dyn_p, pi_inf, gamma, rho, qv, rhoyks(:), pres, t, &
1664 & dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t)
1667 if (hypoelasticity)
then
1668 tau_e(1) =
q_cons_vf(eqn_idx%stress%end)%sf(
j - 2,
k,
l)/rho
1671 if (bubbles_euler)
then
1673 if (num_fluids == 3)
then
1686 nr3 = nr3 +
weight(s)*(nr(s)**3._wp)
1689 nbub = sqrt((4._wp*
pi/3._wp)*nr3/alf)
1692 print *,
'In probe, nbub: ', nbub
1706 varr = m20 - m10**2._wp
1707 varv = m02 - m01**2._wp
1710 rdot(:) = nrdot(:)/nbub
1713 ptot = pres - ptilde
1720 accel = accel_mag(
j - 2,
k,
l)
1722 else if (p == 0)
then
1724 do d = 1, num_species
1725 rhoyks(d) =
q_cons_vf(eqn_idx%species%beg + d - 1)%sf(
j - 2,
k - 2,
l)
1729 if ((probe(i)%x >=
x_cb(-1)) .and. (probe(i)%x <=
x_cb(m)))
then
1730 if ((probe(i)%y >=
y_cb(-1)) .and. (probe(i)%y <=
y_cb(n)))
then
1732 distx(s) =
x_cb(s) - probe(i)%x
1733 if (distx(s) < 0._wp) distx(s) = 1000._wp
1736 disty(s) =
y_cb(s) - probe(i)%y
1737 if (disty(s) < 0._wp) disty(s) = 1000._wp
1739 j = minloc(distx, 1)
1740 k = minloc(disty, 1)
1746 call s_convert_to_mixture_variables(
q_cons_vf,
j - 2,
k - 2,
l, rho, gamma, pi_inf, qv, re, g_local, &
1749 vel(s) =
q_cons_vf(eqn_idx%cont%end + s)%sf(
j - 2,
k - 2,
l)/rho
1751 do s = 1, num_fluids
1752 alpha(s) =
q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(
j - 2,
k - 2,
l)
1755 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
1757 if (hypoelasticity)
then
1758 if (cont_damage)
then
1759 damage_state =
q_cons_vf(eqn_idx%damage)%sf(
j - 2,
k - 2,
l)
1760 g_local = g_local*max((1._wp - damage_state), 0._wp)
1764 &
k - 2,
l), dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t, &
1766 &
q_cons_vf(eqn_idx%mom%beg)%sf(
j - 2,
k - 2,
l), g_local)
1769 &
k - 2,
l), dyn_p, pi_inf, gamma, rho, qv, rhoyks, pres, t)
1772 if (hypoelasticity)
then
1774 tau_e(s) =
q_cons_vf(eqn_idx%stress%beg + s - 1)%sf(
j - 2,
k - 2,
l)/rho
1778 if (bubbles_euler)
then
1790 nr3 = nr3 +
weight(s)*(nr(s)**3._wp)
1793 nbub = sqrt((4._wp*
pi/3._wp)*nr3/alf)
1797 rdot(:) = nrdot(:)/nbub
1801 & 0._wp, 0._wp, c, qv)
1805 if ((probe(i)%x >=
x_cb(-1)) .and. (probe(i)%x <=
x_cb(m)))
then
1806 if ((probe(i)%y >=
y_cb(-1)) .and. (probe(i)%y <=
y_cb(n)))
then
1807 if ((probe(i)%z >=
z_cb(-1)) .and. (probe(i)%z <=
z_cb(p)))
then
1809 distx(s) =
x_cb(s) - probe(i)%x
1810 if (distx(s) < 0._wp) distx(s) = 1000._wp
1813 disty(s) =
y_cb(s) - probe(i)%y
1814 if (disty(s) < 0._wp) disty(s) = 1000._wp
1817 distz(s) =
z_cb(s) - probe(i)%z
1818 if (distz(s) < 0._wp) distz(s) = 1000._wp
1820 j = minloc(distx, 1)
1821 k = minloc(disty, 1)
1822 l = minloc(distz, 1)
1829 & g_local, fluid_pp(:)%G)
1831 vel(s) =
q_cons_vf(eqn_idx%cont%end + s)%sf(
j - 2,
k - 2,
l - 2)/rho
1833 do s = 1, num_fluids
1834 alpha(s) =
q_cons_vf(eqn_idx%adv%beg + s - 1)%sf(
j - 2,
k - 2,
l - 2)
1837 dyn_p = 0.5_wp*rho*dot_product(vel, vel)
1840 do d = 1, num_species
1841 rhoyks(d) =
q_cons_vf(eqn_idx%species%beg + d - 1)%sf(
j - 2,
k - 2,
l - 2)
1845 if (hypoelasticity)
then
1846 if (cont_damage)
then
1847 damage_state =
q_cons_vf(eqn_idx%damage)%sf(
j - 2,
k - 2,
l - 2)
1848 g_local = g_local*max((1._wp - damage_state), 0._wp)
1852 &
q_cons_vf(eqn_idx%alf)%sf(
j - 2,
k - 2,
l - 2), dyn_p, pi_inf, gamma, &
1853 & rho, qv, rhoyks, pres, t,
q_cons_vf(eqn_idx%stress%beg)%sf(
j - 2, &
1854 &
k - 2,
l - 2),
q_cons_vf(eqn_idx%mom%beg)%sf(
j - 2,
k - 2,
l - 2), &
1858 &
q_cons_vf(eqn_idx%alf)%sf(
j - 2,
k - 2,
l - 2), dyn_p, pi_inf, gamma, &
1859 & rho, qv, rhoyks, pres, t)
1864 & 0._wp, 0._wp, c, qv)
1866 accel = accel_mag(
j - 2,
k - 2,
l - 2)
1872# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1874 call s_mpi_allreduce_sum(tmp, rho)
1875# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1877 call s_mpi_allreduce_sum(tmp, pres)
1878# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1880 call s_mpi_allreduce_sum(tmp, gamma)
1881# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1883 call s_mpi_allreduce_sum(tmp, pi_inf)
1884# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1886 call s_mpi_allreduce_sum(tmp, qv)
1887# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1889 call s_mpi_allreduce_sum(tmp, c)
1890# 1458 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1892 call s_mpi_allreduce_sum(tmp, accel)
1893# 1461 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1897 call s_mpi_allreduce_sum(tmp, vel(s))
1900 if (bubbles_euler)
then
1901# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1903 call s_mpi_allreduce_sum(tmp, alf)
1904# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1906 call s_mpi_allreduce_sum(tmp, alfgr)
1907# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1909 call s_mpi_allreduce_sum(tmp, nbub)
1910# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1912 call s_mpi_allreduce_sum(tmp, nr(1))
1913# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1915 call s_mpi_allreduce_sum(tmp, nrdot(1))
1916# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1918 call s_mpi_allreduce_sum(tmp, m00)
1919# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1921 call s_mpi_allreduce_sum(tmp, r(1))
1922# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1924 call s_mpi_allreduce_sum(tmp, rdot(1))
1925# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1927 call s_mpi_allreduce_sum(tmp, ptilde)
1928# 1469 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1930 call s_mpi_allreduce_sum(tmp, ptot)
1931# 1472 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1934# 1475 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1936 call s_mpi_allreduce_sum(tmp, varr)
1937# 1475 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1939 call s_mpi_allreduce_sum(tmp, varv)
1940# 1475 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1942 call s_mpi_allreduce_sum(tmp, m10)
1943# 1475 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1945 call s_mpi_allreduce_sum(tmp, m01)
1946# 1475 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1948 call s_mpi_allreduce_sum(tmp, m20)
1949# 1475 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1951 call s_mpi_allreduce_sum(tmp, m02)
1952# 1478 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1956 if (hypoelasticity)
then
1957 do s = 1, (num_dims*(num_dims + 1))/2
1959 call s_mpi_allreduce_sum(tmp, tau_e(s))
1963 if (cont_damage)
then
1965 call s_mpi_allreduce_sum(tmp, damage_state)
1970 if (bubbles_euler .and. (num_fluids <= 2))
then
1972 write (i + 30,
'(6x,f12.6,14f28.16)') nondim_time, rho, vel(1), pres, alf, r(1), rdot(1), nr(1), &
1973 & nrdot(1), varr, varv, m10, m01, m20, m02
1975 write (i + 30,
'(6x,f12.6,8f24.8)') nondim_time, rho, vel(1), pres, alf, r(1), rdot(1), nr(1), nrdot(1)
1978 else if (bubbles_euler .and. (num_fluids == 3))
then
1980 &
'(6x,f12.6,f24.8,f24.8,f24.8,f24.8,f24.8,' //
'f24.8,f24.8,f24.8,f24.8,f24.8, f24.8)') &
1981 & nondim_time, rho, vel(1), pres, alf, alfgr, nr(1), nrdot(1), r(1), rdot(1), ptilde, ptot
1982 else if (bubbles_euler .and. num_fluids == 4)
then
1984 &
'(6x,f12.6,f24.8,f24.8,f24.8,f24.8,' //
'f24.8,f24.8,f24.8,f24.8,f24.8,f24.8,f24.8,f24.8,f24.8)') &
1986 & 0, 0),
q_cons_vf(4)%sf(
j - 2, 0, 0),
q_cons_vf(5)%sf(
j - 2, 0, 0),
q_cons_vf(6)%sf(
j - 2, 0, 0), &
1988 &
q_cons_vf(10)%sf(
j - 2, 0, 0), nbub, r(1), rdot(1)
1990 write (i + 30,
'(6X,F12.6,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), pres
1992 else if (p == 0)
then
1993 if (bubbles_euler)
then
1994# 1520 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1995 write (i + 30,
'(6X,10F24.8)') nondim_time, rho, vel(1), vel(2), pres, alf, nr(1), nrdot(1), r(1), &
1997# 1523 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
1998 else if (hypoelasticity)
then
1999# 1525 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2000 write (i + 30,
'(6X,F12.6,F24.8,F24.8,F24.8,F24.8,' //
'F24.8,F24.8,F24.8)') nondim_time, rho, &
2001 & vel(1), vel(2), pres, tau_e(1), tau_e(2), tau_e(3)
2002# 1528 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2004 write (i + 30,
'(6X,F12.6,F24.8,F24.8,F24.8)') nondim_time, rho, vel(1), pres
2005 print *,
'time =', nondim_time,
'rho =', rho,
'pres =', pres
2008# 1534 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2010 &
'(6X,F12.6,F24.8,F24.8,F24.8,F24.8,' //
'F24.8,F24.8,F24.8,F24.8,F24.8,' //
'F24.8)') &
2011 & nondim_time, rho, vel(1), vel(2), vel(3), pres, gamma, pi_inf, qv, c, accel
2012# 1538 "/home/runner/work/MFC/MFC/src/simulation/m_data_output.fpp"
2017 if (integral_wrt .and. bubbles_euler)
then
2019 do i = 1, num_integrals
2035 if ((integral(i)%xmin <=
x_cb(
j)) .and. (integral(i)%xmax >=
x_cb(
j)))
then
2039 vel(s) =
q_cons_vf(eqn_idx%cont%end + s)%sf(
j,
k,
l)/rho
2043 &
l)**2._wp)/rho)/(1._wp -
q_cons_vf(eqn_idx%alf)%sf(
j,
k,
l)) - pi_inf - qv)/gamma
2044 int_pres = int_pres + (pres - 1._wp)**2._wp
2047 int_pres = sqrt(int_pres/(1._wp*npts))
2051 call s_mpi_allreduce_sum(tmp, int_pres)
2055 if (bubbles_euler .and. (num_fluids <= 2))
then
2056 write (i + 70,
'(6x,f12.6,f24.8)') nondim_time, int_pres
2060 else if (p == 0)
then
2061 if (num_integrals /= 3)
then
2062 call s_mpi_abort(
'Incorrect number of integrals')
2065 rad = integral(1)%xmax
2066 thickness = integral(1)%xmin
2068 do i = 1, num_integrals
2078 if (sqrt(
x_cb(
j)**2._wp +
y_cb(
k)**2._wp) < (rad - 0.5_wp*thickness)) trigger = .true.
2079 else if (i == 2)
then
2081 if (sqrt(
x_cb(
j)**2._wp +
y_cb(
k)**2._wp) > (rad - 0.5_wp*thickness) .and. sqrt(
x_cb(
j)**2._wp &
2082 & +
y_cb(
k)**2._wp) < (rad + 0.5_wp*thickness)) trigger = .true.
2083 else if (i == 3)
then
2085 if (sqrt(
x_cb(
j)**2._wp +
y_cb(
k)**2._wp) > (rad + 0.5_wp*thickness)) trigger = .true.
2102 vel(s) =
q_cons_vf(eqn_idx%cont%end + s)%sf(
j,
k,
l)/rho
2106 &
l)**2._wp)/rho)/(1._wp -
q_cons_vf(eqn_idx%alf)%sf(
j,
k,
l)) - pi_inf - qv)/gamma
2107 int_pres = int_pres + abs(pres - 1._wp)
2108 max_pres = max(max_pres, abs(pres - 1._wp))
2114 int_pres = int_pres/(1._wp*npts)
2121 call s_mpi_allreduce_sum(tmp, int_pres)
2124 call s_mpi_allreduce_max(tmp, max_pres)
2128 if (bubbles_euler .and. (num_fluids <= 2))
then
2129 write (i + 70,
'(6x,f12.6,f24.8,f24.8)') nondim_time, int_pres, max_pres