OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
alefvm_sfint3.F File Reference
#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "vect01_c.inc"
#include "inter22.inc"
#include "param_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine alefvm_sfint3 (ixs, nv46, ale_connect, ialefvm_flg, ipm, iparg, ng, x, iad22, nel)

Function/Subroutine Documentation

◆ alefvm_sfint3()

subroutine alefvm_sfint3 ( integer, dimension(nixs,*) ixs,
integer nv46,
type(t_ale_connectivity), intent(in) ale_connect,
integer ialefvm_flg,
integer, dimension(npropmi,*) ipm,
integer, dimension(nparg,*) iparg,
integer ng,
x,
iad22,
integer nel )

Definition at line 33 of file alefvm_sfint3.F.

37C-----------------------------------------------
38C D e s c r i p t i o n
39C-----------------------------------------------
40C 'alefvm' is related to a collocated scheme (built from FVM and based on Godunov scheme)
41C which was temporarily introduced for experimental option /INTER/TYPE22 (FSI coupling with cut cell method)
42C This cut cell method is not completed, abandoned, and is not an official option.
43C There is no other use for this scheme which is automatically enabled when /INTER/TYPE22 is defined (INT22>0 => IALEFVM=1).
44C
45C This subroutine is treating an uncut cell.
46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE alefvm_mod
50 USE i22bufbric_mod !debug
51 USE i22tri_mod
52 USE ale_connectivity_mod
53C-----------------------------------------------
54C I m p l i c i t T y p e s
55C-----------------------------------------------
56#include "implicit_f.inc"
57C-----------------------------------------------
58C G l o b a l P a r a m e t e r s
59C-----------------------------------------------
60#include "mvsiz_p.inc"
61C-----------------------------------------------
62C C o m m o n B l o c k s
63C-----------------------------------------------
64#include "vect01_c.inc"
65#include "inter22.inc"
66#include "param_c.inc"
67C-----------------------------------------------
68C D e s c r i p t i o n
69C-----------------------------------------------
70C This subroutines computes internal forces for
71C finite volume scheme (IALEFVM==1)
72C
73C If option is not detected in input file then
74C subroutine is unplugged
75C-----------------------------------------------
76C D u m m y A r g u m e n t s
77C-----------------------------------------------
78 INTEGER NEL
79 INTEGER :: IXS(NIXS,*),NV46,IALEFVM_FLG, IPM(NPROPMI,*),IPARG(NPARG,*),NG
80 my_real :: x(3,*)
81 my_real :: iad22(*)
82 TYPE(t_ale_connectivity), INTENT(IN) :: ALE_CONNECT
83C-----------------------------------------------
84C L o c a l V a r i a b l e s
85C-----------------------------------------------
86 INTEGER :: IAD2, IAD3
87 INTEGER :: I, II, IV, J, JV, IMAT, ILAW,IFLG_ALE,IFLG_EUL
88 INTEGER :: NIN, IB, IPRES_MOM
89 INTEGER :: NC1,NC2,NC3,NC4,NC5,NC6,NC7,NC8
90 my_real :: f0(3,mvsiz), fface(3,nv46,mvsiz)
91 my_real :: nx(6,mvsiz), ny(6,mvsiz), nz(6,mvsiz), p1,p2,denom,pf
92 my_real :: theta, m1, m2, mf
93 my_real ::
94 . x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz), x5(mvsiz), x6(mvsiz), x7(mvsiz), x8(mvsiz),
95 . y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz), y5(mvsiz), y6(mvsiz), y7(mvsiz), y8(mvsiz),
96 . z1(mvsiz), z2(mvsiz), z3(mvsiz), z4(mvsiz), z5(mvsiz), z6(mvsiz), z7(mvsiz), z8(mvsiz),
97 . z_1, z_2, u1n1, u2n1
98
99C-----------------------------------------------
100C P r e - C o n d i t i o n s
101C-----------------------------------------------
102 IF(alefvm_param%IEnabled==0) RETURN
103 imat = ixs(1,1+nft)
104 ilaw = ipm(2,imat)
105 ialefvm_flg = ipm(251,imat)
106 IF(ialefvm_flg <= 1)RETURN
107 IF(ilaw==11) RETURN
108 iflg_ale = iparg(7,ng)
109 iflg_eul = iparg(11,ng)
110C-----------------------------------------------
111C S o u r c e L i n e s
112C-----------------------------------------------
113
114 nin = 1
115
116 !-------------------------------------------------------------!
117 ! NORMAL VECTOR FOR ALE !
118 !-------------------------------------------------------------!
119! IF(JALE==1 .OR. JEUL==0)THEN !lagrangian bricks also
120 DO i=1,nel
121 ii = i + nft
122 !---8 local node numbers NC1 TO NC8 for solid element I ---!
123 nc1=ixs(2,ii)
124 nc2=ixs(3,ii)
125 nc3=ixs(4,ii)
126 nc4=ixs(5,ii)
127 nc5=ixs(6,ii)
128 nc6=ixs(7,ii)
129 nc7=ixs(8,ii)
130 nc8=ixs(9,ii)
131 !
132 !---Coordinates of the 8 nodes
133 x1(i)=x(1,nc1)
134 y1(i)=x(2,nc1)
135 z1(i)=x(3,nc1)
136 !
137 x2(i)=x(1,nc2)
138 y2(i)=x(2,nc2)
139 z2(i)=x(3,nc2)
140 !
141 x3(i)=x(1,nc3)
142 y3(i)=x(2,nc3)
143 z3(i)=x(3,nc3)
144 !
145 x4(i)=x(1,nc4)
146 y4(i)=x(2,nc4)
147 z4(i)=x(3,nc4)
148 !
149 x5(i)=x(1,nc5)
150 y5(i)=x(2,nc5)
151 z5(i)=x(3,nc5)
152 !
153 x6(i)=x(1,nc6)
154 y6(i)=x(2,nc6)
155 z6(i)=x(3,nc6)
156 !
157 x7(i)=x(1,nc7)
158 y7(i)=x(2,nc7)
159 z7(i)=x(3,nc7)
160 !
161 x8(i)=x(1,nc8)
162 y8(i)=x(2,nc8)
163 z8(i)=x(3,nc8)
164 ENDDO
165 DO i=1,nel
166 ! Face-1
167 nx(1,i)=(y3(i)-y1(i))*(z2(i)-z4(i)) - (z3(i)-z1(i))*(y2(i)-y4(i))
168 ny(1,i)=(z3(i)-z1(i))*(x2(i)-x4(i)) - (x3(i)-x1(i))*(z2(i)-z4(i))
169 nz(1,i)=(x3(i)-x1(i))*(y2(i)-y4(i)) - (y3(i)-y1(i))*(x2(i)-x4(i))
170 ! Face-2
171 nx(2,i)=(y7(i)-y4(i))*(z3(i)-z8(i)) - (z7(i)-z4(i))*(y3(i)-y8(i))
172 ny(2,i)=(z7(i)-z4(i))*(x3(i)-x8(i)) - (x7(i)-x4(i))*(z3(i)-z8(i))
173 nz(2,i)=(x7(i)-x4(i))*(y3(i)-y8(i)) - (y7(i)-y4(i))*(x3(i)-x8(i))
174 ! Face-3
175 nx(3,i)=(y6(i)-y8(i))*(z7(i)-z5(i)) - (z6(i)-z8(i))*(y7(i)-y5(i))
176 ny(3,i)=(z6(i)-z8(i))*(x7(i)-x5(i)) - (x6(i)-x8(i))*(z7(i)-z5(i))
177 nz(3,i)=(x6(i)-x8(i))*(y7(i)-y5(i)) - (y6(i)-y8(i))*(x7(i)-x5(i))
178 ! Face-4
179 nx(4,i)=(y2(i)-y5(i))*(z6(i)-z1(i)) - (z2(i)-z5(i))*(y6(i)-y1(i))
180 ny(4,i)=(z2(i)-z5(i))*(x6(i)-x1(i)) - (x2(i)-x5(i))*(z6(i)-z1(i))
181 nz(4,i)=(x2(i)-x5(i))*(y6(i)-y1(i)) - (y2(i)-y5(i))*(x6(i)-x1(i))
182 ! Face-5
183 nx(5,i)=(y7(i)-y2(i))*(z6(i)-z3(i)) - (z7(i)-z2(i))*(y6(i)-y3(i))
184 ny(5,i)=(z7(i)-z2(i))*(x6(i)-x3(i)) - (x7(i)-x2(i))*(z6(i)-z3(i))
185 nz(5,i)=(x7(i)-x2(i))*(y6(i)-y3(i)) - (y7(i)-y2(i))*(x6(i)-x3(i))
186 ! Face-6
187 nx(6,i)=(y8(i)-y1(i))*(z4(i)-z5(i)) - (z8(i)-z1(i))*(y4(i)-y5(i))
188 ny(6,i)=(z8(i)-z1(i))*(x4(i)-x5(i)) - (x8(i)-x1(i))*(z4(i)-z5(i))
189 nz(6,i)=(x8(i)-x1(i))*(y4(i)-y5(i)) - (y8(i)-y1(i))*(x4(i)-x5(i))
190 ENDDO
191! ENDIF
192
193 !-------------------------------------------------------------!
194 ! Assembling Forces !
195 !-------------------------------------------------------------!
196
197 !IPRES_MOM = 0,1,2,3,4 : (P1+P2)/2
198 !IPRES_MOM = 5 : (rho1c1P1+rho2c2P2)/(rho1c1+rho2c2) + (rho1c1*rho2c2/(rho1c1+rho2c2)) <uel-uadj,nel> acoustic solver
199
200 ipres_mom = alefvm_param%ISOLVER
201
202! print *, "IPRES_MOM, SFINT3", IPRES_MOM
203
204 SELECT CASE (ipres_mom)
205 CASE(5) ! acoustic solver
206 DO i=1,nel
207 ii = i + nft
208 iad2 = ale_connect%ee_connect%iad_connect(ii)
209 p1 = alefvm_buffer%F_FACE(1,4,ii)
210 m1 = alefvm_buffer%F_FACE(1,5,ii)
211 z_1 = alefvm_buffer%F_FACE(1,3,ii) ! impedance
212 DO j=1,nv46
213 iv = ale_connect%ee_connect%connected(iad2 + j - 1)
214 IF(iv > 0)THEN
215 !--ADJACENT ELEM DOES EXIST
216 iad3 = ale_connect%ee_connect%iad_connect(iv)
217 DO jv=1,nv46
218 IF(ale_connect%ee_connect%connected(iad3 + jv - 1)==ii)EXIT
219 ENDDO
220 u1n1 = + alefvm_buffer%F_FACE(3, j,ii)
221 u2n1 = - alefvm_buffer%F_FACE(3,jv,iv) !UN2 = <U2,N2> = <U2,-N1>
222 z_2 = alefvm_buffer%F_FACE(1, 3,iv)
223 p2 = alefvm_buffer%F_FACE(1, 4,iv)
224 m2 = alefvm_buffer%F_FACE(1, 5,iv)
225 denom = z_1 + z_2
226 mf = min(m1,m2)
227 theta = min(one,mf) !see dellacherie,omnes,raviart : fixing low mach issue
228 pf = (z_1*p2 + z_2*p1)/denom + theta*(z_1*z_2*(u1n1-u2n1)/denom)
229 ELSE
230 !--ADJACENT ELEM DOES NOT EXIST
231 !SLIDING RWALL BC
232 u1n1 = + alefvm_buffer%F_FACE(3, j,ii)
233 !Pf = P1
234 mf = m1
235 theta = min(one,mf) !see dellacherie,omnes,raviart : fixing low mach issue
236 pf = p1 + theta*half*z_1*u1n1
237 !Pf = P1
238 ENDIF
239 fface(1,j,i) = -half*pf*nx(j,i) !N= 2.SURF.n , where abs(n)=1 => abs(N) = 2S , this lines leads to FFACE(1J,I) = Pf.S.n
240 fface(2,j,i) = -half*pf*ny(j,i) ! negative sign because int(div.sigma,DV) => int(-P,dS)
241 fface(3,j,i) = -half*pf*nz(j,i)
242 enddo!next J
243 enddo!next I
244 CASE DEFAULT
245 DO i=1,nel
246 ii = i + nft
247 iad2 = ale_connect%ee_connect%iad_connect(ii)
248 p1 = alefvm_buffer%F_FACE(1,4,ii)
249 DO j=1,nv46
250 iv = ale_connect%ee_connect%connected(iad2 + j - 1)
251 IF(iv > 0)THEN
252 !--ADJACENT ELEM DOES EXIST
253 iad3 = ale_connect%ee_connect%iad_connect(iv)
254 DO jv=1,nv46
255 IF(ale_connect%ee_connect%connected(iad3 + jv - 1)==ii)EXIT
256 ENDDO
257 p2 = alefvm_buffer%F_FACE(1,4,iv)
258 pf = half*(p1+p2)
259 ELSE
260 !--ADJACENT ELEM DOES NOT EXIST
261 !SLIDING RWALL BC
262 pf = p1
263 ENDIF
264 fface(1,j,i) = -half*pf*nx(j,i) ! negative sign because int(div.sigma,DV) => int(-P,dS)
265 fface(2,j,i) = -half*pf*ny(j,i)
266 fface(3,j,i) = -half*pf*nz(j,i)
267 !print *, "P1,P2,Pf", P1,P2,Pf
268 !print *, " normal ", NX(J,I), NY(J,I), NZ(J,I)
269c write ( *, FMT='(A,I10,A,I3,A,3F20.12,A,F20.12)')
270c ." id=", ixs(11,ii)," ncycle=", NCYCLE, " F=", FFACE(1:3,J,I) , "Pf=", Pf
271 enddo!next J
272 enddo!next I
273 END SELECT
274
275 !-------------------------------------------------------------!
276 ! INTEGRAL ON EACH FACE from Integral(DIV(SIGMA),Volume) !
277 !-------------------------------------------------------------!
278 DO i=1,nel
279 ii = i + nft
280 IF(int22/=0)THEN
281 ib = nint(iad22(i))
282 IF(ib>0) cycle
283 ENDIF
284 f0(1,i) = sum(fface(1,1:nv46,i))
285 f0(2,i) = sum(fface(2,1:nv46,i))
286 f0(3,i) = sum(fface(3,1:nv46,i))
287
288 ! write (*,FMT='(A,I6,A,3F30.16)') " brick ID=", ixs(11,II)," Fint=", F0(1:3,I)
289
290
291 alefvm_buffer%FINT_CELL(1,ii) = f0(1,i)
292 alefvm_buffer%FINT_CELL(2,ii) = f0(2,i)
293 alefvm_buffer%FINT_CELL(3,ii) = f0(3,i)
294
295 !Fcell already contains gravity force, so do not erase but add
296 alefvm_buffer%FCELL(1,ii) = alefvm_buffer%FCELL(1,ii) + f0(1,i) + alefvm_buffer%FEXT_CELL(1,ii)
297 alefvm_buffer%FCELL(2,ii) = alefvm_buffer%FCELL(2,ii) + f0(2,i) + alefvm_buffer%FEXT_CELL(2,ii)
298 alefvm_buffer%FCELL(3,ii) = alefvm_buffer%FCELL(3,ii) + f0(3,i) + alefvm_buffer%FEXT_CELL(3,ii)
299 enddo!next I
300
301
302
303 !-------------------------------------------------------------!
304 ! DBUG OUTPUT !
305 !-------------------------------------------------------------!
306! if(IALEFVM_OUTP_FINT /= 0)then
307!
308! if(itask==0)then
309! print *, " |----alefvm_sfint3.F-----|"
310! print *, " | THREAD INFORMATION |"
311! print *, " |------------------------|"
312! print *, " NCYCLE =", NCYCLE
313!
314! do i=1,nb
315! ii = nft + i
316! IF(INT22/=0)THEN
317! IB = NINT(IAD22(I))
318! IF(IB>0) CYCLE
319! ENDIF
320! print *, " brique=", ixs(11,nft+i)
321! write(*,FMT='(A34,6A26)') " ",
322! . "#--- internal force -----#"
323! write (*,FMT='(A,8E26.14)') " force-X =", FCELL(1,II)
324! write (*,FMT='(A,8E26.14)') " force-Y =", FCELL(2,II)
325! write (*,FMT='(A,8E26.14)') " force-Z =", FCELL(3,II)
326! print *, " "
327! enddo
328!
329! endif!if(itask)
330!
331! endif!if(IALEFVM_OUTP_FINT /= 0)
332 !-----------------------------------------!
333
334 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
min
#define min(a, b)
Definition macros.h:20
ale_connectivity_mod
Definition ale_connectivity_mod.F:223
alefvm_mod
Definition alefvm_mod.F:88
alefvm_mod::alefvm_buffer
type(alefvm_buffer_), target alefvm_buffer
Definition alefvm_mod.F:120
alefvm_mod::alefvm_param
type(alefvm_param_), target alefvm_param
Definition alefvm_mod.F:121
i22bufbric_mod
Definition cut-cell-search_mod.F:95
i22tri_mod
Definition cut-cell-search_mod.F:212