40
41
42
43 USE python_funct_mod
44 use python_call_funct_cload_mod
45 use nodal_arrays_mod
46 USE sensor_mod
47
48
49
50#include "implicit_f.inc"
51
52
53
54#include "com01_c.inc"
55#include "com04_c.inc"
56#include "com08_c.inc"
57#include "param_c.inc"
58#include "sms_c.inc"
59#include "units_c.inc"
60
61
62
63
64 TYPE(python_), INTENT(IN) :: PYTHON
65 TYPE(nodal_arrays_) :: NODES
66 INTEGER ,INTENT(IN) :: NSENSOR
67 INTEGER NPC(*),CPTREAC,NODREAC(*)
68 INTEGER IBFV(NIFV,*), NODNX_SMS(*)
69 my_real :: tf(*), vel(lfxvelr,*), fthreac(6,*)
70
71
72
73 INTEGER N, , J, K, L, I1, N1, N2, ISENS,
74 . ILENC, IPOSC, IADC, ICOOR, ISMOOTH
75 my_real ax, axi, vv, a0, aa, fac, facx, startt, stopt, ts,
76 . dydx, dw, dw2, dd,
77 . yc, tsc, dydxc, rw_sms,
78 . skew1, skew2, skew3,
79 . xi, yi, zi, xf, yf, zf, xa, ya, za,
80 . aold0(3),
81 . vx1, vy1, vz1, vx2, vy2, vz2
82 my_real,
DIMENSION(:),
ALLOCATABLE :: mass, vx, vy, vz
83 TYPE (SENSOR_STR_) ,DIMENSION(NSENSOR) ,INTENT(IN) :: SENSOR_TAB
84 DOUBLE PRECISION, INTENT(INOUT) :: WFEXT
85
86
87
89 EXTERNAL finter2,finter2_smooth
90
91 IF(n2d==1)THEN
92 ax=one
93 ELSE
94 ax=zero
95 ENDIF
96
97
98
99 ALLOCATE(mass(numnod),vx(numnod),vy(numnod),vz(numnod))
100 mass(1:numnod)=zero
101 vx(1:numnod)=zero
102 vy(1:numnod)=zero
103 vz(1:numnod)=zero
104 DO n=1,nfxvel
105 IF (ibfv(13,n) /= 2) cycle
106 stopt = vel(3,n)
107 IF(tt>stopt) THEN
108 n1 = iabs(ibfv(1,n))
109 n2 = ibfv(14,n)
110
111 vx1= nodes%V(1,n1)+nodes%A(1,n1)*dt12
112 vy1= nodes%V(2,n1)+nodes%A(2,n1)*dt12
113 vz1= nodes%V(3,n1)+nodes%A(3,n1)*dt12
114 vx2= nodes%V(1,n2)+nodes%A(1,n2)*dt12
115 vy2= nodes%V(2,n2)+nodes%A(2,n2)*dt12
116 vz2= nodes%V(3,n2)+nodes%A(3,n2)*dt12
117
118 vx(n2) = vx(n2) + nodes%MS(n1)*vx1+ nodes%MS(n2)*vx2/ibfv(16,n)
119 vy(n2) = vy(n2) + nodes%MS(n1)*vy1+ nodes%MS(n2)*vy2/ibfv(16,n)
120 vz(n2) = vz(n2) + nodes%MS(n1)*vz1+ nodes%MS(n2)*vz2/ibfv(16,n)
121 mass(n2) = mass(n2) + nodes%MS(n1) + nodes%MS(n2)/ibfv(16,n)
122 ENDIF
123 ENDDO
124 DO n=1,numnod
125 IF(mass(n) == zero) cycle
126 vx(n)=vx(n)/mass(n)
127 vy(n)=vy(n)/mass(n)
128 vz(n)=vz(n)/mass(n)
129 ENDDO
130
131 DO n=1,nfxvel
132 IF (ibfv(13,n) /= 2) cycle
133 stopt = vel(3,n)
134 IF(tt>stopt) THEN
135 n1 = iabs(ibfv(1,n))
136 n2 = ibfv(14,n)
137 nodes%A(1,n1)=(vx(n2)-nodes%V(1,n1))/dt12
138 nodes%A(2,n1)=(vy(n2)-nodes%V(2,n1))/dt12
139 nodes%A(3,n1)=(vz(n2)-nodes%V(3,n1))/dt12
140 nodes%A(1,n2)=(vx(n2)-nodes%V(1,n2))/dt12
141 nodes%A(2,n2)=(vy(n2)-nodes%V(2,n2))/dt12
142 nodes%A(3,n2)=(vz(n2)-nodes%V(3,n2))/dt12
143 ENDIF
144 ENDDO
145 DEALLOCATE(mass, vx, vy, vz)
146
147
148
149 DO n=1,nfxvel
150 IF (ibfv(13,n) == 0) cycle
151 startt = vel(2,n)
152 stopt = vel(3,n)
153 facx = vel(5,n)
154 i=iabs(ibfv(1,n))
155 IF(tt<startt) cycle
156 IF(tt>stopt) cycle
157 IF(nsensor>0) THEN
158 isens = ibfv(4,n)
159 IF(isens==0)THEN
160 ts=tt
161 ELSE
162 ts = tt - sensor_tab(isens)%TSTART
163 IF(ts<0.0) cycle
164 ENDIF
165 ELSE
166 ts=tt
167 ENDIF
168 tsc = (ts+dt2)*facx
169 IF(idtmins==0.AND.idtmins_int==0)THEN
170 rw_sms = one
171 ELSE
172 IF(nodnx_sms(i)==0)THEN
173 rw_sms = one
174 ELSE
175 rw_sms = zero
176 ENDIF
177 ENDIF
178
179 IF(cptreac/=0) THEN
180 aold0(1)=nodes%A(1,i)
181 aold0(2)=nodes%A(2,i)
182 aold0(3)=nodes%A(3,i)
183 ENDIF
184
185 l = ibfv(3,n)
186 IF(ncycle==0)THEN
187 iposc = 0
188 ELSE
189 iposc = ibfv(5,n)
190 ENDIF
191 iadc = npc(l) / 2 + 1
192 ilenc = npc(l+1) / 2 - iadc - iposc
193
194 ismooth = 0
195 IF (l > 0) ismooth = npc(2*nfunct+l+1)
196 IF (ismooth == 0) THEN
197 yc = finter2(tf,iadc,iposc,ilenc,tsc,dydxc)
198 ELSE IF(ismooth > 0) THEN
199 yc = finter2_smooth(tf,iadc,iposc,ilenc,tsc,dydxc)
200 ELSE
201 ismooth = -ismooth
202
203 CALL python_call_funct_cload(python, ismooth,tsc, yc,i,nodes)
204 ENDIF
205 ibfv(5,n) = iposc
206
207 IF(ibfv(13,n) == 1) THEN
208
209 fac = vel(1,n)
210 skew1= vel(7,n)
211 skew2= vel(8,n)
212 skew3= vel(9,n)
213 ELSEIF(ibfv(13,n) == 2) THEN
214
215 xa = nodes%X(1,i)
216 ya = nodes%X(2,i)
217 za = nodes%X(3,i)
218 i1 = ibfv(14,n)
219 xf = nodes%X(1,i1)
220 yf = nodes%X(2,i1)
221 zf = nodes%X(3,i1)
222 fac= sqrt((xf-xa)**2+(yf-ya)**2+(zf-za)**2)
223 IF(fac < vel(7,n)) THEN
224 vel(3,n) = tt
225 yc = zero
226 WRITE(iout,'(A,I10,1X,I10,A,1PE12.5)')
227 . ' RIGID LINK ON NODES',nodes%ITAB(ibfv(1,n)),ibfv(14,n),
228 . ' ACTIVATED AT TIME',tt
229 ENDIF
230 skew1= (xf-xa)/fac
231 skew2= (yf-ya)/fac
232 skew3= (zf-za)/fac
233 fac = vel(1,n)
234 ENDIF
235 vel(6,n) = yc
236
237
238
239 dw = vel(4,n)
240 wfext = wfext + rw_sms*dw
241 vel(4,n)= (one-rw_sms)*vel(4,n)
242
243 yc = yc * fac
244 axi=one-ax+ax*nodes%X(2,i)
245
246 dd = skew1*nodes%D(1,i) + skew2*nodes%D(2,i) + skew3*nodes%D(3,i)
247 vv = skew1*nodes%V(1,i) + skew2*nodes%V(2,i) + skew3*nodes%V(3,i)
248 a0 = skew1*nodes%A(1,i) + skew2*nodes%A(2,i) + skew3
249
250 IF(ibfv(13,n) == 1) yc =(yc-dd)/dt2
251 yc =(yc-vv)/dt12
252 aa = yc - a0
253
254 nodes%A(1,i)=skew1*yc
255 nodes%A(2,i)=skew2*yc
256 nodes%A(3,i)=skew3*yc
257
258 dw = fourth*nodes%MS(i)*(yc*dt12+two*vv)*aa*axi*nodes%WEIGHT(i)
259 wfext = wfext + dt1*dw
260
261
262
263 vel(4,n) = vel(4,n)+dt2*dw
264
265 IF (cptreac/=0) THEN
266 i=iabs(ibfv(1,n))
267 IF (nodreac(i)/=0) THEN
268 fthreac(1,nodreac(i)) = fthreac(1,nodreac(i)) + (nodes%A(1,i) -
269 & aold0(1))*nodes%MS(i)*dt12
270 fthreac(2,nodreac(i)) = fthreac(2,nodreac(i)) + (nodes%A(2,i) -
271 & aold0(2))*nodes%MS(i)*dt12
272 fthreac(3,nodreac(i)) = fthreac(3,nodreac(i)) + (nodes%A(3,i) -
273 & aold0(3))*nodes%MS(i)*dt12
274 ENDIF
275 ENDIF
276
277 ENDDO
278
279 RETURN