OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
ebcs6_inip.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| ebcs6_inip_mod ../engine/source/boundary_conditions/ebcs/ebcs6_inip.F
25!||--- called by ------------------------------------------------------
26!|| ebcs_main ../engine/source/boundary_conditions/ebcs/ebcs_main.F
27!||====================================================================
29 IMPLICIT NONE
30 CONTAINS
31!||====================================================================
32!|| ebcs6_inip ../engine/source/boundary_conditions/ebcs/ebcs6_inip.F
33!||--- called by ------------------------------------------------------
34!|| ebcs_main ../engine/source/boundary_conditions/ebcs/ebcs_main.F
35!||--- uses -----------------------------------------------------
36!|| ebcs_mod ../common_source/modules/boundary_conditions/ebcs_mod.F90
37!|| output_mod ../common_source/modules/output/output_mod.f90
38!|| segvar_mod ../engine/share/modules/segvar_mod.F
39!||====================================================================
40 SUBROUTINE ebcs6_inip(NSEG,ISEG,SEGVAR,
41 . A,V,X,
42 . LISTE,NOD,IRECT,
43 . RO0,EN0,P0,
44 . VO,PO,LA,
45 . MS,STIFN,EBCS,OUTPUT,DT1,TIME)
46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE ebcs_mod
50 USE segvar_mod
51 USE output_mod , ONLY : output_
52C-----------------------------------------------
53C Pression arret initiale imposee avec impedance
54C-----------------------------------------------
55C I m p l i c i t T y p e s
56C-----------------------------------------------
57#include "implicit_f.inc"
58C-----------------------------------------------
59C C o m m o n B l o c k s
60C-----------------------------------------------
61#include "param_c.inc"
62C-----------------------------------------------
63C D u m m y A r g u m e n t s
64C-----------------------------------------------
65 INTEGER NSEG,NOD,ISEG(NSEG),LISTE(NOD),IRECT(4,NSEG)
66 my_real :: A(3,*),X(3,*),V(3,*),LA(3,NOD),RO0(NSEG),EN0(NSEG),P0(NOD),VO(NOD),PO(NOD),MS(*),STIFN(*)
67 TYPE(t_ebcs_inip), INTENT(IN) :: EBCS
68 TYPE(t_segvar) :: SEGVAR
69 my_real,INTENT(IN) :: DT1 !< time step
70 my_real,INTENT(IN) :: time !< current time
71 TYPE(output_), INTENT(INOUT) :: OUTPUT !< output structure
72C-----------------------------------------------
73C L o c a l V a r i a b l e s
74C-----------------------------------------------
75 INTEGER :: I,IS,KSEG,N1,N2,N3,N4,NG1,NG2,NG3,NG4,N
76 my_real :: orient,rho,c,lcar,roc,alp,fac,
77 . x13,y13,z13,x24,y24,z24,nx,ny,nz,s,
78 . roou,enou,vmx,vmy,vmz,fluxi,fluxo,vn,pn,du,dp,p,dpdv
79 my_real :: de, de_in, de_out, dm_in, dm_out
80C-----------------------------------------------
81 c=ebcs%c
82 rho=ebcs%rho
83 roc=rho*c
84 lcar=ebcs%lcar
85 alp=zero
86 de_in = zero
87 de_out = zero
88 dm_in = zero
89 dm_out = zero
90
91 IF (lcar > zero) alp=c*dt1/lcar
92C Save densite and initial energy
93 IF(time == zero)THEN
94 DO is=1,nseg
95 kseg=abs(iseg(is))
96 ro0(is) = segvar%RHO(kseg)
97 en0(is) = segvar%EINT(kseg)
98 ENDDO
99 ENDIF
100C SURFACE NORMALE NODALES
101 DO i=1,nod
102 la(1,i)=zero
103 la(2,i)=zero
104 la(3,i)=zero
105 ENDDO
106
107 DO is=1,nseg
108 kseg=abs(iseg(is))
109 orient=float(iseg(is)/kseg)
110 n1=irect(1,is)
111 n2=irect(2,is)
112 n3=irect(3,is)
113 n4=irect(4,is)
114 IF(n4 == 0 .OR. n4 == n3) THEN
115 fac=one_over_6*orient
116 n4=n3
117 ELSE
118 fac=one_over_8*orient
119 ENDIF
120c
121 ng1=liste(n1)
122 ng2=liste(n2)
123 ng3=liste(n3)
124 ng4=liste(n4)
125 x13=x(1,ng3)-x(1,ng1)
126 y13=x(2,ng3)-x(2,ng1)
127 z13=x(3,ng3)-x(3,ng1)
128 x24=x(1,ng4)-x(1,ng2)
129 y24=x(2,ng4)-x(2,ng2)
130 z24=x(3,ng4)-x(3,ng2)
131c
132 nx=(y13*z24-z13*y24)*fac
133 ny=(z13*x24-x13*z24)*fac
134 nz=(x13*y24-y13*x24)*fac
135c
136 la(1,n1)=la(1,n1)+nx
137 la(2,n1)=la(2,n1)+ny
138 la(3,n1)=la(3,n1)+nz
139 la(1,n2)=la(1,n2)+nx
140 la(2,n2)=la(2,n2)+ny
141 la(3,n2)=la(3,n2)+nz
142 la(1,n3)=la(1,n3)+nx
143 la(2,n3)=la(2,n3)+ny
144 la(3,n3)=la(3,n3)+nz
145C
146 vmx=v(1,ng1)+v(1,ng2)+v(1,ng3)
147 vmy=v(2,ng1)+v(2,ng2)+v(2,ng3)
148 vmz=v(3,ng1)+v(3,ng2)+v(3,ng3)
149 IF(n4/=n3) THEN
150 la(1,n4)=la(1,n4)+nx
151 la(2,n4)=la(2,n4)+ny
152 la(3,n4)=la(3,n4)+nz
153 vmx=vmx+v(1,ng4)
154 vmy=vmy+v(2,ng4)
155 vmz=vmz+v(3,ng4)
156 ENDIF
157C
158c mass and total energy balance
159c
160 roou = segvar%RHO(kseg)
161 enou = segvar%EINT(kseg)
162C
163 ! post-treatment
164 fluxo=(vmx*nx+vmy*ny+vmz*nz)*dt1
165 fluxi=min(fluxo,zero)
166 fluxo=max(fluxo,zero)
167 dm_out=dm_out-fluxo*roou
168 dm_in=dm_in-fluxi*ro0(is)
169 de_out=de_out-fluxo*enou
170 de_in=de_in-fluxi*en0(is)
171
172C
173C storage of density and incoming energy in facet buffer
174C
175 segvar%RHO(kseg)=ro0(is)
176 segvar%EINT(kseg)=en0(is)
177
178 ENDDO
179
180C
181C compute initial nodal pressure
182C
183 IF (time == zero) THEN
184 DO i=1,nod
185 n=liste(i)
186 s=la(1,i)**2+la(2,i)**2+la(3,i)**2
187 vn=(v(1,n)*la(1,i)+v(2,n)*la(2,i)+v(3,n)*la(3,i))/sqrt(s)
188c We memorize the stop pressure
189 p0(i)=p0(i)/s
190 IF(vn<zero)p0(i)=p0(i)+half*rho*vn**2
191c write(6,*)'nodal P', P0(I)
192 vo(i)=vn
193 po(i)=p0(i)
194 ENDDO
195 ENDIF
196C
197 DO i=1,nod
198 n=liste(i)
199 s=sqrt(la(1,i)**2+la(2,i)**2+la(3,i)**2)
200 vn=(v(1,n)*la(1,i)+v(2,n)*la(2,i)+v(3,n)*la(3,i))/s
201 dpdv=roc
202c stopping condition
203 pn=p0(i)
204 IF(vn<0)THEN
205 pn=p0(i)-half*rho*vn**2
206 dpdv=dpdv-rho*vn
207 ENDIF
208c silent boundary
209 du=roc*(vn-vo(i))
210 dp=alp*(pn-po(i))
211 p=po(i)+dp+du
212c write(6,*)'compute', P,P0(I)
213c
214 a(1,n)=a(1,n)-p*la(1,i)
215 a(2,n)=a(2,n)-p*la(2,i)
216 a(3,n)=a(3,n)-p*la(3,i)
217 stifn(n)=stifn(n)+(two*(s*dpdv)**2)/ms(n)
218C
219 de=-half*(po(i)+p)*dt1*vn*s
220 de_in = de_in + max(de, zero)
221 de_out = de_out + min(de, zero)
222C
223 vo(i)=vn
224 po(i)=p
225 ENDDO
226
227!$OMP CRITICAL
228 output%DATA%INOUT%DM_IN = output%DATA%INOUT%DM_IN + dm_in
229 output%DATA%INOUT%DM_OUT = output%DATA%INOUT%DM_OUT + dm_out
230 output%DATA%INOUT%DE_IN = output%DATA%INOUT%DE_IN + de_in
231 output%DATA%INOUT%DE_OUT = output%DATA%INOUT%DE_OUT + de_out
232!$OMP END CRITICAL
233
234c-----------
235 RETURN
236
237 END SUBROUTINE ebcs6_inip
238 END MODULE ebcs6_inip_mod
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine ebcs6_inip(nseg, iseg, segvar, a, v, x, liste, nod, irect, ro0, en0, p0, vo, po, la, ms, stifn, ebcs, output, dt1, time)
Definition ebcs6_inip.F:46