1!
2!     CalculiX - A 3-dimensional finite element program
3!              Copyright (C) 1998-2021 Guido Dhondt
4!
5!     This program is free software; you can redistribute it and/or
6!     modify it under the terms of the GNU General Public License as
7!     published by the Free Software Foundation(version 2);
8!
9!
10!     This program is distributed in the hope that it will be useful,
11!     but WITHOUT ANY WARRANTY; without even the implied warranty of
12!     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13!     GNU General Public License for more details.
14!
15!     You should have received a copy of the GNU General Public License
16!     along with this program; if not, write to the Free Software
17!     Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
18!
19      subroutine getdesiinfobou(ndesibou,nodedesibou,nodedesiinv,
20     &  lakon,ipkon,kon,ipoface,nodface,nodedesiinvbou,ndesi,
21     &  nodedesi,nk)
22!
23!     storing the design variables on the boundary in nodedesibou
24!
25!     create the field nodedesiinvbou which contains
26!     the following information:
27!     entry > 0: node contains to the set ndesibou
28!     number of entry: specifies the position of the node
29!            in the set nodedesiinvbou
30!
31      implicit none
32!
33      character*8 lakon(*)
34!
35      character*81 setname
36!
37      integer ndesibou,node,nodedesibou(*),i,k,m,index,
38     &  nodedesiinv(*),nelem,nope,nopedesi,ipkon(*),nnodes,
39     &  kon(*),konl(26),iaux,kflag,ipoface(*),nodface(5,*),jfacem,
40     &  nopesurf(9),ifaceq(8,6),ifacet(6,4),ifacew1(4,5),
41     &  ifacew2(8,5),nopem,nodedesiinvbou(*),ndesi,nodedesi(*),nk
42!
43      setname(1:1)=' '
44      ndesibou=0
45!
46!     nodes per face for hex elements
47!
48      data ifaceq /4,3,2,1,11,10,9,12,
49     &            5,6,7,8,13,14,15,16,
50     &            1,2,6,5,9,18,13,17,
51     &            2,3,7,6,10,19,14,18,
52     &            3,4,8,7,11,20,15,19,
53     &            4,1,5,8,12,17,16,20/
54!
55!     nodes per face for tet elements
56!
57      data ifacet /1,3,2,7,6,5,
58     &             1,2,4,5,9,8,
59     &             2,3,4,6,10,9,
60     &             1,4,3,8,10,7/
61!
62!     nodes per face for linear wedge elements
63!
64      data ifacew1 /1,3,2,0,
65     &             4,5,6,0,
66     &             1,2,5,4,
67     &             2,3,6,5,
68     &             3,1,4,6/
69!
70!     nodes per face for quadratic wedge elements
71!
72      data ifacew2 /1,3,2,9,8,7,0,0,
73     &             4,5,6,10,11,12,0,0,
74     &             1,2,5,4,7,14,10,13,
75     &             2,3,6,5,8,15,11,14,
76     &             3,1,4,6,9,13,12,15/
77!
78!     A design node is considered to be at the border if the external
79!     does not contain as many designvariables as nodes on the surface
80!
81      do i=1,nk
82         node=i
83         if(ipoface(node).eq.0) cycle
84         index=ipoface(node)
85         do
86            nelem=nodface(3,index)
87            jfacem=nodface(4,index)
88!
89            if(lakon(nelem)(4:4).eq.'8') then
90               nope=8
91               nopedesi=4
92               nopem=4
93            elseif(lakon(nelem)(4:5).eq.'20') then
94               nope=20
95               nopedesi=8
96               nopem=8
97            elseif(lakon(nelem)(4:5).eq.'10') then
98               nope=10
99               nopedesi=4
100               nopem=6
101            elseif(lakon(nelem)(4:4).eq.'4') then
102               nope=4
103               nopedesi=3
104               nopem=3
105            elseif(lakon(nelem)(4:4).eq.'6') then
106               nope=6
107               if(jfacem.le.2) then
108                  nopem=3
109                  nopedesi=3
110               else
111                  nopem=4
112                  nopedesi=4
113               endif
114            elseif(lakon(nelem)(4:5).eq.'15') then
115               nope=15
116               if(jfacem.le.2) then
117                  nopem=6
118                  nopedesi=6
119               else
120                  nopem=8
121                  nopedesi=8
122               endif
123            endif
124!
125!     actual position of the nodes belonging to the
126!     master surface
127!
128            do k=1,nope
129               konl(k)=kon(ipkon(nelem)+k)
130            enddo
131!
132            if((nope.eq.20).or.(nope.eq.8)) then
133               do m=1,nopem
134                  nopesurf(m)=konl(ifaceq(m,jfacem))
135               enddo
136            elseif((nope.eq.10).or.(nope.eq.4)) then
137               do m=1,nopem
138                  nopesurf(m)=konl(ifacet(m,jfacem))
139               enddo
140            elseif(nope.eq.15) then
141               do m=1,nopem
142                  nopesurf(m)=konl(ifacew2(m,jfacem))
143               enddo
144            else
145               do m=1,nopem
146                  nopesurf(m)=konl(ifacew1(m,jfacem))
147               enddo
148            endif
149!
150!     sum up how many designvariables are on that surface
151!
152            nnodes=0
153            do m=1,nopem
154               if(nodedesiinv(nopesurf(m)).ne.0) then
155                  nnodes=nnodes+1
156               endif
157            enddo
158!
159            if(nnodes.lt.nopedesi) then
160               do m=1,nopem
161                  if(nodedesiinv(nopesurf(m)).eq.1) then
162                     nodedesiinv(nopesurf(m))=-1
163                     ndesibou=ndesibou+1
164                     nodedesibou(ndesibou)=nopesurf(m)
165                     nodedesiinvbou(nopesurf(m))=1
166                  endif
167               enddo
168            endif
169            index=nodface(5,index)
170            if(index.eq.0) exit
171         enddo
172      enddo
173!
174!     sort the boundary nodes
175!
176      kflag=1
177      call isortii(nodedesibou,iaux,ndesibou,kflag)
178!
179      do i=1,ndesibou
180         index=nodedesibou(i)
181         nodedesiinvbou(index)=i
182         if(nodedesiinv(index).eq.-1) then
183            nodedesiinv(index)=1
184         endif
185      enddo
186!
187      return
188      end
189
190