MFC:Pre_process  v1.0
m_global_parameters.f90
Go to the documentation of this file.
1 !! __ _______________
2 !! / |/ / ____/ ____/
3 !! / /|_/ / /_ / /
4 !! / / / / __/ / /___
5 !! /_/ /_/_/ \____/
6 !!
7 !! This file is part of MFC.
8 !!
9 !! MFC is the legal property of its developers, whose names
10 !! are listed in the copyright file included with this source
11 !! distribution.
12 !!
13 !! MFC is free software: you can redistribute it and/or modify
14 !! it under the terms of the GNU General Public License as published
15 !! by the Free Software Foundation, either version 3 of the license
16 !! or any later version.
17 !!
18 !! MFC is distributed in the hope that it will be useful,
19 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
20 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 !! GNU General Public License for more details.
22 !!
23 !! You should have received a copy of the GNU General Public License
24 !! along with MFC (LICENSE).
25 !! If not, see <http://www.gnu.org/licenses/>.
26 
33 
38 
39 
40  ! Dependencies =============================================================
41  USE mpi ! Message passing interface (MPI) module
42 
43  USE m_derived_types ! Definitions of the derived types
44  ! ==========================================================================
45 
46 
47  IMPLICIT NONE
48 
49 
50  ! Logistics ================================================================
51  INTEGER :: num_procs
52  INTEGER , PARAMETER :: num_stcls_min = 5
53  INTEGER , PARAMETER :: path_len = 400
54  INTEGER , PARAMETER :: name_len = 50
55  REAL(KIND(0d0)), PARAMETER :: dflt_real = -1d6
56  INTEGER , PARAMETER :: dflt_int = -100
57  CHARACTER(LEN = path_len) :: case_dir
58  LOGICAL :: old_grid
59  LOGICAL :: old_ic
60  INTEGER :: t_step_old
61  ! ==========================================================================
62 
63 
64  ! Computational Domain Parameters ==========================================
65 
66  INTEGER :: proc_rank
67 
68  INTEGER :: m
69  INTEGER :: n
70  INTEGER :: p
72 
73 
74  INTEGER :: m_glb, n_glb, p_glb
76 
77  INTEGER :: num_dims
78 
79 
80  LOGICAL :: cyl_coord
81  INTEGER :: grid_geometry
83 
84 
85  REAL(KIND(0d0)), ALLOCATABLE, DIMENSION(:) :: x_cc, y_cc, z_cc
87 
88 
89  REAL(KIND(0d0)), ALLOCATABLE, DIMENSION(:) :: x_cb, y_cb, z_cb
91 
92 
93  REAL(KIND(0d0)) :: dx, dy, dz
95 
96 
99 
100 
103 
104  ! Parameters of the grid stretching function for the x-, y- and z-coordinate
105  ! directions. The "a" parameters are a measure of the rate at which the grid
106  ! is stretched while the remaining parameters are indicative of the location
107  ! on the grid at which the stretching begins.
108  REAL(KIND(0d0)) :: a_x, a_y, a_z
109  INTEGER :: loops_x, loops_y, loops_z
110  REAL(KIND(0d0)) :: x_a, y_a, z_a
111  REAL(KIND(0d0)) :: x_b, y_b, z_b
112 
113  ! ==========================================================================
114 
115 
116  ! Simulation Algorithm Parameters ==========================================
117  INTEGER :: model_eqns
118  INTEGER :: num_fluids
119  LOGICAL :: adv_alphan
120  LOGICAL :: mpp_lim
121  INTEGER :: sys_size
122  INTEGER :: weno_order
123 
124  ! Annotations of the structure, i.e. the organization, of the state vectors
127  INTEGER :: e_idx
128  INTEGER :: alf_idx
132  INTEGER :: gamma_idx
133  INTEGER :: pi_inf_idx
134 
135 
138 
139 
140  LOGICAL :: parallel_io
141  INTEGER :: precision
142 
143  ! Perturb density of surrounding air so as to break symmetry of grid
144  LOGICAL :: perturb_flow
146  LOGICAL :: perturb_sph
147  INTEGER :: perturb_sph_fluid
148  REAL(KIND(0d0)), DIMENSION(num_fluids_max) :: fluid_rho
149 
150 
151  INTEGER, ALLOCATABLE, DIMENSION(:) :: proc_coords
153 
154  INTEGER, ALLOCATABLE, DIMENSION(:) :: start_idx
156 
157  TYPE(mpi_io_var), PUBLIC :: mpi_io_data
158 
159 
160  CHARACTER(LEN = name_len) :: mpiiofs
161  INTEGER :: mpi_info_int
163 
164  INTEGER, PRIVATE :: ierr
165  ! ==========================================================================
166 
167 
168  ! Initial Condition Parameters =============================================
169  INTEGER :: num_patches
170 
171  TYPE(ic_patch_parameters), DIMENSION(num_patches_max) :: patch_icpp
176  ! ==========================================================================
177 
178 
179  ! Fluids Physical Parameters ===============================================
180  TYPE(physical_parameters), DIMENSION(num_fluids_max) :: fluid_pp
184 
185  ! ==========================================================================
186 
187 
188  REAL(KIND(0d0)) :: rhoref, pref
189 
192  INTEGER :: nb
193  REAL(KIND(0d0)) :: r0ref
194  REAL(KIND(0d0)) :: ca, web, re_inv
195  REAL(KIND(0d0)), DIMENSION(:), ALLOCATABLE :: weight, r0, v0
196  LOGICAL :: bubbles
198 
201  LOGICAL :: polytropic
202  LOGICAL :: polydisperse
203  INTEGER :: thermal !1 = adiabatic, 2 = isotherm, 3 = transfer
204  REAL(KIND(0d0)) :: r_n, r_v, phi_vn, phi_nv, pe_c, tw
205  REAL(KIND(0d0)), DIMENSION(:), ALLOCATABLE :: k_n, k_v, pb0, mass_n0, mass_v0, pe_t
206  REAL(KIND(0d0)), DIMENSION(:), ALLOCATABLE :: re_trans_t, re_trans_c, im_trans_t, im_trans_c, omegan
207  REAL(KIND(0d0)) :: poly_sigma
209 
210  INTEGER, ALLOCATABLE, DIMENSION(:,:,:) :: logic_grid
211 
212  ! Mathematical and Physical Constants ======================================
213  REAL(KIND(0d0)), PARAMETER :: pi = 3.141592653589793d0
214  ! ==========================================================================
215 
216  CONTAINS
217 
218 
222  SUBROUTINE s_assign_default_values_to_user_inputs() ! ------------------
223 
224  INTEGER :: i
225 
226  ! Logistics
227  case_dir = ' '
228  old_grid = .false.
229  old_ic = .false.
231 
232 
233  ! Computational domain parameters
234  m = dflt_int
235  n = dflt_int
236  p = dflt_int
237 
238  cyl_coord = .false.
239 
240  x_domain%beg = dflt_real
241  x_domain%end = dflt_real
242  y_domain%beg = dflt_real
243  y_domain%end = dflt_real
244  z_domain%beg = dflt_real
245  z_domain%end = dflt_real
246 
247  stretch_x = .false.
248  stretch_y = .false.
249  stretch_z = .false.
250 
251  a_x = dflt_real
252  a_y = dflt_real
253  a_z = dflt_real
254  loops_x = 1
255  loops_y = 1
256  loops_z = 1
257  x_a = dflt_real
258  x_b = dflt_real
259  y_a = dflt_real
260  y_b = dflt_real
261  z_a = dflt_real
262  z_b = dflt_real
263 
264  ! Simulation algorithm parameters
267  adv_alphan = .false.
269 
270  bc_x%beg = dflt_real
271  bc_x%end = dflt_real
272  bc_y%beg = dflt_real
273  bc_y%end = dflt_real
274  bc_z%beg = dflt_real
275  bc_z%end = dflt_real
276 
277  parallel_io = .false.
278  precision = 2
279  perturb_flow = .false.
281  perturb_sph = .false.
284 
285 
286  ! Initial condition parameters
288 
289  DO i = 1, num_patches_max
290  patch_icpp(i)%geometry = dflt_int
291  patch_icpp(i)%x_centroid = dflt_real
292  patch_icpp(i)%y_centroid = dflt_real
293  patch_icpp(i)%z_centroid = dflt_real
294  patch_icpp(i)%length_x = dflt_real
295  patch_icpp(i)%length_y = dflt_real
296  patch_icpp(i)%length_z = dflt_real
297  patch_icpp(i)%radius = dflt_real
298  patch_icpp(i)%epsilon = dflt_real
299  patch_icpp(i)%beta = dflt_real
300  patch_icpp(i)%normal = dflt_real
301  patch_icpp(i)%radii = dflt_real
302  patch_icpp(i)%alter_patch = .false.
303  patch_icpp(i)%alter_patch(0) = .true.
304  patch_icpp(i)%smoothen = .false.
305  patch_icpp(i)%smooth_patch_id = i
306  patch_icpp(i)%smooth_coeff = dflt_real
307  patch_icpp(i)%alpha_rho = dflt_real
308  patch_icpp(i)%rho = dflt_real
309  patch_icpp(i)%vel = dflt_real
310  patch_icpp(i)%pres = dflt_real
311  patch_icpp(i)%alpha = dflt_real
312  patch_icpp(i)%gamma = dflt_real
313  patch_icpp(i)%pi_inf = dflt_real
314  !should get all of r0's and v0's
315  patch_icpp(i)%r0 = dflt_real
316  patch_icpp(i)%v0 = dflt_real
317  END DO
318 
319  ! Tait EOS
320  rhoref = dflt_real
321  pref = dflt_real
322 
323  ! Bubble modeling
324  bubbles = .false.
325  polytropic = .true.
326  polydisperse= .false.
327 
328  thermal = dflt_int
329  r0ref = dflt_real
330  nb = dflt_int
331 
332  ca = dflt_real
333  re_inv = dflt_real
334  web = dflt_real
336 
337  r_n = dflt_real
338  r_v = dflt_real
339  phi_vn = dflt_real
340  phi_nv = dflt_real
341  pe_c = dflt_real
342  tw = dflt_real
343 
344  ! Fluids physical parameters
345  DO i = 1, num_fluids_max
346  fluid_pp(i)%gamma = dflt_real
347  fluid_pp(i)%pi_inf = dflt_real
348  fluid_pp(i)%mul0 = dflt_real
349  fluid_pp(i)%ss = dflt_real
350  fluid_pp(i)%pv = dflt_real
351  fluid_pp(i)%gamma_v = dflt_real
352  fluid_pp(i)%M_v = dflt_real
353  fluid_pp(i)%mu_v = dflt_real
354  fluid_pp(i)%k_v = dflt_real
355  END DO
356 
357 
358  END SUBROUTINE s_assign_default_values_to_user_inputs ! ----------------
359 
362  SUBROUTINE s_initialize_global_parameters_module() ! ----------------------
363 
364  INTEGER :: i, fac
365 
366  ! Determining the layout of the state vectors and overall size of
367  ! the system of equations, given the dimensionality and choice of
368  ! the equations of motion
369 
370  ! Gamma/Pi_inf Model ===============================================
371  IF(model_eqns == 1) THEN
372 
373  ! Setting number of fluids
374  num_fluids = 1
375 
376  ! Annotating structure of the state and flux vectors belonging
377  ! to the system of equations defined by the selected number of
378  ! spatial dimensions and the gamma/pi_inf model
379  cont_idx%beg = 1
380  cont_idx%end = cont_idx%beg
381  mom_idx%beg = cont_idx%end + 1
382  mom_idx%end = cont_idx%end + num_dims
383  e_idx = mom_idx%end + 1
384  adv_idx%beg = e_idx + 1
385  adv_idx%end = adv_idx%beg + 1
386  gamma_idx = adv_idx%beg
387  pi_inf_idx = adv_idx%end
388  sys_size = adv_idx%end
389 
390  ! ==================================================================
391 
392 
393  ! Volume Fraction Model (5-equation model) =========================
394  ELSE IF(model_eqns == 2) THEN
395 
396  ! Annotating structure of the state and flux vectors belonging
397  ! to the system of equations defined by the selected number of
398  ! spatial dimensions and the volume fraction model
399  cont_idx%beg = 1
400  cont_idx%end = num_fluids
401  mom_idx%beg = cont_idx%end + 1
402  mom_idx%end = cont_idx%end + num_dims
403  e_idx = mom_idx%end + 1
404  adv_idx%beg = e_idx + 1
405  adv_idx%end = e_idx + num_fluids
406 
407  IF( (adv_alphan .NEQV. .true.) .and. &
408  (num_fluids > 1)) adv_idx%end = adv_idx%end - 1
409 
410  sys_size = adv_idx%end
411 
412  IF (bubbles) THEN
413  alf_idx = adv_idx%end
414  ELSE
415  alf_idx = 0
416  END IF
417 
418  IF (bubbles) THEN
419  bub_idx%beg = sys_size+1
420  bub_idx%end = sys_size+2*nb
421  IF (polytropic .NEQV. .true.) THEN
422  bub_idx%end = sys_size+4*nb
423  END IF
424  sys_size = bub_idx%end
425 
426  ALLOCATE( bub_idx%rs(nb), bub_idx%vs(nb) )
427  ALLOCATE( bub_idx%ps(nb), bub_idx%ms(nb) )
428  ALLOCATE( weight(nb),r0(nb),v0(nb) )
429 
430  DO i = 1, nb
431  IF (polytropic .NEQV. .true.) THEN
432  fac = 4
433  ELSE
434  fac = 2
435  END IF
436 
437  bub_idx%rs(i) = bub_idx%beg+(i-1)*fac
438  bub_idx%vs(i) = bub_idx%rs(i)+1
439 
440  IF (polytropic .NEQV. .true.) THEN
441  bub_idx%ps(i) = bub_idx%vs(i)+1
442  bub_idx%ms(i) = bub_idx%ps(i)+1
443  END IF
444  END DO
445 
446  IF (nb == 1) THEN
447  weight(:) = 1d0
448  r0(:) = 1d0
449  v0(:) = 0d0
450  ELSE IF (nb > 1) THEN
451  CALL s_simpson(nb)
452  v0(:) = 0d0
453  ELSE
454  stop 'Invalid value of nb'
455  END IF
456 
457  IF (polytropic .NEQV. .true.) THEN
459  ELSE
460  rhoref = 1.d0
461  pref = 1.d0
462  END IF
463  END IF
464  ! ==================================================================
465 
466 
467  ! Volume Fraction Model (6-equation model) =========================
468  ELSE IF(model_eqns == 3) THEN
469 
470  ! Annotating structure of the state and flux vectors belonging
471  ! to the system of equations defined by the selected number of
472  ! spatial dimensions and the volume fraction model
473  cont_idx%beg = 1
474  cont_idx%end = num_fluids
475  mom_idx%beg = cont_idx%end + 1
476  mom_idx%end = cont_idx%end + num_dims
477  e_idx = mom_idx%end + 1
478  adv_idx%beg = e_idx + 1
479  adv_idx%end = e_idx + num_fluids
480  IF(adv_alphan .NEQV. .true.) adv_idx%end = adv_idx%end - 1
481  internalenergies_idx%beg = adv_idx%end + 1
484  !========================
485  ELSE IF(model_eqns == 4) THEN
486  ! 4 equation model with subgrid bubbles
487  cont_idx%beg = 1 ! one continuity equation
488  cont_idx%end = 1 ! num_fluids
489  mom_idx%beg = cont_idx%end + 1 ! one momentum equation in each direction
490  mom_idx%end = cont_idx%end + num_dims
491  e_idx = mom_idx%end + 1 ! one energy equation
492  adv_idx%beg = e_idx + 1
493  adv_idx%end = adv_idx%beg !one volume advection equation
494  alf_idx = adv_idx%end
495  sys_size = alf_idx !adv_idx%end
496 
497  IF (bubbles) THEN
498  bub_idx%beg = sys_size+1
499  bub_idx%end = sys_size+2*nb
500  IF (polytropic .NEQV. .true.) THEN
501  bub_idx%end = sys_size+4*nb
502  END IF
503  sys_size = bub_idx%end
504 
505  ALLOCATE( bub_idx%rs(nb), bub_idx%vs(nb) )
506  ALLOCATE( bub_idx%ps(nb), bub_idx%ms(nb) )
507  ALLOCATE( weight(nb),r0(nb),v0(nb) )
508 
509  DO i = 1, nb
510  IF (polytropic .NEQV. .true.) THEN
511  fac = 4
512  ELSE
513  fac = 2
514  END IF
515 
516  bub_idx%rs(i) = bub_idx%beg+(i-1)*fac
517  bub_idx%vs(i) = bub_idx%rs(i)+1
518 
519  IF (polytropic .NEQV. .true.) THEN
520  bub_idx%ps(i) = bub_idx%vs(i)+1
521  bub_idx%ms(i) = bub_idx%ps(i)+1
522  END IF
523  END DO
524 
525  IF (nb == 1) THEN
526  weight(:) = 1d0
527  r0(:) = 1d0
528  v0(:) = 0d0
529  ELSE IF (nb > 1) THEN
530  CALL s_simpson(nb)
531  v0(:) = 0d0
532  ELSE
533  stop 'Invalid value of nb'
534  END IF
535 
536  IF (polytropic .NEQV. .true.) THEN
538  ELSE
539  rhoref = 1.d0
540  pref = 1.d0
541  END IF
542  END IF
543  END IF
544  ! ==================================================================
545 
546  ALLOCATE(mpi_io_data%view(1:sys_size))
547  ALLOCATE(mpi_io_data%var(1:sys_size))
548 
549  DO i = 1, sys_size
550  ALLOCATE(mpi_io_data%var(i)%sf(0:m,0:n,0:p))
551  mpi_io_data%var(i)%sf => null()
552  END DO
553 
554  ! Allocating grid variables for the x-direction
555  ALLOCATE(x_cc(0:m), x_cb(-1:m))
556  ! Allocating grid variables for the y- and z-directions
557  IF(n > 0) THEN
558  ALLOCATE(y_cc(0:n), y_cb(-1:n))
559  IF(p > 0) THEN
560  ALLOCATE(z_cc(0:p), z_cb(-1:p))
561  END IF
562  END IF
563 
564  IF (cyl_coord .NEQV. .true.) THEN ! Cartesian grid
565  grid_geometry = 1
566  ELSEIF (cyl_coord .AND. p == 0) THEN ! Axisymmetric cylindrical grid
567  grid_geometry = 2
568  ELSE ! Fully 3D cylindrical grid
569  grid_geometry = 3
570  END IF
571 
572 
573  ALLOCATE( logic_grid(0:m,0:n,0:p) )
574 
575  END SUBROUTINE s_initialize_global_parameters_module ! --------------------
576 
577 
580  SUBROUTINE s_initialize_nonpoly
581  INTEGER :: ir
582  REAL(KIND(0.D0)) :: rhol0
583  REAL(KIND(0.D0)) :: pl0
584  REAL(KIND(0.D0)) :: uu
585  REAL(KIND(0.D0)) :: D_m
586  REAL(KIND(0.D0)) :: temp
587  REAL(KIND(0.D0)) :: omega_ref
588  REAL(KIND(0.D0)), DIMENSION(Nb) :: chi_vw0
589  REAL(KIND(0.D0)), DIMENSION(Nb) :: cp_m0
590  REAL(KIND(0.D0)), DIMENSION(Nb) :: k_m0
591  REAL(KIND(0.D0)), DIMENSION(Nb) :: rho_m0
592  REAL(KIND(0.D0)), DIMENSION(Nb) :: x_vw
593  ! polytropic index used to compute isothermal natural frequency
594  REAL(KIND(0.D0)), PARAMETER :: k_poly = 1.d0
595  ! universal gas constant
596  REAL(KIND(0.D0)), PARAMETER :: Ru = 8314.d0
597 
598 
599  ! liquid physical properties
600  REAL(KIND(0.D0)) :: mul0, ss, pv, gamma_v, M_v, mu_v
601 
602  ! gas physical properties
603  REAL(KIND(0.D0)) :: gamma_m, gamma_n, M_n, mu_n
604 
605  rhol0 = rhoref
606  pl0 = pref
607 
608  ALLOCATE( pb0(nb), mass_n0(nb), mass_v0(nb), pe_t(nb) )
609  ALLOCATE( k_n(nb), k_v(nb), omegan(nb) )
610  ALLOCATE( re_trans_t(nb), re_trans_c(nb), im_trans_t(nb), im_trans_c(nb) )
611 
612  pb0(:) = dflt_real
613  mass_n0(:) = dflt_real
614  mass_v0(:) = dflt_real
615  pe_t(:) = dflt_real
616  omegan(:) = dflt_real
617 
618  mul0 = fluid_pp(1)%mul0
619  ss = fluid_pp(1)%ss
620  pv = fluid_pp(1)%pv
621  gamma_v = fluid_pp(1)%gamma_v
622  m_v = fluid_pp(1)%M_v
623  mu_v = fluid_pp(1)%mu_v
624  k_v(:) = fluid_pp(1)%k_v
625 
626  gamma_n = fluid_pp(2)%gamma_v
627  m_n = fluid_pp(2)%M_v
628  mu_n = fluid_pp(2)%mu_v
629  k_n(:) = fluid_pp(2)%k_v
630 
631  gamma_m = gamma_n
632  IF (thermal==2 ) gamma_m = 1.d0 !isothermal
633 
634  temp = 293.15d0
635  d_m = 0.242d-4
636  uu = dsqrt( pl0/rhol0 )
637 
638  omega_ref = 3.d0*k_poly*ca + 2.d0*( 3.d0*k_poly-1.d0 )/web
639 
640  !!! thermal properties !!!
641  ! gas constants
642  r_n = ru/m_n
643  r_v = ru/m_v
644  ! phi_vn & phi_nv (phi_nn = phi_vv = 1)
645  phi_vn = ( 1.d0+dsqrt(mu_v/mu_n)*(m_n/m_v)**(0.25d0) )**2 &
646  / ( dsqrt(8.d0)*dsqrt(1.d0+m_v/m_n) )
647  phi_nv = ( 1.d0+dsqrt(mu_n/mu_v)*(m_v/m_n)**(0.25d0) )**2 &
648  / ( dsqrt(8.d0)*dsqrt(1.d0+m_n/m_v) )
649  ! internal bubble pressure
650  pb0 = pl0 + 2.d0*ss/( r0ref*r0 )
651 
652  ! mass fraction of vapor
653  chi_vw0 = 1.d0/( 1.d0+r_v/r_n*(pb0/pv-1.d0) )
654  ! specific heat for gas/vapor mixture
655  cp_m0 = chi_vw0*r_v*gamma_v/( gamma_v-1.d0 ) &
656  + ( 1.d0-chi_vw0 )*r_n*gamma_n/( gamma_n-1.d0 )
657  ! mole fraction of vapor
658  x_vw = m_n*chi_vw0/( m_v+(m_n-m_v)*chi_vw0 )
659  ! thermal conductivity for gas/vapor mixture
660  k_m0 = x_vw*k_v/( x_vw+(1.d0-x_vw)*phi_vn ) &
661  + ( 1.d0-x_vw )*k_n/( x_vw*phi_nv+1.d0-x_vw )
662  ! mixture density
663  rho_m0 = pv/( chi_vw0*r_v*temp )
664 
665  ! mass of gas/vapor computed using dimensional quantities
666  mass_n0 = 4.d0*( pb0-pv )*pi/( 3.d0*r_n*temp*rhol0 )*r0**3
667  mass_v0 = 4.d0*pv*pi/( 3.d0*r_v*temp*rhol0 )*r0**3
668  ! Peclet numbers
669  pe_t = rho_m0*cp_m0*uu*r0ref/k_m0
670  pe_c = uu*r0ref/d_m
671  ! nondimensional properties
672  r_n = rhol0*r_n*temp/pl0
673  r_v = rhol0*r_v*temp/pl0
674  k_n = k_n/k_m0
675  k_v = k_v/k_m0
676  pb0 = pb0/pl0
677  pv = pv/pl0
678 
679  ! bubble wall temperature, normalized by T0, in the liquid
680  ! keeps a constant (cold liquid assumption)
681  tw = 1.d0
682  ! natural frequencies
683  omegan = dsqrt( 3.d0*k_poly*ca+2.d0*(3.d0*k_poly-1.d0)/(web*r0) )/r0
684 
685  pl0 = 1.d0
686  DO ir = 1,nb
687  CALL s_transcoeff( omegan(ir)*r0(ir),pe_t(ir)*r0(ir), &
688  re_trans_t(ir),im_trans_t(ir) )
689  CALL s_transcoeff( omegan(ir)*r0(ir),pe_c*r0(ir), &
690  re_trans_c(ir),im_trans_c(ir) )
691  END DO
692  im_trans_t = 0d0
693  im_trans_c = 0d0
694 
695  rhoref = 1.d0
696  pref = 1.d0
697  END SUBROUTINE s_initialize_nonpoly
698 
704  SUBROUTINE s_transcoeff( omega,peclet,Re_trans,Im_trans )
706  REAL(KIND(0.D0)), INTENT(IN) :: omega
707  REAL(KIND(0.D0)), INTENT(IN) :: peclet
708  REAL(KIND(0.D0)), INTENT(OUT) :: Re_trans
709  REAL(KIND(0.D0)), INTENT(OUT) :: Im_trans
710  COMPLEX :: trans, c1, c2, c3
711  COMPLEX :: imag = ( 0.,1. )
712  REAL(KIND(0.D0)) :: f_transcoeff
713 
714  c1 = imag*omega*peclet
715  c2 = csqrt( c1 )
716  c3 = ( cexp(c2)-cexp(-c2) )/( cexp(c2)+cexp(-c2) ) ! TANH(c2)
717  trans = ( (c2/c3-1.d0)**(-1)-3.d0/c1 )**( -1 ) ! transfer function
718 
719  re_trans = dble( trans )
720  im_trans = aimag( trans )
721 
722  END SUBROUTINE s_transcoeff
723 
724 
725  SUBROUTINE s_initialize_parallel_io() ! --------------------------------
727  num_dims = 1 + min(1,n) + min(1,p)
728 
729  ALLOCATE(proc_coords(1:num_dims))
730 
731  IF (parallel_io .NEQV. .true.) RETURN
732 
733  ! Option for Lustre file system (Darter/Comet/Stampede)
734  WRITE(mpiiofs, '(A)') '/lustre_'
735  mpiiofs = trim(mpiiofs)
736  CALL mpi_info_create(mpi_info_int, ierr)
737  CALL mpi_info_set(mpi_info_int, 'romio_ds_write', 'disable', ierr)
738 
739  ! Option for UNIX file system (Hooke/Thomson)
740  ! WRITE(mpiiofs, '(A)') '/ufs_'
741  ! mpiiofs = TRIM(mpiiofs)
742  ! mpi_info_int = MPI_INFO_NULL
743 
744  ALLOCATE(start_idx(1:num_dims))
745 
746  END SUBROUTINE s_initialize_parallel_io ! ------------------------------
747 
748 
749  SUBROUTINE s_finalize_global_parameters_module() ! ------------------------
751  INTEGER :: i
752 
753  ! Deallocating grid variables for the x-direction
754  DEALLOCATE(x_cc, x_cb)
755  ! Deallocating grid variables for the y- and z-directions
756  IF(n > 0) THEN
757  DEALLOCATE(y_cc, y_cb)
758  IF(p > 0) THEN
759  DEALLOCATE(z_cc, z_cb)
760  END IF
761  END IF
762 
763  DEALLOCATE(proc_coords)
764  IF (parallel_io) THEN
765  DEALLOCATE(start_idx)
766  DO i = 1, sys_size
767  mpi_io_data%var(i)%sf => null()
768  END DO
769 
770  DEALLOCATE(mpi_io_data%var)
771  DEALLOCATE(mpi_io_data%view)
772  END IF
773 
774  END SUBROUTINE s_finalize_global_parameters_module ! ----------------------
775 
780  SUBROUTINE s_comp_n_from_cons( vftmp,nRtmp,ntmp )
781  REAL(KIND(0.D0)), INTENT(IN) :: vftmp
782  REAL(KIND(0.D0)), DIMENSION(nb), INTENT(IN) :: nRtmp
783  REAL(KIND(0.D0)), INTENT(OUT) :: ntmp
784  REAL(KIND(0.D0)) :: nR3
785 
786  CALL s_quad( nrtmp**3,nr3 ) !returns itself if NR0 = 1
787  ntmp = dsqrt( (4.d0*pi/3.d0)*nr3/vftmp )
788 
789  END SUBROUTINE s_comp_n_from_cons
790 
795  SUBROUTINE s_comp_n_from_prim( vftmp,Rtmp,ntmp )
796  REAL(KIND(0.D0)), INTENT(IN) :: vftmp
797  REAL(KIND(0.D0)), DIMENSION(nb), INTENT(IN) :: Rtmp
798  REAL(KIND(0.D0)), INTENT(OUT) :: ntmp
799  REAL(KIND(0.D0)) :: R3
800 
801  CALL s_quad( rtmp**3,r3 ) !returns itself if NR0 = 1
802  ntmp = (3.d0/(4.d0*pi)) * vftmp/r3
803 
804  END SUBROUTINE s_comp_n_from_prim
805 
809  SUBROUTINE s_quad( func,mom )
811  REAL(KIND(0.D0)), DIMENSION(nb), INTENT(IN) :: func
812  REAL(KIND(0.D0)), INTENT(OUT) :: mom
813 
814  mom = dot_product( weight,func )
815 
816  END SUBROUTINE s_quad
817 
820  SUBROUTINE s_simpson( Npt )
822  INTEGER, INTENT(IN) :: Npt
823  INTEGER :: ir
824  REAL(KIND(0.D0)) :: R0mn
825  REAL(KIND(0.D0)) :: R0mx
826  REAL(KIND(0.D0)) :: dphi
827  REAL(KIND(0.D0)) :: tmp
828  REAL(KIND(0.D0)) :: sd
829  REAL(KIND(0.D0)), DIMENSION(Npt) :: phi
830 
831  ! nondiml. min. & max. initial radii for numerical quadrature
832  !sd = 0.05D0
833  !R0mn = 0.75D0
834  !R0mx = 1.3D0
835 
836  !sd = 0.3D0
837  !R0mn = 0.3D0
838  !R0mx = 6.D0
839 
840  !sd = 0.7D0
841  !R0mn = 0.12D0
842  !R0mx = 150.D0
843 
844  sd = poly_sigma
845  r0mn = 0.8d0*dexp(-2.8d0 * sd)
846  r0mx = 0.2d0*dexp( 9.5d0 * sd) + 1.d0
847 
848  ! phi = ln( R0 ) & return R0
849  DO ir = 1,npt
850  phi(ir) = dlog( r0mn ) &
851  + dble( ir-1 )*dlog( r0mx/r0mn )/dble( npt-1 )
852  r0(ir) = dexp( phi(ir) )
853  END DO
854  dphi = phi(2) - phi(1)
855 
856  ! weights for quadrature using Simpson's rule
857  DO ir = 2,npt-1
858  ! Gaussian
859  tmp = dexp( -0.5d0*(phi(ir)/sd)**2 )/dsqrt( 2.d0*pi )/sd
860  IF ( mod(ir,2)==0 ) THEN
861  weight(ir) = tmp*4.d0*dphi/3.d0
862  ELSE
863  weight(ir) = tmp*2.d0*dphi/3.d0
864  END IF
865  END DO
866 
867  tmp = dexp( -0.5d0*(phi(1)/sd)**2 )/dsqrt( 2.d0*pi )/sd
868  weight(1) = tmp*dphi/3.d0
869  tmp = dexp( -0.5d0*(phi(npt)/sd)**2 )/dsqrt( 2.d0*pi )/sd
870  weight(npt) = tmp*dphi/3.d0
871 
872  END SUBROUTINE s_simpson
873 
874 
875 
876 END MODULE m_global_parameters
Derived type adding initial condition (ic) patch parameters as attributes NOTE: The requirements for ...
real(kind(0d0)), dimension(num_fluids_max) fluid_rho
integer, parameter num_patches_max
Maximum number of patches allowed.
real(kind(0d0)), dimension(:), allocatable v0
integer, dimension(:), allocatable proc_coords
Processor coordinates in MPI_CART_COMM.
integer t_step_old
Existing IC/grid folder.
subroutine s_finalize_global_parameters_module()
Integer boounds for variables.
real(kind(0d0)), dimension(:), allocatable re_trans_t
subroutine s_initialize_parallel_io()
integer perturb_sph_fluid
Fluid to be perturbed with perturb_sph flag.
subroutine s_initialize_nonpoly
Initializes and computes bubble properties for non-polytropic processes.
type(physical_parameters), dimension(num_fluids_max) fluid_pp
Database of the physical parameters of each of the fluids that is present in the flow. These include the stiffened gas equation of state parameters, the Reynolds numbers and the Weber numbers.
integer pi_inf_idx
Index of liquid stiffness func. eqn.
integer mpi_info_int
MPI info for parallel IO with Lustre file systems.
real(kind(0d0)), dimension(:), allocatable y_cc
This module contains all of the parameters characterizing the computational domain, simulation algorithm, initial condition and the stiffened equation of state.
integer, parameter num_stcls_min
Mininum # of stencils.
subroutine s_initialize_global_parameters_module()
Computation of parameters, allocation procedures, and/or any other tasks needed to properly setup the...
logical old_grid
Use existing grid data.
type(bounds_info) z_domain
Locations of the domain bounds in the x-, y- and z-coordinate directions.
type(bounds_info) x_domain
integer, parameter path_len
Maximum path length.
bounds for the bubble dynamic variables
real(kind(0d0)), dimension(:), allocatable y_cb
real(kind(0d0)), dimension(:), allocatable k_v
real(kind(0d0)), dimension(:), allocatable re_trans_c
real(kind(0d0)), dimension(:), allocatable x_cb
type(int_bounds_info) cont_idx
Indexes of first & last continuity eqns.
real(kind(0d0)), dimension(:), allocatable pe_t
integer p_glb
Global number of cells in each direction.
This file contains the definitions of all of the custom-defined types used in the pre-process code...
subroutine s_assign_default_values_to_user_inputs()
Assigns default values to user inputs prior to reading them in. This allows for an easier consistency...
Derived type annexing the physical parameters (PP) of the fluids. These include the specific heat rat...
type(int_bounds_info) adv_idx
Indexes of first & last advection eqns.
real(kind(0d0)), dimension(:), allocatable z_cb
Locations of cell-boundaries (cb) in x-, y- and z-directions, respectively.
integer p
Number of cells in the x-, y- and z-coordinate directions.
real(kind(0d0)), parameter dflt_real
Default real value.
integer alf_idx
Index of void fraction.
integer num_fluids
Number of different fluids present in the flow.
real(kind(0d0)), dimension(:), allocatable mass_v0
real(kind(0d0)), dimension(:), allocatable im_trans_c
real(kind(0d0)), dimension(:), allocatable x_cc
real(kind(0d0)), dimension(:), allocatable omegan
logical old_ic
Use existing IC data.
integer, dimension(:), allocatable start_idx
Starting cell-center index of local processor in global grid.
integer sys_size
Number of unknowns in the system of equations.
logical mpp_lim
Alpha limiter.
subroutine s_simpson(Npt)
Computes the Simpson weights for quadrature.
integer model_eqns
Multicomponent flow model.
subroutine s_comp_n_from_prim(vftmp, Rtmp, ntmp)
Computes the bubble number density n from the primitive variables.
type(bub_bounds_info) bub_idx
Indexes of first & last bubble variable eqns.
integer gamma_idx
Index of specific heat ratio func. eqn.
integer precision
Precision of output files.
type(mpi_io_var), public mpi_io_data
integer grid_geometry
Cylindrical coordinates (either axisymmetric or full 3D)
real(kind(0d0)), dimension(:), allocatable im_trans_t
character(len=path_len) case_dir
Case folder location.
integer, dimension(:,:,:), allocatable logic_grid
subroutine s_quad(func, mom)
Computes the quadrature for polydisperse bubble populations.
character(len=name_len) mpiiofs
integer e_idx
Index of total energy equation.
real(kind(0d0)) pref
Reference parameters for Tait EOS.
type(int_bounds_info) mom_idx
Indexes of first & last momentum eqns.
integer, parameter num_fluids_max
Maximum number of fluids allowed.
integer num_procs
Number of processors.
real(kind(0d0)), dimension(:), allocatable mass_n0
type(bounds_info) y_domain
real(kind(0d0)), dimension(:), allocatable pb0
integer num_patches
Number of patches composing initial condition.
type(ic_patch_parameters), dimension(num_patches_max) patch_icpp
Database of the initial condition patch parameters (icpp) for each of the patches employed in the con...
integer weno_order
Order of accuracy for the WENO reconstruction.
real(kind(0d0)), dimension(:), allocatable r0
real(kind(0d0)), dimension(:), allocatable z_cc
Locations of cell-centers (cc) in x-, y- and z-directions, respectively.
real(kind(0d0)) dz
Minimum cell-widths in the x-, y- and z-coordinate directions.
logical parallel_io
Format of the data files.
subroutine s_comp_n_from_cons(vftmp, nRtmp, ntmp)
Computes the bubble number density n from the conservative variables.
type(int_bounds_info) internalenergies_idx
Indexes of first & last internal energy eqns.
integer num_dims
Number of spatial dimensions.
integer, parameter name_len
Maximum name length.
real(kind(0d0)), dimension(:), allocatable k_n
integer, parameter dflt_int
Default integer value.
real(kind(0d0)), parameter pi
logical adv_alphan
Advection of the last volume fraction.
Derived type adding beginning (beg) and end bounds info as attributes.
type(bounds_info) bc_z
Boundary conditions in the x-, y- and z-coordinate directions.
subroutine s_transcoeff(omega, peclet, Re_trans, Im_trans)
Computes the transfer coefficient for the non-polytropic bubble compression process.
integer proc_rank
Rank of the local processor.
logical stretch_z
Grid stretching flags for the x-, y- and z-coordinate directions.
integer perturb_flow_fluid
Fluid to be perturbed with perturb_flow flag.
real(kind(0d0)), dimension(:), allocatable weight