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!     construction of the B matrix
20!
21!     author: Yannick Muller
22!
23      subroutine flowoutput(itg,ieg,ntg,nteq,bc,lakon,ntmat_,
24     &     v,shcon,nshcon,ipkon,kon,co,nflow, dtime,ttime,time,
25     &     ielmat,prop,ielprop,nactdog,nacteq,iin,physcon,
26     &     camt,camf,camp,rhcon,nrhcon,vold,jobnamef,set,istartset,
27     &     iendset,ialset,nset,mi,iaxial,istep,iit)
28!
29      implicit none
30!
31      logical identity
32!
33      character*8 lakon(*)
34      character*81 set(*)
35      character*132 jobnamef(*),fnnet
36!
37      integer mi(*),itg(*),ieg(*),ntg,nflow,ielmat(mi(3),*),i,
38     &     nrhcon(*),node,iaxial,ider,idirf(8),ieq,imat,kflag,
39     &     ntmat_,nteq,nshcon(*),nelem,index,ipkon(*),kon(*),iin,
40     &     nactdog(0:3,*),nacteq(0:3,*),ielprop(*),node1,nodem,node2,
41     &     istartset(*),iendset(*),ialset(*),nset,nodef(8),numf,
42     &     istep,iit,iplausi
43!
44      real*8 physcon(*),v(0:mi(2),*),shcon(0:3,ntmat_,*),co(3,*),
45     &     prop(*),dtime,ttime,time,xflow,camp(*),camt(*),camf(*),
46     &     rhcon(0:1,ntmat_,*),vold(0:mi(2),*),eta,
47     &     bc(*),cp,dvi,df(8),gastemp,f,g(3),r,rho,ts1,ts2
48!
49!
50      do i=1,132
51         if(jobnamef(1)(i:i).eq.' ') exit
52      enddo
53      i=i-1
54      fnnet=jobnamef(1)(1:i)//'.net'
55      open(1,file=fnnet,status='unknown')
56!
57      kflag=3
58!
59      do i=1,nflow
60         nelem=ieg(i)
61!
62!        output for gas networks
63!
64         if((lakon(nelem)(2:5).ne.'LIPI').and.
65     &      (lakon(nelem)(2:5).ne.'LICH')) then
66!
67            index=ipkon(nelem)
68            node1=kon(index+1)
69            nodem=kon(index+2)
70            node2=kon(index+3)
71!
72            xflow=v(1,nodem)
73!
74            if(lakon(nelem)(2:3).ne.'LP') then
75!
76!              compressible
77!
78               if(node1.eq.0) then
79                  ts1=v(3,node2)
80                  ts2=ts1
81               elseif(node2.eq.0) then
82                  ts1=v(3,node1)
83                  ts2=ts1
84               else
85                  ts1=v(3,node1)
86                  ts2=v(3,node2)
87               endif
88               gastemp=(ts1+ts2)/2.d0
89            else
90!
91!              incompressible
92!
93               if(xflow.gt.0) then
94                  gastemp=v(3,node1)
95               else
96                  gastemp=v(3,node2)
97               endif
98            endif
99!
100            imat=ielmat(1,nelem)
101!
102            call materialdata_tg(imat,ntmat_,gastemp,shcon,nshcon,cp,r,
103     &         dvi,rhcon,nrhcon,rho)
104!
105            if(nacteq(2,nodem).ne.0) then
106               ieq=nacteq(2,nodem)
107!
108!              dummy set number
109!
110               numf=1
111!
112               call flux(node1,node2,nodem,nelem,lakon,kon,ipkon,
113     &              nactdog,identity,
114     &              ielprop,prop,kflag,v,xflow,f,nodef,idirf,df,
115     &              cp,r,rho,physcon,g,co,dvi,numf,vold,set,shcon,
116     &              nshcon,rhcon,nrhcon,ntmat_,mi,ider,ttime,time,
117     &              iaxial,iplausi)
118            endif
119         endif
120!
121         if(lakon(ieg(i))(2:5).eq.'LICH') then
122            if((lakon(ieg(i))(6:7).eq.'SG').or.
123     &           (lakon(ieg(i))(6:7).eq.'WE').or.
124     &           (lakon(ieg(i))(6:7).eq.'DS')) then
125               index=ipkon(ieg(i))
126               node=kon(index+2)
127               if(nactdog(3,node).eq.0) cycle
128               index=ielprop(ieg(i))
129               if(lakon(ieg(i))(6:7).eq.'SG') then
130                  eta=prop(index+4)
131                  nelem=int(prop(index+7))
132               elseif(lakon(ieg(i))(6:7).eq.'WE') then
133                  eta=prop(index+4)
134                  nelem=int(prop(index+7))
135               elseif(lakon(ieg(i))(6:7).eq.'DS') then
136                  eta=prop(index+7)
137                  nelem=int(prop(index+9))
138               endif
139               if(nelem.ne.0) then
140                  write(*,*) '     *INFO in flowoutput: hydraulic jump'
141                  write(*,*) '           in element ',nelem,'.'
142                  write(*,*) '           relative location:',eta
143                  write(*,*)
144               endif
145            endif
146         endif
147      enddo
148!
149      close(1)
150!
151      return
152      end
153