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

Go to the source code of this file.

Functions/Subroutines

subroutine solidsections (inpc, textpart, set, istartset, iendset, ialset, nset, ielmat, matname, nmat, ielorien, orname, norien, lakon, thicke, kon, ipkon, irstrt, istep, istat, n, iline, ipol, inl, ipoinp, inp, cs, mcs, iaxial, ipoinpc, mi, co, ixfree, xnor, iponor, nelcon)
 

Function/Subroutine Documentation

◆ solidsections()

subroutine solidsections ( character*1, dimension(*)  inpc,
character*132, dimension(16)  textpart,
character*81, dimension(*)  set,
integer, dimension(*)  istartset,
integer, dimension(*)  iendset,
integer, dimension(*)  ialset,
integer  nset,
integer, dimension(mi(3),*)  ielmat,
character*80, dimension(*)  matname,
integer  nmat,
integer, dimension(mi(3),*)  ielorien,
character*80, dimension(*)  orname,
integer  norien,
character*8, dimension(*)  lakon,
real*8, dimension(mi(3),*)  thicke,
integer, dimension(*)  kon,
integer, dimension(*)  ipkon,
integer  irstrt,
integer  istep,
integer  istat,
integer  n,
integer  iline,
integer  ipol,
integer  inl,
integer, dimension(2,*)  ipoinp,
integer, dimension(3,*)  inp,
real*8, dimension(17,*)  cs,
integer  mcs,
integer  iaxial,
integer, dimension(0:*)  ipoinpc,
integer, dimension(*)  mi,
real*8, dimension(3,*)  co,
integer  ixfree,
real*8, dimension(*)  xnor,
integer, dimension(2,*)  iponor,
integer, dimension(2,*)  nelcon 
)
24 !
25 ! reading the input deck: *SOLID SECTION
26 !
27  implicit none
28 !
29  character*1 inpc(*)
30  character*8 lakon(*)
31  character*80 matname(*),orname(*),material,orientation
32  character*81 set(*),elset
33  character*132 textpart(16)
34 !
35  integer mi(*),istartset(*),iendset(*),ialset(*),ielmat(mi(3),*),
36  & ielorien(mi(3),*),kon(*),ipkon(*),indexe,irstrt,nset,nmat,
37  & norien,ielem,node1,node2,m,indexx,ixfree,iponor(2,*),
38  & istep,istat,n,key,i,j,k,l,imaterial,iorientation,ipos,
39  & iline,ipol,inl,ipoinp(2,*),inp(3,*),mcs,iaxial,ipoinpc(0:*),
40  & nelcon(2,*)
41 !
42  real*8 thicke(mi(3),*),thickness,pi,cs(17,*),xn(3),co(3,*),p(3),
43  & dd,xnor(*)
44 !
45  if((istep.gt.0).and.(irstrt.ge.0)) then
46  write(*,*) '*ERROR reading *SOLID SECTION: *SOLID SECTION'
47  write(*,*)' should be placed before all step definitions'
48  call exit(201)
49  endif
50 !
51  pi=4.d0*datan(1.d0)
52 !
53  orientation='
54  & '
55  elset='
56  & '
57  ipos=0
58 !
59  do i=2,n
60  if(textpart(i)(1:9).eq.'MATERIAL=') then
61  material=textpart(i)(10:89)
62  elseif(textpart(i)(1:12).eq.'ORIENTATION=') then
63  orientation=textpart(i)(13:92)
64  elseif(textpart(i)(1:6).eq.'ELSET=') then
65  elset=textpart(i)(7:86)
66  elset(81:81)=' '
67  ipos=index(elset,' ')
68  elset(ipos:ipos)='E'
69  else
70  write(*,*)
71  & '*WARNING reading *SOLID SECTION: parameter not recognized:'
72  write(*,*) ' ',
73  & textpart(i)(1:index(textpart(i),' ')-1)
74  call inputwarning(inpc,ipoinpc,iline,
75  &"*SOLID SECTION%")
76  endif
77  enddo
78 !
79 ! check for the existence of the set,the material and orientation
80 !
81  do i=1,nmat
82  if(matname(i).eq.material) exit
83  enddo
84  if(i.gt.nmat) then
85  do i=1,nmat
86  if(matname(i)(1:11).eq.'ANISO_CREEP') then
87  if(matname(i)(12:20).eq.material(1:9)) exit
88  elseif(matname(i)(1:10).eq.'ANISO_PLAS') then
89  if(matname(i)(11:20).eq.material(1:10)) exit
90  endif
91  enddo
92  endif
93  if(i.gt.nmat) then
94  write(*,*)
95  & '*ERROR reading *SOLID SECTION: nonexistent material'
96  write(*,*) ' '
97  call inputerror(inpc,ipoinpc,iline,
98  &"*SOLID SECTION%")
99  call exit(201)
100  endif
101  imaterial=i
102 !
103  if(orientation.eq.' ') then
104  iorientation=0
105 c
106 c next lines were removed since some people use local systems
107 c to the the stresses stored in that system (e.g. cylindrical
108 c system)
109 c
110 c elseif(nelcon(1,i).eq.2) then
111 c write(*,*) '*INFO reading *SOLID SECTION: an orientation'
112 c write(*,*) ' is for isotropic materials irrelevant'
113 c call inputinfo(inpc,ipoinpc,iline,
114 c &"*SOLID SECTION%")
115 c iorientation=0
116  else
117  do i=1,norien
118  if(orname(i).eq.orientation) exit
119  enddo
120  if(i.gt.norien) then
121  write(*,*)
122  & '*ERROR reading *SOLID SECTION: nonexistent orientation'
123  write(*,*) ' '
124  call inputerror(inpc,ipoinpc,iline,
125  &"*SOLID SECTION%")
126  call exit(201)
127  endif
128  iorientation=i
129  endif
130 !
131  if(ipos.eq.0) then
132  write(*,*) '*ERROR reading *SOLID SECTION: no element set ',
133  & elset
134  write(*,*) ' was been defined. '
135  call inputerror(inpc,ipoinpc,iline,
136  &"*SOLID SECTION%")
137  call exit(201)
138  endif
139  do i=1,nset
140  if(set(i).eq.elset) exit
141  enddo
142  if(i.gt.nset) then
143  elset(ipos:ipos)=' '
144  write(*,*) '*ERROR reading *SOLID SECTION: element set ',elset
145  write(*,*) ' has not yet been defined. '
146  call inputerror(inpc,ipoinpc,iline,
147  &"*SOLID SECTION%")
148  call exit(201)
149  endif
150 !
151 ! assigning the elements of the set the appropriate material
152 ! and orientation number
153 !
154  do j=istartset(i),iendset(i)
155  if(ialset(j).gt.0) then
156  if((lakon(ialset(j))(1:1).eq.'B').or.
157  & (lakon(ialset(j))(1:1).eq.'S')) then
158  write(*,*)
159  & '*ERROR reading *SOLID SECTION: *SOLID SECTION can'
160  write(*,*) ' not be used for beam or shell elements
161  &'
162  write(*,*) ' Faulty element: ',ialset(j)
163  call exit(201)
164  endif
165  ielmat(1,ialset(j))=imaterial
166  ielorien(1,ialset(j))=iorientation
167  else
168  k=ialset(j-2)
169  do
170  k=k-ialset(j)
171  if(k.ge.ialset(j-1)) exit
172  if((lakon(k)(1:1).eq.'B').or.
173  & (lakon(k)(1:1).eq.'S')) then
174  write(*,*) '*ERROR reading *SOLID SECTION: *SOLID SECT
175  &ION can'
176  write(*,*) ' not be used for beam or shell eleme
177  &nts'
178  write(*,*) ' Faulty element: ',k
179  call exit(201)
180  endif
181  ielmat(1,k)=imaterial
182  ielorien(1,k)=iorientation
183  enddo
184  endif
185  enddo
186 !
187  call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl,
188  & ipoinp,inp,ipoinpc)
189 !
190 ! assigning a thickness to plane stress elements and an angle to
191 ! axisymmetric elements
192 !
193 cccc if((key.eq.0).or.(lakon(ialset(istartset(i)))(1:2).eq.'CA')) then
194  if(key.eq.0) then
195  read(textpart(1)(1:20),'(f20.0)',iostat=istat) thickness
196  if(istat.gt.0) call inputerror(inpc,ipoinpc,iline,
197  &"*SOLID SECTION%")
198 !
199 ! for axial symmetric structures:
200 ! thickness for axial symmetric elements: 2 degrees
201 ! thickness for plane stress elements: reduced by 180
202 ! thickness for plane strain elements: reduced by 180
203 !
204  if(iaxial.eq.180) then
205  if(lakon(ialset(istartset(i)))(1:2).eq.'CA') then
206  thickness=datan(1.d0)*8.d0/iaxial
207  elseif(lakon(ialset(istartset(i)))(1:3).eq.'CPS') then
208  thickness=thickness/iaxial
209  elseif(lakon(ialset(istartset(i)))(1:3).eq.'CPE') then
210  thickness=thickness/iaxial
211  endif
212  endif
213 cccc else
214 cccc thickness=datan(1.d0)*8.d0/iaxial
215 cccc endif
216 !
217 ! assigning the thickness to each node of the corresponding
218 ! elements (thickness specified)
219 !
220  do j=istartset(i),iendset(i)
221  if(ialset(j).gt.0) then
222  if((lakon(ialset(j))(1:2).eq.'CP').or.
223  & (lakon(ialset(j))(1:2).eq.'CA')) then
224 !
225 ! plane stress/strain or axisymmetric elements
226 !
227  indexe=ipkon(ialset(j))
228  do l=1,8
229  thicke(1,indexe+l)=thickness
230  enddo
231  elseif(lakon(ialset(j))(1:1).eq.'T') then
232  ielem=ialset(j)
233 !
234 ! default cross section for trusses is the
235 ! rectangular cross section
236 !
237  lakon(ielem)(8:8)='R'
238  indexe=ipkon(ielem)
239  node1=kon(indexe+1)
240  if(lakon(ielem)(4:4).eq.'2') then
241  node2=kon(indexe+2)
242  else
243  node2=kon(indexe+3)
244  endif
245 !
246 ! determining a vector orthogonal to the truss
247 ! element
248 !
249  do l=1,3
250  xn(l)=co(l,node2)-co(l,node1)
251  enddo
252  if(dabs(xn(1)).gt.0.d0) then
253  p(1)=-xn(3)
254  p(2)=0.d0
255  p(3)=xn(1)
256  elseif(dabs(xn(2)).gt.0.d0) then
257  p(1)=xn(2)
258  p(2)=-xn(1)
259  p(3)=0.d0
260  else
261  p(1)=0.d0
262  p(2)=xn(3)
263  p(3)=-xn(2)
264  endif
265  dd=dsqrt(p(1)*p(1)+p(2)*p(2)+p(3)*p(3))
266  if(dd.lt.1.d-10) then
267  write(*,*)
268  & '*ERROR reading *SOLID SECTION: normal'
269  write(*,*) ' in direction 1 has zero size'
270  call exit(201)
271  endif
272  do l=1,3
273  p(l)=p(l)/dd
274  enddo
275  do l=1,3
276  thicke(1,indexe+l)=dsqrt(thickness)
277  thicke(2,indexe+l)=dsqrt(thickness)
278  indexx=ixfree
279  do m=1,3
280  xnor(indexx+m)=p(m)
281  enddo
282  ixfree=ixfree+6
283  iponor(1,indexe+l)=indexx
284  enddo
285  endif
286  else
287  k=ialset(j-2)
288  do
289  k=k-ialset(j)
290  if(k.ge.ialset(j-1)) exit
291  if((lakon(k)(1:2).eq.'CP').or.
292  & (lakon(k)(1:2).eq.'CA')) then
293  indexe=ipkon(k)
294  do l=1,8
295  thicke(1,indexe+l)=thickness
296  enddo
297  elseif(lakon(k)(1:1).eq.'T') then
298 !
299 ! default cross section for trusses is the
300 ! rectangular cross section
301 !
302  lakon(k)(8:8)='R'
303  indexe=ipkon(k)
304  node1=kon(indexe+1)
305  if(lakon(k)(4:4).eq.'2') then
306  node2=kon(indexe+2)
307  else
308  node2=kon(indexe+3)
309  endif
310 !
311 ! determining a vector orthogonal to the truss
312 ! element
313 !
314  do l=1,3
315  xn(l)=co(l,node2)-co(l,node1)
316  enddo
317  if(dabs(xn(1)).gt.0.d0) then
318  p(1)=-xn(3)
319  p(2)=0.d0
320  p(3)=xn(1)
321  elseif(dabs(xn(2)).gt.0.d0) then
322  p(1)=xn(2)
323  p(2)=-xn(1)
324  p(3)=0.d0
325  else
326  p(1)=0.d0
327  p(2)=xn(3)
328  p(3)=-xn(2)
329  endif
330  dd=dsqrt(p(1)*p(1)+p(2)*p(2)+p(3)*p(3))
331  if(dd.lt.1.d-10) then
332  write(*,*)
333  & '*ERROR reading *SOLID SECTION: normal'
334  write(*,*) ' in direction 1 has zero size'
335  call exit(201)
336  endif
337  do l=1,3
338  p(l)=p(l)/dd
339  enddo
340  do l=1,3
341  thicke(1,indexe+l)=dsqrt(thickness)
342  thicke(2,indexe+l)=dsqrt(thickness)
343  indexx=ixfree
344  do m=1,3
345  xnor(indexx+m)=p(m)
346  enddo
347  ixfree=ixfree+6
348  iponor(1,indexe+l)=indexx
349  enddo
350  endif
351  enddo
352  endif
353  enddo
354  else
355 !
356 ! assigning the thickness to each node of the corresponding
357 ! elements (thickness not specified: only axisymmetric elements)
358 !
359  thickness=datan(1.d0)*8.d0/iaxial
360  do j=istartset(i),iendset(i)
361  if(ialset(j).gt.0) then
362  if(lakon(ialset(j))(1:2).eq.'CA') then
363 c write(*,*) 'solidsections ',
364 c & ialset(j),lakon(ialset(j)),thickness
365 !
366 ! axisymmetric elements
367 !
368  indexe=ipkon(ialset(j))
369  do l=1,8
370  thicke(1,indexe+l)=thickness
371  enddo
372  endif
373  else
374  k=ialset(j-2)
375  do
376  k=k-ialset(j)
377  if(k.ge.ialset(j-1)) exit
378  if(lakon(k)(1:2).eq.'CA') then
379 !
380 ! axisymmetric elements
381 !
382  indexe=ipkon(k)
383  do l=1,8
384  thicke(1,indexe+l)=thickness
385  enddo
386  endif
387  enddo
388  endif
389  enddo
390  endif
391 !
392 ! defining cyclic symmetric conditions for axisymmetric
393 ! elements (needed for cavity radiation)
394 !
395  do j=istartset(i),iendset(i)
396  if(ialset(j).gt.0) then
397  if(lakon(ialset(j))(1:2).eq.'CA') then
398  if(mcs.gt.1) then
399  write(*,*) '*ERROR reading *SOLID SECTION: '
400  write(*,*) ' axisymmetric elements cannot be
401  &combined with cyclic symmetry'
402  call exit(201)
403  elseif(mcs.eq.1) then
404  if(int(cs(1,1)).ne.int(2.d0*pi/thickness+0.5d0))
405  & then
406  write(*,*) '*ERROR reading *SOLID SECTION: '
407  write(*,*) ' it is not allowed to define t
408  &wo different'
409  write(*,*) ' angles for an axisymmetric st
410  &ructure'
411  call exit(201)
412  else
413  exit
414  endif
415  endif
416  mcs=1
417  cs(1,1)=2.d0*pi/thickness+0.5d0
418  cs(2,1)=-0.5d0
419  cs(3,1)=-0.5d0
420  cs(5,1)=1.5d0
421  do k=6,9
422  cs(k,1)=0.d0
423  enddo
424  cs(10,1)=1.d0
425  cs(11,1)=0.d0
426  cs(12,1)=-1.d0
427  cs(14,1)=0.5
428  cs(15,1)=dcos(thickness)
429  cs(16,1)=dsin(thickness)
430  exit
431  endif
432  else
433  k=ialset(j-2)
434  do
435  k=k-ialset(j)
436  if(k.ge.ialset(j-1)) exit
437 c if(lakon(ialset(j))(1:2).eq.'CA') then
438  if(lakon(k)(1:2).eq.'CA') then
439  if(mcs.gt.1) then
440  write(*,*) '*ERROR reading *SOLID SECTION: '
441  write(*,*) ' axisymmetric elements cannot
442  &be combined with cyclic symmetry'
443  call exit(201)
444  elseif(mcs.eq.1) then
445  if(int(cs(1,1)).ne.int(2.d0*pi/thickness+0.5d0))
446  & then
447  write(*,*) '*ERROR reading *SOLID SECTION: '
448  write(*,*) ' it is not allowed to defin
449  &e two different'
450  write(*,*) ' angles for an axisymmetric
451  & structure'
452  call exit(201)
453  else
454  exit
455  endif
456  endif
457  mcs=1
458  cs(1,1)=2.d0*pi/thickness+0.5d0
459  cs(2,1)=-0.5d0
460  cs(3,1)=-0.5d0
461  cs(5,1)=1.5d0
462  do k=6,9
463  cs(k,1)=0.d0
464  enddo
465  cs(10,1)=1.d0
466  cs(11,1)=0.d0
467  cs(12,1)=-1.d0
468  cs(14,1)=0.5
469  cs(15,1)=dcos(thickness)
470  cs(16,1)=dsin(thickness)
471  exit
472  endif
473  enddo
474  endif
475  enddo
476 !
477  if(key.eq.0) then
478  call getnewline(inpc,textpart,istat,n,key,iline,ipol,inl,
479  & ipoinp,inp,ipoinpc)
480  endif
481 cccc endif
482 !
483  return
inputwarning
subroutine inputwarning(inpc, ipoinpc, iline, text)
Definition: inputwarning.f:20
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
thickness
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)