1!$Id:$
2      subroutine pblend2a(iblend,iside,isd)
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
11!     Purpose:  Construct two dimensional interpolation using blending
12
13!     Inputs:
14!        iblend(*) - Blending functions parameters/sides
15!        isd       - Dimension for sides array
16
17!     Outputs:
18!        iside(4)  - Side connection numbers for face
19!-----[--.----+----.----+----.-----------------------------------------]
20      implicit  none
21
22      include  'cblend.h'
23      include  'pointer.h'
24      include  'comblk.h'
25
26      logical   setvar,palloc, flsd
27      integer   isd,i,i1,i2
28
29      integer   iblend(*),iside(*), lblend(4)
30
31      save
32
33!     Save values of vertex nodes in case pointer to iblend changes
34
35      do i = 1,4
36        lblend(i) = iblend(10+i)
37      end do
38
39!     Check each side
40
41      do i = 1,4
42        i1 = lblend(i)
43        i2 = lblend(mod(i,4)+1)
44        call pblenda1(i,i1,i2,mr(np(162)),iside,isd, flsd)
45
46        if(flsd) then
47          numsd = numsd + 1
48          setvar = palloc(162,'BSIDE',numsd*isd,1)
49          call pblenda2(i,i1,i2,mr(np(162)),iside,isd)
50        endif
51      end do ! i
52
53      end
54
55      subroutine pblenda1(i,i1,i2,is,iside,isd,flsd)
56
57!-----[--.----+----.----+----.-----------------------------------------]
58!     Purpose:  Find side number for 2-d blend.
59
60!     Inputs:
61!        is(isd,*) - Blending side supernode lists
62!        iblend(*) - Blending functions parameters/sides
63!        isd       - Dimension for sides array
64
65!     Outputs:
66!        iside(4)  - Side connection numbers for face
67!        flsd      - false if no new side needed
68!-----[--.----+----.----+----.-----------------------------------------]
69      implicit  none
70
71      include  'cblend.h'
72
73      logical   flsd
74      integer   i,i1,i2,i3,i4, j,k, isd
75      integer   is(isd,*),iside(*)
76
77      save
78
79!     Check each side
80
81      do j = 1,numsd
82        i3 = is(1,j)
83        if(i3.eq.0 .or. i3.eq.1 .or.i3.eq.3) then
84          k = 3
85        elseif(i3.eq.2) then
86          do i4 = 3,isd
87            if(is(i4,j).ne.0) then
88              k = i4
89            endif
90          end do ! i4
91        endif
92        if((i1.eq.is(2,j) .and. i2.eq.is(k,j)) .or.
93     &     (i1.eq.is(k,j) .and. i2.eq.is(2,j))) then
94           iside(i) = j
95           flsd     = .false.
96           return
97        endif
98      end do ! j
99      flsd = .true.
100
101      end
102
103      subroutine pblenda2(i,i1,i2,is,iside,isd)
104
105!-----[--.----+----.----+----.-----------------------------------------]
106!     Purpose:  Set new side number for 2-d blend.
107
108!     Inputs:
109!        is(isd,*) - Blending side supernode lists
110!        iblend(*) - Blending functions parameters/sides
111!        isd       - Dimension for sides array
112
113!     Outputs:
114!        iside(4)  - Side connection numbers for face
115!        flsd      - false if no new side needed
116!-----[--.----+----.----+----.-----------------------------------------]
117      implicit  none
118
119      include  'cblend.h'
120
121      integer   i,i1,i2,i3, isd
122      integer   is(isd,*),iside(*)
123
124      save
125
126      do i3 = 1,isd
127        is(i3,numsd) = 0
128      end do ! i3
129      is(2,numsd) = i1
130      is(3,numsd) = i2
131      iside(i)    = numsd
132
133      end
134
135      subroutine pblend2b(n,xs,is,trb,iblend,ilr,x,ix,
136     &                    iside,isd,ndm,nen1,prt,prth,eflag,nflag)
137
138!-----[--.----+----.----+----.-----------------------------------------]
139
140!     Purpose:  Construct two dimensional interpolation using blending
141
142!     Inputs:
143!        n         - Block number
144!        xs(3,*)   - Blending supernode connections
145!        is(isd,*) - Blending side supernode lists
146!        trb       - Transformation for blending coordinates
147!        iblend(*) - Blending functions parameters/sides
148!        ilr(*)    - Material quantities for blends
149!        iside(*)  - Side connections for face
150!        isd       - Dimension for sides array
151!        ndm       - Spatial dimension of mesh
152!        nen1      - Dimension of ix array
153
154!     Outputs:
155!        x(ndm,*)  - Nodal coordinates for blended patch
156!        ix(nen1,*)- Element connections
157
158!-----[--.----+----.----+----.-----------------------------------------]
159
160      implicit   none
161
162      include   'cdata.h'
163      include   'iofile.h'
164      include   'pointer.h'
165      include   'region.h'
166      include   'trdata.h'
167      include   'comblk.h'
168
169      logical    prt,prth,eflag,nflag, setvar, palloc
170      character  ctype*15, etype*5, pelabl*5
171      integer    i,ii,in, j,jj, k, ma,m1,m2
172      integer    n,ne,nf,ni,nm,nn,nr,ns,nodinc,ntyp, styp, dlayer
173      real*8     trdeto
174
175      integer    isd,ndm,nen1
176      integer    is(isd,*),iblend(*), ix(nen1,*), ilr(*)
177
178      real*8     xs(3,*),trb(3,4),x(ndm,*)
179
180      integer    nsn(4), iside(4)
181
182      save
183
184!     Set up face
185
186      do j = 1,4
187        iblend(j+10) = iside(j)
188      end do
189
190!     Set for 4 sides
191
192      in = 4
193
194!     Check signs on sides for blend
195
196      call mkside(n,iblend(11),is,isd)
197
198!     Check for matching end nodes for blending functions
199
200      do i = 1,in
201        ii = iblend(i+10)
202        if(ii.gt.0) then
203          if(is(1,abs(ii)).eq.2) then
204            do k = 3,isd
205              if(is(k,abs(ii)).ne.0) m1 = k
206            end do
207          else
208            m1 = 3
209          endif
210        elseif(ii.lt.0) then
211          m1 = 2
212        else
213          write(*,*) ' ERROR - Zero side number specified'
214        endif
215
216        j  = mod(i,in) + 1
217        jj = iblend(j+10)
218        if(jj.gt.0) then
219          m2 = 2
220        elseif(jj.lt.0) then
221          if(is(1,abs(jj)).eq.2) then
222            do k = 3,isd
223              if(is(k,abs(jj)).ne.0) m2 = k
224            end do
225          else
226            m2 = 3
227          endif
228        else
229          write(*,*) ' ERROR - Zero side number specified'
230        endif
231
232        if(is(m1,abs(ii)).ne.is(m2,abs(jj))) then
233          if(ior.lt.0) write(*,2000) n,i,j,ii,jj
234          write(iow,2000) n,i,j,ii,jj
235        endif
236
237      end do ! i
238
239!     Quadrilateral region blends
240
241      if(in.eq.4) then
242
243!       Determine the number of nodes for each side
244        do i = 1,in
245          ii = abs(iblend(i+10))
246          do j = isd,2,-1
247            if( is(j,ii).ne.0 ) go to 110
248          end do ! j
249          write(*,*) ' ERROR - No side nodes found for side',i
250          call plstop(.true.)
251110       nsn(i) = j-1
252        end do ! i
253
254!       Get edge interpolations
255
256        nr   = iblend(1)
257        ns   = iblend(2)
258        ntyp = iblend(6)
259
260        setvar = palloc (111, 'TEMP1',(nr+1)*ndm  ,2)
261        setvar = palloc (112, 'TEMP2',(ns+1)*ndm  ,2)
262        setvar = palloc (113, 'TEMP3',(nr+1)*ndm  ,2)
263        setvar = palloc (114, 'TEMP4',(ns+1)*ndm  ,2)
264        setvar = palloc (115, 'TEMP5',max(nr,ns)+1,2)
265        setvar = palloc (116, 'TEMP6',(nr+1)*3    ,2)
266
267        nreg = iblend(10)
268        ii   = iblend(11)
269        jj   = abs(ii)
270        styp = is(1,jj)
271
272        call pside1(nr, xs, trb, ii,is(2,jj), nsn(1),ndm,
273     &              hr(np(115)),hr(np(116)), hr(np(111)),styp)
274
275        ii   = iblend(12)
276        jj   = abs(ii)
277        styp = is(1,jj)
278
279        call pside1(ns, xs, trb, ii,is(2,jj), nsn(2),ndm,
280     &              hr(np(115)),hr(np(116)), hr(np(112)),styp)
281
282        ii   =-iblend(13)
283        jj   = abs(ii)
284        styp = is(1,jj)
285
286        call pside1(nr, xs, trb, ii,is(2,jj), nsn(3),ndm,
287     &              hr(np(115)),hr(np(116)), hr(np(113)),styp)
288
289        ii   =-iblend(14)
290        jj   = abs(ii)
291        styp = is(1,jj)
292        call pside1(ns, xs, trb, ii,is(2,jj), nsn(4),ndm,
293     &              hr(np(115)),hr(np(116)), hr(np(114)),styp)
294
295        ni = iblend(3)
296        call pblendx(nn,nr,ns,ni,ntyp,ndm, hr(np(111)),hr(np(112)),
297     &               hr(np(113)),hr(np(114)),mr(np(190)),x,
298     &               nflag,prt,prth)
299
300        setvar = palloc (116, 'TEMP6',0 ,2)
301        setvar = palloc (115, 'TEMP5',0 ,2)
302        setvar = palloc (114, 'TEMP4',0 ,2)
303        setvar = palloc (113, 'TEMP3',0 ,2)
304        setvar = palloc (112, 'TEMP2',0 ,2)
305        setvar = palloc (111, 'TEMP1',0 ,2)
306
307        if(eflag) then
308          ne     = iblend(4)
309          ma     = iblend(5)
310          nm     = 4
311          nodinc = 0
312          ctype  = 'blen'
313
314          nr     = nr + 1
315          ns     = ns + 1
316
317          trdeto = trdet
318          trdet  = trb(1,1)*(trb(2,2)*trb(3,3) - trb(2,3)*trb(3,2))
319     &           + trb(1,2)*(trb(2,3)*trb(3,1) - trb(2,1)*trb(3,3))
320     &           + trb(1,3)*(trb(2,1)*trb(3,2) - trb(2,2)*trb(3,1))
321
322          if(ma.lt.0) then
323            dlayer = -ma
324          else
325            dlayer =  0
326          endif
327
328          call sblke(nr,ns,x,ix,ni,ne,nn,ndm,nen1,nodinc,ntyp,nm,ma,
329     &               dlayer,ilr,ctype)
330
331          trdet  = trdeto
332
333          nf     = nn
334        endif
335
336      endif
337
338!     Set region numbers
339
340      if(eflag) then
341        do nn = ne,nf
342          ix(nen1-1,nn) = nreg
343        enddo
344
345!       Print lists if wanted
346
347        if(prt.and.ne.gt.0) then
348          do nn = ne,nf,50
349            call prtitl(prth)
350            write(iow,2003) (i,i=1,nen)
351            if(ior.lt.0) then
352              write(  *,2003) (i,i=1,nen)
353            endif
354            j = min(nf,nn+49)
355            do i = nn,j
356              ma = ix(nen1,i)
357              etype = pelabl(ix(nen+7,i))
358              write(iow,2004) i,ma,nreg,etype,(ix(k,i),k=1,nen)
359              if(ior.lt.0) then
360                write(*,2004) i,ma,nreg,etype,(ix(k,i),k=1,nen)
361              endif
362            end do ! i
363          end do ! nn
364        endif
365      endif
366
367!     Formats
368
3692000  format(' ERROR - Blending function',i3/
370     &       '         End node 2 for side',i2,' not same as'/
371     &       '         End node 1 for side',i2/
372     &       '         Node-2 =',i8,' Node-1 =',i8)
373
3742003  format('   E l e m e n t   C o n n e c t i o n s'//
375     &   '   Elmt Mat Reg  Type',7(i3,' node'):/(21x,7(i3,' node')))
376
3772004  format(i7,2i4,1x,a5,i8:/(21x,7i8))
378
379      end
380