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

Go to the source code of this file.

Functions/Subroutines

subroutine elements (inpc, textpart, kon, ipkon, lakon, nkon, ne, ne_, set, istartset, iendset, ialset, nset, nset_, nalset, nalset_, mi, ixfree, iponor, xnor, istep, istat, n, iline, ipol, inl, ipoinp, inp, iaxial, ipoinpc, solid, cfd, network, filab, nlabel, out3d, iuel, nuel_)
 

Function/Subroutine Documentation

◆ elements()

subroutine elements ( character*1, dimension(*)  inpc,
character*132, dimension(16)  textpart,
integer, dimension(*)  kon,
integer, dimension(*)  ipkon,
character*8, dimension(*)  lakon,
integer  nkon,
integer  ne,
integer  ne_,
character*81, dimension(*)  set,
integer, dimension(*)  istartset,
integer, dimension(*)  iendset,
integer, dimension(*)  ialset,
integer  nset,
integer  nset_,
integer  nalset,
integer  nalset_,
integer, dimension(*)  mi,
integer  ixfree,
integer, dimension(2,*)  iponor,
real*8, dimension(*)  xnor,
integer  istep,
integer  istat,
integer  n,
integer  iline,
integer  ipol,
integer  inl,
integer, dimension(2,*)  ipoinp,
integer, dimension(3,*)  inp,
integer  iaxial,
integer, dimension(0:*)  ipoinpc,
logical  solid,
integer  cfd,
integer  network,
character*87, dimension(*)  filab,
integer  nlabel,
logical  out3d,
integer, dimension(4,*)  iuel,
integer  nuel_ 
)
24 !
25 ! reading the input deck: *ELEMENT
26 !
27  implicit none
28 !
29  logical solid,beamshell,out3d
30 !
31  character*1 inpc(*)
32  character*8 lakon(*),label
33  character*81 set(*),elset
34  character*87 filab(*)
35  character*132 textpart(16)
36 !
37  integer kon(*),istartset(*),iendset(*),ialset(*),ne,ne_,nset,
38  & nset_,nalset,nalset_,istep,istat,n,key,i,ielset,js,k,nn,
39  & nteller,j,ipkon(*),nkon,nope,indexe,mi(*),ipos,indexy,ixfree,
40  & iponor(2,*),nopeexp,iline,ipol,inl,ipoinp(2,*),inp(3,*),
41  & iaxial,ipoinpc(0:*),cfd,nlabel,network,iuel(4,*),nuel_,
42  & id,four,number,ndof,intpoints
43 !
44  real*8 xnor(*)
45 !
46  if(istep.gt.0) then
47  write(*,*) '*ERROR reading *ELEMENT: *ELEMENT should be placed'
48  write(*,*) ' before all step definitions'
49  call exit(201)
50  endif
51 !
52  indexy=-1
53  ielset=0
54  beamshell=.false.
55  four=4
56 !
57  label=' '
58 !
59 ! checking for set definition
60 !
61  loop: do i=2,n
62  if(textpart(i)(1:6).eq.'ELSET=') then
63  elset=textpart(i)(7:86)
64  if(textpart(i)(87:87).ne.' ') then
65  write(*,*) '*ERROR reading *ELEMENT: set name too long'
66  write(*,*) ' (more than 80 characters)'
67  write(*,*) ' set name:',textpart(i)(1:132)
68  call exit(201)
69  endif
70  elset(81:81)=' '
71  ipos=index(elset,' ')
72  elset(ipos:ipos)='E'
73  ielset=1
74  do js=1,nset
75  if(set(js).eq.elset) then
76 !
77 ! existent set
78 !
79  if(iendset(js).eq.nalset) then
80  cycle loop
81  else
82 !
83 ! rearranging set information towards the end
84 !
85  nn=iendset(js)-istartset(js)+1
86  if(nalset+nn.gt.nalset_) then
87  write(*,*)
88  & '*ERROR reading *ELEMENT: increase nalset_'
89  call exit(201)
90  endif
91  do k=1,nn
92  ialset(nalset+k)=ialset(istartset(js)+k-1)
93  enddo
94  do k=istartset(js),nalset
95  ialset(k)=ialset(k+nn)
96  enddo
97  do k=1,nset
98  if(istartset(k).gt.iendset(js)) then
99  istartset(k)=istartset(k)-nn
100  iendset(k)=iendset(k)-nn
101  endif
102  enddo
103  istartset(js)=nalset-nn+1
104  iendset(js)=nalset
105  cycle loop
106  endif
107  endif
108  enddo
109 !
110 ! new set
111 !
112  nset=nset+1
113  if(nset.gt.nset_) then
114  write(*,*) '*ERROR reading *ELEMENT: increase nset_'
115  call exit(201)
116  endif
117  js=nset
118  istartset(js)=nalset+1
119  iendset(js)=nalset
120  set(js)=elset
121  cycle
122  elseif(textpart(i)(1:5).eq.'TYPE=') then
123  read(textpart(i)(6:13),'(a8)') label
124 !
125 ! removing the ABAQUS label for heat transfer elements
126 !
127  if((label(1:2).eq.'DC').and.(label(1:7).ne.'DCOUP3D')) then
128  label(1:7)=label(2:8)
129  label(8:8)=' '
130 !
131 ! mapping 2D beam and truss elements onto 3D elements
132 !
133  elseif(label(1:3).eq.'B21') then
134  label='B31 '
135  elseif(label(1:4).eq.'T2D2') then
136  label='T3D2 '
137  endif
138 !
139 ! full integration quadratic hexahedral element
140 ! (including such which are expanded into one)
141 !
142  if((label.eq.'C3D20 ').or.
143  & (label.eq.'CPE8 ').or.
144  & (label.eq.'CPS8 ').or.
145  & (label.eq.'CAX8 ').or.
146  & (label.eq.'M3D8 ').or.
147  & (label.eq.'S8 ').or.
148  & (label.eq.'T3D3 ').or.
149  & (label.eq.'B32 ').or.
150 !
151 ! reduced integration quadratic hexahedral element
152 ! (including such which are expanded into one)
153 !
154  & (label.eq.'C3D20R ').or.
155  & (label.eq.'CPE8R ').or.
156  & (label.eq.'CPS8R ').or.
157  & (label.eq.'CAX8R ').or.
158  & (label.eq.'M3D8R ').or.
159  & (label.eq.'S8R ').or.
160  & (label.eq.'B32R ').or.
161 !
162 ! full integration linear hexahedral element
163 ! (including such which are expanded into one)
164 !
165  & (label.eq.'C3D8 ').or.
166  & (label.eq.'CPE4 ').or.
167  & (label.eq.'CPS4 ').or.
168  & (label.eq.'CAX4 ').or.
169  & (label.eq.'M3D4 ').or.
170  & (label.eq.'S4 ').or.
171  & (label.eq.'T3D2 ').or.
172  & (label.eq.'B31 ').or.
173 !
174 ! reduced integration linear hexahedral element
175 ! (including such which are expanded into one)
176 !
177  & (label.eq.'C3D8R ').or.
178  & (label.eq.'CPE4R ').or.
179  & (label.eq.'CPS4R ').or.
180  & (label.eq.'CAX4R ').or.
181  & (label.eq.'M3D4R ').or.
182  & (label.eq.'S4R ').or.
183  & (label.eq.'B31R ').or.
184 c Bernhardi start
185 c
186 c incompatible modes element
187 c
188  & (label.eq.'C3D8I ').or.
189 c Bernhardi end
190 !
191 ! quadratic tetrahedral element
192 !
193  & (label.eq.'C3D10 ').or.
194  & (label.eq.'C3D10T ').or.
195 !
196 ! linear tetrahedral element
197 !
198  & (label.eq.'C3D4 ').or.
199 !
200 ! quadratic wedge
201 ! (including such which are expanded into one)
202 !
203  & (label.eq.'C3D15 ').or.
204  & (label.eq.'CPE6 ').or.
205  & (label.eq.'CPS6 ').or.
206  & (label.eq.'CAX6 ').or.
207  & (label.eq.'M3D6 ').or.
208  & (label.eq.'S6 ').or.
209 !
210 ! linear wedge
211 ! (including such which are expanded into one)
212 !
213  & (label.eq.'C3D6 ').or.
214  & (label.eq.'CPE3 ').or.
215  & (label.eq.'CPS3 ').or.
216  & (label.eq.'CAX3 ').or.
217  & (label.eq.'M3D3 ').or.
218  & (label.eq.'S3 ').or.
219 !
220 ! gap element
221 !
222  & (label.eq.'GAPUNI ').or.
223 !
224 ! spring element
225 !
226  & (label.eq.'SPRINGA ').or.
227  & (label.eq.'SPRING1 ').or.
228  & (label.eq.'SPRING2 ').or.
229 !
230 ! dashpot element
231 !
232  & (label.eq.'DASHPOTA'))
233  & then
234  solid=.true.
235 !
236 ! 3D fluid element
237 !
238  elseif((label.eq.'F3D8 ').or.
239  & (label.eq.'F3D4 ').or.
240  & (label.eq.'F3D6 ')) then
241  cfd=1
242 !
243 ! distributing coupling element
244 !
245  elseif(label(1:7).eq.'DCOUP3D') then
246 !
247 ! network element
248 !
249  elseif(label(1:1).eq.'D') then
250  network=3
251 !
252 ! mass element
253 !
254  elseif(label(1:4).eq.'MASS') then
255 !
256 ! user element
257 !
258  elseif(label(1:1).eq.'U') then
259 !
260 ! unknown element type
261 !
262  else
263  write(*,*) '*ERROR reading *ELEMENT:'
264  write(*,*) label,' is an unknown element type'
265  call exit(201)
266  endif
267 !
268  if(label(1:3).eq.'CAX') then
269  iaxial=180
270  elseif((label(1:2).eq.'S3').or.
271  & (label(1:2).eq.'S4').or.
272  & (label(1:2).eq.'S6').or.
273  & (label(1:2).eq.'S8').or.
274  & (label(1:1).eq.'B')) then
275  beamshell=.true.
276  endif
277 !
278  else
279  write(*,*)
280  & '*WARNING reading *ELEMENT: parameter not recognized:'
281  write(*,*) ' ',
282  & textpart(i)(1:index(textpart(i),' ')-1)
283  call inputwarning(inpc,ipoinpc,iline,
284  &"*ELEMENT%")
285  endif
286  enddo loop
287 !
288  if(label.eq.' ') then
289  write(*,*) '*ERROR reading *ELEMENT: element type is lacking'
290  write(*,*) ' '
291  call inputerror(inpc,ipoinpc,iline,
292  &"*ELEMENT%")
293  call exit(201)
294  endif
295 !
296 ! nope is the number of nodes per element as defined in the input
297 ! deck, nopeexp is the number of nodes per element after expansion
298 ! (only applicable to 1D and 2D elements such as beams, shells..)
299 !
300 c Bernhardi start
301  if(label(1:5).eq.'C3D8I') then
302  nope=8
303  nopeexp=11
304  elseif(label(4:5).eq.'20') then
305 c Bernhardi end
306  nope=20
307  nopeexp=20
308  elseif((label(1:4).eq.'CPE8').or.(label(1:4).eq.'CPS8').or.
309  & (label(1:4).eq.'CAX8').or.(label(1:2).eq.'S8').or.
310  & (label(1:4).eq.'M3D8')) then
311  nope=8
312  nopeexp=28
313  elseif((label(1:4).eq.'CPE6').or.(label(1:4).eq.'CPS6').or.
314  & (label(1:4).eq.'CAX6').or.(label(1:2).eq.'S6').or.
315  & (label(1:4).eq.'M3D6')) then
316  nope=6
317  nopeexp=21
318  elseif((label(1:3).eq.'B32').or.
319  & (label(1:4).eq.'T3D3')) then
320  nope=3
321  nopeexp=23
322  elseif((label(1:4).eq.'B31 ').or.
323  & (label(1:4).eq.'T3D2')) then
324  nope=2
325 ! expanded into C3D8I: 11 nodes per element
326  nopeexp=13
327  elseif(label(1:4).eq.'B31R') then
328  nope=2
329  nopeexp=10
330  elseif(label(4:4).eq.'8') then
331  nope=8
332  nopeexp=8
333 c elseif((label(1:5).eq.'CPE4 ').or.(label(1:5).eq.'CPS4 ').or.
334 c & (label(1:5).eq.'CAX4 ').or.(label(1:3).eq.'S4 ')) then
335  elseif(label(1:3).eq.'S4 ') then
336  nope=4
337 ! expanded into C3D8I: 11 nodes per element
338  nopeexp=15
339  elseif((label(1:4).eq.'CPE4').or.(label(1:4).eq.'CPS4').or.
340  & (label(1:4).eq.'CAX4').or.(label(1:2).eq.'S4').or.
341  & (label(1:4).eq.'M3D4')) then
342  nope=4
343  nopeexp=12
344  elseif(label(4:5).eq.'10') then
345  nope=10
346  nopeexp=10
347  elseif(label(4:4).eq.'4') then
348  nope=4
349  nopeexp=4
350  elseif(label(4:5).eq.'15') then
351  nope=15
352  nopeexp=15
353  elseif(label(4:4).eq.'6') then
354  nope=6
355  nopeexp=6
356  elseif((label(1:4).eq.'CPE3').or.(label(1:4).eq.'CPS3').or.
357  & (label(1:4).eq.'CAX3').or.(label(1:2).eq.'S3').or.
358  & (label(1:4).eq.'M3D3')) then
359  nope=3
360  nopeexp=9
361  elseif(label(1:8).eq.'DASHPOTA') then
362  label='EDSHPTA1'
363  nope=2
364  nopeexp=2
365  elseif(label(1:7).eq.'DCOUP3D') then
366  nope=1
367  nopeexp=1
368  elseif(label(1:1).eq.'D') then
369  nope=3
370  nopeexp=3
371  elseif(label(1:1).eq.'G') then
372  label='ESPGAPA1'
373  nope=2
374  nopeexp=2
375  elseif(label(1:7).eq.'SPRINGA') then
376  label='ESPRNGA1'
377  nope=2
378  nopeexp=2
379  elseif(label(1:7).eq.'SPRING1') then
380  label='ESPRNG10'
381  nope=1
382  nopeexp=1
383  elseif(label(1:7).eq.'SPRING2') then
384  label='ESPRNG21'
385  nope=2
386  nopeexp=2
387  elseif(label(1:4).eq.'MASS') then
388  nope=1
389  nopeexp=1
390  elseif(label(1:1).eq.'U') then
391  number=ichar(label(2:2))*256**3+
392  & ichar(label(3:3))*256**2+
393  & ichar(label(4:4))*256+
394  & ichar(label(5:5))
395 c read(label(2:3),'(i2)',iostat=istat) number
396  call nidentk(iuel,number,nuel_,id,four)
397  intpoints=iuel(2,id)
398  ndof=iuel(3,id)
399  nope=iuel(4,id)
400  nopeexp=nope
401  write(label(6:6),'(a1)') char(intpoints)
402  write(label(7:7),'(a1)') char(ndof)
403  write(label(8:8),'(a1)') char(nope)
404  endif
405 !
406  do
407  call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl,
408  & ipoinp,inp,ipoinpc)
409  if((istat.lt.0).or.(key.eq.1)) then
410 !
411 ! default for frd output of beams and shell is the 3d-expansion
412 !
413  if(beamshell) then
414  do j=1,nlabel
415  filab(j)(5:5)='E'
416  enddo
417  out3d=.true.
418  endif
419  return
420  endif
421  read(textpart(1)(1:10),'(i10)',iostat=istat) i
422  if(istat.gt.0) call inputerror(inpc,ipoinpc,iline,
423  &"*ELEMENT%")
424  if(i.gt.ne_) then
425  write(*,*) '*ERROR reading *ELEMENT: increase ne_'
426  call exit(201)
427  endif
428 !
429 ! check whether element was already defined
430 !
431  if(ipkon(i).ne.-1) then
432  write(*,*) '*ERROR reading *ELEMENT: element',i
433  write(*,*) ' is already defined'
434  write(*,*) ' '
435  call inputerror(inpc,ipoinpc,iline,
436  &"*ELEMENT%")
437  endif
438 !
439 ! new element
440 !
441  ipkon(i)=nkon
442  lakon(i)=label
443  indexe=nkon
444 !
445  nkon=nkon+nopeexp
446 !
447  do j=2,min(n,nope+1)
448  read(textpart(j)(1:10),'(i10)',iostat=istat) kon(indexe+j-1)
449  if(istat.gt.0) call inputerror(inpc,ipoinpc,iline,
450  &"*ELEMENT%")
451  enddo
452  nteller=n-1
453  if(nteller.lt.nope) then
454  do
455  call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl,
456  & ipoinp,inp,ipoinpc)
457  if((istat.lt.0).or.(key.eq.1)) then
458  write(*,*)
459  & '*ERROR reading *ELEMENT: element definition'
460  write(*,*) ' incomplete for element ',i
461  call exit(201)
462  endif
463  if(nteller+n.gt.nope) n=nope-nteller
464  do j=1,n
465  read(textpart(j)(1:10),'(i10)',iostat=istat)
466  & kon(indexe+nteller+j)
467  if(istat.gt.0) call inputerror(inpc,ipoinpc,iline,
468  &"*ELEMENT%")
469  enddo
470  nteller=nteller+n
471  if(nteller.eq.nope) exit
472  enddo
473  endif
474  ne=max(ne,i)
475 !
476 ! assigning element to set
477 !
478  if(ielset.eq.1) then
479  if(nalset+1.gt.nalset_) then
480  write(*,*) '*ERROR reading *ELEMENT: increase nalset_'
481  call exit(201)
482  endif
483  nalset=nalset+1
484  ialset(nalset)=i
485  iendset(js)=nalset
486  endif
487 !
488 ! for plane stress, plane strain and axisymmetric elements:
489 ! define the normal
490 !
491  if((label(1:2).eq.'CP').or.(label(1:2).eq.'CA')) then
492  if(indexy.eq.-1) then
493  indexy=ixfree
494  xnor(indexy+1)=0.d0
495  xnor(indexy+2)=0.d0
496  xnor(indexy+3)=1.d0
497  ixfree=ixfree+3
498  endif
499  do j=1,nope
500  iponor(1,indexe+j)=indexy
501  enddo
502  endif
503 !
504  enddo
505 !
506  return
max
#define max(a, b)
Definition: cascade.c:32
inputwarning
subroutine inputwarning(inpc, ipoinpc, iline, text)
Definition: inputwarning.f:20
min
#define min(a, b)
Definition: cascade.c:31
getnewline
subroutine getnewline(inpc, textpart, istat, n, key, iline, ipol, inl, ipoinp, inp, ipoinpc)
Definition: getnewline.f:21
inputerror
subroutine inputerror(inpc, ipoinpc, iline, text)
Definition: inputerror.f:20
nidentk
subroutine nidentk(x, px, n, id, k)
Definition: nidentk.f:27
Hosted by OpenAircraft.com, (Michigan UAV, LLC)