MFC:Pre_process  v1.0
m_variables_conversion.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 
39 
40 
41  ! Dependencies =============================================================
42  USE m_derived_types !< definitions of the derived types
43 
44  USE m_global_parameters !< global parameters for the code
45  ! ==========================================================================
46 
47 
48  IMPLICIT NONE
49 
50 
53  abstract INTERFACE
54 
64  SUBROUTINE s_convert_xxxxx_to_mixture_variables( q_vf, i,j,k, &
65  rho,gamma,pi_inf )
66 
67  ! Importing the derived type scalar_field from m_derived_types.f90
68  ! and global variable sys_size, from m_global_variables.f90, as
69  ! the abstract interface does not inherently have access to them
70  IMPORT :: scalar_field, sys_size
71 
72  TYPE(scalar_field), DIMENSION(sys_size), INTENT(IN) :: q_vf
73 
74  INTEGER, INTENT(IN) :: i,j,k
75 
76  REAL(KIND(0d0)), INTENT(OUT) :: rho
77  REAL(KIND(0d0)), INTENT(OUT) :: gamma
78  REAL(KIND(0d0)), INTENT(OUT) :: pi_inf
79 
81 
82  END INTERFACE
83 
84  ! NOTE: These abstract interfaces allow for the declaration of a pointer to
85  ! a procedure such that the choice of the model equations does not have to
86  ! be queried every time the mixture quantities are needed.
87 
88 
89 
91  POINTER :: s_convert_to_mixture_variables => null()
94 
95  CONTAINS
96 
109  SUBROUTINE s_convert_mixture_to_mixture_variables( q_vf, i,j,k, &
110  rho,gamma,pi_inf )
112  TYPE(scalar_field), DIMENSION(sys_size), INTENT(IN) :: q_vf
113  INTEGER, INTENT(IN) :: i,j,k
114 
115  REAL(KIND(0d0)), INTENT(OUT) :: rho
116  REAL(KIND(0d0)), INTENT(OUT) :: gamma
117  REAL(KIND(0d0)), INTENT(OUT) :: pi_inf
118 
119 
120  ! Transfering the density, the specific heat ratio function and the
121  ! liquid stiffness function, respectively
122  rho = q_vf(1)%sf(i,j,k)
123  gamma = q_vf(gamma_idx)%sf(i,j,k)
124  pi_inf = q_vf(pi_inf_idx)%sf(i,j,k)
125 
126 
127  END SUBROUTINE s_convert_mixture_to_mixture_variables ! ----------------
128 
143  j,k,l,&
144  rho_K,gamma_K, pi_inf_K &
145  )
147  TYPE(scalar_field), DIMENSION(sys_size), INTENT(IN) :: qK_vf
148 
149  REAL(KIND(0d0)), INTENT(OUT) :: rho_K, gamma_K, pi_inf_K
150 
151  REAL(KIND(0d0)), DIMENSION(num_fluids) :: alpha_rho_K, alpha_K
153 
154  INTEGER, INTENT(IN) :: j,k,l
155  INTEGER :: i
156 
157  ! Constraining the partial densities and the volume fractions within
158  ! their physical bounds to make sure that any mixture variables that
159  ! are derived from them result within the limits that are set by the
160  ! fluids physical parameters that make up the mixture
161  ! alpha_rho_K(1) = qK_vf(i)%sf(i,j,k)
162  ! alpha_K(1) = qK_vf(E_idx+i)%sf(i,j,k)
163 
164  ! Performing the transfer of the density, the specific heat ratio
165  ! function as well as the liquid stiffness function, respectively
166 
167  IF (model_eqns == 4) THEN
168  rho_k = qk_vf(1)%sf(j,k,l)
169  gamma_k = fluid_pp(1)%gamma !qK_vf(gamma_idx)%sf(i,j,k)
170  pi_inf_k = fluid_pp(1)%pi_inf !qK_vf(pi_inf_idx)%sf(i,j,k)
171  ELSE IF ((model_eqns == 2) .AND. bubbles .AND. adv_alphan) THEN
172  rho_k = 0d0; gamma_k = 0d0; pi_inf_k = 0d0
173 
174  IF (mpp_lim .AND. (num_fluids > 2)) THEN
175  DO i = 1,num_fluids
176  rho_k = rho_k + qk_vf(i)%sf(j,k,l)
177  gamma_k = gamma_k + qk_vf(i+e_idx)%sf(j,k,l)*fluid_pp(i)%gamma
178  pi_inf_k = pi_inf_k + qk_vf(i+e_idx)%sf(j,k,l)*fluid_pp(i)%pi_inf
179  END DO
180  ELSE IF (num_fluids == 2) THEN
181  rho_k = qk_vf(1)%sf(j,k,l)
182  gamma_k = fluid_pp(1)%gamma
183  pi_inf_k = fluid_pp(1)%pi_inf
184  ELSE IF (num_fluids > 2) THEN
185  DO i = 1, num_fluids-1 !leave out bubble part of mixture
186  rho_k = rho_k + qk_vf(i)%sf(j,k,l)
187  gamma_k = gamma_k + qk_vf(i+e_idx)%sf(j,k,l)*fluid_pp(i)%gamma
188  pi_inf_k = pi_inf_k + qk_vf(i+e_idx)%sf(j,k,l)*fluid_pp(i)%pi_inf
189  END DO
190  !rho_K = qK_vf(1)%sf(j,k,l)
191  !gamma_K = fluid_pp(1)%gamma
192  !pi_inf_K = fluid_pp(1)%pi_inf
193  ELSE
194  rho_k = qk_vf(1)%sf(j,k,l)
195  gamma_k = fluid_pp(1)%gamma
196  pi_inf_k = fluid_pp(1)%pi_inf
197  END IF
198  END IF
199 
200  END SUBROUTINE s_convert_species_to_mixture_variables_bubbles ! ----------------
201 
202 
215  SUBROUTINE s_convert_species_to_mixture_variables( q_vf, j,k,l, &
216  rho,gamma,pi_inf )
217 
218 
219  TYPE(scalar_field), DIMENSION(sys_size), INTENT(IN) :: q_vf
220 
221 
222  INTEGER, INTENT(IN) :: j,k,l
223  ! Indices of the cell for which to compute the mixture variables
224 
225  REAL(KIND(0d0)), INTENT(OUT) :: rho
226  REAL(KIND(0d0)), INTENT(OUT) :: gamma
227  REAL(KIND(0d0)), INTENT(OUT) :: pi_inf
228 
229  INTEGER :: i
230 
231 
232  ! Computing the density, the specific heat ratio function and the
233  ! liquid stiffness function, respectively
234  IF(adv_alphan) THEN
235 
236  rho = 0d0; gamma = 0d0; pi_inf = 0d0
237 
238  DO i = 1, num_fluids
239  rho = rho + q_vf(i)%sf(j,k,l)
240  gamma = gamma + q_vf(i+e_idx)%sf(j,k,l)*fluid_pp(i)%gamma
241  pi_inf = pi_inf + q_vf(i+e_idx)%sf(j,k,l)*fluid_pp(i)%pi_inf
242  END DO
243 
244  ELSE
245 
246  rho = q_vf(num_fluids)%sf(j,k,l)
247  gamma = fluid_pp(num_fluids)%gamma
248  pi_inf = fluid_pp(num_fluids)%pi_inf
249 
250  DO i = 1, num_fluids-1
251  rho = rho + q_vf(i)%sf(j,k,l)
252  gamma = gamma + q_vf(i+e_idx)%sf(j,k,l) &
253  * ( fluid_pp( i )%gamma &
254  - fluid_pp(num_fluids)%gamma )
255  pi_inf = pi_inf + q_vf(i+e_idx)%sf(j,k,l) &
256  * ( fluid_pp( i )%pi_inf &
257  - fluid_pp(num_fluids)%pi_inf )
258  END DO
259 
260  END IF
261 
262 
263  END SUBROUTINE s_convert_species_to_mixture_variables ! ----------------
264 
265 
268  SUBROUTINE s_initialize_variables_conversion_module() ! -------------------
270  ! Depending on the model selection for the equations of motion, the
271  ! appropriate procedure for the conversion to the mixture variables
272  ! is targeted by the procedure pointer
273 
274  IF (model_eqns == 1) THEN ! Gamma/pi_inf model
277 
278  ELSE IF (bubbles) THEN
281  ELSE
282  ! Volume fraction model
285  END IF
286 
287 
288  END SUBROUTINE s_initialize_variables_conversion_module ! -----------------
289 
290 
291 
292 
296  SUBROUTINE s_convert_conservative_to_primitive_variables( q_cons_vf, &
297  q_prim_vf )
299  TYPE(scalar_field), &
300  DIMENSION(sys_size), &
301  INTENT(IN) :: q_cons_vf
302 
303  TYPE(scalar_field), &
304  DIMENSION(sys_size), &
305  INTENT(INOUT) :: q_prim_vf
306 
307  ! Density, specific heat ratio function, liquid stiffness function
308  ! and dynamic pressure, as defined in the incompressible flow sense,
309  ! respectively
310  REAL(KIND(0d0)) :: rho
311  REAL(KIND(0d0)) :: gamma
312  REAL(KIND(0d0)) :: pi_inf
313  REAL(KIND(0d0)) :: dyn_pres
314  REAL(KIND(0d0)) :: nbub
315  REAL(KIND(0d0)), DIMENSION(nb) :: nRtmp
316 
317  ! Generic loop iterators
318  INTEGER :: i,j,k,l
319 
320 
321  ! Converting the conservative variables to the primitive variables
322  DO l = 0, p
323  DO k = 0, n
324  DO j = 0, m
325 
326  ! Obtaining the density, specific heat ratio function
327  ! and the liquid stiffness function, respectively
328  IF (model_eqns .ne. 4 ) THEN
329  CALL s_convert_to_mixture_variables( q_cons_vf,j,k,l, &
330  rho,gamma,pi_inf )
331  END IF
332 
333  ! Transferring the continuity equation(s) variable(s)
334  DO i = 1, cont_idx%end
335  q_prim_vf(i)%sf(j,k,l) = q_cons_vf(i)%sf(j,k,l)
336  END DO
337 
338  ! Zeroing out the dynamic pressure since it is computed
339  ! iteratively by cycling through the momentum equations
340  dyn_pres = 0d0
341 
342  ! Computing velocity and dynamic pressure from momenta
343  DO i = mom_idx%beg, mom_idx%end
344  IF (model_eqns .ne. 4 ) THEN
345  q_prim_vf(i)%sf(j,k,l) = q_cons_vf(i)%sf(j,k,l)/rho
346  dyn_pres = dyn_pres + q_cons_vf(i)%sf(j,k,l) * &
347  q_prim_vf(i)%sf(j,k,l) / 2d0
348  ELSE
349  q_prim_vf(i)%sf(j,k,l) = q_cons_vf(i)%sf(j,k,l) &
350  / q_cons_vf(1)%sf(j,k,l)
351  END IF
352  END DO
353 
354  IF ( (model_eqns .ne. 4) .AND. (bubbles .NEQV. .true.) ) THEN
355  ! Computing the pressure from the energy
356  q_prim_vf(e_idx)%sf(j,k,l) = &
357  (q_cons_vf(e_idx)%sf(j,k,l)-dyn_pres-pi_inf)/gamma
358  ELSE IF ( (model_eqns .ne. 4) .AND. bubbles) THEN
359  print*, 'getting model_eqns 2 with bubbles. cons to prim'
360  ! p = ( E/(1-alf) - 0.5 rho u u/(1-alf) - pi_inf_k )/gamma_k
361  q_prim_vf(e_idx)%sf(j,k,l) = &
362  ((q_cons_vf(e_idx)%sf(j,k,l)-dyn_pres)/(1.d0 - q_cons_vf(alf_idx)%sf(j,k,l)) &
363  - pi_inf )/gamma
364  ELSE
365  ! Tait EOS
366  ! p = (pl0 + pi_infty)(rho/(rho_l0(1-alf)))^gamma - pi_infty
367  q_prim_vf(e_idx)%sf(j,k,l) = &
368  (pref+fluid_pp(1)%pi_inf) * &
369  ( &
370  q_prim_vf(1)%sf(j,k,l)/ &
371  (rhoref*(1-q_prim_vf(alf_idx)%sf(j,k,l))) &
372  ) ** (1/fluid_pp(1)%gamma + 1) - fluid_pp(1)%pi_inf
373  END IF
374 
375  ! Set partial pressures to mixture pressure
376  IF(model_eqns == 3) THEN
378  q_prim_vf(i)%sf(j,k,l) = q_prim_vf(e_idx)%sf(j,k,l)
379  END DO
380  END IF
381 
382  ! Transfer the advection equation(s) variable(s)
383  DO i = adv_idx%beg, adv_idx%end
384  q_prim_vf(i)%sf(j,k,l) = q_cons_vf(i)%sf(j,k,l)
385  END DO
386 
387  IF (bubbles) THEN
388  DO i = 1,nb
389  nrtmp(i) = q_cons_vf(bub_idx%rs(i))%sf(j,k,l)
390  END DO
391  CALL s_comp_n_from_cons( q_cons_vf(alf_idx)%sf(j,k,l), nrtmp, nbub)
392  DO i = bub_idx%beg, sys_size
393  q_prim_vf(i)%sf(j,k,l) = q_cons_vf(i)%sf(j,k,l)/nbub
394  END DO
395  END IF
396  END DO
397  END DO
398  END DO
399 
400 
401  END SUBROUTINE s_convert_conservative_to_primitive_variables ! ---------
402 
407  SUBROUTINE s_convert_primitive_to_conservative_variables( q_prim_vf, &
408  q_cons_vf )
410  TYPE(scalar_field), &
411  DIMENSION(sys_size), &
412  INTENT(IN) :: q_prim_vf
413 
414  TYPE(scalar_field), &
415  DIMENSION(sys_size), &
416  INTENT(INOUT) :: q_cons_vf
417 
418  ! Density, specific heat ratio function, liquid stiffness function
419  ! and dynamic pressure, as defined in the incompressible flow sense,
420  ! respectively
421  REAL(KIND(0d0)) :: rho
422  REAL(KIND(0d0)) :: gamma
423  REAL(KIND(0d0)) :: pi_inf
424  REAL(KIND(0d0)) :: dyn_pres
425  REAL(KIND(0d0)) :: nbub
426  REAL(KIND(0d0)), DIMENSION(nb) :: Rtmp
427 
428  INTEGER :: i,j,k,l
429 
430  ! Converting the primitive variables to the conservative variables
431  DO l = 0, p
432  DO k = 0, n
433  DO j = 0, m
434 
435  ! Obtaining the density, specific heat ratio function
436  ! and the liquid stiffness function, respectively
437  CALL s_convert_to_mixture_variables( q_prim_vf,j,k,l, &
438  rho,gamma,pi_inf )
439 
440  ! Transferring the continuity equation(s) variable(s)
441  DO i = 1, cont_idx%end
442  q_cons_vf(i)%sf(j,k,l) = q_prim_vf(i)%sf(j,k,l)
443  END DO
444 
445  ! Zeroing out the dynamic pressure since it is computed
446  ! iteratively by cycling through the velocity equations
447  dyn_pres = 0d0
448 
449  ! Computing momenta and dynamic pressure from velocity
450  DO i = mom_idx%beg, mom_idx%end
451  q_cons_vf(i)%sf(j,k,l) = rho*q_prim_vf(i)%sf(j,k,l)
452  dyn_pres = dyn_pres + q_cons_vf(i)%sf(j,k,l) * &
453  q_prim_vf(i)%sf(j,k,l) / 2d0
454  END DO
455 
456  ! Computing the energy from the pressure
457  IF ( (model_eqns .ne. 4) .AND. (bubbles .NEQV. .true.) ) THEN
458  ! E = Gamma*P + \rho u u /2 + \pi_inf
459  q_cons_vf(e_idx)%sf(j,k,l) = &
460  gamma*q_prim_vf(e_idx)%sf(j,k,l)+dyn_pres+pi_inf
461  ELSE IF ( (model_eqns .ne. 4) .AND. (bubbles) ) THEN
462  ! \tilde{E} = dyn_pres + (1-\alf)(\Gamma p_l + \Pi_inf)
463  q_cons_vf(e_idx)%sf(j,k,l) = dyn_pres + &
464  (1.d0 - q_prim_vf(alf_idx)%sf(j,k,l)) * &
465  (gamma*q_prim_vf(e_idx)%sf(j,k,l) + pi_inf)
466  ELSE
467  !Tait EOS, no conserved energy variable
468  q_cons_vf(e_idx)%sf(j,k,l) = 0.
469  END IF
470 
471  ! Computing the internal energies from the pressure and continuities
472  IF(model_eqns == 3) THEN
474  q_cons_vf(i)%sf(j,k,l) = q_cons_vf(i-adv_idx%end)%sf(j,k,l) * &
475  fluid_pp(i-adv_idx%end)%gamma*q_prim_vf(e_idx)%sf(j,k,l)+fluid_pp(i-adv_idx%end)%pi_inf
476  END DO
477  END IF
478 
479  ! Transferring the advection equation(s) variable(s)
480  DO i = adv_idx%beg, adv_idx%end
481  q_cons_vf(i)%sf(j,k,l) = q_prim_vf(i)%sf(j,k,l)
482  END DO
483 
484  IF (bubbles) THEN
485  DO i = 1,nb
486  rtmp(i) = q_prim_vf(bub_idx%rs(i))%sf(j,k,l)
487  END DO
488  CALL s_comp_n_from_prim( q_prim_vf(alf_idx)%sf(j,k,l), rtmp, nbub)
489 
490  DO i = bub_idx%beg, sys_size
491  q_cons_vf(i)%sf(j,k,l) = q_prim_vf(i)%sf(j,k,l)*nbub
492  END DO
493  END IF
494 
495  END DO
496  END DO
497  END DO
498 
499  END SUBROUTINE s_convert_primitive_to_conservative_variables ! ---------
500 
501 
503  SUBROUTINE s_finalize_variables_conversion_module() ! ----------------
505  ! Nullifying the procedure pointer to the subroutine transfering/
506  ! computing the mixture/species variables to the mixture variables
508 
509 
510  END SUBROUTINE s_finalize_variables_conversion_module ! --------------
511 
512 
513 END MODULE m_variables_conversion
subroutine s_initialize_variables_conversion_module()
Computation of parameters, allocation procedures, and/or any other tasks needed to properly setup the...
Abstract interface to two subroutines designed for the transfer/conversion of the mixture/species var...
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.
This module contains all of the parameters characterizing the computational domain, simulation algorithm, initial condition and the stiffened equation of state.
Derived type annexing a scalar field (SF)
subroutine s_convert_species_to_mixture_variables(q_vf, j, k, l, rho, gamma, pi_inf)
This subroutine is designed for the volume fraction model and provided a set of either conservative o...
procedure(s_convert_xxxxx_to_mixture_variables), pointer s_convert_to_mixture_variables
type(int_bounds_info) cont_idx
Indexes of first & last continuity eqns.
This file contains the definitions of all of the custom-defined types used in the pre-process code...
subroutine s_convert_primitive_to_conservative_variables(q_prim_vf, q_cons_vf)
Converts the primitive variables to the conservative ones. Used when initializing patches...
type(int_bounds_info) adv_idx
Indexes of first & last advection eqns.
integer p
Number of cells in the x-, y- and z-coordinate directions.
integer alf_idx
Index of void fraction.
integer num_fluids
Number of different fluids present in the flow.
subroutine s_finalize_variables_conversion_module()
Deallocation procedures for the module.
integer sys_size
Number of unknowns in the system of equations.
logical mpp_lim
Alpha limiter.
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.
subroutine s_convert_mixture_to_mixture_variables(q_vf, i, j, k, rho, gamma, pi_inf)
This subroutine is designed for the gamma/pi_inf model and provided a set of either conservative or p...
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.
subroutine s_convert_conservative_to_primitive_variables(q_cons_vf, q_prim_vf)
Converts the conservative variables to the primitive ones.
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.
subroutine s_convert_species_to_mixture_variables_bubbles(qK_vf, j, k, l, rho_K, gamma_K, pi_inf_K)
This procedure is used alongside with the gamma/pi_inf model to transfer the density, the specific heat ratio function and liquid stiffness function from the vector of conservative or primitive variables to their scalar counterparts. Specifially designed for when subgrid bubbles must be included.
This module consists of subroutines used in the conversion of the conservative variables into the pri...
logical adv_alphan
Advection of the last volume fraction.