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

Go to the source code of this file.

Functions/Subroutines

subroutine blockanalysis (set, nset, istartset, iendset, ialset, nblk, ipkon, kon, ielfa, nodface, neiel, neij, neifa, ipoface, ipnei, konf, istartblk, iendblk, nactdoh, nblket, nblkze, neielsize, ielblk, nk, nactdohinv)
 

Function/Subroutine Documentation

◆ blockanalysis()

subroutine blockanalysis ( character*81, dimension(*)  set,
integer  nset,
integer, dimension(*)  istartset,
integer, dimension(*)  iendset,
integer, dimension(*)  ialset,
integer  nblk,
integer, dimension(*)  ipkon,
integer, dimension(*)  kon,
integer, dimension(4,*)  ielfa,
integer, dimension(5,*)  nodface,
integer, dimension(*)  neiel,
integer, dimension(*)  neij,
integer, dimension(*)  neifa,
integer, dimension(*)  ipoface,
integer, dimension(*)  ipnei,
integer, dimension(*)  konf,
integer, dimension(*)  istartblk,
integer, dimension(*)  iendblk,
integer, dimension(*)  nactdoh,
integer, dimension(*)  nblket,
integer, dimension(*)  nblkze,
integer  neielsize,
integer, dimension(*)  ielblk,
integer  nk,
integer, dimension(*)  nactdohinv 
)
23 !
24 ! renumbering the blocks in a CFD-analysis
25 !
26  implicit none
27 !
28  logical iexternal(6)
29 !
30  character*81 set(*)
31 !
32  integer nset,istartset(*),iendset(*),ialset(*),
33  & nblk,i,j,k,neighbor(8,8,8),iel,ii,indexe,ifreenei,ifour,
34  & iaux,iset,j2,index,indexold,ifree,nodes(4),ipkon(*),kon(*),
35  & kflag,ielfa(4,*),nodface(5,*),neiel(*),neij(*),neifa(*),
36  & ipoface(*),ipnei(*),ifaceq(8,6),ifreenei2,iel2,ifa,
37  & kon2(8),kon3(8),kon4(8),kon5(8),kon6(8),nactdohinv(*),
38  & kon7(8),kon8(8),ielsav,ifasav,jopposite8(6),
39  & konf(*),nef,indexet,indexeet,indexze,indexeze,ifaf,
40  & indexe2,istartblk(*),iendblk(*),nactdoh(*),nblket(*),
41  & nblkze(*),nblketl,nblkzel,neielsize,ielblk(*),nk
42 !
43 ! nodes belonging to the element faces
44 !
45  data ifaceq /4,3,2,1,11,10,9,12,
46  & 5,6,7,8,13,14,15,16,
47  & 1,2,6,5,9,18,13,17,
48  & 2,3,7,6,10,19,14,18,
49  & 3,4,8,7,11,20,15,19,
50  & 4,1,5,8,12,17,16,20/
51  data jopposite8 /2,1,5,6,3,4/
52  data kon2 /2,3,4,1,6,7,8,5/
53  data kon3 /3,4,1,2,7,8,5,6/
54  data kon4 /4,1,2,3,8,5,6,7/
55  data kon5 /5,8,7,6,1,4,3,2/
56  data kon6 /6,5,8,7,2,1,4,3/
57  data kon7 /7,6,5,8,3,2,1,4/
58  data kon8 /8,7,6,5,4,3,2,1/
59 !
60  nblk=0
61  do iset=1,nset
62  if(set(iset)(1:10).eq.'FLUIDBLOCK') then
63  nblk=nblk+1
64  endif
65  enddo
66 !
67  if(nblk.eq.0) return
68 !
69 ! initialize neighbor
70 !
71  do i=1,8
72  do j=1,8
73  do k=1,8
74  neighbor(i,j,k)=0
75  enddo
76  enddo
77  enddo
78  neighbor(5,1,2)=4
79  neighbor(2,1,5)=4
80  neighbor(1,2,6)=3
81  neighbor(6,2,1)=3
82  neighbor(2,6,5)=7
83  neighbor(5,6,2)=7
84  neighbor(6,5,1)=8
85  neighbor(1,5,6)=8
86  neighbor(2,1,4)=5
87  neighbor(4,1,2)=5
88  neighbor(1,4,3)=8
89  neighbor(3,4,1)=8
90  neighbor(4,3,2)=7
91  neighbor(2,3,4)=7
92  neighbor(3,2,1)=6
93  neighbor(1,2,3)=6
94  neighbor(5,6,7)=2
95  neighbor(7,6,5)=2
96  neighbor(6,7,8)=3
97  neighbor(8,7,6)=3
98  neighbor(7,8,5)=4
99  neighbor(5,8,7)=4
100  neighbor(8,5,6)=1
101  neighbor(6,5,8)=1
102  neighbor(2,3,7)=4
103  neighbor(7,3,2)=4
104  neighbor(3,7,6)=8
105  neighbor(6,7,3)=8
106  neighbor(7,6,2)=5
107  neighbor(2,6,7)=5
108  neighbor(6,2,3)=1
109  neighbor(3,2,6)=1
110  neighbor(8,7,3)=6
111  neighbor(3,7,8)=6
112  neighbor(7,3,4)=2
113  neighbor(4,3,7)=2
114  neighbor(3,4,8)=1
115  neighbor(8,4,3)=1
116  neighbor(4,8,7)=5
117  neighbor(7,8,4)=5
118  neighbor(5,8,4)=7
119  neighbor(4,8,5)=7
120  neighbor(8,4,1)=3
121  neighbor(1,4,8)=3
122  neighbor(4,1,5)=2
123  neighbor(5,1,4)=2
124  neighbor(1,5,8)=6
125  neighbor(8,5,1)=6
126 !
127  ifour=4
128  kflag=1
129  nef=0
130 !
131 ! loop over all blocks
132 !
133  nblk=0
134  do iset=1,nset
135  if(set(iset)(1:10).ne.'FLUIDBLOCK') cycle
136  nblk=nblk+1
137  iel=ialset(istartset(iset))
138 !
139 !
140 ! determining the external element faces of the fluid mesh
141 ! the faces are catalogued by the three lowes nodes numbers
142 ! in ascending order. ipoface(i) points to a face for which
143 ! node i is the lowest node and nodface(1,ipoface(i)) and
144 ! nodface(2,ipoface(i)) are the next lower ones.
145 ! nodface(3,ipoface(i)) contains the element number,
146 ! nodface(4,ipoface(i)) the face number and nodface(5,ipoface(i))
147 ! is a pointer to the next surface for which node i is the
148 ! lowest node; if there are no more such surfaces the pointer
149 ! has the value zero
150 ! An external element face is one which belongs to one element
151 ! only
152 !
153  ifree=1
154  ifreenei=0
155  do i=1,nk
156  ipoface(i)=0
157  enddo
158  do i=1,6*neielsize
159  neiel(i)=0
160  enddo
161 !
162  do ii=istartset(iset),iendset(iset)
163  i=ialset(ii)
164  indexe=ipkon(i)
165 !
166 ! hex element
167 !
168  ipnei(i)=ifreenei
169  do j=1,6
170  do k=1,4
171  nodes(k)=kon(indexe+ifaceq(k,j))
172  enddo
173  call isortii(nodes,iaux,ifour,kflag)
174  indexold=0
175  index=ipoface(nodes(1))
176  do
177 !
178 ! adding a surface which has not been
179 ! catalogued so far
180 !
181  if(index.eq.0) then
182  nodface(1,ifree)=nodes(2)
183  nodface(2,ifree)=nodes(3)
184  nodface(3,ifree)=i
185  nodface(4,ifree)=j
186  nodface(5,ifree)=ipoface(nodes(1))
187  ipoface(nodes(1))=ifree
188  ifreenei=ifreenei+1
189  neiel(ifreenei)=0
190  neifa(ifreenei)=ifree
191  ielfa(1,ifree)=i
192  ielfa(4,ifree)=j
193  ifree=ifree+1
194  exit
195  endif
196 !
197 ! removing a surface which has already
198 ! been catalogued
199 !
200  if((nodface(1,index).eq.nodes(2)).and.
201  & (nodface(2,index).eq.nodes(3))) then
202  ifreenei=ifreenei+1
203 !
204 ! completing the facial info in neifa
205 !
206  neifa(ifreenei)=index
207  ielfa(2,index)=i
208 !
209 ! the neighboring elements to the face are i and iel2
210 ! with local face number for the face of j and j2
211 !
212  iel2=ielfa(1,index)
213  j2=ielfa(4,index)
214 !
215 ! completing the neighboring info for (element i,side j)
216 !
217  neiel(ifreenei)=iel2
218  neij(ifreenei)=j2
219 !
220 ! completing the neighboring info for (element iel2,side j2)
221 !
222  ifreenei2=ipnei(iel2)+j2
223  neiel(ifreenei2)=i
224  neij(ifreenei2)=j
225  exit
226  endif
227  indexold=index
228  index=nodface(5,index)
229  enddo
230  enddo
231  enddo
232 !
233 ! looking for a corner element of the block
234 !
235  index=ipnei(iel)
236  ifa=1
237 !
238 ! marching in one direction up to the wall
239 !
240  do
241  if(neiel(index+ifa).eq.0) then
242  if((ifa.eq.1).or.(ifa.eq.2)) then
243  ifa=3
244  elseif((ifa.eq.3).or.(ifa.eq.5)) then
245  ifa=4
246  else
247  ifa=1
248  endif
249  exit
250  endif
251  iel=neiel(index+ifa)
252  ifa=neij(index+ifa)
253  ifa=jopposite8(ifa)
254  index=ipnei(iel)
255  enddo
256 !
257 ! marching in the next direction up to the wall
258 !
259  do
260  if(neiel(index+ifa).eq.0) then
261  if((ifa.eq.1).or.(ifa.eq.2)) then
262  ifa=3
263  elseif((ifa.eq.3).or.(ifa.eq.5)) then
264  ifa=4
265  else
266  ifa=1
267  endif
268  exit
269  endif
270  iel=neiel(index+ifa)
271  ifa=neij(index+ifa)
272  ifa=jopposite8(ifa)
273  index=ipnei(iel)
274  enddo
275 !
276 ! marching in the third direction up to the wall
277 !
278  ielsav=iel
279  ifasav=ifa
280 !
281  do
282  if(neiel(index+ifa).eq.0) then
283  exit
284  endif
285  iel=neiel(index+ifa)
286  ifa=neij(index+ifa)
287  ifa=jopposite8(ifa)
288  index=ipnei(iel)
289  enddo
290 !
291 ! check whether the element is a corner element and
292 ! modifying the topology in such
293 ! a way that the local axes are pointing inwards into the
294 ! block
295 !
296  index=ipnei(iel)
297  do j=1,6
298  if(neiel(index+j).eq.0) then
299  iexternal(j)=.true.
300  else
301  iexternal(j)=.false.
302  endif
303  enddo
304  indexe=ipkon(iel)
305  if(iexternal(1).and.iexternal(6).and.iexternal(3)) then
306  do j=1,8
307  konf(indexe+j)=kon(indexe+j)
308  enddo
309  elseif(iexternal(1).and.iexternal(3).and.iexternal(4)) then
310  do j=1,8
311  konf(indexe+j)=kon(indexe+kon2(j))
312  enddo
313  elseif(iexternal(1).and.iexternal(4).and.iexternal(5)) then
314  do j=1,8
315  konf(indexe+j)=kon(indexe+kon3(j))
316  enddo
317  elseif(iexternal(1).and.iexternal(5).and.iexternal(6)) then
318  do j=1,8
319  konf(indexe+j)=kon(indexe+kon4(j))
320  enddo
321  elseif(iexternal(2).and.iexternal(3).and.iexternal(6)) then
322  do j=1,8
323  konf(indexe+j)=kon(indexe+kon5(j))
324  enddo
325  elseif(iexternal(2).and.iexternal(3).and.iexternal(6)) then
326  do j=1,8
327  konf(indexe+j)=kon(indexe+kon5(j))
328  enddo
329  elseif(iexternal(2).and.iexternal(3).and.iexternal(4)) then
330  do j=1,8
331  konf(indexe+j)=kon(indexe+kon6(j))
332  enddo
333  elseif(iexternal(2).and.iexternal(4).and.iexternal(5)) then
334  do j=1,8
335  konf(indexe+j)=kon(indexe+kon7(j))
336  enddo
337  elseif(iexternal(2).and.iexternal(5).and.iexternal(6)) then
338  do j=1,8
339  konf(indexe+j)=kon(indexe+kon8(j))
340  enddo
341  else
342 !
343 ! if no corner: look in the other direction
344 !
345  if((ifasav.eq.1).or.(ifasav.eq.2)) then
346  ifa=4
347  elseif((ifasav.eq.3).or.(ifasav.eq.5)) then
348  ifa=1
349  else
350  ifa=3
351  endif
352  iel=ielsav
353 !
354  do
355  if(neiel(index+ifa).eq.0) then
356  exit
357  endif
358  iel=neiel(index+ifa)
359  ifa=neij(index+ifa)
360  ifa=jopposite8(ifa)
361  index=ipnei(iel)
362  enddo
363 !
364 ! iel is now a corner element
365 !
366 ! modifying the topology of the corner element in such
367 ! a way that the local axes are pointing inwards into the
368 ! block
369 !
370  index=ipnei(iel)
371  do j=1,6
372  if(neiel(index+j).eq.0) then
373  iexternal(j)=.true.
374  else
375  iexternal(j)=.false.
376  endif
377  enddo
378  indexe=ipkon(iel)
379  if(iexternal(1).and.iexternal(6).and.iexternal(3)) then
380  do j=1,8
381  konf(indexe+j)=kon(indexe+j)
382  enddo
383  elseif(iexternal(1).and.iexternal(3).and.iexternal(4)) then
384  do j=1,8
385  konf(indexe+j)=kon(indexe+kon2(j))
386  enddo
387  elseif(iexternal(1).and.iexternal(4).and.iexternal(5)) then
388  do j=1,8
389  konf(indexe+j)=kon(indexe+kon3(j))
390  enddo
391  elseif(iexternal(1).and.iexternal(5).and.iexternal(6)) then
392  do j=1,8
393  konf(indexe+j)=kon(indexe+kon4(j))
394  enddo
395  elseif(iexternal(2).and.iexternal(3).and.iexternal(6)) then
396  do j=1,8
397  konf(indexe+j)=kon(indexe+kon5(j))
398  enddo
399  elseif(iexternal(2).and.iexternal(3).and.iexternal(6)) then
400  do j=1,8
401  konf(indexe+j)=kon(indexe+kon5(j))
402  enddo
403  elseif(iexternal(2).and.iexternal(3).and.iexternal(4)) then
404  do j=1,8
405  konf(indexe+j)=kon(indexe+kon6(j))
406  enddo
407  elseif(iexternal(2).and.iexternal(4).and.iexternal(5)) then
408  do j=1,8
409  konf(indexe+j)=kon(indexe+kon7(j))
410  enddo
411  elseif(iexternal(2).and.iexternal(5).and.iexternal(6)) then
412  do j=1,8
413  konf(indexe+j)=kon(indexe+kon8(j))
414  enddo
415  else
416  write(*,*) '*ERROR in blockanalysis'
417  write(*,*) ' no block structure'
418  write(*,*) ' no corner element found'
419  call exit(201)
420  endif
421  endif
422 !
423  istartblk(nblk)=nef+1
424  index=ipnei(iel)
425  nblketl=0
426  nblkzel=0
427 !
428 ! renumbering the elements
429 !
430  do
431  indexze=index
432  indexeze=indexe
433 !
434 ! check whether all layers have the same number of
435 ! elements
436 !
437  if(nblkzel.ne.0) then
438  if(nblkze(nblk).eq.0) then
439  nblkze(nblk)=nblkzel
440  elseif(nblkze(nblk).ne.nblkzel) then
441  write(*,*) '*ERROR in blockanalysis'
442  write(*,*) ' block ',set(iset)
443  write(*,*) ' is no block'
444  call exit(201)
445  endif
446  nblkzel=0
447  endif
448 !
449  do
450  indexet=index
451  indexeet=indexe
452 !
453  if(nblketl.ne.0) then
454  if(nblket(nblk).eq.0) then
455  nblket(nblk)=nblketl
456  elseif(nblket(nblk).ne.nblketl) then
457  write(*,*) '*ERROR in blockanalysis'
458  write(*,*) ' block ',set(iset)
459  write(*,*) ' is no block'
460  call exit(201)
461  endif
462  nblketl=0
463  endif
464 !
465  nblketl=nblketl+1
466  nblkzel=nblkzel+1
467 !
468  nef=nef+1
469  nactdoh(iel)=nef
470  nactdohinv(nef)=iel
471  ielblk(nef)=nblk
472  ifaf=4
473  call identifyface(konf(indexe+1),kon(indexe+1),ifaf,ifa)
474 !
475  do
476  iel=neiel(index+ifa)
477  if(iel.eq.0) exit
478  ifa=neij(index+ifa)
479 !
480  nblketl=nblketl+1
481  nblkzel=nblkzel+1
482 !
483  nef=nef+1
484  nactdoh(iel)=nef
485  nactdohinv(nef)=iel
486  ielblk(nef)=nblk
487  ifa=jopposite8(ifa)
488 !
489  indexe2=ipkon(iel)
490  call adaptconnectivity(konf(indexe+1),ifaf,
491  & kon(indexe2+1),konf(indexe2+1),neighbor)
492  indexe=indexe2
493  index=ipnei(iel)
494  enddo
495 !
496  ifaf=5
497  call identifyface(konf(indexeet+1),kon(indexeet+1),
498  & ifaf,ifa)
499  iel=neiel(indexet+ifa)
500  if(iel.eq.0) exit
501 !
502  indexe2=ipkon(iel)
503  call adaptconnectivity(konf(indexeet+1),ifaf,
504  & kon(indexe2+1),konf(indexe2+1),neighbor)
505  indexe=indexe2
506  index=ipnei(iel)
507  enddo
508 !
509  ifaf=2
510  call identifyface(konf(indexeze+1),kon(indexeze+1),
511  & ifaf,ifa)
512  iel=neiel(indexze+ifa)
513  if(iel.eq.0) exit
514 !
515  indexe2=ipkon(iel)
516  call adaptconnectivity(konf(indexeze+1),ifaf,
517  & kon(indexe2+1),konf(indexe2+1),neighbor)
518  indexe=indexe2
519  index=ipnei(iel)
520  enddo
521  iendblk(nblk)=nef
522 !
523  enddo
524 !
525 c open(91,file='kon',status='unknown')
526 c do i=1,nef
527 c iel=nactdohinv(i)
528 c indexe=ipkon(iel)
529 c write(91,100) i
530 c 100 format(' -1',i10,' 1 0WATER')
531 c write(91,101) (konf(indexe+j),j=1,8)
532 c 101 format(' -2',8i10)
533 c enddo
534 c close(91)
535 !
536  return
isortii
subroutine isortii(ix, iy, n, kflag)
Definition: isortii.f:6
nodes
subroutine nodes(inpc, textpart, co, nk, nk_, set, istartset, iendset, ialset, nset, nset_, nalset, nalset_, istep, istat, n, iline, ipol, inl, ipoinp, inp, ipoinpc)
Definition: nodes.f:22
identifyface
subroutine identifyface(konl1, konl2, iface1, iface2)
Definition: identifyface.f:20
adaptconnectivity
subroutine adaptconnectivity(konl1, iface, konl2, konl2f, neighbor)
Definition: adaptconnectivity.f:21
Hosted by OpenAircraft.com, (Michigan UAV, LLC)