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

Go to the source code of this file.

Functions/Subroutines

subroutine mafillem (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, fext, 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, ncmat_, mass, stiffness, buckling, rhsi, intscheme, physcon, shcon, nshcon, cocon, ncocon, ttime, time, istep, iinc, coriolis, ibody, xloadold, reltime, veold, springarea, nstate_, xstateini, xstate, thicke, integerglob, doubleglob, tieset, istartset, iendset, ialset, ntie, nasym, iactive, h0, pslavsurf, pmastsurf, mortar, clearini, ielprop, prop, iponoel, inoel, network)
 

Function/Subroutine Documentation

◆ mafillem()

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