38
39
40
44
45
46
47#include "implicit_f.inc"
48
49
50
51
52
53 INTEGER IGOTO ,SIZEPN,NPTBI
54 INTEGER NVAR,NDC
56 . cons(*),dcons(*),x(
nvar),uparam(*),fun(
nvar),dfun(
nvar),
57 . xbia(nptbi)
58
59
60
61
62 INTEGER NPC(*), NFUNC, IFUNC(NFUNC)
64 EXTERNAL finter
65
66 INTENT(IN) :: npc,pld,igoto ,xbia,nptbi
67 INTENT(INOUT) :: fun,dfun,cons,dcons
68
69
70
71 INTEGER I,II,J,ID,NP1,NP2,J1,K,K1,NCON, ITER,NITER,
72 . PN1,PN2,ISYM,LENC,LENT
74 . embc,embt,y,fmins,fminf,dembc,embcp,flexp,
75 . flex1,flex2,lc0,lt0,dc0,dt0,hc0,ht0,kc,kt,kfc,kft,
76 . hc,ht,dcc,dc,dt,udc,udt,hdc,hdt,fc,ft,fpc,fpt,kf,
77 . deric,dtt,func,dcp,dtp, t1,t2,a1,a2,dflex,flex2p,flex1p,
78 . fminc,fmint,fminc2,fmins2,fminf2,fminf1,yp,embtp,dembt,yfac(8)
79
81 . ec(lenc),fcu(lenc),lc(lenc),yc(lenc),sigc(lenc),
82 . sigcp(lenc),
83 . xfib(lenc),yfib(lenc),
84 . xcfib(nptbi),ycfib(nptbi),xtfib(nptbi),ytfib(nptbi)
86 . et(lent),ftu(lent),lt(lent),yt(lent),sigt(lent),
87 . sigtp(lent)
88
89
90
91
92 IF (isym == 1)THEN
93 fminf = zero
94
95 fminf1 = zero
96 fminf2 = zero
97 fminc = zero
98 fmint = zero
99 fmins = zero
100 fmins2 = zero
101
102
103 embc = x(1)
104 embt = x(3)
105 flex1= x(2)
106 flex2= x(4)
107 dembc =
max(em03 * embc, em10)
108 dflex =
max(em03 * flex1, em10)
109 dembt =
max(em03 * embt, em10)
110
111 niter= 5
112 DO i=1,8
113 yfac(i)= uparam(8+i)
114 ENDDO
115 dc0 = one+embc
116 hc0 = sqrt(dc0*dc0 - one)
117
118 dt0 = one+embt
119 ht0 = sqrt(dt0*dt0 - one)
120
121
122
123 CALL fct_fiber_2(npc,pld ,ifunc(7),ifunc(8),yfac(7),yfac(8),xbia,nptbi,
124 . flex1,flex2,dc0,hc0,dt0,ht0,xcfib,ycfib,xtfib,ytfib )
125
126 pn1 = npc(ifunc(1))
127 pn2 = npc(ifunc(1)+1)
128 lenc = (pn2 - pn1) / 2
129
130 DO i = 1,lenc
131 ii = 2*(i-1)
132
133 ec(i) = pld(npc(ifunc(1))+ ii )
134 fcu(i)= pld(npc(ifunc(1))+ ii+1 )*yfac(1)
136 . yfac(1),flex1,flex2,embc ,sigc(i))
137 a1 = (fcu(i) - sigc(i))/
max(em20,fcu(i))
138 fminc = fminc + a1**2
139 ENDDO
140
141 pn1 = npc(ifunc(2))
142 pn2 = npc(ifunc(2)+1)
143 lent = half*(pn2 - pn1)
144 DO i = 1,lent
145 ii = 2*(i-1)
146 et(i) = pld(npc(ifunc(2))+ ii )
147 ftu(i)= pld(npc(ifunc(2))+ ii+1 )*yfac(2)
149 . yfac(2),flex1,flex2,embt,sigt(i))
150 a2 = (ftu(i) - sigt(i))/
max(em20,ftu(i))
151 fmint = fmint + a2**2
152 ENDDO
153
154 IF (igoto == 5 .or. igoto == 8) fun = (fminc + fmint )/two
155
156
157
158
159 IF (igoto == 3 .or. igoto == 8) THEN
160
161 embcp = x(1) + dembc
162 flex1 = x(2)
163 dc0 = one + embcp
164 hc0 = sqrt(dc0*dc0 - one)
165 DO i = 1,lenc
166 ii = 2*(i-1)
167 ec(i) = pld(npc(ifunc(1))+ ii )
168 fcu(i)= pld(npc(ifunc(1))+ ii +1)*yfac(1)
170 . yfac(1),flex1,flex2,embcp,sigc(i))
171
172 a1 = (fcu(i) - sigc(i))/
max(em20,fcu(i))
173 fmins = fmins + a1**2
174 ENDDO
175
176 dfun(1) = (fmins - fminc) / dembc
177
178
179
180
181 embc = x(1)
182 embt = x(3)
183 flexp = x(2) + dflex
184 flex2 = x(4)
185 dc0 = one + embc
186 hc0 = sqrt(dc0*dc0 - one)
187 DO i = 1,lenc
188 ii = 2*(i-1)
189 ec(i) = pld(npc(ifunc(1))+ ii )
190 fcu(i)= pld(npc(ifunc(1))+ ii +1)*yfac(1)
192 . yfac(1),flexp,flex2,embc,sigc(i))
193
194 a1 = (fcu(i) - sigc(i))/
max(em20,fcu(i))
195
196 fminf1 = fminf1 + a1**2
197 ENDDO
198
199 dt0 = one + embt
200 ht0 = sqrt(dt0*dt0 - one)
201 DO i = 1,lent
202 ii = 2*(i-1)
203 et(i) = pld(npc(ifunc(2))+ ii )
204 ftu(i)= pld(npc(ifunc(2))+ ii +1)*yfac(2)
206 . yfac(2),flexp,flex2,embt,sigt(i))
207
208 a1 = (ftu(i) - sigt(i))/
max(em20,ftu(i))
209 fminf2 = fminf2 + a1**2
210 ENDDO
211
212 dfun(2) = (fminf1 + fminf2 - fminc - fmint) / dflex
213
214
215
216 fminf1 = zero
217 fminf2 = zero
218 embc = x(1)
219 flex1= x(2)
220
221 embtp = x(3) + dembt
222 flex2 = x(4)
223 dt0 = one + embtp
224 ht0 = sqrt(dt0*dt0 - one)
225 DO i = 1,lent
226 ii = 2*(i-1)
227 et(i) = pld(npc(ifunc(2))+ ii )
228 ftu(i)= pld(npc(ifunc(2))+ ii +1)*yfac(2)
230 . yfac(2),flex1,flex2,embtp,sigt(i))
231
232 a1 = (ftu(i) - sigt(i))/
max(em20,ftu(i))
233 fmins2 = fmins2 + a1**2
234 ENDDO
235
236 dfun(3) = (fmins2 - fmint) / dembt
237
238
239
240 embc = x(1)
241 embt = x(3)
242 flex1= x(2)
243 flexp = x(4) + dflex
244
245 dc0 = one + embc
246 hc0 = sqrt(dc0*dc0 - one)
247 DO i = 1,lenc
248 ii = 2*(i-1)
249 ec(i) = pld(npc(ifunc(1))+ ii )
250 fcu(i)= pld(npc(ifunc(1))+ ii +1)*yfac(1)
252 . yfac(1),flex1,flexp,embc,sigc(i
253
254 a1 = (fcu(i) - sigc(i))/
max(em20,fcu(i
255
256 fminf1 = fminf1 + a1**2
257 ENDDO
258 dt0 = one + embt
259 ht0 = sqrt(dt0*dt0 - one)
260 DO i = 1,lent
261 ii = 2*(i-1)
262 et(i) = pld(npc(ifunc(2))+ ii )
263 ftu(i)= pld(npc(ifunc(2))+ ii +1)*yfac(2)
265 . yfac(2),flex1,flexp,embt,sigt(i))
266
267 a1 = (ftu(i) - sigt(i))/
max(em20,ftu(i))
268
269 fminf2 = fminf2 + a1**2
270 ENDDO
271 dfun(4) = (fminf1 + fminf2 - fminc - fmint) / dflex
272 ENDIF
273
274
275
276 ELSE
277 fminc = zero
278 fmins = zero
279 fminf = zero
280
281 embc = x(1)
282 embt = x(1)
283 flex1= x(2)
284 flex2= x(2)
285 dembc =
max(em03 * embc, em10)
286 dflex =
max(em03 * flex1, em10)
287 niter= 5
288 DO i=1,8
289 yfac(i)= uparam(8+i)
290 ENDDO
291 dc0 = one+embc
292 hc0 = sqrt(dc0*dc0 - one)
293
294 dt0 = one+embt
295 ht0 = sqrt(dt0*dt0 - one)
296
297 pn1 = npc(ifunc(7))
298 pn2 = npc(ifunc(7)+1)
299 sizepn = half*(pn2 - pn1)
300
301
302
303 CALL fct_fiber(npc,pld ,ifunc(7),sizepn,dc0,hc0,xfib,yfib)
304
305
306
307 DO i = 1,sizepn
308 ii = 2*(i-1)
309
310 ec(i) = pld(npc(ifunc(1))+ ii )
311 fcu(i)= pld(npc(ifunc(1))+ ii+1 )*yfac(1)
312
313 CALL calc_uniax(ec(i),xfib,yfib,sizepn,dc0,hc0,
314 . yfac(1),flex1,embc,sigc(i))
315 a1 = (fcu(i) - sigc(i))/
max(em20,fcu(i))
316 fminc = fminc + a1**2
317 ENDDO
318
319 IF (igoto == 5 .or. igoto == 8) fun = fminc
320
321
322
323
324 IF (igoto == 3 .or. igoto == 8) THEN
325 embcp = x(1) + dembc
326 flex1 = x(2)
327 dc0 = one + embcp
328 hc0 = sqrt(dc0*dc0 - one)
329 dt0 = dc0
330 ht0 = hc0
331 DO i = 1,sizepn
332 ii = 2*(i-1)
333 ec(i) = pld(npc(ifunc(1))+ ii )
334 fcu(i)= pld(npc(ifunc(1))+ ii +1)*yfac(1)
335 CALL calc_uniax(ec(i),xfib,yfib,sizepn,dc0,hc0,
336 . yfac(1),flex1,embcp,sigc(i))
337 a1 = (fcu(i) - sigc(i))/
max(em20,fcu
338 fmins = fmins + a1**2
339 ENDDO
340
341 dfun(1) = (fmins - fminc) / dembc
342
343
344 embc = x(1)
345 flexp = x(2) + dflex
346 dc0 = one + embc
347 hc0 = sqrt(dc0*dc0 - one)
348 dt0 = dc0
349 ht0 = hc0
350 DO i = 1,sizepn
351 ii = 2*(i-1)
352 ec(i) = pld
353 fcu(i)= pld(npc(ifunc(1))+ ii +1)*yfac(1)
354 CALL calc_uniax(ec(i),xfib,yfib,sizepn,dc0,hc0,
355 . yfac(1),flexp,embc,sigc(i))
356
357 a1 = (fcu(i) - sigc(i))/
max(em20,fcu(i))
358
359 fminf = fminf + a1**2
360 ENDDO
361
362 dfun(2) = (fminf - fminc) / dflex
363
364 ENDIF
365 ENDIF
366
367
368
369 RETURN
subroutine fct_fiber(npc, pld, idn, sizepn, dc0, hc0, xfib, yfib)
subroutine calc_uniax_2(ec, xfib, yfib, sizepn, dc0, hc0, yfac, flex1, flex2, embc, sigc)
subroutine calc_uniax(ec, xfib, yfib, sizepn, dc0, hc0, yfac, flex, stretch, sigc)
subroutine fct_fiber_2(npc, pld, idn1, idn2, yfac1, yfac2, xbia, nptbi, kfc, kft, dc0, hc0, dt0, ht0, xcfib, ycfib, xtfib, ytfib)
integer, parameter nchartitle
integer function nvar(text)