CalculiX  2.13
A Free Software Three-Dimensional Structural Finite Element Program
resultsmech_se.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine resultsmech_se (co, kon, ipkon, lakon, ne, v, stx, elcon, nelcon, rhcon, nrhcon, alcon, nalcon, alzero, ielmat, ielorien, norien, orab, ntmat_, t0, t1, ithermal, prestr, iprestr, eme, iperturb, fn, iout, vold, nmethod, veold, dtime, time, ttime, plicon, nplicon, plkcon, nplkcon, xstateini, xstiff, xstate, npmat_, matname, mi, ielas, icmd, ncmat_, nstate_, stiini, vini, ener, eei, enerini, istep, iinc, springarea, reltime, calcul_fn, calcul_cauchy, nener, ikin, ne0, thicke, emeini, pslavsurf, pmastsurf, mortar, clearini, nea, neb, ielprop, prop, dfn, idesvar, nodedesi, fn0, sti, icoordinate, dxstiff, ialdesi, xdesi)
 

Function/Subroutine Documentation

◆ resultsmech_se()

subroutine resultsmech_se ( real*8, dimension(3,*)  co,
integer, dimension(*)  kon,
integer, dimension(*)  ipkon,
character*8, dimension(*)  lakon,
integer  ne,
real*8, dimension(0:mi(2),*)  v,
real*8, dimension(6,mi(1),*)  stx,
real*8, dimension(0:ncmat_,ntmat_,*)  elcon,
integer, dimension(2,*)  nelcon,
real*8, dimension(0:1,ntmat_,*)  rhcon,
integer, dimension(*)  nrhcon,
real*8, dimension(0:6,ntmat_,*)  alcon,
integer, dimension(2,*)  nalcon,
real*8, dimension(*)  alzero,
integer, dimension(mi(3),*)  ielmat,
integer, dimension(mi(3),*)  ielorien,
integer  norien,
real*8, dimension(7,*)  orab,
integer  ntmat_,
real*8, dimension(*)  t0,
real*8, dimension(*)  t1,
integer, dimension(2)  ithermal,
real*8, dimension(6,mi(1),*)  prestr,
integer  iprestr,
real*8, dimension(6,mi(1),*)  eme,
integer, dimension(*)  iperturb,
real*8, dimension(0:mi(2),*)  fn,
integer  iout,
real*8, dimension(0:mi(2),*)  vold,
integer  nmethod,
real*8, dimension(0:mi(2),*)  veold,
real*8  dtime,
real*8  time,
real*8  ttime,
real*8, dimension(0:2*npmat_,ntmat_,*)  plicon,
integer, dimension(0:ntmat_,*)  nplicon,
real*8, dimension(0:2*npmat_,ntmat_,*)  plkcon,
integer, dimension(0:ntmat_,*)  nplkcon,
real*8, dimension(nstate_,mi(1),*)  xstateini,
real*8, dimension(27,mi(1),*)  xstiff,
real*8, dimension(nstate_,mi(1),*)  xstate,
integer  npmat_,
character*80, dimension(*)  matname,
integer, dimension(*)  mi,
integer  ielas,
integer  icmd,
integer  ncmat_,
integer  nstate_,
real*8, dimension(6,mi(1),*)  stiini,
real*8, dimension(0:mi(2),*)  vini,
real*8, dimension(mi(1),*)  ener,
real*8, dimension(6,mi(1),*)  eei,
real*8, dimension(mi(1),*)  enerini,
integer  istep,
integer  iinc,
real*8, dimension(2,*)  springarea,
real*8  reltime,
integer  calcul_fn,
integer  calcul_cauchy,
integer  nener,
integer  ikin,
integer  ne0,
real*8, dimension(mi(3),*)  thicke,
real*8, dimension(6,mi(1),*)  emeini,
real*8, dimension(3,*)  pslavsurf,
real*8, dimension(6,*)  pmastsurf,
integer  mortar,
real*8, dimension(3,9,*)  clearini,
integer  nea,
integer  neb,
integer, dimension(*)  ielprop,
real*8, dimension(*)  prop,
real*8, dimension(0:mi(2),*)  dfn,
integer  idesvar,
integer, dimension(*)  nodedesi,
real*8, dimension(0:mi(2),*)  fn0,
real*8, dimension(6,mi(1),*)  sti,
integer  icoordinate,
real*8, dimension(27,mi(1),ne,*)  dxstiff,
integer, dimension(*)  ialdesi,
real*8, dimension(3,*)  xdesi 
)
31 !
32 ! calculates stresses and the material tangent at the integration
33 ! points and the internal forces at the nodes
34 !
35  implicit none
36 !
37  integer cauchy
38 !
39  character*8 lakon(*),lakonl
40  character*80 amat,matname(*)
41 !
42  integer kon(*),konl(26),nea,neb,mi(*),mint2d,nopes,
43  & nelcon(2,*),nrhcon(*),nalcon(2,*),ielmat(mi(3),*),
44  & ielorien(mi(3),*),ntmat_,ipkon(*),ne0,iflag,null,
45  & istep,iinc,mt,ne,mattyp,ithermal(2),iprestr,i,j,k,m1,m2,jj,
46  & i1,m3,m4,kk,nener,indexe,nope,norien,iperturb(*),iout,
47  & icmd,ihyper,nmethod,kode,imat,mint3d,iorien,ielas,
48  & istiff,ncmat_,nstate_,ikin,ilayer,nlayer,ki,kl,ielprop(*),
49  & nplicon(0:ntmat_,*),nplkcon(0:ntmat_,*),npmat_,calcul_fn,
50  & calcul_cauchy,nopered,mortar,jfaces,iloc,igauss,
51  & idesvar,node,nodedesi(*),kscale,idir,
52  & iactive,icoordinate,ialdesi(*),ii
53 !
54  real*8 co(3,*),v(0:mi(2),*),shp(4,26),stiini(6,mi(1),*),
55  & stx(6,mi(1),*),xl(3,26),vl(0:mi(2),26),stre(6),prop(*),
56  & elcon(0:ncmat_,ntmat_,*),rhcon(0:1,ntmat_,*),xs2(3,7),
57  & alcon(0:6,ntmat_,*),vini(0:mi(2),*),thickness,
58  & alzero(*),orab(7,*),elas(21),rho,fn(0:mi(2),*),
59  & fnl(3,10),skl(3,3),beta(6),xl2(3,8),qa(4),
60  & vkl(0:3,3),t0(*),t1(*),prestr(6,mi(1),*),eme(6,mi(1),*),
61  & ckl(3,3),vold(0:mi(2),*),eloc(9),veold(0:mi(2),*),
62  & springarea(2,*),elconloc(21),eth(6),xkl(3,3),voldl(0:mi(2),26),
63  & xikl(3,3),ener(mi(1),*),emec(6),eei(6,mi(1),*),enerini(mi(1),*),
64  & emec0(6),veoldl(0:mi(2),26),xsj2(3),shp2(7,8),
65  & e,un,al,um,am1,xi,et,ze,tt,exx,eyy,ezz,exy,exz,eyz,
66  & xsj,vj,t0l,t1l,dtime,weight,pgauss(3),vij,time,ttime,
67  & plicon(0:2*npmat_,ntmat_,*),plkcon(0:2*npmat_,ntmat_,*),
68  & xstiff(27,mi(1),*),xstate(nstate_,mi(1),*),plconloc(802),
69  & vokl(3,3),xstateini(nstate_,mi(1),*),vikl(3,3),
70  & gs(8,4),a,reltime,tlayer(4),dlayer(4),xlayer(mi(3),4),
71  & thicke(mi(3),*),emeini(6,mi(1),*),clearini(3,9,*),
72  & pslavsurf(3,*),pmastsurf(6,*),sti(6,mi(1),*),
73  & fn0(0:mi(2),*),dfn(0:mi(2),*),hglf(3,4),ahr,
74  & dxstiff(27,mi(1),ne,*),xdesi(3,*)
75 !
76  include "gauss.f"
77 !
78  iflag=3
79  null=0
80  qa(3)=-1.d0
81  qa(4)=0.d0
82 !
83  mt=mi(2)+1
84 !
85 ! -------------------------------------------------------------
86 ! Initialisation of the loop for the variation of
87 ! the internal forces
88 ! -------------------------------------------------------------
89 !
90 !
91 ! Loop over all elements in thread
92  do ii=nea,neb
93  if(idesvar.gt.0) then
94  i=ialdesi(ii)
95  else
96  i=ii
97  endif
98 !
99  lakonl=lakon(i)
100 !
101  if(ipkon(i).lt.0) cycle
102 !
103 ! no 3D-fluid elements
104 !
105  if(lakonl(1:1).eq.'F') cycle
106  if(lakonl(1:7).eq.'DCOUP3D') cycle
107 !
108  if(lakonl(7:8).ne.'LC') then
109 !
110  imat=ielmat(1,i)
111  amat=matname(imat)
112  if(norien.gt.0) then
113  iorien=ielorien(1,i)
114  else
115  iorien=0
116  endif
117 !
118  if(nelcon(1,imat).lt.0) then
119  ihyper=1
120  else
121  ihyper=0
122  endif
123  elseif(lakonl(4:5).eq.'20') then
124 !
125 ! composite materials
126 !
127  mint2d=4
128  nopes=8
129 ! determining the number of layers
130 !
131  nlayer=0
132  do k=1,mi(3)
133  if(ielmat(k,i).ne.0) then
134  nlayer=nlayer+1
135  endif
136  enddo
137 !
138 ! determining the layer thickness and global thickness
139 ! at the shell integration points
140 !
141  iflag=1
142  indexe=ipkon(i)
143  do kk=1,mint2d
144  xi=gauss3d2(1,kk)
145  et=gauss3d2(2,kk)
146  call shape8q(xi,et,xl2,xsj2,xs2,shp2,iflag)
147  tlayer(kk)=0.d0
148  do k=1,nlayer
149  thickness=0.d0
150  do j=1,nopes
151  thickness=thickness+thicke(k,indexe+j)*shp2(4,j)
152  enddo
153  tlayer(kk)=tlayer(kk)+thickness
154  xlayer(k,kk)=thickness
155  enddo
156  enddo
157  iflag=3
158 !
159  ilayer=0
160  do k=1,4
161  dlayer(k)=0.d0
162  enddo
163  elseif(lakonl(4:5).eq.'15') then
164 !
165 ! composite materials
166 !
167  mint2d=3
168  nopes=6
169 ! determining the number of layers
170 !
171  nlayer=0
172  do k=1,mi(3)
173  if(ielmat(k,i).ne.0) then
174  nlayer=nlayer+1
175  endif
176  enddo
177 !
178 ! determining the layer thickness and global thickness
179 ! at the shell integration points
180 !
181  iflag=1
182  indexe=ipkon(i)
183  do kk=1,mint2d
184  xi=gauss3d10(1,kk)
185  et=gauss3d10(2,kk)
186  call shape6tri(xi,et,xl2,xsj2,xs2,shp2,iflag)
187  tlayer(kk)=0.d0
188  do k=1,nlayer
189  thickness=0.d0
190  do j=1,nopes
191  thickness=thickness+thicke(k,indexe+j)*shp2(4,j)
192  enddo
193  tlayer(kk)=tlayer(kk)+thickness
194  xlayer(k,kk)=thickness
195  enddo
196  enddo
197  iflag=3
198 !
199  ilayer=0
200  do k=1,3
201  dlayer(k)=0.d0
202  enddo
203 !
204  endif
205 !
206  indexe=ipkon(i)
207 c Bernhardi start
208  if(lakonl(1:5).eq.'C3D8I') then
209  nope=11
210  elseif(lakonl(4:5).eq.'20') then
211 c Bernhardi end
212  nope=20
213  elseif(lakonl(4:4).eq.'8') then
214  nope=8
215  elseif(lakonl(4:5).eq.'10') then
216  nope=10
217  elseif(lakonl(4:4).eq.'4') then
218  nope=4
219  elseif(lakonl(4:5).eq.'15') then
220  nope=15
221  elseif(lakonl(4:4).eq.'6') then
222  nope=6
223  elseif((lakonl(1:1).eq.'E').and.(lakonl(7:7).ne.'F')) then
224 !
225 ! spring elements, contact spring elements and
226 ! dashpot elements
227 !
228  if(lakonl(7:7).eq.'C') then
229 !
230 ! contact spring elements
231 !
232  if(mortar.eq.1) then
233 !
234 ! face-to-face penalty
235 !
236  nope=kon(ipkon(i))
237  elseif(mortar.eq.0) then
238 !
239 ! node-to-face penalty
240 !
241  nope=ichar(lakonl(8:8))-47
242  konl(nope+1)=kon(indexe+nope+1)
243  endif
244  else
245 !
246 ! genuine spring elements and dashpot elements
247 !
248  nope=ichar(lakonl(8:8))-47
249  endif
250  else
251  cycle
252  endif
253 !
254 ! computation of the coordinates of the local nodes w/ variation
255 ! if the designnode belongs to the considered element
256 !
257  if(idesvar.gt.0) then
258 !
259  if(icoordinate.eq.1) then
260 !
261 ! the coordinates of the nodes are the design variables
262 !
263  do j=1,nope
264  node=kon(indexe+j)
265  if(node.eq.nodedesi(idesvar)) then
266  iactive=j
267  exit
268  endif
269  enddo
270  endif
271 !
272  endif
273 !
274  if(lakonl(4:5).eq.'8R') then
275  mint3d=1
276  elseif(lakonl(4:7).eq.'20RB') then
277  if((lakonl(8:8).eq.'R').or.(lakonl(8:8).eq.'C')) then
278  mint3d=50
279  else
280  call beamintscheme(lakonl,mint3d,ielprop(i),prop,
281  & null,xi,et,ze,weight)
282  endif
283  elseif((lakonl(4:4).eq.'8').or.
284  & (lakonl(4:6).eq.'20R')) then
285  if(lakonl(7:8).eq.'LC') then
286  mint3d=8*nlayer
287  else
288  mint3d=8
289  endif
290  elseif(lakonl(4:4).eq.'2') then
291  mint3d=27
292  elseif(lakonl(4:5).eq.'10') then
293  mint3d=4
294  elseif(lakonl(4:4).eq.'4') then
295  mint3d=1
296  elseif(lakonl(4:5).eq.'15') then
297  if(lakonl(7:8).eq.'LC') then
298  mint3d=6*nlayer
299  else
300  mint3d=9
301  endif
302  elseif(lakonl(4:4).eq.'6') then
303  mint3d=2
304  elseif(lakonl(1:1).eq.'E') then
305  mint3d=0
306  endif
307 
308  do j=1,nope
309  konl(j)=kon(indexe+j)
310  do k=1,3
311  xl(k,j)=co(k,konl(j))
312  vl(k,j)=v(k,konl(j))
313  voldl(k,j)=vold(k,konl(j))
314  enddo
315  enddo
316 !
317 ! computation of the coordinates of the local nodes w/ variation
318 ! if the designnode belongs to the considered element
319 !
320  if((idesvar.gt.0).and.(icoordinate.eq.1)) then
321  do j=1,3
322  xl(j,iactive)=xl(j,iactive)+xdesi(j,idesvar)
323  enddo
324  endif
325 !
326 ! calculating the forces for the contact elements
327 !
328  if(mint3d.eq.0) then
329 !
330 ! "normal" spring and dashpot elements
331 !
332  kode=nelcon(1,imat)
333  if(lakonl(7:7).eq.'A') then
334  t0l=0.d0
335  t1l=0.d0
336  if(ithermal(1).eq.1) then
337  t0l=(t0(konl(1))+t0(konl(2)))/2.d0
338  t1l=(t1(konl(1))+t1(konl(2)))/2.d0
339  elseif(ithermal(1).ge.2) then
340  t0l=(t0(konl(1))+t0(konl(2)))/2.d0
341  t1l=(vold(0,konl(1))+vold(0,konl(2)))/2.d0
342  endif
343  endif
344 !
345 ! spring elements (including contact springs)
346 !
347  if(lakonl(2:2).eq.'S') then
348  if((lakonl(7:7).eq.'A').or.(lakonl(7:7).eq.'1').or.
349  & (lakonl(7:7).eq.'2').or.((mortar.eq.0).and.
350  & ((nmethod.ne.1).or.(iperturb(1).ge.2).or.(iout.ne.-1))))
351  & then
352  call springforc_n2f(xl,konl,vl,imat,elcon,nelcon,elas,
353  & fnl,ncmat_,ntmat_,nope,lakonl,t1l,kode,elconloc,
354  & plicon,nplicon,npmat_,ener(1,i),nener,
355  & stx(1,1,i),mi,springarea(1,konl(nope+1)),nmethod,
356  & ne0,nstate_,xstateini,xstate,reltime,
357  & ielas,ener(1,i+ne),ielorien,orab,norien,i)
358  elseif((mortar.eq.1).and.
359  & ((nmethod.ne.1).or.(iperturb(1).ge.2).or.(iout.ne.-1)))
360  & then
361  iloc=kon(indexe+nope+1)
362  jfaces=kon(indexe+nope+2)
363  igauss=kon(indexe+nope+1)
364  call springforc_f2f(xl,vl,imat,elcon,nelcon,elas,
365  & fnl,ncmat_,ntmat_,nope,lakonl,t1l,kode,elconloc,
366  & plicon,nplicon,npmat_,ener(1,i),nener,
367  & stx(1,1,i),mi,springarea(1,iloc),nmethod,
368  & ne0,nstate_,xstateini,xstate,reltime,
369  & ielas,iloc,jfaces,igauss,pslavsurf,pmastsurf,
370  & clearini,ener(1,i+ne),kscale,konl,iout,i)
371  endif
372 !
373 ! next lines are not executed in linstatic.c before the
374 ! setup of the stiffness matrix (i.e. nmethod=1 and
375 ! iperturb(1)<1 and iout=-1).
376 !
377  if((lakonl(7:7).eq.'A').or.
378  & ((nmethod.ne.1).or.(iperturb(1).ge.2).or.(iout.ne.-1)))
379  & then
380  if(idesvar.eq.0) then
381  do j=1,nope
382  do k=1,3
383  fn0(k,indexe+j)=fn0(k,indexe+j)+fnl(k,j)
384  enddo
385  enddo
386  else
387  do j=1,nope
388  do k=1,3
389  dfn(k,konl(j))=dfn(k,konl(j))+fnl(k,j)
390  enddo
391  enddo
392  endif
393  endif
394 !
395 ! dashpot elements (including contact dashpots)
396 !
397  elseif((nmethod.eq.4).or.(nmethod.eq.5).or.
398  & ((abs(nmethod).eq.1).and.(iperturb(1).ge.2))) then
399  endif
400  elseif(ikin.eq.1) then
401  do j=1,nope
402  do k=1,3
403  veoldl(k,j)=veold(k,konl(j))
404  enddo
405  enddo
406  endif
407 !
408  do jj=1,mint3d
409  if(lakonl(4:5).eq.'8R') then
410  xi=gauss3d1(1,jj)
411  et=gauss3d1(2,jj)
412  ze=gauss3d1(3,jj)
413  weight=weight3d1(jj)
414  elseif(lakonl(4:7).eq.'20RB') then
415  if((lakonl(8:8).eq.'R').or.(lakonl(8:8).eq.'C')) then
416  xi=gauss3d13(1,jj)
417  et=gauss3d13(2,jj)
418  ze=gauss3d13(3,jj)
419  weight=weight3d13(jj)
420  else
421  call beamintscheme(lakonl,mint3d,ielprop(i),prop,
422  & jj,xi,et,ze,weight)
423  endif
424  elseif((lakonl(4:4).eq.'8').or.
425  & (lakonl(4:6).eq.'20R'))
426  & then
427  if(lakonl(7:8).ne.'LC') then
428  xi=gauss3d2(1,jj)
429  et=gauss3d2(2,jj)
430  ze=gauss3d2(3,jj)
431  weight=weight3d2(jj)
432  else
433  kl=mod(jj,8)
434  if(kl.eq.0) kl=8
435 !
436  xi=gauss3d2(1,kl)
437  et=gauss3d2(2,kl)
438  ze=gauss3d2(3,kl)
439  weight=weight3d2(kl)
440 !
441  ki=mod(jj,4)
442  if(ki.eq.0) ki=4
443 !
444  if(kl.eq.1) then
445  ilayer=ilayer+1
446  if(ilayer.gt.1) then
447  do k=1,4
448  dlayer(k)=dlayer(k)+xlayer(ilayer-1,k)
449  enddo
450  endif
451  endif
452  ze=2.d0*(dlayer(ki)+(ze+1.d0)/2.d0*xlayer(ilayer,ki))/
453  & tlayer(ki)-1.d0
454  weight=weight*xlayer(ilayer,ki)/tlayer(ki)
455 !
456 ! material and orientation
457 !
458  imat=ielmat(ilayer,i)
459  amat=matname(imat)
460  if(norien.gt.0) then
461  iorien=ielorien(ilayer,i)
462  else
463  iorien=0
464  endif
465 !
466  if(nelcon(1,imat).lt.0) then
467  ihyper=1
468  else
469  ihyper=0
470  endif
471  endif
472  elseif(lakonl(4:4).eq.'2') then
473  xi=gauss3d3(1,jj)
474  et=gauss3d3(2,jj)
475  ze=gauss3d3(3,jj)
476  weight=weight3d3(jj)
477  elseif(lakonl(4:5).eq.'10') then
478  xi=gauss3d5(1,jj)
479  et=gauss3d5(2,jj)
480  ze=gauss3d5(3,jj)
481  weight=weight3d5(jj)
482  elseif(lakonl(4:4).eq.'4') then
483  xi=gauss3d4(1,jj)
484  et=gauss3d4(2,jj)
485  ze=gauss3d4(3,jj)
486  weight=weight3d4(jj)
487  elseif(lakonl(4:5).eq.'15') then
488  if(lakonl(7:8).ne.'LC') then
489  xi=gauss3d8(1,jj)
490  et=gauss3d8(2,jj)
491  ze=gauss3d8(3,jj)
492  weight=weight3d8(jj)
493  else
494  kl=mod(jj,6)
495  if(kl.eq.0) kl=6
496 !
497  xi=gauss3d10(1,kl)
498  et=gauss3d10(2,kl)
499  ze=gauss3d10(3,kl)
500  weight=weight3d10(kl)
501 !
502  ki=mod(jj,3)
503  if(ki.eq.0) ki=3
504 !
505  if(kl.eq.1) then
506  ilayer=ilayer+1
507  if(ilayer.gt.1) then
508  do k=1,3
509  dlayer(k)=dlayer(k)+xlayer(ilayer-1,k)
510  enddo
511  endif
512  endif
513  ze=2.d0*(dlayer(ki)+(ze+1.d0)/2.d0*xlayer(ilayer,ki))/
514  & tlayer(ki)-1.d0
515  weight=weight*xlayer(ilayer,ki)/tlayer(ki)
516 !
517 ! material and orientation
518 !
519  imat=ielmat(ilayer,i)
520  amat=matname(imat)
521  if(norien.gt.0) then
522  iorien=ielorien(ilayer,i)
523  else
524  iorien=0
525  endif
526 !
527  if(nelcon(1,imat).lt.0) then
528  ihyper=1
529  else
530  ihyper=0
531  endif
532  endif
533  elseif(lakonl(4:4).eq.'6') then
534  xi=gauss3d7(1,jj)
535  et=gauss3d7(2,jj)
536  ze=gauss3d7(3,jj)
537  weight=weight3d7(jj)
538  endif
539 !
540 c Bernhardi start
541  if(lakonl(1:5).eq.'C3D8R') then
542  call shape8hr(xl,xsj,shp,gs,a)
543  elseif(lakonl(1:5).eq.'C3D8I') then
544  call shape8hu(xi,et,ze,xl,xsj,shp,iflag)
545  elseif(nope.eq.20) then
546 c Bernhardi end
547  if(lakonl(7:7).eq.'A') then
548  call shape20h_ax(xi,et,ze,xl,xsj,shp,iflag)
549  elseif((lakonl(7:7).eq.'E').or.
550  & (lakonl(7:7).eq.'S')) then
551  call shape20h_pl(xi,et,ze,xl,xsj,shp,iflag)
552  else
553  call shape20h(xi,et,ze,xl,xsj,shp,iflag)
554  endif
555  elseif(nope.eq.8) then
556  call shape8h(xi,et,ze,xl,xsj,shp,iflag)
557  elseif(nope.eq.10) then
558  call shape10tet(xi,et,ze,xl,xsj,shp,iflag)
559  elseif(nope.eq.4) then
560  call shape4tet(xi,et,ze,xl,xsj,shp,iflag)
561  elseif(nope.eq.15) then
562  call shape15w(xi,et,ze,xl,xsj,shp,iflag)
563  else
564  call shape6w(xi,et,ze,xl,xsj,shp,iflag)
565  endif
566 !
567 ! vkl(m2,m3) contains the derivative of the m2-
568 ! component of the displacement with respect to
569 ! direction m3
570 !
571  do m2=1,3
572  do m3=1,3
573  vkl(m2,m3)=0.d0
574  enddo
575  enddo
576 !
577  do m1=1,nope
578  do m2=1,3
579  do m3=1,3
580  vkl(m2,m3)=vkl(m2,m3)+shp(m3,m1)*vl(m2,m1)
581  enddo
582 c write(*,*) 'vnoeie',i,konl(m1),(vkl(m2,k),k=1,3)
583  enddo
584  enddo
585 !
586 ! for frequency analysis or buckling with preload the
587 ! strains are calculated with respect to the deformed
588 ! configuration
589 ! for a linear iteration within a nonlinear increment:
590 ! the tangent matrix is calculated at strain at the end
591 ! of the previous increment
592 !
593  if((iperturb(1).eq.1).or.(iperturb(1).eq.-1))then
594  do m2=1,3
595  do m3=1,3
596  vokl(m2,m3)=0.d0
597  enddo
598  enddo
599 !
600  do m1=1,nope
601  do m2=1,3
602  do m3=1,3
603  vokl(m2,m3)=vokl(m2,m3)+
604  & shp(m3,m1)*voldl(m2,m1)
605  enddo
606  enddo
607  enddo
608  endif
609 !
610  kode=nelcon(1,imat)
611 !
612 ! calculating the strain
613 !
614 ! attention! exy,exz and eyz are engineering strains!
615 !
616  exx=vkl(1,1)
617  eyy=vkl(2,2)
618  ezz=vkl(3,3)
619  exy=vkl(1,2)+vkl(2,1)
620  exz=vkl(1,3)+vkl(3,1)
621  eyz=vkl(2,3)+vkl(3,2)
622 !
623  if(iperturb(2).eq.1) then
624 !
625 ! Lagrangian strain
626 !
627  exx=exx+(vkl(1,1)**2+vkl(2,1)**2+vkl(3,1)**2)/2.d0
628  eyy=eyy+(vkl(1,2)**2+vkl(2,2)**2+vkl(3,2)**2)/2.d0
629  ezz=ezz+(vkl(1,3)**2+vkl(2,3)**2+vkl(3,3)**2)/2.d0
630  exy=exy+vkl(1,1)*vkl(1,2)+vkl(2,1)*vkl(2,2)+
631  & vkl(3,1)*vkl(3,2)
632  exz=exz+vkl(1,1)*vkl(1,3)+vkl(2,1)*vkl(2,3)+
633  & vkl(3,1)*vkl(3,3)
634  eyz=eyz+vkl(1,2)*vkl(1,3)+vkl(2,2)*vkl(2,3)+
635  & vkl(3,2)*vkl(3,3)
636 !
637 ! for frequency analysis or buckling with preload the
638 ! strains are calculated with respect to the deformed
639 ! configuration
640 !
641  elseif(iperturb(1).eq.1) then
642  exx=exx+vokl(1,1)*vkl(1,1)+vokl(2,1)*vkl(2,1)+
643  & vokl(3,1)*vkl(3,1)
644  eyy=eyy+vokl(1,2)*vkl(1,2)+vokl(2,2)*vkl(2,2)+
645  & vokl(3,2)*vkl(3,2)
646  ezz=ezz+vokl(1,3)*vkl(1,3)+vokl(2,3)*vkl(2,3)+
647  & vokl(3,3)*vkl(3,3)
648  exy=exy+vokl(1,1)*vkl(1,2)+vokl(1,2)*vkl(1,1)+
649  & vokl(2,1)*vkl(2,2)+vokl(2,2)*vkl(2,1)+
650  & vokl(3,1)*vkl(3,2)+vokl(3,2)*vkl(3,1)
651  exz=exz+vokl(1,1)*vkl(1,3)+vokl(1,3)*vkl(1,1)+
652  & vokl(2,1)*vkl(2,3)+vokl(2,3)*vkl(2,1)+
653  & vokl(3,1)*vkl(3,3)+vokl(3,3)*vkl(3,1)
654  eyz=eyz+vokl(1,2)*vkl(1,3)+vokl(1,3)*vkl(1,2)+
655  & vokl(2,2)*vkl(2,3)+vokl(2,3)*vkl(2,2)+
656  & vokl(3,2)*vkl(3,3)+vokl(3,3)*vkl(3,2)
657  endif
658 !
659 c if(iperturb(2).eq.1) then
660 c!
661 c! Lagrangian strain
662 c!
663 c exx=exx+(vkl(1,1)**2+vkl(2,1)**2+vkl(3,1)**2)/2.d0
664 c eyy=eyy+(vkl(1,2)**2+vkl(2,2)**2+vkl(3,2)**2)/2.d0
665 c ezz=ezz+(vkl(1,3)**2+vkl(2,3)**2+vkl(3,3)**2)/2.d0
666 c exy=exy+vkl(1,1)*vkl(1,2)+vkl(2,1)*vkl(2,2)+
667 c & vkl(3,1)*vkl(3,2)
668 c exz=exz+vkl(1,1)*vkl(1,3)+vkl(2,1)*vkl(2,3)+
669 c & vkl(3,1)*vkl(3,3)
670 c eyz=eyz+vkl(1,2)*vkl(1,3)+vkl(2,2)*vkl(2,3)+
671 c & vkl(3,2)*vkl(3,3)
672 c!
673 c endif
674 !
675 ! storing the local strains
676 !
677  if(iperturb(1).ne.-1) then
678  eloc(1)=exx
679  eloc(2)=eyy
680  eloc(3)=ezz
681  eloc(4)=exy/2.d0
682  eloc(5)=exz/2.d0
683  eloc(6)=eyz/2.d0
684  else
685 !
686 ! linear iteration within a nonlinear increment:
687 !
688  eloc(1)=vokl(1,1)+
689  & (vokl(1,1)**2+vokl(2,1)**2+vokl(3,1)**2)/2.d0
690  eloc(2)=vokl(2,2)+
691  & (vokl(1,2)**2+vokl(2,2)**2+vokl(3,2)**2)/2.d0
692  eloc(3)=vokl(3,3)+
693  & (vokl(1,3)**2+vokl(2,3)**2+vokl(3,3)**2)/2.d0
694  eloc(4)=(vokl(1,2)+vokl(2,1)+vokl(1,1)*vokl(1,2)+
695  & vokl(2,1)*vokl(2,2)+vokl(3,1)*vokl(3,2))/2.d0
696  eloc(5)=(vokl(1,3)+vokl(3,1)+vokl(1,1)*vokl(1,3)+
697  & vokl(2,1)*vokl(2,3)+vokl(3,1)*vokl(3,3))/2.d0
698  eloc(6)=(vokl(2,3)+vokl(3,2)+vokl(1,2)*vokl(1,3)+
699  & vokl(2,2)*vokl(2,3)+vokl(3,2)*vokl(3,3))/2.d0
700  endif
701 !
702 ! calculating the deformation gradient (needed to
703 ! convert the element stiffness matrix from spatial
704 ! coordinates to material coordinates
705 ! deformation plasticity)
706 !
707  if((kode.eq.-50).or.(kode.le.-100)) then
708 !
709 ! calculating the deformation gradient
710 !
711 c Bernhardi start
712  xkl(1,1)=vkl(1,1)+1.0d0
713  xkl(2,2)=vkl(2,2)+1.0d0
714  xkl(3,3)=vkl(3,3)+1.0d0
715 c Bernhardi end
716  xkl(1,2)=vkl(1,2)
717  xkl(1,3)=vkl(1,3)
718  xkl(2,3)=vkl(2,3)
719  xkl(2,1)=vkl(2,1)
720  xkl(3,1)=vkl(3,1)
721  xkl(3,2)=vkl(3,2)
722 !
723 ! calculating the Jacobian
724 !
725  vj=xkl(1,1)*(xkl(2,2)*xkl(3,3)-xkl(2,3)*xkl(3,2))
726  & -xkl(1,2)*(xkl(2,1)*xkl(3,3)-xkl(2,3)*xkl(3,1))
727  & +xkl(1,3)*(xkl(2,1)*xkl(3,2)-xkl(2,2)*xkl(3,1))
728 !
729 ! inversion of the deformation gradient (only for
730 ! deformation plasticity)
731 !
732  if(kode.eq.-50) then
733 !
734  ckl(1,1)=(xkl(2,2)*xkl(3,3)-xkl(2,3)*xkl(3,2))/vj
735  ckl(2,2)=(xkl(1,1)*xkl(3,3)-xkl(1,3)*xkl(3,1))/vj
736  ckl(3,3)=(xkl(1,1)*xkl(2,2)-xkl(1,2)*xkl(2,1))/vj
737  ckl(1,2)=(xkl(1,3)*xkl(3,2)-xkl(1,2)*xkl(3,3))/vj
738  ckl(1,3)=(xkl(1,2)*xkl(2,3)-xkl(2,2)*xkl(1,3))/vj
739  ckl(2,3)=(xkl(2,1)*xkl(1,3)-xkl(1,1)*xkl(2,3))/vj
740  ckl(2,1)=(xkl(3,1)*xkl(2,3)-xkl(2,1)*xkl(3,3))/vj
741  ckl(3,1)=(xkl(2,1)*xkl(3,2)-xkl(2,2)*xkl(3,1))/vj
742  ckl(3,2)=(xkl(3,1)*xkl(1,2)-xkl(1,1)*xkl(3,2))/vj
743 !
744 ! converting the Lagrangian strain into Eulerian
745 ! strain (only for deformation plasticity)
746 !
747  cauchy=0
748  call str2mat(eloc,ckl,vj,cauchy)
749  endif
750 !
751  endif
752 !
753 ! calculating fields for incremental plasticity
754 !
755  if(kode.le.-100) then
756 !
757 ! calculating the deformation gradient at the
758 ! start of the increment
759 !
760 ! calculating the displacement gradient at the
761 ! start of the increment
762 !
763  do m2=1,3
764  do m3=1,3
765  vikl(m2,m3)=0.d0
766  enddo
767  enddo
768 !
769  do m1=1,nope
770  do m2=1,3
771  do m3=1,3
772  vikl(m2,m3)=vikl(m2,m3)
773  & +shp(m3,m1)*vini(m2,konl(m1))
774  enddo
775  enddo
776  enddo
777 !
778 ! calculating the deformation gradient of the old
779 ! fields
780 !
781  xikl(1,1)=vikl(1,1)+1
782  xikl(2,2)=vikl(2,2)+1.
783  xikl(3,3)=vikl(3,3)+1.
784  xikl(1,2)=vikl(1,2)
785  xikl(1,3)=vikl(1,3)
786  xikl(2,3)=vikl(2,3)
787  xikl(2,1)=vikl(2,1)
788  xikl(3,1)=vikl(3,1)
789  xikl(3,2)=vikl(3,2)
790 !
791 ! calculating the Jacobian
792 !
793  vij=xikl(1,1)*(xikl(2,2)*xikl(3,3)
794  & -xikl(2,3)*xikl(3,2))
795  & -xikl(1,2)*(xikl(2,1)*xikl(3,3)
796  & -xikl(2,3)*xikl(3,1))
797  & +xikl(1,3)*(xikl(2,1)*xikl(3,2)
798  & -xikl(2,2)*xikl(3,1))
799 !
800 ! stresses at the start of the increment
801 !
802  do m1=1,6
803  stre(m1)=sti(m1,jj,i)
804  enddo
805 !
806  endif
807 !
808 ! prestress values
809 !
810  if(iprestr.eq.1) then
811  do kk=1,6
812  beta(kk)=-prestr(kk,jj,i)
813  enddo
814  else
815  do kk=1,6
816  beta(kk)=0.d0
817  enddo
818  endif
819 !
820  if(ithermal(1).ge.1) then
821 !
822 ! calculating the temperature difference in
823 ! the integration point
824 !
825  t0l=0.d0
826  t1l=0.d0
827  if(ithermal(1).eq.1) then
828  if((lakonl(4:5).eq.'8 ').or.
829  & (lakonl(4:5).eq.'8I')) then
830  do i1=1,8
831  t0l=t0l+t0(konl(i1))/8.d0
832  t1l=t1l+t1(konl(i1))/8.d0
833  enddo
834  elseif(lakonl(4:6).eq.'20 ') then
835  nopered=20
836  call lintemp(t0,t1,konl,nopered,jj,t0l,t1l)
837  elseif(lakonl(4:6).eq.'10T') then
838  call linscal10(t0,konl,t0l,null,shp)
839  call linscal10(t1,konl,t1l,null,shp)
840  else
841  do i1=1,nope
842  t0l=t0l+shp(4,i1)*t0(konl(i1))
843  t1l=t1l+shp(4,i1)*t1(konl(i1))
844  enddo
845  endif
846  elseif(ithermal(1).ge.2) then
847  if((lakonl(4:5).eq.'8 ').or.
848  & (lakonl(4:5).eq.'8I')) then
849  do i1=1,8
850  t0l=t0l+t0(konl(i1))/8.d0
851  t1l=t1l+vold(0,konl(i1))/8.d0
852  enddo
853  elseif(lakonl(4:6).eq.'20 ') then
854  nopered=20
855  call lintemp_th(t0,vold,konl,nopered,jj,t0l,t1l,mi)
856  elseif(lakonl(4:6).eq.'10T') then
857  call linscal10(t0,konl,t0l,null,shp)
858  call linscal10(vold,konl,t1l,mi(2),shp)
859  else
860  do i1=1,nope
861  t0l=t0l+shp(4,i1)*t0(konl(i1))
862  t1l=t1l+shp(4,i1)*vold(0,konl(i1))
863  enddo
864  endif
865  endif
866  tt=t1l-t0l
867  endif
868 !
869 ! calculating the coordinates of the integration point
870 ! for material orientation purposes (for cylindrical
871 ! coordinate systems)
872 !
873  if((iorien.gt.0).or.(kode.le.-100)) then
874  do j=1,3
875  pgauss(j)=0.d0
876  do i1=1,nope
877  pgauss(j)=pgauss(j)+shp(4,i1)*co(j,konl(i1))
878  enddo
879  enddo
880  endif
881 !
882 ! material data; for linear elastic materials
883 ! this includes the calculation of the stiffness
884 ! matrix
885 !
886  istiff=0
887 !
888  call materialdata_me(elcon,nelcon,rhcon,nrhcon,alcon,
889  & nalcon,imat,amat,iorien,pgauss,orab,ntmat_,
890  & elas,rho,i,ithermal,alzero,mattyp,t0l,t1l,ihyper,
891  & istiff,elconloc,eth,kode,plicon,nplicon,
892  & plkcon,nplkcon,npmat_,plconloc,mi(1),dtime,jj,
893  & xstiff,ncmat_)
894 !
895 ! determining the mechanical strain
896 !
897  if(ithermal(1).ne.0) then
898  do m1=1,6
899  emec(m1)=eloc(m1)-eth(m1)
900 c emec0(m1)=emeini(m1,jj,i)
901  enddo
902  else
903  do m1=1,6
904  emec(m1)=eloc(m1)
905 c emec0(m1)=emeini(m1,jj,i)
906  enddo
907  endif
908  if(kode.le.-100) then
909  do m1=1,6
910  emec0(m1)=emeini(m1,jj,i)
911  enddo
912  endif
913 !
914 ! subtracting the plastic initial strains
915 !
916  if(iprestr.eq.2) then
917  do m1=1,6
918  emec(m1)=emec(m1)-prestr(m1,jj,i)
919  enddo
920  endif
921 !
922 ! calculating the local stiffness and stress
923 !
924  call mechmodel(elconloc,elas,emec,kode,emec0,ithermal,
925  & icmd,beta,stre,xkl,ckl,vj,xikl,vij,
926  & plconloc,xstate,xstateini,ielas,
927  & amat,t1l,dtime,time,ttime,i,jj,nstate_,mi(1),
928  & iorien,pgauss,orab,eloc,mattyp,qa(3),istep,iinc,
929  & ipkon,nmethod,iperturb,qa(4))
930 !
931  if(((nmethod.ne.4).or.(iperturb(1).ne.0)).and.
932  & (nmethod.ne.5)) then
933  if(idesvar.eq.0) then
934  do m1=1,21
935  xstiff(m1,jj,i)=elas(m1)
936  enddo
937  elseif(icoordinate.ne.1) then
938 !
939 ! for orientation design variables: store the modified
940 ! stiffness matrix for use in mafillsmse
941 !
942  idir=idesvar-3*((idesvar-1)/3)
943  do m1=1,21
944  dxstiff(m1,jj,i,idir)=elas(m1)
945  enddo
946  endif
947  endif
948 !
949  if(iperturb(1).eq.-1) then
950 !
951 ! if the forced displacements were changed at
952 ! the start of a nonlinear step, the nodal
953 ! forces due do this displacements are
954 ! calculated in a purely linear way, and
955 ! the first iteration is purely linear in order
956 ! to allow the displacements to redistribute
957 ! in a quasi-static way (only applies to
958 ! quasi-static analyses (*STATIC))
959 !
960  eloc(1)=exx-vokl(1,1)
961  eloc(2)=eyy-vokl(2,2)
962  eloc(3)=ezz-vokl(3,3)
963  eloc(4)=exy-(vokl(1,2)+vokl(2,1))
964  eloc(5)=exz-(vokl(1,3)+vokl(3,1))
965  eloc(6)=eyz-(vokl(2,3)+vokl(3,2))
966 !
967  if(mattyp.eq.1) then
968  e=elas(1)
969  un=elas(2)
970  um=e/(1.d0+un)
971  al=un*um/(1.d0-2.d0*un)
972  um=um/2.d0
973  am1=al*(eloc(1)+eloc(2)+eloc(3))
974  stre(1)=am1+2.d0*um*eloc(1)
975  stre(2)=am1+2.d0*um*eloc(2)
976  stre(3)=am1+2.d0*um*eloc(3)
977  stre(4)=um*eloc(4)
978  stre(5)=um*eloc(5)
979  stre(6)=um*eloc(6)
980  elseif(mattyp.eq.2) then
981  stre(1)=eloc(1)*elas(1)+eloc(2)*elas(2)
982  & +eloc(3)*elas(4)
983  stre(2)=eloc(1)*elas(2)+eloc(2)*elas(3)
984  & +eloc(3)*elas(5)
985  stre(3)=eloc(1)*elas(4)+eloc(2)*elas(5)
986  & +eloc(3)*elas(6)
987  stre(4)=eloc(4)*elas(7)
988  stre(5)=eloc(5)*elas(8)
989  stre(6)=eloc(6)*elas(9)
990  elseif(mattyp.eq.3) then
991  stre(1)=eloc(1)*elas(1)+eloc(2)*elas(2)+
992  & eloc(3)*elas(4)+eloc(4)*elas(7)+
993  & eloc(5)*elas(11)+eloc(6)*elas(16)
994  stre(2)=eloc(1)*elas(2)+eloc(2)*elas(3)+
995  & eloc(3)*elas(5)+eloc(4)*elas(8)+
996  & eloc(5)*elas(12)+eloc(6)*elas(17)
997  stre(3)=eloc(1)*elas(4)+eloc(2)*elas(5)+
998  & eloc(3)*elas(6)+eloc(4)*elas(9)+
999  & eloc(5)*elas(13)+eloc(6)*elas(18)
1000  stre(4)=eloc(1)*elas(7)+eloc(2)*elas(8)+
1001  & eloc(3)*elas(9)+eloc(4)*elas(10)+
1002  & eloc(5)*elas(14)+eloc(6)*elas(19)
1003  stre(5)=eloc(1)*elas(11)+eloc(2)*elas(12)+
1004  & eloc(3)*elas(13)+eloc(4)*elas(14)+
1005  & eloc(5)*elas(15)+eloc(6)*elas(20)
1006  stre(6)=eloc(1)*elas(16)+eloc(2)*elas(17)+
1007  & eloc(3)*elas(18)+eloc(4)*elas(19)+
1008  & eloc(5)*elas(20)+eloc(6)*elas(21)
1009  endif
1010  endif
1011 !
1012  skl(1,1)=stre(1)
1013  skl(2,2)=stre(2)
1014  skl(3,3)=stre(3)
1015  skl(2,1)=stre(4)
1016  skl(3,1)=stre(5)
1017  skl(3,2)=stre(6)
1018 !
1019  skl(1,2)=skl(2,1)
1020  skl(1,3)=skl(3,1)
1021  skl(2,3)=skl(3,2)
1022 !
1023 ! calculation of the nodal forces
1024 !
1025  if(calcul_fn.eq.1)then
1026  if(idesvar.eq.0) then
1027 !
1028 ! unperturbed state
1029 !
1030  do m1=1,nope
1031  do m2=1,3
1032 !
1033 ! linear elastic part
1034 !
1035  do m3=1,3
1036  fn0(m2,indexe+m1)=fn0(m2,indexe+m1)+
1037  & xsj*skl(m2,m3)*shp(m3,m1)*weight
1038  enddo
1039 !
1040 ! nonlinear geometric part
1041 !
1042  if(iperturb(2).eq.1) then
1043  do m3=1,3
1044  do m4=1,3
1045  fn0(m2,indexe+m1)=fn0(m2,indexe+m1)+
1046  & xsj*skl(m4,m3)*weight*
1047  & (vkl(m2,m4)*shp(m3,m1)+
1048  & vkl(m2,m3)*shp(m4,m1))/2.d0
1049  enddo
1050  enddo
1051  endif
1052 !
1053  enddo
1054  enddo
1055 c Bernhardi start
1056  if(lakonl(1:5).eq.'C3D8R') then
1057  ahr=elas(1)*a
1058 !
1059  do i1=1,3
1060  do k=1,4
1061  hglf(i1,k)=0.0d0
1062  do j=1,8
1063  hglf(i1,k)=hglf(i1,k)+gs(j,k)*vl(i1,j)
1064  enddo
1065  hglf(i1,k)=hglf(i1,k)*ahr
1066  enddo
1067  enddo
1068  do i1=1,3
1069  do j=1,8
1070  do k=1,4
1071  fn0(i1,indexe+j)=fn0(i1,indexe+j)+
1072  & hglf(i1,k)*gs(j,k)
1073  enddo
1074  enddo
1075  enddo
1076  endif
1077 c Bernhardi end
1078  else
1079 !
1080 ! perturbed state
1081 !
1082  do m1=1,nope
1083  do m2=1,3
1084 !
1085 ! linear elastic part
1086 !
1087  do m3=1,3
1088  dfn(m2,konl(m1))=dfn(m2,konl(m1))+
1089  & xsj*skl(m2,m3)*shp(m3,m1)*weight
1090  enddo
1091 !
1092 ! nonlinear geometric part
1093 !
1094  if(iperturb(2).eq.1) then
1095  do m3=1,3
1096  do m4=1,3
1097  dfn(m2,konl(m1))=dfn(m2,konl(m1))+
1098  & xsj*skl(m4,m3)*weight*
1099  & (vkl(m2,m4)*shp(m3,m1)+
1100  & vkl(m2,m3)*shp(m4,m1))/2.d0
1101  enddo
1102  enddo
1103  endif
1104 !
1105  enddo
1106  enddo
1107 c Bernhardi start
1108  if(lakonl(1:5).eq.'C3D8R') then
1109 !
1110  ahr=elas(1)*a
1111 !
1112  do i1=1,3
1113  do k=1,4
1114  hglf(i1,k)=0.0d0
1115  do j=1,8
1116  hglf(i1,k)=hglf(i1,k)+gs(j,k)*vl(i1,j)
1117  enddo
1118  hglf(i1,k)=hglf(i1,k)*ahr
1119  enddo
1120  enddo
1121  do i1=1,3
1122  do j=1,8
1123  do k=1,4
1124  dfn(i1,konl(j))=dfn(i1,konl(j))
1125  & +hglf(i1,k)*gs(j,k)
1126  enddo
1127  enddo
1128  enddo
1129  endif
1130 c Bernhardi end
1131  endif
1132  endif
1133 !
1134  enddo ! enddo loop over the integration points
1135 !
1136 ! subtracting the unperturbed state of the element at stake
1137 !
1138  if(idesvar.gt.0) then
1139  do m1=1,nope
1140  do m2=1,3
1141  dfn(m2,konl(m1))=dfn(m2,konl(m1))
1142  & -fn0(m2,indexe+m1)
1143  enddo
1144  enddo
1145  endif
1146 !
1147 ! end of loop over all elements in thread
1148  enddo
1149 !
1150  return
subroutine shape20h_pl(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape20h_pl.f:20
subroutine shape6w(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape6w.f:20
subroutine linscal10(scal, konl, scall, idim, shp)
Definition: linscal10.f:20
subroutine shape8q(xi, et, xl, xsj, xs, shp, iflag)
Definition: shape8q.f:20
subroutine springforc_n2f(xl, konl, vl, imat, elcon, nelcon, elas, fnl, ncmat_, ntmat_, nope, lakonl, t1l, kode, elconloc, plicon, nplicon, npmat_, senergy, nener, cstr, mi, springarea, nmethod, ne0, nstate_, xstateini, xstate, reltime, ielas, venergy, ielorien, orab, norien, nelem)
Definition: springforc_n2f.f:24
subroutine mechmodel(elconloc, elas, emec, kode, emec0, ithermal, icmd, beta, stre, xkl, ckl, vj, xikl, vij, plconloc, xstate, xstateini, ielas, amat, t1l, dtime, time, ttime, iel, iint, nstate_, mi, iorien, pgauss, orab, eloc, mattyp, pnewdt, istep, iinc, ipkon, nmethod, iperturb, depvisc)
Definition: mechmodel.f:24
subroutine shape10tet(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape10tet.f:20
subroutine shape8h(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape8h.f:20
subroutine shape15w(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape15w.f:20
subroutine shape20h_ax(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape20h_ax.f:20
subroutine lintemp_th(t0, vold, konl, nope, jj, t0l, t1l, mi)
Definition: lintemp_th.f:20
subroutine materialdata_me(elcon, nelcon, rhcon, nrhcon, alcon, nalcon, imat, amat, iorien, pgauss, orab, ntmat_, elas, rho, iel, ithermal, alzero, mattyp, t0l, t1l, ihyper, istiff, elconloc, eth, kode, plicon, nplicon, plkcon, nplkcon, npmat_, plconloc, mi, dtime, iint, xstiff, ncmat_)
Definition: materialdata_me.f:23
subroutine lintemp(t0, t1, konl, nope, jj, t0l, t1l)
Definition: lintemp.f:20
subroutine shape20h(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape20h.f:20
subroutine str2mat(str, ckl, vj, cauchy)
Definition: str2mat.f:20
subroutine shape8hr(xl, xsj, shp, gs, a)
Definition: shape8hr.f:20
subroutine beamintscheme(lakonl, mint3d, npropstart, prop, kk, xi, et, ze, weight)
Definition: beamintscheme.f:21
subroutine shape4tet(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape4tet.f:20
subroutine cauchy(n, x, l, u, nbd, g, iorder, iwhere, t, d, xcp, m, wy, ws, sy, wt, theta, col, head, p, c, wbp, v, nseg, iprint, sbgnrm, info, epsmch)
Definition: lbfgsb.f:1225
subroutine shape8hu(xi, et, ze, xl, xsj, shp, iflag)
Definition: shape8hu.f:20
subroutine shape6tri(xi, et, xl, xsj, xs, shp, iflag)
Definition: shape6tri.f:20
subroutine springforc_f2f(xl, vl, imat, elcon, nelcon, elas, fnl, ncmat_, ntmat_, nope, lakonl, t1l, kode, elconloc, plicon, nplicon, npmat_, senergy, nener, cstr, mi, springarea, nmethod, ne0, nstate_, xstateini, xstate, reltime, ielas, iloc, jfaces, igauss, pslavsurf, pmastsurf, clearini, venergy, kscale, konl, iout, nelem)
Definition: springforc_f2f.f:26
subroutine thickness(dgdx, nobject, nodedesiboun, ndesiboun, objectset, xo, yo, zo, x, y, z, nx, ny, nz, co, ifree, ndesia, ndesib, iobject, ndesi, dgdxglob, nk)
Definition: thickness.f:22
Hosted by OpenAircraft.com, (Michigan UAV, LLC)