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 e_c3d_v2rhs(konl,lakonl,bb,nelem,v,dtimef,mi,
20     &     ipvar,var,nk)
21!
22!     computation of the velocity element matrix and rhs for the element with
23!     element with the topology in konl: step 3 (correction **)
24!
25!     bb: rhs
26!
27      implicit none
28!
29      character*8 lakonl
30!
31      integer konl(8),nelem,i,j,k,nope,mint3d,mi(*),ipvar(*),index,nk
32!
33      real*8 shp(4,8),bb(3,8),xsjmod,vl(0:mi(2),8),v(nk,0:mi(2)),
34     &     var(*),dtimef,xsj,ddpress(3)
35!
36      if(lakonl(4:4).eq.'4') then
37        nope=4
38        mint3d=1
39      elseif(lakonl(4:4).eq.'6') then
40        nope=6
41        mint3d=2
42      elseif(lakonl(4:5).eq.'8R') then
43        nope=8
44        mint3d=1
45      elseif(lakonl(4:4).eq.'8') then
46        nope=8
47        mint3d=8
48      endif
49!
50!     initialisation of the rhs
51!
52      do i=1,nope
53        do j=1,3
54          bb(j,i)=0.d0
55        enddo
56      enddo
57!
58!     change in pressure
59!
60      do i=1,nope
61        vl(4,i)=v(konl(i),4)
62      enddo
63!
64!     computation of the matrix: loop over the Gauss points
65!
66      index=ipvar(nelem)
67      do k=1,mint3d
68!
69!     copying the shape functions, their derivatives and the
70!     Jacobian determinant from field var
71!
72        do j=1,nope
73          do i=1,4
74            index=index+1
75            shp(i,j)=var(index)
76          enddo
77        enddo
78        index=index+1
79        xsj=var(index)
80!
81        xsjmod=dtimef*xsj
82!
83        index=index+1
84!
85!     only for the semi-implicit procedure: calculate ddpress;
86!
87        do j=1,3
88          ddpress(j)=0.d0
89        enddo
90        do i=1,nope
91          do j=1,3
92            ddpress(j)=ddpress(j)+shp(j,i)*vl(4,i)
93          enddo
94        enddo
95!
96!     3. contribution to the second part of the
97!     momentum equation
98!
99        do j=1,nope
100          bb(1,j)=bb(1,j)-xsjmod*shp(4,j)*ddpress(1)
101          bb(2,j)=bb(2,j)-xsjmod*shp(4,j)*ddpress(2)
102          bb(3,j)=bb(3,j)-xsjmod*shp(4,j)*ddpress(3)
103        enddo
104!
105      enddo
106!
107c      write(*,*) 'e_c3d_v2rhs'
108c      do k=1,8
109c        write(*,*) nelem,k,(bb(j,k),j=1,3)
110c      enddo
111!
112      return
113      end
114