OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
iniebcs_dp.F File Reference
#include "implicit_f.inc"
#include "param_c.inc"
#include "com04_c.inc"
#include "com01_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine iniebcs_dp (nseg, nod, iseg, ielem, irect, liste, iparg, elbuf_str, x, ixs, ixq, ixtg, dp0, iparts, ipartq, iparttg)

Function/Subroutine Documentation

◆ iniebcs_dp()

subroutine iniebcs_dp ( integer nseg,
integer nod,
integer, dimension(nseg) iseg,
integer, dimension(nseg) ielem,
integer, dimension(4,nseg) irect,
integer, dimension(nod) liste,
integer, dimension(nparg,ngroup) iparg,
type (elbuf_struct_), dimension(ngroup), target elbuf_str,
x,
integer, dimension(nixs,numels), intent(in) ixs,
integer, dimension(nixq,numelq), intent(in) ixq,
integer, dimension(nixtg,numeltg), intent(in) ixtg,
dp0,
integer, dimension(numels), intent(in) iparts,
integer, dimension(numelq), intent(in) ipartq,
integer, dimension(numeltg), intent(in) iparttg )

Definition at line 30 of file iniebcs_dp.F.

33C-----------------------------------------------
34C D e s c r i p t i o n
35C-----------------------------------------------
36C INPUT : SEGMENT defined from /EBCS/NRF
37C OUTPUT : EBCS%DP0 : static hydropressure increment due to gravity loading if defined (/INIGRAV)
38C pre-condition : /EBCS/NRF defined (already check from parent subroutine)
39C NSEG >0 : LOOP is from 1 to NSEG ; NSEG=0 : subroutine does nothing
40C
41C Comments: DP0 = Rho0 * Grav0 * Dist
42C rho0 is retrieved from adjacent cell, initialized possibly by /INIGRAV, so it is in GBUF%RHO and not PM( )
43C grav0 is stored during INIGRAV procedure
44C dist is computed here, it is 2*dh where h=distance(cell_centroid, cell_face). It corresponds to distance to the centroid of a ghost cell
45C-----------------------------------------------
46C M o d u l e s
47C-----------------------------------------------
48 USE elbufdef_mod
49 USE inigrav
50 use element_mod , only : nixs,nixq,nixtg
51C-----------------------------------------------
52C I m p l i c i t T y p e s
53C-----------------------------------------------
54#include "implicit_f.inc"
55C-----------------------------------------------
56C C o m m o n B l o c k s
57C-----------------------------------------------
58#include "param_c.inc"
59#include "com04_c.inc"
60#include "com01_c.inc"
61C-----------------------------------------------
62C D u m m y A r g u m e n t s
63C-----------------------------------------------
64 INTEGER :: NSEG,NOD,ISEG(NSEG),IRECT(4,NSEG),LISTE(NOD),IPARG(NPARG,NGROUP), IELEM(NSEG)
65 INTEGER,INTENT(IN) :: IXS(NIXS,NUMELS),IXQ(NIXQ,NUMELQ), IXTG(NIXTG,NUMELTG)
66 INTEGER, INTENT(IN) :: IPARTS(NUMELS), IPARTQ(NUMELQ), IPARTTG(NUMELTG)
67 my_real :: x(3,numnod), dp0(nseg)
68 TYPE (ELBUF_STRUCT_),DIMENSION(NGROUP), TARGET :: ELBUF_STR
69C-----------------------------------------------
70C L o c a l V a r i a b l e s
71C-----------------------------------------------
72 INTEGER :: NodeLOC(4),NG,IS,KSEG,NodeG(8),ESEG,EAD,KTY,KLT,MFT,ISOLNOD,ITY,IP
73 my_real :: orient, fac,zf(3),zc(3),vec(3),rho ,dist,grav0
74 my_real :: ngx,ngy,ngz
75 TYPE(G_BUFEL_) ,POINTER :: GBUF
76 LOGICAL IS_TRIA
77C=======================================================================
78
79 DO is=1,nseg
80
81 kseg = abs(iseg(is))
82 orient = zero
83 IF(kseg /= 0)orient=float(iseg(is)/kseg)
84
85 !local ids in IRECT array
86 nodeloc(1:4) = irect(1:4,is)
87
88 !global internal node ids
89 nodeg(1:4)=liste(nodeloc(1:4))
90
91 is_tria = .false.
92 IF(nodeloc(4) == 0 .OR. nodeloc(3) == nodeloc(4))is_tria=.true.
93
94 !centroid at face
95 IF(is_tria) THEN
96 fac=third*orient
97 !NodeLOC(3)=NodeLOC(4)
98 zf(1) = fac*sum( x(1,nodeg(1:3)) )
99 zf(2) = fac*sum( x(2,nodeg(1:3)) )
100 zf(3) = fac*sum( x(3,nodeg(1:3)) )
101 ELSE
102 fac=fourth*orient
103 zf(1) = fac*sum( x(1,nodeg(1:4)) )
104 zf(2) = fac*sum( x(2,nodeg(1:4)) )
105 zf(3) = fac*sum( x(3,nodeg(1:4)) )
106 ENDIF
107
108 !centroid at cell
109 eseg=ielem(is)
110 !get density
111 mft = 0
112 ity = 0
113 isolnod = 0
114 klt = 0
115 DO ng=1,ngroup
116 kty = iparg(5,ng)
117 klt = iparg(2,ng)
118 mft = iparg(3,ng)
119 ity = iparg(5,ng)
120 isolnod = iparg(28,ng)
121 IF (eseg<=klt+mft) EXIT
122 IF(n2d==0)THEN
123 IF(ity/=1)cycle
124 IF(isolnod == 0)THEN
125 print *,"**ERROR /EBCS/NRF : #2205"
126 cycle
127 ENDIF
128 ELSE
129 IF(ity /= 2 .AND. ity /= 7)cycle
130 ENDIF
131 ENDDO
132 ead = eseg-mft ! at this step it is ensured that EAD \in [1,LLT]
133 IF(ead<=0)cycle
134 gbuf => elbuf_str(ng)%GBUF
135 rho = gbuf%RHO(ead) !retrieve rho brom GBUF and not PM since it could be initialized by INIGRAV
136 !cell centroid
137 IF(ity==1)THEN
138 nodeg(1:8)=ixs(2:9,eseg)
139 zc(1)= sum(x(1,nodeg(1:isolnod)))/isolnod
140 zc(2)= sum(x(2,nodeg(1:isolnod)))/isolnod
141 zc(3)= sum(x(3,nodeg(1:isolnod)))/isolnod
142 ip = iparts(eseg)
143 ELSEIF(ity == 2)THEN
144 nodeg(1:4)=ixq(2:5,eseg)
145 zc(1)= fourth*sum(x(1,nodeg(1:4)))
146 zc(2)= fourth*sum(x(2,nodeg(1:4)))
147 zc(3)= fourth*sum(x(3,nodeg(1:4)))
148 ip = ipartq(eseg)
149 ELSEIF(ity == 7)THEN
150 nodeg(1:3)=ixtg(2:4,eseg)
151 zc(1)= third*sum(x(1,nodeg(1:3)))
152 zc(2)= third*sum(x(2,nodeg(1:3)))
153 zc(3)= third*sum(x(3,nodeg(1:3)))
154 ip = iparttg(eseg)
155 ELSE
156 !not supposed to happen
157 print *, "**ERROR /EBCS/NRF : ONE SEGMENT IS LOCATED TO AN UNEXPECTED TYPE OF ELEMENTS"
158 cycle !next IS (segment)
159 ENDIF
160
161 !-- distance from cell centroid to ghost cell is twice the distance from centroid to face centroid
162 vec(1) = -zc(1)+zf(1)
163 vec(2) = -zc(2)+zf(2)
164 vec(3) = -zc(3)+zf(3)
165 !DIST = VEC(1)*VEC(1) + VEC(2)*VEC(2) + VEC(3)*VEC(3)
166 !DIST = TWO*SQRT(DIST)
167
168 !-- increment of static pressure due to initial gravity loading
169 ! stored in DP0 => EBCS_TAB%..%DP0
170 IF(inigrav_parts%IS_ALLOCATED)THEN
171 IF(inigrav_parts%TAGPART(ip) == 1)THEN
172 grav0 = inigrav_parts%GRAV0(ip) !0.0 if not related to inigrav option
173 ELSE
174 grav0 = zero
175 ENDIF
176 ELSE
177 grav0 = zero
178 ENDIF
179 IF(grav0 /= zero)THEN
180
181 ngx=inigrav_parts%NG(1,ip)
182 ngy=inigrav_parts%NG(2,ip)
183 ngz=inigrav_parts%NG(3,ip)
184
185 dist = vec(1)*ngx + vec(2)*ngy +vec(3)*ngz
186 dist = two * dist
187
188 dp0(is) = rho*abs(grav0)*dist
189 ELSE
190 !no increment to hydrostatic pressure
191 dp0(is)=zero
192 ENDIF
193
194
195 ENDDO
196
197
my_real
#define my_real
Definition cppsort.cpp:32
inigrav::inigrav_parts
type(t_inigrav_parts) inigrav_parts
Definition inigrav_mod.F:52