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_plhs(lakonl,sm,nelem,ipvar,var)
20!
21!     computation of the pressure element matrix for the element with
22!     the topology in konl
23!
24      implicit none
25!
26      character*8 lakonl
27!
28      integer nelem,i,j,k,nope,mint3d,ipvar(*),index
29!
30      real*8 shp(4,20),sm(8,8),weight,var(*)
31!
32      if(lakonl(4:4).eq.'8') then
33        nope=8
34      elseif(lakonl(4:4).eq.'4') then
35        nope=4
36      elseif(lakonl(4:4).eq.'6') then
37        nope=6
38      endif
39!
40      if(lakonl(4:5).eq.'8R') then
41        mint3d=1
42      elseif(lakonl(4:4).eq.'8') then
43        mint3d=8
44      elseif(lakonl(4:4).eq.'4') then
45        mint3d=1
46      elseif(lakonl(4:5).eq.'6 ') then
47        mint3d=2
48      endif
49!
50!     initialisation of sm
51!
52      do i=1,nope
53        do j=1,nope
54          sm(i,j)=0.d0
55        enddo
56      enddo
57!
58!     computation of the matrix: loop over the Gauss points
59!
60      index=ipvar(nelem)
61      do k=1,mint3d
62!
63!     copying the shape functions, their derivatives and the
64!     Jacobian determinant from field var
65!
66        do j=1,nope
67          do i=1,4
68            index=index+1
69            shp(i,j)=var(index)
70          enddo
71        enddo
72!
73!     Jacobian determinant
74!
75        index=index+1
76        weight=var(index)
77        index=index+1
78        do j=1,nope
79!
80          do i=1,j
81!
82!     lhs pressure matrix (only for semi-implicit incompressible
83!     calculations)
84!
85            sm(i,j)=sm(i,j)
86     &           +(shp(1,i)*shp(1,j)+
87     &           shp(2,i)*shp(2,j)+
88     &           shp(3,i)*shp(3,j))*weight
89          enddo
90        enddo
91      enddo
92!
93      return
94      end
95