1!$Id:$
2      subroutine blkgen(ndm,nen1,x,ix,prt,prth)
3
4!      * * F E A P * * A Finite Element Analysis Program
5
6!....  Copyright (c) 1984-2017: Regents of the University of California
7!                               All rights reserved
8
9!-----[--.----+----.----+----.-----------------------------------------]
10!      Purpose: Generate a block of elements and nodes for mesh
11!               descriptions.
12
13!      Inputs:
14!         ndm        - Dimension of 'x' array
15!         nen1       - Dimension for 'ix' array
16!         prt        - Print generated data if true
17!         prth       - Print headers if true
18
19!      Outputs:
20!         x(ndm,*)   - Block of nodal coordinates
21!         ix(nen1,*) - Block of elements
22
23!-----[--.----+----.----+----.-----------------------------------------]
24      implicit  none
25
26      include  'cblktr.h'
27      include  'cdata.h'
28      include  'crotas.h'
29      include  'iofile.h'
30      include  'pointer.h'
31      include  'region.h'
32      include  'comblk.h'
33
34      logical   prt,prth,pcomp,errck,pinput,tinput,palloc
35      logical   eltype, pblktyp
36      character xh*6,ctype*15,layer*15, etype*5, pelabl*5
37      integer   i,j,k,l,n,nm,nn,nr,ns,nt,nf,ng,ni,ntyp,nodinc
38      integer   ndm,nen1,ne,ma,mab, dlayer,nlay
39      real*8    dr,ds,dt
40
41      integer   ix(nen1,*)
42      integer   ixl(27)
43      real*8    x(ndm,*),xl(3,27),shp(3,9),tb(7),tc(4),td(15)
44
45      save
46
47      data      xh/' coord'/
48
49!     Block mesh generation routine
50
51100   if(ior.lt.0) then
52        if(ndm.le.2) write(*,5000)
53        if(ndm.ge.3) write(*,5001)
54      endif
55      errck = tinput(ctype,1,tb,7)
56      if(errck) go to 100
57
58!     Check for form of coordinate system
59
60      if(.not.pcomp(ctype,'pola',4) .and.
61     &   .not.pcomp(ctype,'cyli',4) .and.
62     &   .not.pcomp(ctype,'sphe',4) .and.
63     &   .not.pcomp(ctype,'surf',4)) then
64        ctype = 'cartesian'
65      endif
66
67!     Set parameters for block
68
69      nr     = nint(tb(1))
70      ns     = nint(tb(2))
71      nt     = nint(tb(3))
72      ntyp   = nint(tb(7))
73      do n = 1,27
74        ixl(n) = 0
75        do j = 1,ndm
76          xl(j,n) = 0.0
77        end do ! j
78      end do ! n
79      nn     = 0   ! Number of block nodes
80      nm     = 0   ! Number of block nodes
81      nlay   = 0   ! Number of material layers in block
82      dlayer = 0   ! Direction of layering
83      netyp  = 0   ! Default undefined element type
84
85!     Check for layered material properties
86
87110   errck = tinput(layer,1,td,3)
88      if(errck) go to 110
89
90      if(pcomp(layer,'laye',4)) then
91        dlayer = nint(td(1))
92        if(dlayer.eq.1) then
93          nlay = nr
94        elseif(dlayer.eq.2) then
95          nlay = ns
96        elseif(dlayer.eq.3) then
97          nlay = nt
98        endif
99        errck = palloc(111,'TEMP1',nlay,1)
100
101        j = 1
102        do while(j.le.nlay)
103120       errck = tinput(layer,0,td,15)
104          if(errck) go to 120
105          do i = j,min(j+15,nlay)
106            mr(np(111)-1+i) = nint(td(i-j+1))
107          end do ! i
108          j = j + 16
109        end do ! while
110        go to 110
111
112!     Data is a block coordinate
113
114      else
115        mab    = 0
116        eltype = .false.
117        do while (.not.eltype)
118          eltype = pblktyp(layer,td, ntyp,ns,mab)
119          if(.not.eltype) then
120            errck = tinput(layer,1,td,5)
121          endif
122        end do !
123        if(eltype) then
124          call setval(layer,15,tc(1))
125          l  = nint(tc(1))
126          if(l.eq.0 ) go to 22
127          if(l.gt.27) go to 402
128          nm = max(nm,l)
129          nn = nn + 1
130          ixl(l)  = l
131          xl(1,l) = td(1)
132          xl(2,l) = td(2)
133          xl(3,l) = td(3)
134        else
135          go to 110  ! Input another record
136        endif
137      end if ! no layers
138
139!     Input block coordinates
140
14121    if(ior.lt.0) write(*,5002)
142        errck = pinput(tc,4)
143        if(errck) go to 21
144        l = nint(tc(1))
145        if(l.eq.0 ) go to 22
146        if(l.gt.27) go to 402
147        nm = max(nm,l)
148        nn = nn + 1
149        ixl(l)  = l
150        xl(1,l) = tc(2)
151        xl(2,l) = tc(3)
152        xl(3,l) = tc(4)
153      go to 21
154
155!     Set a default value for unspecified block type
156
15722    if(ntyp.eq.0 .and. nt.gt.0 .and. ndm.eq.3 .and. ixl(8).ne.0) then
158        ntyp = 10
159      endif
160
161      if(pcomp(ctype,'surf',4) .and. ntyp.eq.10) then
162        ntyp = 0
163      endif
164
165!     2-d generations
166
167      if(ntyp.lt.10) then
168        nt     = 1
169        ni     = nint(tb(3))
170        ne     = nint(tb(4))
171        ma     = nint(tb(5))
172        nodinc = nint(tb(6))
173
174!     3-d generations
175
176      elseif(ntyp.lt.20) then
177        nt     = nint(tb(3))
178        ni     = nint(tb(4))
179        ne     = nint(tb(5))
180        ma     = nint(tb(6))
181        nodinc = 0
182
183!     Shell generations
184
185      elseif(ntyp.lt.30) then
186
187        nt     = 1
188        ni     = nint(tb(3))
189        ne     = nint(tb(4))
190        ma     = nint(tb(5))
191        nodinc = nint(tb(6))
192
193!     User generations
194
195      else
196        if(ndm.le.2) then
197          nt     = 1
198          ni     = nint(tb(3))
199          ne     = nint(tb(4))
200          ma     = nint(tb(5))
201          nodinc = nint(tb(6))
202        else
203          nt     = nint(tb(3))
204          ni     = nint(tb(4))
205          ne     = nint(tb(5))
206          ma     = nint(tb(6))
207          nodinc = 0
208        endif
209      endif
210
211!     Set the final material number
212
213      if(mab.gt.0) then
214        ma = mab
215      endif
216
217!     Reset to default values if necessary
218
219      if(ni.eq.0) ni = nio + 1
220      if(ne.eq.0) ne = neo + 1
221      if(ma.eq.0) ma = mao
222
223      ma     = max(ma,1)
224      nodinc = max(nodinc,0)
225      nr     = max(nr,1)
226      ns     = max(ns,1)
227      nt     = max(nt,1)
228      ni     = max(ni,1)
229
230!     Output block data
231
232      if(prt) then
233        call prtitl(prth)
234        write(iow,2000) nr,ns,nt,ni,ne,ma,nodinc,ntyp
235        if(ne.le.0) write(iow,3000)
236        if(ior.lt.0) then
237          write(*,2000) nr,ns,nt,ni,ne,ma,nodinc,ntyp
238          if(ne.le.0) write(*,3000)
239        endif
240        if(nlay.gt.0) then
241          write(iow,2005) (j,j=1,min(5,nlay))
242          write(iow,2006) (j,mr(np(111)-1+j),j=1,nlay)
243          if(ior.lt.0) then
244            write(*,2005) (j,j=1,min(5,nlay))
245            write(*,2006) (j,mr(np(111)-1+j),j=1,nlay)
246          endif
247        end if ! nlay > 0
248        write(iow,2002) ctype,(i,xh,i=1,ndm)
249        if(ior.lt.0) then
250          write(*,2002) ctype,(i,xh,i=1,ndm)
251        endif
252        do l = 1,27
253          if(ixl(l).gt.0) then
254            write(iow,2001) l,(xl(i,l),i=1,ndm)
255            if(ior.lt.0) then
256              write(*,2001) l,(xl(i,l),i=1,ndm)
257            endif
258          endif
259        end do ! l
260      endif
261
262!     Set generation increments of natural coordinates
263
264      dr = 2.d0/nr
265      ds = 2.d0/ns
266
267!     Determine last element number to be generated
268
269      if(ntyp.lt.10) then
270        if(nn.le.3) then
271          ng = ni + nr
272          if(ns.eq.1) then
273            nf = ne + nr - 1
274          else
275            nf = ne + nr/2 - 1
276          endif
277          nr = nr + 1
278        else
279          if (ntyp.eq.0) then
280            nf = ne + nr*ns - 1
281          elseif (abs(ntyp).eq.7) then
282            nf = ne + (nr*ns)/2 - 1
283          elseif (ntyp.ge.8) then
284            nf = ne + (nr*ns)/4 - 1
285          elseif (ntyp.lt.0) then
286            nf = ne + 4*nr*ns - 1
287          else
288            nf = ne + 2*nr*ns - 1
289          endif
290
291!         Determine last node number to be generated
292
293          nr = nr + 1
294          ns = ns + 1
295          if(ndm.eq.1) ns = 1
296          ng = nr*ns + nodinc*(ns-1) + ni -1
297          if(ntyp.eq. -7) then
298            ng = ng + ((nr-1)*(ns-1))/2
299          elseif(ntyp .eq. -1) then
300            ng = ng + (nr-1)*(ns-1)
301          elseif(ntyp .eq.  8) then
302            ng = ng - ((nr-1)*(ns-1))/4
303          endif
304        endif
305        if(nf.gt.numel.and.ne.gt.0) go to 401
306        if(ng.gt.numnp) go to 400
307
308!       Generate nodes
309
310        call sblkn(nr,ns,xl,ixl,shp,x,dr,ds,ni,n,ndm,nodinc,ntyp,
311     &             nm,ctype,prt,prth)
312
313!       Generate elements
314
315        call sblke(nr,ns,x,ix,ni,ne,n,ndm,nen1,nodinc,ntyp,nm,ma,
316     &             dlayer,mr(np(111)),ctype)
317
318!     3-d generations
319
320      elseif(ntyp.lt.20) then
321        dt = 2.d0/nt
322        if(ntyp.eq.10) then
323          nf = ne + nr*ns*nt - 1
324        elseif(ntyp.eq.11) then
325          nf = ne + 6*nr*ns*nt - 1
326        else
327          write(iow,4003) ntyp
328          call plstop(.true.)
329        endif
330        if(nf.gt.numel.and.ne.gt.0) go to 401
331        nr = nr + 1
332        ns = ns + 1
333        nt = nt + 1
334        ng = nr*ns*nt + ni -1
335        if(ng.gt.numnp) go to 400
336          call vblkn(nr,ns,nt,xl,x,ixl,dr,ds,dt,
337     &               ni,ndm,ctype,prt,prth)
338
339          if(ne.gt.0) then
340            call vblke(nr,ns,nt,ix,ni,ne,nf,nen1,ma,ntyp,
341     &                 dlayer,mr(np(111)))
342          endif
343
344!     User generations
345
346      else
347        dt = 2.d0/nt
348        nf = ne
349        ng = ni
350        call ublk(td,nn,nr,ns,nt,xl,x,ixl,ix,dr,ds,dt,
351     &            ng,nf,ndm,nen1,ma,ctype,prt, 2)
352        if(nf.gt.numel.and.ne.gt.0) go to 401
353        if(ng.gt.numnp) go to 400
354      endif
355
356!     Set old numbers
357
358      if(ne.gt.0) neo = nf
359      nio = ng
360      mao = ma
361
362!     Set node type to active
363
364      do n = ni,ng
365        mr(np(190)-1+n) = 0
366      end do ! n
367
368!     Set region number
369
370      do n = ne,nf
371        ix(nen1-1,n) = nreg
372      end do ! n
373
374      if(nlay.gt.0) then
375        errck = palloc(111,'TEMP1',0,1)
376      endif
377
378!     Print lists if wanted
379
380      if(prt.and.ne.gt.0) then
381        do n = ne,nf,50
382          call prtitl(prth)
383          write(iow,2003) (i,i=1,nen)
384          if(ior.lt.0) then
385              write(  *,2003) (i,i=1,nen)
386          endif
387          j = min(nf,n+49)
388          do i = n,j
389            ma = ix(nen1,i)
390            etype = pelabl(ix(nen+7,i))
391            write(iow,2004) i,ma,nreg,etype,(ix(k,i),k=1,nen)
392            if(ior.lt.0) then
393              write(*,2004) i,ma,nreg,etype,(ix(k,i),k=1,nen)
394            endif
395          end do ! i
396        end do ! n
397      endif
398
399      return
400
401!     Error messages
402
403400   write(iow,4000) ng,numnp
404      if(ior.lt.0) write(  *,4000) ng,numnp
405      return
406
407401   write(iow,4001) nf,numel
408      if(ior.lt.0) write(  *,4001) nf,numel
409      return
410
411402   write(iow,4002) l
412      if(ior.lt.0) write(*,4002) l
413
414!     Formats
415
4162000  format('   N o d e   G e n e r a t i o n s'//
417     &   10x,'Number of xi_1-increments ',i5/
418     &   10x,'Number of xi_2-increments ',i5/
419     &   10x,'Number of xi_3-increments ',i5/
420     &   10x,'First node number         ',i5/
421     &   10x,'First element number      ',i5/
422     &   10x,'Element material number   ',i5/
423     &   10x,'Node line increment       ',i5/
424     &   10x,'Block type (0-10)         ',i5/1x)
425
4262001  format(i9,1p,3e12.3)
427
4282002  format(/5x,'Block coordinates specified as: ',a,//,
429     &       5x,'node',3(i6,a6))
430
4312003  format('   E l e m e n t   C o n n e c t i o n s'//
432     &   '   Elmt Mat Reg  Type',7(i3,'-node'):/(21x,7(i3,'-node')))
433
4342004  format(i7,2i4,1x,a5,7i8:/(21x,7i8))
435
4362005  format(/5x,'Layered Material Properties'/
437     &      /(7x,4(i2,'-Layer Matl')))
438
4392006  format(7x,4(i8,i5))
440
4413000  format(' *WARNING* No elements are generated ')
442
4434000  format(' *ERROR* Insufficient storage for nodes'/
444     &   10x,'final node =',i5,5x,'numnp =',i5)
445
4464001  format(' *ERROR* Insufficient storage for elements'/
447     &   10x,'final element =',i5,5x,'numel =',i5)
448
4494002  format(' *ERROR* Block node has number > 27.  No. =',i8)
450
4514003  format(' *ERROR* Block type incorrect: Ntype =',i8)
452
4535000  format(' Input: type,nr,ns,ni,ne,ma,nodinc,ntyp'/3x,'>',$)
454
4555001  format(' Input: type,nr,ns,nt,ni,ne,ma,ntyp'/3x,'>',$)
456
4575002  format(' Input: node, x-1, x-2, x-3'/3x,'>',$)
458
459      end
460