File: /nfs/home/0/users/jenkins/mfix.git/model/qmomk/qmomk_initial_conditions.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: QMOMK_INITIAL_CONDITIONS                               C
4     !  Author: Alberto Passalacqua                        Date: 30-Jul-08  C
5     !  Reviewer:                                          Date:            C
6     !                                                                      C
7     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
8     
9           SUBROUTINE QMOMK_INITIAL_CONDITIONS
10     
11     !-----------------------------------------------
12     !     M o d u l e s
13     !-----------------------------------------------
14           USE param
15           USE param1
16           USE constant
17           USE physprop
18           USE fldvar
19           USE geometry
20           USE compar
21           USE indices
22           USE qmom_kinetic_equation
23           USE qmomk_quadrature
24           USE qmomk_parameters
25           USE ic
26           USE functions
27     
28           IMPLICIT NONE
29     
30           DOUBLE PRECISION InitVal
31           INTEGER :: I, M, IJK
32     
33           DO IJK = IJKSTART3, IJKEND3
34     
35             IF(.NOT.FLUID_AT(IJK)) cycle
36             DO M = 1, MMAX
37               ! A.P.
38               ! Initializing weights as volume_fraction/number_of_nodes
39               ! Note that MFIX doesn't initialize the volume fraction in
40               ! boundaries, and QMOM need a positive value there.
41               ! Here the boundaries are not considered, and they will
42               ! be initialized separetely in the correspoding module
43     
44               DO I = 1, QMOMK_NN
45                  QMOMK_N0 (I, IJK, M) = ROP_S(IJK, M)/(QMOMK_NN * RO_s(IJK,M))
46               END DO
47     
48               ! A.P. Granular temperature minimum value is bounded
49               !      This is required by QMOM
50               InitVal = MAX(SQRT(THETA_M(IJK,M)), MINIMUM_THETA)
51     
52               QMOMK_U0(1, IJK, M) = -InitVal + U_S(IJK, M)
53               QMOMK_U0(2, IJK, M) = +InitVal + U_S(IJK, M)
54               QMOMK_U0(3, IJK, M) = -InitVal + U_S(IJK, M)
55               QMOMK_U0(4, IJK, M) = +InitVal + U_S(IJK, M)
56               QMOMK_U0(5, IJK, M) = -InitVal + U_S(IJK, M)
57               QMOMK_U0(6, IJK, M) = +InitVal + U_S(IJK, M)
58               QMOMK_U0(7, IJK, M) = -InitVal + U_S(IJK, M)
59               QMOMK_U0(8, IJK, M) = +InitVal + U_S(IJK, M)
60     
61               QMOMK_V0(1, IJK, M) = -InitVal + V_S(IJK, M)
62               QMOMK_V0(2, IJK, M) = -InitVal + V_S(IJK, M)
63               QMOMK_V0(3, IJK, M) = +InitVal + V_S(IJK, M)
64               QMOMK_V0(4, IJK, M) = +InitVal + V_S(IJK, M)
65               QMOMK_V0(5, IJK, M) = -InitVal + V_S(IJK, M)
66               QMOMK_V0(6, IJK, M) = -InitVal + V_S(IJK, M)
67               QMOMK_V0(7, IJK, M) = +InitVal + V_S(IJK, M)
68               QMOMK_V0(8, IJK, M) = +InitVal + V_S(IJK, M)
69     
70               QMOMK_W0(1, IJK, M) = -InitVal + W_S(IJK, M)
71               QMOMK_W0(2, IJK, M) = -InitVal + W_S(IJK, M)
72               QMOMK_W0(3, IJK, M) = -InitVal + W_S(IJK, M)
73               QMOMK_W0(4, IJK, M) = -InitVal + W_S(IJK, M)
74               QMOMK_W0(5, IJK, M) = +InitVal + W_S(IJK, M)
75               QMOMK_W0(6, IJK, M) = +InitVal + W_S(IJK, M)
76               QMOMK_W0(7, IJK, M) = +InitVal + W_S(IJK, M)
77               QMOMK_W0(8, IJK, M) = +InitVal + W_S(IJK, M)
78             END DO
79           END DO
80     
81           DO IJK = IJKSTART3, IJKEND3
82            DO M = 1, MMAX
83              IF (FLUID_AT(IJK)) THEN
84                CALL MOMENTS_TWENTY_EIGHT_NODES (QMOMK_N0(:,IJK,M), &
85                     QMOMK_U0(:,IJK,M), QMOMK_V0(:,IJK,M), QMOMK_W0(:,IJK,M), QMOMK_M0(:,IJK,M))
86     
87                CALL EIGHT_NODE_3D (QMOMK_M0(:,IJK,M), QMOMK_N0(:,IJK,M), &
88                     QMOMK_U0(:,IJK,M), QMOMK_V0(:,IJK,M), QMOMK_W0(:,IJK,M))
89     
90                CALL MOMENTS_TWENTY_EIGHT_NODES (QMOMK_N0(:,IJK,M), &
91                     QMOMK_U0(:,IJK,M), QMOMK_V0(:,IJK,M), QMOMK_W0(:,IJK,M), QMOMK_M0(:,IJK,M))
92              END IF
93            END DO
94           END DO
95     
96           END SUBROUTINE QMOMK_INITIAL_CONDITIONS
97