CalculiX  2.13
A Free Software Three-Dimensional Structural Finite Element Program
Functions/Subroutines
mafillsmcs.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine mafillsmcs (co, nk, kon, ipkon, lakon, ne, nodeboun, ndirboun, xboun, nboun, ipompc, nodempc, coefmpc, nmpc, nodeforc, ndirforc, xforc, nforc, nelemload, sideload, xload, nload, xbody, ipobody, nbody, cgr, ad, au, bb, nactdof, icol, jq, irow, neq, nzl, nmethod, ikmpc, ilmpc, ikboun, ilboun, elcon, nelcon, rhcon, nrhcon, alcon, nalcon, alzero, ielmat, ielorien, norien, orab, ntmat_, t0, t1, ithermal, prestr, iprestr, vold, iperturb, sti, nzs, stx, adb, aub, iexpl, plicon, nplicon, plkcon, nplkcon, xstiff, npmat_, dtime, matname, mi, ics, cs, nm, ncmat_, labmpc, mass, stiffness, buckling, rhsi, intscheme, mcs, coriolis, ibody, xloadold, reltime, ielcs, veold, springarea, thicke, integerglob, doubleglob, tieset, istartset, iendset, ialset, ntie, nasym, pslavsurf, pmastsurf, mortar, clearini, ielprop, prop, ne0, kscale, xstateini, xstate, nstate_)
 

Function/Subroutine Documentation

◆ mafillsmcs()

subroutine mafillsmcs ( real*8, dimension(3,*)  co,
integer  nk,
integer, dimension(*)  kon,
integer, dimension(*)  ipkon,
character*8, dimension(*)  lakon,
integer  ne,
integer, dimension(*)  nodeboun,
integer, dimension(*)  ndirboun,
real*8, dimension(*)  xboun,
integer  nboun,
integer, dimension(*)  ipompc,
integer, dimension(3,*)  nodempc,
real*8, dimension(*)  coefmpc,
integer  nmpc,
integer, dimension(2,*)  nodeforc,
integer, dimension(*)  ndirforc,
real*8, dimension(*)  xforc,
integer  nforc,
integer, dimension(2,*)  nelemload,
character*20, dimension(*)  sideload,
real*8, dimension(2,*)  xload,
integer  nload,
real*8, dimension(7,*)  xbody,
integer, dimension(2,*)  ipobody,
integer  nbody,
real*8, dimension(4,*)  cgr,
real*8, dimension(*)  ad,
real*8, dimension(*)  au,
real*8, dimension(*)  bb,
integer, dimension(0:mi(2),*)  nactdof,
integer, dimension(*)  icol,
integer, dimension(*)  jq,
integer, dimension(*)  irow,
integer  neq,
integer  nzl,
integer  nmethod,
integer, dimension(*)  ikmpc,
integer, dimension(*)  ilmpc,
integer, dimension(*)  ikboun,
integer, dimension(*)  ilboun,
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  ithermal,
real*8, dimension(6,mi(1),*)  prestr,
integer  iprestr,
real*8, dimension(0:mi(2),*)  vold,
integer, dimension(*)  iperturb,
real*8, dimension(6,mi(1),*)  sti,
integer  nzs,
real*8, dimension(6,mi(1),*)  stx,
real*8, dimension(*)  adb,
real*8, dimension(*)  aub,
integer  iexpl,
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(27,mi(1),*)  xstiff,
integer  npmat_,
real*8  dtime,
character*80, dimension(*)  matname,
integer, dimension(*)  mi,
integer, dimension(*)  ics,
real*8, dimension(17,*)  cs,
integer  nm,
integer  ncmat_,
character*20, dimension(*)  labmpc,
logical  mass,
logical  stiffness,
logical  buckling,
logical  rhsi,
integer  intscheme,
integer  mcs,
logical  coriolis,
integer, dimension(3,*)  ibody,
real*8, dimension(2,*)  xloadold,
real*8  reltime,
integer, dimension(*)  ielcs,
real*8, dimension(0:mi(2),*)  veold,
real*8, dimension(2,*)  springarea,
real*8, dimension(mi(3),*)  thicke,
integer, dimension(*)  integerglob,
real*8, dimension(*)  doubleglob,
character*81, dimension(3,*)  tieset,
integer, dimension(*)  istartset,
integer, dimension(*)  iendset,
integer, dimension(*)  ialset,
integer  ntie,
integer  nasym,
real*8, dimension(3,*)  pslavsurf,
real*8, dimension(6,*)  pmastsurf,
integer  mortar,
real*8, dimension(3,9,*)  clearini,
integer, dimension(*)  ielprop,
real*8, dimension(*)  prop,
integer  ne0,
integer  kscale,
real*8, dimension(nstate_,mi(1),*)  xstateini,
real*8, dimension(nstate_,mi(1),*)  xstate,
integer  nstate_ 
)
35 !
36 ! filling the stiffness matrix in spare matrix format (sm)
37 ! for cyclic symmetry calculations
38 !
39  implicit none
40 !
41  logical mass,stiffness,buckling,rhsi,coriolis
42 !
43  character*8 lakon(*)
44  character*20 labmpc(*),sideload(*)
45  character*80 matname(*)
46  character*81 tieset(3,*)
47 !
48  integer kon(*),nodeboun(*),ndirboun(*),ipompc(*),nodempc(3,*),
49  & nodeforc(2,*),ndirforc(*),nelemload(2,*),icol(*),jq(*),ikmpc(*),
50  & ilmpc(*),ikboun(*),ilboun(*),mi(*),nstate_,ne0,ielprop(*),
51  & nactdof(0:mi(2),*),irow(*),istartset(*),iendset(*),kscale,
52  & nelcon(2,*),nrhcon(*),nalcon(2,*),ielmat(mi(3),*),
53  & ielorien(mi(3),*),integerglob(*),ialset(*),ntie,
54  & ipkon(*),ics(*),ij,ilength,lprev,ipobody(2,*),nbody,
55  & ibody(3,*),nk,ne,nboun,nmpc,nforc,nload,neq,nzl,nmethod,
56  & ithermal,iprestr,iperturb(*),nzs,i,j,k,l,m,idist,jj,
57  & ll,id,id1,id2,ist,ist1,ist2,index,jdof1,jdof2,idof1,idof2,
58  & mpc1,mpc2,index1,index2,node1,node2,kflag,nasym,mortar,
59  & ntmat_,indexe,nope,norien,iexpl,i0,nm,inode,icomplex,
60  & inode1,icomplex1,inode2,icomplex2,ner,ncmat_,intscheme,istep,
61  & iinc,mcs,ielcs(*),nplicon(0:ntmat_,*),nplkcon(0:ntmat_,*),npmat_
62 !
63  real*8 co(3,*),xboun(*),coefmpc(*),xforc(*),xload(2,*),p1(3),
64  & p2(3),ad(*),au(*),bodyf(3),bb(*),xbody(7,*),cgr(4,*),prop(*),
65  & t0(*),t1(*),prestr(6,mi(1),*),vold(0:mi(2),*),s(60,60),
66  & xstate(nstate_,mi(1),*),xstateini(nstate_,mi(1),*),ff(60),
67  & sti(6,mi(1),*),sm(60,60),stx(6,mi(1),*),adb(*),aub(*),
68  & elcon(0:ncmat_,ntmat_,*),rhcon(0:1,ntmat_,*),xloadold(2,*),
69  & alcon(0:6,ntmat_,*),cs(17,*),alzero(*),orab(7,*),reltime,
70  & springarea(2,*),plicon(0:2*npmat_,ntmat_,*),
71  & plkcon(0:2*npmat_,ntmat_,*),thicke(mi(3),*),doubleglob(*),
72  & xstiff(27,mi(1),*),pi,theta,ti,tr,veold(0:mi(2),*),om,valu2,
73  & value,dtime,walue,walu2,time,ttime,clearini(3,9,*),
74  & pslavsurf(3,*),pmastsurf(6,*)
75 !
76 !
77 ! calculating the scaling factors for the cyclic symmetry calculation
78 !
79  pi=4.d0*datan(1.d0)
80 !
81  do i=1,mcs
82  write(*,*) '*INFO in mafillsmcs: calculating nodal diameter',nm
83  write(*,*) ' for',cs(1,i),'sectors'
84  write(*,*) ' (cyclic symmetry definition number',i,')'
85  write(*,*)
86  theta=nm*2.d0*pi/cs(1,i)
87  cs(15,i)=dcos(theta)
88  cs(16,i)=dsin(theta)
89  enddo
90 !
91  kflag=2
92  i0=0
93 !
94 ! determining nzl
95 !
96  nzl=0
97  do i=neq,1,-1
98  if(icol(i).gt.0) then
99  nzl=i
100  exit
101  endif
102  enddo
103 !
104 ! initializing the matrices
105 !
106  do i=1,neq
107  ad(i)=0.d0
108  enddo
109  do i=1,nzs
110  au(i)=0.d0
111  enddo
112 !
113  do i=1,neq
114  adb(i)=0.d0
115  enddo
116  do i=1,nzs
117  aub(i)=0.d0
118  enddo
119 !
120  ner=neq/2
121 !
122 ! loop over all elements
123 !
124 ! initialisation of the error parameter
125 !
126 c ne0=0
127  do i=1,ne
128 !
129  if(ipkon(i).lt.0) cycle
130  indexe=ipkon(i)
131 c Bernhardi start
132  if(lakon(i)(4:5).eq.'8I') then
133  nope=11
134 c Bernhardi end
135  elseif(lakon(i)(4:5).eq.'20') then
136  nope=20
137  elseif(lakon(i)(4:4).eq.'2') then
138  nope=26
139  elseif(lakon(i)(4:4).eq.'8') then
140  nope=8
141  elseif(lakon(i)(4:5).eq.'10') then
142  nope=10
143  elseif(lakon(i)(4:4).eq.'4') then
144  nope=4
145  elseif(lakon(i)(4:5).eq.'15') then
146  nope=15
147  elseif(lakon(i)(4:4).eq.'6') then
148  nope=6
149  elseif(lakon(i)(1:2).eq.'ES') then
150  read(lakon(i)(8:8),'(i1)') nope
151  nope=nope+1
152 !
153 ! local contact spring number
154 !
155  if(lakon(i)(7:7).eq.'C') then
156  if(nasym.eq.1) cycle
157  endif
158  else
159  cycle
160  endif
161 !
162  om=0.d0
163 !
164  if((nbody.gt.0).and.(lakon(i)(1:1).ne.'E')) then
165 !
166 ! assigning centrifugal forces
167 !
168  index=i
169  do
170  j=ipobody(1,index)
171  if(j.eq.0) exit
172  if(ibody(1,j).eq.1) then
173  om=xbody(1,j)
174  p1(1)=xbody(2,j)
175  p1(2)=xbody(3,j)
176  p1(3)=xbody(4,j)
177  p2(1)=xbody(5,j)
178  p2(2)=xbody(6,j)
179  p2(3)=xbody(7,j)
180  endif
181  index=ipobody(2,index)
182  if(index.eq.0) exit
183  enddo
184  endif
185 !
186  call e_c3d(co,kon,lakon(i),p1,p2,om,bodyf,nbody,s,sm,ff,i,
187  & nmethod,elcon,nelcon,rhcon,nrhcon,alcon,nalcon,
188  & alzero,ielmat,ielorien,norien,orab,ntmat_,
189  & t0,t1,ithermal,vold,iperturb,nelemload,sideload,xload,
190  & nload,idist,sti,stx,iexpl,plicon,
191  & nplicon,plkcon,nplkcon,xstiff,npmat_,
192  & dtime,matname,mi(1),ncmat_,mass,stiffness,buckling,rhsi,
193  & intscheme,ttime,time,istep,iinc,coriolis,xloadold,
194  & reltime,ipompc,nodempc,coefmpc,nmpc,ikmpc,ilmpc,veold,
195  & springarea,nstate_,xstateini,xstate,ne0,ipkon,thicke,
196  & integerglob,doubleglob,tieset,istartset,
197  & iendset,ialset,ntie,nasym,pslavsurf,pmastsurf,mortar,
198  & clearini,ielprop,prop,kscale)
199 !
200  do jj=1,3*nope
201 !
202  j=(jj-1)/3+1
203  k=jj-3*(j-1)
204 !
205  node1=kon(indexe+j)
206  jdof1=nactdof(k,node1)
207 !
208  do ll=jj,3*nope
209  if (mcs.gt.1)then
210  if(ielcs(i).gt.0) then
211  s(jj,ll)=(cs(1,(ielcs(i)+1))/cs(1,1))*s(jj,ll)
212  sm(jj,ll)=(cs(1,(ielcs(i)+1))/cs(1,1))*sm(jj,ll)
213  endif
214  endif
215 !
216  l=(ll-1)/3+1
217  m=ll-3*(l-1)
218 !
219  node2=kon(indexe+l)
220  jdof2=nactdof(m,node2)
221 !
222 ! check whether one of the DOF belongs to a SPC or MPC
223 !
224  if((jdof1.gt.0).and.(jdof2.gt.0)) then
225  call add_sm_ei(au,ad,aub,adb,jq,irow,jdof1,jdof2,
226  & s(jj,ll),sm(jj,ll),jj,ll)
227  call add_sm_ei(au,ad,aub,adb,jq,irow,jdof1+ner,jdof2+ner,
228  & s(jj,ll),sm(jj,ll),jj,ll)
229  elseif((jdof1.gt.0).or.(jdof2.gt.0)) then
230 !
231 ! idof1: genuine DOF
232 ! idof2: nominal DOF of the SPC/MPC
233 !
234  if(jdof1.le.0) then
235  idof1=jdof2
236 c idof2=(node1-1)*8+k
237  idof2=jdof1
238  else
239  idof1=jdof1
240 c idof2=(node2-1)*8+m
241  idof2=jdof2
242  endif
243 !
244  if(nmpc.gt.0) then
245 c call nident(ikmpc,idof2,nmpc,id)
246 c if((id.gt.0).and.(ikmpc(id).eq.idof2)) then
247  if(idof2.ne.2*(idof2/2)) then
248 !
249 ! regular DOF / MPC
250 !
251 c id1=ilmpc(id)
252  id1=(-idof2+1)/2
253  ist=ipompc(id1)
254  index=nodempc(3,ist)
255  if(index.eq.0) cycle
256  do
257  inode=nodempc(1,index)
258  icomplex=0
259 c write(*,*) id1,labmpc(id1)(1:9)
260  if(labmpc(id1)(1:6).eq.'CYCLIC') then
261  read(labmpc(id1)(7:20),'(i14)') icomplex
262  elseif(labmpc(id1)(1:9).eq.'SUBCYCLIC') then
263  do ij=1,mcs
264  ilength=int(cs(4,ij))
265  lprev=int(cs(14,ij))
266  call nident(ics(lprev+1),inode,ilength,id)
267  if(id.gt.0) then
268  if(ics(lprev+id).eq.inode) then
269  icomplex=ij
270  exit
271  endif
272  endif
273  enddo
274  endif
275  idof2=nactdof(nodempc(2,index),inode)
276  if(idof2.gt.0) then
277  value=-coefmpc(index)*s(jj,ll)/coefmpc(ist)
278  valu2=-coefmpc(index)*sm(jj,ll)/
279  & coefmpc(ist)
280  if(idof1.eq.idof2) then
281  value=2.d0*value
282  valu2=2.d0*valu2
283  endif
284  if(icomplex.eq.0) then
285  call add_sm_ei(au,ad,aub,adb,jq,irow,
286  & idof1,idof2,value,valu2,i0,i0)
287  call add_sm_ei(au,ad,aub,adb,jq,irow,
288  & idof1+ner,idof2+ner,value,valu2,i0,i0)
289  else
290  walue=value*cs(15,icomplex)
291  walu2=valu2*cs(15,icomplex)
292  call add_sm_ei(au,ad,aub,adb,jq,irow,
293  & idof1,idof2,walue,walu2,i0,i0)
294  call add_sm_ei(au,ad,aub,adb,jq,irow,
295  & idof1+ner,idof2+ner,walue,walu2,i0,i0)
296  if(idof1.ne.idof2) then
297  walue=value*cs(16,icomplex)
298  walu2=valu2*cs(16,icomplex)
299  call add_sm_ei(au,ad,aub,adb,jq,irow,
300  & idof1,idof2+ner,walue,walu2,i0,i0)
301  walue=-walue
302  walu2=-walu2
303  call add_sm_ei(au,ad,aub,adb,jq,irow,
304  & idof1+ner,idof2,walue,walu2,i0,i0)
305  endif
306  endif
307  endif
308  index=nodempc(3,index)
309  if(index.eq.0) exit
310  enddo
311  cycle
312  endif
313  endif
314 !
315  else
316 c idof1=(node1-1)*8+k
317 c idof2=(node2-1)*8+m
318  idof1=jdof1
319  idof2=jdof2
320 !
321  mpc1=0
322  mpc2=0
323  if(nmpc.gt.0) then
324 c call nident(ikmpc,idof1,nmpc,id1)
325 c if((id1.gt.0).and.(ikmpc(id1).eq.idof1)) mpc1=1
326 c call nident(ikmpc,idof2,nmpc,id2)
327 c if((id2.gt.0).and.(ikmpc(id2).eq.idof2)) mpc2=1
328  if(idof1.ne.2*(idof1/2)) mpc1=1
329  if(idof2.ne.2*(idof2/2)) mpc2=1
330  endif
331  if((mpc1.eq.1).and.(mpc2.eq.1)) then
332 c id1=ilmpc(id1)
333 c id2=ilmpc(id2)
334  id1=(-idof1+1)/2
335  id2=(-idof2+1)/2
336  if(id1.eq.id2) then
337 !
338 ! MPC id1 / MPC id1
339 !
340  ist=ipompc(id1)
341  index1=nodempc(3,ist)
342  if(index1.eq.0) cycle
343  do
344  inode1=nodempc(1,index1)
345  icomplex1=0
346  if(labmpc(id1)(1:6).eq.'CYCLIC') then
347  read(labmpc(id1)(7:20),'(i14)') icomplex1
348  elseif(labmpc(id1)(1:9).eq.'SUBCYCLIC') then
349  do ij=1,mcs
350  ilength=int(cs(4,ij))
351  lprev=int(cs(14,ij))
352  call nident(ics(lprev+1),inode1,
353  & ilength,id)
354  if(id.gt.0) then
355  if(ics(lprev+id).eq.inode1) then
356  icomplex1=ij
357  exit
358  endif
359  endif
360  enddo
361  endif
362  idof1=nactdof(nodempc(2,index1),inode1)
363  index2=index1
364  do
365  inode2=nodempc(1,index2)
366  icomplex2=0
367  if(labmpc(id1)(1:6).eq.'CYCLIC') then
368  read(labmpc(id1)(7:20),'(i14)') icomplex2
369  elseif(labmpc(id1)(1:9).eq.'SUBCYCLIC') then
370  do ij=1,mcs
371  ilength=int(cs(4,ij))
372  lprev=int(cs(14,ij))
373  call nident(ics(lprev+1),inode2,
374  & ilength,id)
375  if(id.gt.0) then
376  if(ics(lprev+id).eq.inode2) then
377  icomplex2=ij
378  exit
379  endif
380  endif
381  enddo
382  endif
383  idof2=nactdof(nodempc(2,index2),inode2)
384  if((idof1.gt.0).and.(idof2.gt.0)) then
385  value=coefmpc(index1)*coefmpc(index2)*
386  & s(jj,ll)/coefmpc(ist)/coefmpc(ist)
387  valu2=coefmpc(index1)*coefmpc(index2)*
388  & sm(jj,ll)/coefmpc(ist)/coefmpc(ist)
389  if((icomplex1.eq.0).and.
390  & (icomplex2.eq.0)) then
391  call add_sm_ei(au,ad,aub,adb,jq,
392  & irow,idof1,idof2,value,valu2,i0,i0)
393  call add_sm_ei(au,ad,aub,adb,jq,
394  & irow,idof1+ner,idof2+ner,value,
395  & valu2,i0,i0)
396  elseif((icomplex1.ne.0).and.
397  & (icomplex2.ne.0)) then
398  if(icomplex1.eq.icomplex2) then
399  call add_sm_ei(au,ad,aub,adb,jq,
400  & irow,idof1,idof2,value,valu2,i0,i0)
401  call add_sm_ei(au,ad,aub,adb,jq,
402  & irow,idof1+ner,idof2+ner,value,
403  & valu2,i0,i0)
404  else
405  tr=cs(15,icomplex1)*cs(15,icomplex2)
406  & +cs(16,icomplex1)*cs(16,icomplex2)
407 c write(*,*) 'tr= ',tr
408  walue=value*tr
409  walu2=valu2*tr
410  call add_sm_ei(au,ad,aub,adb,jq,
411  & irow,idof1,idof2,walue,walu2,i0,i0)
412  call add_sm_ei(au,ad,aub,adb,jq,
413  & irow,idof1+ner,idof2+ner,walue,
414  & walu2,i0,i0)
415  ti=cs(15,icomplex1)*cs(16,icomplex2)
416  & -cs(15,icomplex2)*cs(16,icomplex1)
417 c write(*,*) icomplex1,icomplex2,
418 c & cs(15,icomplex1),cs(16,icomplex1),
419 c & cs(15,icomplex2),cs(16,icomplex2)
420 c write(*,*) 'ti= ',ti
421  walue=value*ti
422  walu2=valu2*ti
423 c write(*,'(2i8,2(1x,e11.4))')
424 c & idof1,idof2+ner,
425 c & walue,walu2
426  call add_sm_ei(au,ad,aub,adb,jq,irow
427  & ,idof1,idof2+ner,walue,walu2,i0,i0)
428  walue=-walue
429  walu2=-walu2
430  call add_sm_ei(au,ad,aub,adb,jq,irow
431  & ,idof1+ner,idof2,walue,walu2,i0,i0)
432  endif
433  elseif((icomplex1.eq.0).or.
434  & (icomplex2.eq.0)) then
435  if(icomplex2.ne.0) then
436  walue=value*cs(15,icomplex2)
437  walu2=valu2*cs(15,icomplex2)
438  else
439  walue=value*cs(15,icomplex1)
440  walu2=valu2*cs(15,icomplex1)
441  endif
442  call add_sm_ei(au,ad,aub,adb,jq,irow,
443  & idof1,idof2,walue,walu2,i0,i0)
444  call add_sm_ei(au,ad,aub,adb,jq,irow,
445  & idof1+ner,idof2+ner,walue,walu2,i0,i0)
446  if(icomplex2.ne.0) then
447  walue=value*cs(16,icomplex2)
448  walu2=valu2*cs(16,icomplex2)
449  else
450  walue=-value*cs(16,icomplex1)
451  walu2=-valu2*cs(16,icomplex1)
452  endif
453 c walue=value*st
454 c walu2=valu2*st
455  call add_sm_ei(au,ad,aub,adb,jq,irow,
456  & idof1,idof2+ner,walue,walu2,i0,i0)
457  walue=-walue
458  walu2=-walu2
459  call add_sm_ei(au,ad,aub,adb,jq,irow,
460  & idof1+ner,idof2,walue,walu2,i0,i0)
461  endif
462  endif
463  index2=nodempc(3,index2)
464  if(index2.eq.0) exit
465  enddo
466  index1=nodempc(3,index1)
467  if(index1.eq.0) exit
468  enddo
469  else
470 !
471 ! MPC id1 / MPC id2
472 !
473  ist1=ipompc(id1)
474  index1=nodempc(3,ist1)
475  if(index1.eq.0) cycle
476  do
477  inode1=nodempc(1,index1)
478  icomplex1=0
479  if(labmpc(id1)(1:6).eq.'CYCLIC') then
480  read(labmpc(id1)(7:20),'(i14)') icomplex1
481  elseif(labmpc(id1)(1:9).eq.'SUBCYCLIC') then
482  do ij=1,mcs
483  ilength=int(cs(4,ij))
484  lprev=int(cs(14,ij))
485  call nident(ics(lprev+1),inode1,
486  & ilength,id)
487  if(id.gt.0) then
488  if(ics(lprev+id).eq.inode1) then
489  icomplex1=ij
490  exit
491  endif
492  endif
493  enddo
494  endif
495  idof1=nactdof(nodempc(2,index1),inode1)
496  ist2=ipompc(id2)
497  index2=nodempc(3,ist2)
498  if(index2.eq.0) then
499  index1=nodempc(3,index1)
500  if(index1.eq.0) then
501  exit
502  else
503  cycle
504  endif
505  endif
506  do
507  inode2=nodempc(1,index2)
508  icomplex2=0
509  if(labmpc(id2)(1:6).eq.'CYCLIC') then
510  read(labmpc(id2)(7:20),'(i14)') icomplex2
511  elseif(labmpc(id2)(1:9).eq.'SUBCYCLIC') then
512  do ij=1,mcs
513  ilength=int(cs(4,ij))
514  lprev=int(cs(14,ij))
515  call nident(ics(lprev+1),inode2,
516  & ilength,id)
517  if(id.gt.0) then
518  if(ics(lprev+id).eq.inode2) then
519  icomplex2=ij
520  exit
521  endif
522  endif
523  enddo
524  endif
525  idof2=nactdof(nodempc(2,index2),inode2)
526  if((idof1.gt.0).and.(idof2.gt.0)) then
527  value=coefmpc(index1)*coefmpc(index2)*
528  & s(jj,ll)/coefmpc(ist1)/coefmpc(ist2)
529  valu2=coefmpc(index1)*coefmpc(index2)*
530  & sm(jj,ll)/coefmpc(ist1)/coefmpc(ist2)
531  if(idof1.eq.idof2) then
532  value=2.d0*value
533  valu2=2.d0*valu2
534  endif
535  if((icomplex1.eq.0).and.
536  & (icomplex2.eq.0)) then
537  call add_sm_ei(au,ad,aub,adb,jq,
538  & irow,idof1,idof2,value,valu2,i0,i0)
539  call add_sm_ei(au,ad,aub,adb,jq,
540  & irow,idof1+ner,idof2+ner,value,
541  & valu2,i0,i0)
542  elseif((icomplex1.ne.0).and.
543  & (icomplex2.ne.0)) then
544  if(icomplex1.eq.icomplex2) then
545  call add_sm_ei(au,ad,aub,adb,jq,
546  & irow,idof1,idof2,value,valu2,i0,i0)
547  call add_sm_ei(au,ad,aub,adb,jq,
548  & irow,idof1+ner,idof2+ner,value,
549  & valu2,i0,i0)
550  else
551  tr=cs(15,icomplex1)*cs(15,icomplex2)
552  & +cs(16,icomplex1)*cs(16,icomplex2)
553 c write(*,*) 'tr= ',tr
554  walue=value*tr
555  walu2=valu2*tr
556  call add_sm_ei(au,ad,aub,adb,jq,
557  & irow,idof1,idof2,walue,walu2,i0,i0)
558  call add_sm_ei(au,ad,aub,adb,jq,
559  & irow,idof1+ner,idof2+ner,walue,
560  & walu2,i0,i0)
561  ti=cs(15,icomplex1)*cs(16,icomplex2)
562  & -cs(15,icomplex2)*cs(16,icomplex1)
563 c write(*,*) icomplex1,icomplex2,
564 c & cs(15,icomplex1),cs(16,icomplex1),
565 c & cs(15,icomplex2),cs(16,icomplex2)
566 c write(*,*) 'ti= ',ti
567  walue=value*ti
568  walu2=valu2*ti
569 c write(*,'(2i8,2(1x,e11.4))')
570 c & idof1,idof2+ner,
571 c & walue,walu2
572  call add_sm_ei(au,ad,aub,adb,jq,irow
573  & ,idof1,idof2+ner,walue,walu2,i0,i0)
574  walue=-walue
575  walu2=-walu2
576  call add_sm_ei(au,ad,aub,adb,jq,irow
577  & ,idof1+ner,idof2,walue,walu2,i0,i0)
578  endif
579  elseif((icomplex1.eq.0).or.
580  & (icomplex2.eq.0)) then
581  if(icomplex2.ne.0) then
582  walue=value*cs(15,icomplex2)
583  walu2=valu2*cs(15,icomplex2)
584  else
585  walue=value*cs(15,icomplex1)
586  walu2=valu2*cs(15,icomplex1)
587  endif
588  call add_sm_ei(au,ad,aub,adb,jq,irow,
589  & idof1,idof2,walue,walu2,i0,i0)
590  call add_sm_ei(au,ad,aub,adb,jq,irow,
591  & idof1+ner,idof2+ner,walue,walu2,i0,i0)
592  if(idof1.ne.idof2) then
593  if(icomplex2.ne.0) then
594  walue=value*cs(16,icomplex2)
595  walu2=valu2*cs(16,icomplex2)
596  else
597  walue=-value*cs(16,icomplex1)
598  walu2=-valu2*cs(16,icomplex1)
599  endif
600 c walue=value*st
601 c walu2=valu2*st
602  call add_sm_ei(au,ad,aub,adb,jq,
603  & irow,idof1,idof2+ner,walue,
604  & walu2,i0,i0)
605  walue=-walue
606  walu2=-walu2
607  call add_sm_ei(au,ad,aub,adb,jq,
608  & irow,idof1+ner,idof2,walue,
609  & walu2,i0,i0)
610  endif
611  endif
612  endif
613  index2=nodempc(3,index2)
614  if(index2.eq.0) exit
615  enddo
616  index1=nodempc(3,index1)
617  if(index1.eq.0) exit
618  enddo
619  endif
620  endif
621  endif
622  enddo
623 !
624  enddo
625  enddo
626 !
627  return
e_c3d
subroutine e_c3d(co, kon, lakonl, p1, p2, omx, bodyfx, nbody, s, sm, ff, nelem, nmethod, elcon, nelcon, rhcon, nrhcon, alcon, nalcon, alzero, ielmat, ielorien, norien, orab, ntmat_, t0, t1, ithermal, vold, iperturb, nelemload, sideload, xload, nload, idist, sti, stx, iexpl, plicon, nplicon, plkcon, nplkcon, xstiff, npmat_, dtime, matname, mi, ncmat_, mass, stiffness, buckling, rhsi, intscheme, ttime, time, istep, iinc, coriolis, xloadold, reltime, ipompc, nodempc, coefmpc, nmpc, ikmpc, ilmpc, veold, springarea, nstate_, xstateini, xstate, ne0, ipkon, thicke, integerglob, doubleglob, tieset, istartset, iendset, ialset, ntie, nasym, pslavsurf, pmastsurf, mortar, clearini, ielprop, prop, kscale)
Definition: e_c3d.f:31
add_sm_ei
subroutine add_sm_ei(au, ad, aub, adb, jq, irow, i, j, value, valuem, i0, i1)
Definition: add_sm_ei.f:21
idist
static ITG * idist
Definition: radflowload.c:39
nident
subroutine nident(x, px, n, id)
Definition: nident.f:26
Hosted by OpenAircraft.com, (Michigan UAV, LLC)