36
37
38
40 USE elbufdef_mod
41 USE multi_fvm_mod
42 use element_mod , only : nixs,nixq,nixtg
43
44
45
46#include "implicit_f.inc"
47
48
49
50#include "com01_c.inc"
51#include "param_c.inc"
52#include "mvsiz_p.inc"
53
54
55
57 INTEGER, INTENT(IN) :: NG
58 TYPE(ELBUF_STRUCT_), TARGET, DIMENSION(NGROUP) :: ELBUF_TAB
59 INTEGER, INTENT(IN) :: IPARG(NPARG, *)
60 INTEGER, INTENT(IN) :: ITASK
61 INTEGER, INTENT(IN) :: IXS(NIXS, *), IXQ(NIXQ, *), IXTG(NIXTG, *)
62 TYPE(MULTI_FVM_STRUCT), INTENT(INOUT) :: MULTI_FVM
63 my_real,
INTENT(IN) :: pm(npropm, *)
64 INTEGER, INTENT(IN) :: (NPROPMI, *)
65 my_real,
INTENT(IN) :: gravity(4, *)
66 my_real,
INTENT(IN) :: current_time
67
68
69
70 TYPE(G_BUFEL_), POINTER :: GBUF
71 INTEGER :: II, I, J, JJ, KFACE, KFACE2, NB_FACE
72 INTEGER :: IPLA
75 my_real :: volfrac, massfrac, massfracii, massfracjj, eintii, eintjj
77 my_real :: massflux, ss, sumflux(3), tmp2, tmp3, tmp4, normal_vel
78 my_real :: sr, sl, sstar, surf, vii(3), vjj(3), fii(3), fjj(3)
79 my_real :: fiistar(3), fjjstar(3), viistar(3), vjjstar(3)
81 my_real :: normal_velii, normal_veljj
82 my_real :: alphaii, alphajj, rhoii, rhojj, rhoeii, rhoejj, pii, pjj
84 my_real :: pstar, estar, pinf, gam, gravii(3)
85 INTEGER :: LOCAL_MATID, MATLAW, NVERTEX, NODEID, IMAT
86 my_real :: sum_normalvel, pshift
87 LOGICAL ::
88 TYPE(LBUF_PTR) :: LBUFS(MULTI_FVM%NBMAT)
89 INTEGER :: ISOLNOD, NEL, ITY, NFT
90
91
92
93
94 gbuf => elbuf_tab(ng)%GBUF
95 nel = iparg(2, ng)
96 nft = iparg(3, ng)
97 ity = iparg(5, ng)
98 isolnod = iparg(28, ng)
99
100 pshift = multi_fvm%PRES_SHIFT
101
102 IF (multi_fvm%NBMAT > 1) THEN
103
104 DO imat = 1, multi_fvm%NBMAT
105 lbufs(imat)%LBUF => elbuf_tab(ng)%BUFLY(imat)%LBUF(1, 1, 1)
106 ENDDO
107 ENDIF
108 DO imat = 1, multi_fvm%NBMAT
109
110 IF (multi_fvm%SYM == 0) THEN
111 nb_face = 6
112 local_matid = ipm(20 + imat, ixs(1, 1 + nft))
113 ELSEIF (ity == 2) THEN
114
115 nb_face = 4
116 local_matid = ipm(20 + imat, ixq(1, 1 + nft))
117 ELSEIF (ity == 7) THEN
118
119 nb_face = 3
120 local_matid = ipm(20 + imat, ixtg(1, 1 + nft))
121 ENDIF
122 matlaw = ipm(2, local_matid)
123
124
125 DO ii = 1, nel
126 i = ii + nft
127
128 IF (multi_fvm%SYM == 0 .AND. isolnod /= 4) THEN
129
130 sumflux(1) = multi_fvm%SUBVOL_FLUXES(imat, 1, i) + multi_fvm%SUBVOL_FLUXES(imat, 2, i) +
131 . multi_fvm%SUBVOL_FLUXES(imat, 3, i) + multi_fvm%SUBVOL_FLUXES
132 . multi_fvm%SUBVOL_FLUXES(imat, 5, i) + multi_fvm%SUBVOL_FLUXES(imat,
133
134 sumflux(2) = multi_fvm%SUBMASS_FLUXES(imat, 1, i) + multi_fvm%SUBMASS_FLUXES(imat, 2, i) +
135 . multi_fvm%SUBMASS_FLUXES(imat, 3, i) + multi_fvm%SUBMASS_FLUXES(imat, 4, i) +
136 . multi_fvm%SUBMASS_FLUXES(imat, 5, i) + multi_fvm%SUBMASS_FLUXES(imat, 6, i)
137
138 sumflux(3) = multi_fvm%SUBENER_FLUXES(imat, 1, i) + multi_fvm%SUBENER_FLUXES(imat, 2, i) +
139 . multi_fvm%SUBENER_FLUXES(imat, 3, i) + multi_fvm%SUBENER_FLUXES(imat, 4, i) +
140 . multi_fvm%SUBENER_FLUXES(imat, 5, i) + multi_fvm%SUBENER_FLUXES(imat, 6, i)
141 ELSEIF (isolnod == 4) THEN
142
143 sumflux(1) = multi_fvm%SUBVOL_FLUXES(imat, 5, i) + multi_fvm%SUBVOL_FLUXES
144 . multi_fvm%SUBVOL_FLUXES(imat, 2, i) + multi_fvm%SUBVOL_FLUXES(imat, 4, i)
145
146 sumflux(2) = multi_fvm%SUBMASS_FLUXES(imat, 5, i) + multi_fvm%SUBMASS_FLUXES(imat, 6, i) +
147 . multi_fvm%SUBMASS_FLUXES(imat, 2, i) + multi_fvm%SUBMASS_FLUXES(imat, 4, i)
148
149 sumflux(3) = multi_fvm%SUBENER_FLUXES(imat, 5, i) + multi_fvm%SUBENER_FLUXES(imat, 6, i) +
150 . multi_fvm%SUBENER_FLUXES(imat, 2, i) + multi_fvm%SUBENER_FLUXES(imat, 4, i)
151 ELSE
152
153
154 sumflux(1) = multi_fvm%SUBVOL_FLUXES(imat, 1, i) + multi_fvm%SUBVOL_FLUXES(imat, 2, i) +
155 . multi_fvm%SUBVOL_FLUXES(imat, 3, i)
156
157 sumflux(2) = multi_fvm%SUBMASS_FLUXES(imat, 1, i) + multi_fvm%SUBMASS_FLUXES(imat, 2, i) +
158 . multi_fvm%SUBMASS_FLUXES(imat, 3, i)
159
160 sumflux(3) = multi_fvm%SUBENER_FLUXES(imat, 1, i) + multi_fvm%SUBENER_FLUXES(imat, 2, i) +
161 . multi_fvm%SUBENER_FLUXES(imat, 3, i)
162 IF (ity == 2) THEN
163
164
165 sumflux(1) = sumflux(1) + multi_fvm%SUBVOL_FLUXES(imat, 4, i)
166
167 sumflux(2) = sumflux(2) + multi_fvm%SUBMASS_FLUXES(imat, 4, i)
168
169 sumflux(3) = sumflux(3) + multi_fvm%SUBENER_FLUXES(imat, 4, i)
170 ENDIF
171 ENDIF
172
173 sum_normalvel = zero
174 IF (multi_fvm%SYM == 0 .AND. isolnod == 4) THEN
175 kface = 5
176 normal_vel = multi_fvm%FLUXES(6, kface, i)
177 sum_normalvel = sum_normalvel + normal_vel
178 kface = 6
179 normal_vel = multi_fvm%FLUXES(6, kface, i)
180 sum_normalvel = sum_normalvel + normal_vel
181 kface = 2
182 normal_vel = multi_fvm%FLUXES(6, kface, i)
183 sum_normalvel = sum_normalvel + normal_vel
184 kface = 4
185 normal_vel = multi_fvm%FLUXES(6, kface, i)
186 sum_normalvel = sum_normalvel + normal_vel
187 ELSE
188 DO kface = 1, nb_face
189 normal_vel = multi_fvm%FLUXES(6, kface, i)
190 sum_normalvel = sum_normalvel + normal_vel
191 ENDDO
192 ENDIF
193 alphaii = multi_fvm%PHASE_ALPHA(imat, i)
194 pii = multi_fvm%PHASE_PRES(imat, i)
195
196
197 tmp2 = lbufs(imat)%LBUF%VOL(ii) * lbufs(imat)%LBUF%RHO(ii) -
198 . timestep * sumflux(2)
199
200 tmp3 = lbufs(imat)%LBUF%VOL(ii) * lbufs(imat)%LBUF%EINT(ii) -
201 . timestep * (sumflux(3) + alphaii * (pii + pshift) * sum_normalvel)
202
203
204
205 tmp4 = lbufs(imat)%LBUF%VOL(ii) -
206 . timestep * (sumflux(1) - alphaii * sum_normalvel)
207
208
209 multi_fvm%PHASE_ALPHA(imat, i) = tmp4
210
211 multi_fvm%PHASE_RHO(imat, i) = tmp2
212
213 multi_fvm%PHASE_EINT(imat, i) = tmp3
214 ENDDO
215 ENDDO
216
217
218
219
220