42
43
44
47
48
49
50#include "implicit_f.inc"
51#include "comlock.inc"
52
53
54
55#include "com06_c.inc"
56#include "com08_c.inc"
57#include "scr07_c.inc"
58#include "scr14_c.inc"
59#include "scr16_c.inc"
60#include "parit_c.inc"
61#include "scr18_c.inc"
62
63
64
65 INTEGER NELTST,ITYPTST, IFUNC,IFUN2, NOINT, NSN, NMN,HFLAG,ICOR
66 INTEGER MSR(*), NSV(*), NPC(*), ISKY(*)
67
69 . dt2t,ansmx,ansmx0,ff0,fmx,fmy,fmz,xmas,ffac,facx,fac2,stiff,
70 . peni
72 . e(*),es(*),em(*),tf(*),v(*),fsav(*),fskyi(lskyi,nfskyi),
73 . fcont(3,*)
74 TYPE(H3D_DATABASE) :: H3D_DATA
75
76
77
78 INTEGER NPOINT, IL, IG, IG3, IG2, IG1, IL3, IL2, IL1, NISKYL
79
81 . vsmax, vmmax, vmax, ft,fu,ff, xk, dtmi, fac, facdt, dx, finter
82 EXTERNAL finter
83
84 IF (tt == zero) THEN
85 ansmx0 = zero
86 ff0 = zero
87 IF (icor == 1) THEN
88 ansmx0 = -ansmx
89 IF (hflag == 1) ff0 = finter(ifun2,abs(peni)*facx,npc,tf,xk)
90 ENDIF
91 ENDIF
92
93 IF (ansmx > zero)THEN
94 vsmax =zero
95 vmmax =zero
96 ansmx0=zero
97
98 DO il=1,nsn
99 ig=nsv(il)
100 ig3=3*ig
101 ig2=ig3-1
102 ig1=ig2-1
103 vsmax =
max(vsmax,v(ig1)**2+v(ig2)**2+v(ig3)**2)
104 ENDDO
105
106 DO il=1,nmn
107 ig=msr(il)
108 ig3=3*ig
109 ig2=ig3-1
110 ig1=ig2-1
111 vmmax =
max(vmmax,v(ig1)**2+v(ig2)**2+v(ig3)**2)
112 ENDDO
113
114 vmax = sqrt(vsmax)+sqrt(vmmax)+ em15
115 ft = finter(ifunc,zero,npc,tf,xk)
117 dtmi =
max(em01*sqrt(xmas/xk),ansmx/vmax)
118
119 ELSEIF (ansmx == zero) THEN
120 ft = finter(ifunc,ansmx*facx,npc,tf,xk)
121 xk =
max(em15,xk*ffac)
122 ft = ft*ffac
123 dtmi = em01*sqrt(xmas/xk)
124
125 ELSE
126
127 ansmx = -ansmx
128 ft = finter(ifunc,ansmx*facx,npc,tf,xk)
129 xk =
max(em15,xk*ffac)
130 ft = ft*ffac
131 fu = zero
132 IF (hflag == 1) THEN
133 fu = finter(ifun2,ansmx*facx,npc,tf,xk)
134 fu = fu*fac2
135 ENDIF
136
137 IF (hflag > 0) THEN
138 dx = ansmx-ansmx0
139 IF (dx >= zero) THEN
140
141 ft =
min(ft, ff0 + stiff*dx)
142 ELSE
143
144 ft =
max(fu, ff0 + stiff*dx)
145 ENDIF
146 xk = ft - ff0 /
max(em15,dx)
147 ansmx0 = ansmx
148 ff0 = ft
149 ENDIF
150 fac = ft /
max(em15,sqrt(fmx**2+fmy**2+fmz**2))
151 facdt = fac*dt1
152
153 fsav(1)=fsav(1)+fmx*facdt
154 fsav(2)=fsav(2)+fmy*facdt
155 fsav(3)=fsav(3)+fmz*facdt
156 fsav(4)=fsav(4)-fmx*facdt
157 fsav(5)=fsav(5)-fmy*facdt
158 fsav(6)=fsav(6)-fmz*facdt
159
160 IF (iparit == 0) THEN
161 DO 190 il=1,nsn
162 ig=nsv(il)
163 ig3=3*ig
164 ig2=ig3-1
165 ig1=ig2-1
166 il3=3*il
167 il2=il3-1
168 il1=il2-1
169 e(ig1)=es(il1)*fac
170 e(ig2)=es(il2)*fac
171 e(ig3)=es(il3)*fac
172 190 CONTINUE
173
174 DO 200 il=1,nmn
175 ig=msr(il)
176 ig3=3*ig
177 ig2=ig3-1
178 ig1=ig2-1
179 il3=3*il
180 il2=il3-1
181 il1=il2-1
182 e(ig1)=em(il1)*fac
183 e(ig2)=em(il2)*fac
184 e(ig3)=em(il3)*fac
185 fsav(4)=fsav(4)-em(il1)*facdt
186 fsav(5)=fsav(5)-em(il2)*facdt
187 fsav(6)=fsav(6)-em(il3)*facdt
188 200 CONTINUE
189
190 ELSE
191
192#include "lockon.inc"
193 niskyl = nisky
194 nisky = nisky + nsn + nmn
195#include "lockoff.inc"
196 IF (kdtint == 0) THEN
197 DO 220 il=1,nsn
198 il3=3*il
199 il2=il3-1
200 il1=il2-1
201 niskyl = niskyl + 1
202 fskyi(niskyl,1)=es(il1)*fac
203 fskyi(niskyl,2)=es(il2)*fac
204 fskyi(niskyl,3)=es(il3)*fac
205 fskyi(niskyl,4)=zero
206 isky(niskyl) = nsv(il)
207 220 CONTINUE
208
209 DO 240 il=1,nmn
210 il3=3*il
211 il2=il3-1
212 il1=il2-1
213 niskyl = niskyl + 1
214 fskyi(niskyl,1)=em(il1)*fac
215 fskyi(niskyl,2)=em(il2)*fac
216 fskyi(niskyl,3)=em(il3)*fac
217 fskyi(niskyl,4)=zero
218 isky(niskyl) = msr(il)
219 fsav(4)=fsav(4)-em(il1)*facdt
220 fsav(5)=fsav(5)-em(il2)*facdt
221 fsav(6)=fsav(6)-em(il3)*facdt
222 240 CONTINUE
223 ELSE
224 DO il=1,nsn
225 il3=3*il
226 il2=il3-1
227 il1=il2-1
228 niskyl = niskyl + 1
229 fskyi(niskyl,1)=es(il1)*fac
230 fskyi(niskyl,2)=es(il2)*fac
231 fskyi(niskyl,3)=es(il3)*fac
232 fskyi(niskyl,4)=zero
233 fskyi(niskyl,5)=zero
234 isky(niskyl) = nsv(il)
235 ENDDO
236
237 DO il=1,nmn
238 il3=3*il
239 il2=il3-1
240 il1=il2-1
241 niskyl = niskyl + 1
242 fskyi(niskyl,1)=em(il1)*fac
243 fskyi(niskyl,2)=em(il2)*fac
244 fskyi(niskyl,3)=em(il3)*fac
245 fskyi(niskyl,4)=zero
246 fskyi(niskyl,5)=zero
247 isky(niskyl) = msr(il)
248 fsav(4)=fsav(4)-em(il1)*facdt
249 fsav(5)=fsav(5)-em(il2)*facdt
250 fsav(6)=fsav(6)-em(il3)*facdt
251 ENDDO
252 ENDIF
253 ENDIF
254
255 IF(anim_v(4)+outp_v(4)+h3d_data%N_VECT_CONT>0.AND.
256 . ((tt>=tanim .AND. tt<=tanim_stop).OR.tt>=toutp.OR.(tt>=h3d_data%TH3D.AND.tt<=h3d_data%TH3D_STOP).OR.
257 . (manim>=4.AND.manim<=15).OR. h3d_data%MH3D /= 0))THEN
258#include "lockon.inc"
259 DO il=1,nsn
260 il3=3*il
261 il2=il3-1
262 il1=il2-1
263 fcont(1,nsv(il)) =fcont(1,nsv(il)) + es(il1)*fac
264 fcont(2,nsv(il)) =fcont(2,nsv(il)) + es(il2)*fac
265 fcont(3,nsv(il)) =fcont(3,nsv(il)) + es(il3)*fac
266 ENDDO
267
268 DO il=1,nmn
269 il3=3*il
270 il2=il3-1
271 il1=il2-1
272 fcont(1,msr(il)) =fcont(1,msr(il)) + em(il1)*fac
273 fcont(2,msr(il)) =fcont(2,msr(il)) + em(il2)*fac
274 fcont(3,msr(il)) =fcont(3,msr(il)) + em(il3)*fac
275 ENDDO
276#include "lockoff.inc"
277 ENDIF
278
279 xk =
max(xk,ft /
max(em15,ansmx))
280 dtmi = em01*sqrt(xmas/
max(xk,em20))
281 ENDIF
282
283 IF(dtmi<dt2t)THEN
284 dt2t = dtmi
285 neltst = noint
286 ityptst = 10
287 ENDIF
288
289 RETURN