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