xref: /phasta/phSolver/common/fillsparse.f (revision 595995161822a203c8467e0e4a253d7bd7d6df32)
1*59599516SKenneth E. Jansen	subroutine fillsparseI(	iens, xKebe,	lhsK,
2*59599516SKenneth E. Jansen     &                               xGoC,      lhsP,
3*59599516SKenneth E. Jansen     1			             row,	col)
4*59599516SKenneth E. Jansenc
5*59599516SKenneth E. Jansenc
6*59599516SKenneth E. Jansenc
7*59599516SKenneth E. Jansen	include "common.h"
8*59599516SKenneth E. Jansen	real*8	xKebe(npro,9,nshl,nshl), xGoC(npro,4,nshl,nshl)
9*59599516SKenneth E. Jansen	integer	ien(npro,nshl),	col(nshg+1), row(nshg*nnz)
10*59599516SKenneth E. Jansen	real*8	lhsK(9,nnz_tot),	lhsP(4,nnz_tot)
11*59599516SKenneth E. Jansenc
12*59599516SKenneth E. Jansen	integer	aa,	b,	c,	e,	i,	k,	n
13*59599516SKenneth E. Jansenc
14*59599516SKenneth E. Jansen	integer sparseloc
15*59599516SKenneth E. Jansen
16*59599516SKenneth E. Jansen	integer iens(npro,nshl)
17*59599516SKenneth E. Jansenc
18*59599516SKenneth E. Jansenc prefer to show explicit absolute value needed for cubic modes and
19*59599516SKenneth E. Jansenc higher rather than inline abs on pointer as in past versions
20*59599516SKenneth E. Jansenc iens is the signed ien array ien is unsigned
21*59599516SKenneth E. Jansenc
22*59599516SKenneth E. Jansen	ien=abs(iens)
23*59599516SKenneth E. Jansenc
24*59599516SKenneth E. Jansenc.... Accumulate the lhs
25*59599516SKenneth E. Jansenc
26*59599516SKenneth E. Jansen	do e = 1, npro ! loop over the elements
27*59599516SKenneth E. Jansen	    do aa = 1, nshl ! loop over the local equation numbers
28*59599516SKenneth E. Jansen		i = ien(e,aa) ! finds the global equation number or
29*59599516SKenneth E. Jansen			      ! block-row of our matrix
30*59599516SKenneth E. Jansen		c = col(i)    ! starting point to look for the matching column
31*59599516SKenneth E. Jansen		n = col(i+1) - c  !length of the list of entries in rowp
32*59599516SKenneth E. Jansen		do b = 1, nshl ! local variable number tangent respect
33*59599516SKenneth E. Jansen			       ! to
34*59599516SKenneth E. Jansenc function that searches row until it finds the match that gives the
35*59599516SKenneth E. Jansenc		   global equation number
36*59599516SKenneth E. Jansen
37*59599516SKenneth E. Jansen		    k = sparseloc( row(c), n, ien(e,b) ) + c-1
38*59599516SKenneth E. Jansenc
39*59599516SKenneth E. Jansenc                                             *         *
40*59599516SKenneth E. Jansenc                   dimension egmass(npro,ndof,nenl,ndof,nenl)
41*59599516SKenneth E. Jansenc
42*59599516SKenneth E. Jansenc compressible      lhsT(1:5,1:5,k)=lhsT(1:5,1:5,k)+egmass(e,1:5,aa,1:5,b)
43*59599516SKenneth E. Jansenc
44*59599516SKenneth E. Jansen		    lhsK(1,k) = lhsK(1,k) + xKebe(e,1,aa,b)
45*59599516SKenneth E. Jansen		    lhsK(2,k) = lhsK(2,k) + xKebe(e,2,aa,b)
46*59599516SKenneth E. Jansen		    lhsK(3,k) = lhsK(3,k) + xKebe(e,3,aa,b)
47*59599516SKenneth E. Jansen		    lhsK(4,k) = lhsK(4,k) + xKebe(e,4,aa,b)
48*59599516SKenneth E. Jansen		    lhsK(5,k) = lhsK(5,k) + xKebe(e,5,aa,b)
49*59599516SKenneth E. Jansen		    lhsK(6,k) = lhsK(6,k) + xKebe(e,6,aa,b)
50*59599516SKenneth E. Jansen		    lhsK(7,k) = lhsK(7,k) + xKebe(e,7,aa,b)
51*59599516SKenneth E. Jansen		    lhsK(8,k) = lhsK(8,k) + xKebe(e,8,aa,b)
52*59599516SKenneth E. Jansen		    lhsK(9,k) = lhsK(9,k) + xKebe(e,9,aa,b)
53*59599516SKenneth E. Jansenc
54*59599516SKenneth E. Jansen		    lhsP(1,k) = lhsP(1,k) + xGoC(e,1,aa,b)
55*59599516SKenneth E. Jansen		    lhsP(2,k) = lhsP(2,k) + xGoC(e,2,aa,b)
56*59599516SKenneth E. Jansen		    lhsP(3,k) = lhsP(3,k) + xGoC(e,3,aa,b)
57*59599516SKenneth E. Jansen		    lhsP(4,k) = lhsP(4,k) + xGoC(e,4,aa,b)
58*59599516SKenneth E. Jansen		enddo
59*59599516SKenneth E. Jansen	    enddo
60*59599516SKenneth E. Jansen	enddo
61*59599516SKenneth E. Jansenc
62*59599516SKenneth E. Jansenc.... end
63*59599516SKenneth E. Jansenc
64*59599516SKenneth E. Jansen	return
65*59599516SKenneth E. Jansen	end
66*59599516SKenneth E. Jansen
67*59599516SKenneth E. Jansen
68*59599516SKenneth E. Jansen	subroutine fillsparseC(	iens, EGmass,	lhsK,
69*59599516SKenneth E. Jansen     1			             row,	col)
70*59599516SKenneth E. Jansenc
71*59599516SKenneth E. Jansenc-----------------------------------------------------------
72*59599516SKenneth E. Jansenc This routine fills up the spasely stored LHS mass matrix
73*59599516SKenneth E. Jansenc
74*59599516SKenneth E. Jansenc Nahid Razmra, Spring 2000. (Sparse Matrix)
75*59599516SKenneth E. Jansenc-----------------------------------------------------------
76*59599516SKenneth E. Jansenc
77*59599516SKenneth E. Jansenc
78*59599516SKenneth E. Jansen
79*59599516SKenneth E. Jansen	include "common.h"
80*59599516SKenneth E. Jansen
81*59599516SKenneth E. Jansen 	real*8	EGmass(npro,nedof,nedof)
82*59599516SKenneth E. Jansen        integer	ien(npro,nshl),	col(nshg+1), row(nnz*nshg)
83*59599516SKenneth E. Jansen        real*8 lhsK(nflow*nflow,nnz_tot)
84*59599516SKenneth E. Jansen
85*59599516SKenneth E. Jansenc
86*59599516SKenneth E. Jansen        integer aa, b, c, e, i, k, n, n1
87*59599516SKenneth E. Jansen        integer f, g, r, s, t
88*59599516SKenneth E. Jansenc
89*59599516SKenneth E. Jansen	integer sparseloc
90*59599516SKenneth E. Jansen
91*59599516SKenneth E. Jansen	integer iens(npro,nshl)
92*59599516SKenneth E. Jansenc
93*59599516SKenneth E. Jansenc prefer to show explicit absolute value needed for cubic modes and
94*59599516SKenneth E. Jansenc higher rather than inline abs on pointer as in past versions
95*59599516SKenneth E. Jansenc iens is the signed ien array ien is unsigned
96*59599516SKenneth E. Jansenc
97*59599516SKenneth E. Jansen	ien=abs(iens)
98*59599516SKenneth E. Jansenc
99*59599516SKenneth E. Jansenc.... Accumulate the lhsK
100*59599516SKenneth E. Jansenc
101*59599516SKenneth E. Jansen	do e = 1, npro
102*59599516SKenneth E. Jansen	    do aa = 1, nshl
103*59599516SKenneth E. Jansen		i = ien(e,aa)
104*59599516SKenneth E. Jansen		c = col(i)
105*59599516SKenneth E. Jansen		n = col(i+1) - c
106*59599516SKenneth E. Jansen		r = (aa-1)*nflow
107*59599516SKenneth E. Jansen		do b = 1, nshl
108*59599516SKenneth E. Jansen		    s = (b-1)*nflow
109*59599516SKenneth E. Jansen                    k = sparseloc( row(c), n, ien(e,b) ) + c-1
110*59599516SKenneth E. Jansenc
111*59599516SKenneth E. Jansen		    do g = 1, nflow
112*59599516SKenneth E. Jansen		       t = (g-1)*nflow
113*59599516SKenneth E. Jansen		       do f = 1, nflow
114*59599516SKenneth E. Jansen
115*59599516SKenneth E. Jansen			  lhsK(t+f,k) = lhsK(t+f,k) + EGmass(e,r+f,s+g)
116*59599516SKenneth E. Jansenc
117*59599516SKenneth E. Jansen
118*59599516SKenneth E. Jansen                       enddo
119*59599516SKenneth E. Jansen                   enddo
120*59599516SKenneth E. Jansen		enddo
121*59599516SKenneth E. Jansen	    enddo
122*59599516SKenneth E. Jansen	enddo
123*59599516SKenneth E. Jansenc
124*59599516SKenneth E. Jansenc.... end
125*59599516SKenneth E. Jansenc
126*59599516SKenneth E. Jansen	return
127*59599516SKenneth E. Jansen	end
128*59599516SKenneth E. Jansen
129*59599516SKenneth E. Jansen	subroutine fillsparseSclr(   iens,      xSebe,	lhsS,
130*59599516SKenneth E. Jansen     1			             row,	col)
131*59599516SKenneth E. Jansenc
132*59599516SKenneth E. Jansenc
133*59599516SKenneth E. Jansenc
134*59599516SKenneth E. Jansen	include "common.h"
135*59599516SKenneth E. Jansen	real*8	xSebe(npro,nshl,nshl)
136*59599516SKenneth E. Jansen	integer	ien(npro,nshl),	col(nshg+1), row(nshg*nnz)
137*59599516SKenneth E. Jansen	real*8	lhsS(nnz_tot)
138*59599516SKenneth E. Jansenc
139*59599516SKenneth E. Jansen	integer	aa,	b,	c,	e,	i,	k,	n
140*59599516SKenneth E. Jansenc
141*59599516SKenneth E. Jansen	integer sparseloc
142*59599516SKenneth E. Jansen
143*59599516SKenneth E. Jansen	integer iens(npro,nshl)
144*59599516SKenneth E. Jansenc
145*59599516SKenneth E. Jansenc prefer to show explicit absolute value needed for cubic modes and
146*59599516SKenneth E. Jansenc higher rather than inline abs on pointer as in past versions
147*59599516SKenneth E. Jansenc iens is the signed ien array ien is unsigned
148*59599516SKenneth E. Jansenc
149*59599516SKenneth E. Jansen	ien=abs(iens)
150*59599516SKenneth E. Jansenc
151*59599516SKenneth E. Jansenc.... Accumulate the lhs
152*59599516SKenneth E. Jansenc
153*59599516SKenneth E. Jansen	do e = 1, npro
154*59599516SKenneth E. Jansen	    do aa = 1, nshl
155*59599516SKenneth E. Jansen		i = ien(e,aa)
156*59599516SKenneth E. Jansen		c = col(i)
157*59599516SKenneth E. Jansen		n = col(i+1) - c
158*59599516SKenneth E. Jansen		do b = 1, nshl
159*59599516SKenneth E. Jansen		    k = sparseloc( row(c), n, ien(e,b) ) + c-1
160*59599516SKenneth E. Jansenc
161*59599516SKenneth E. Jansen		    lhsS(k) = lhsS(k) + xSebe(e,aa,b)
162*59599516SKenneth E. Jansen		enddo
163*59599516SKenneth E. Jansen	    enddo
164*59599516SKenneth E. Jansen	enddo
165*59599516SKenneth E. Jansenc
166*59599516SKenneth E. Jansenc.... end
167*59599516SKenneth E. Jansenc
168*59599516SKenneth E. Jansen	return
169*59599516SKenneth E. Jansen	end
170*59599516SKenneth E. Jansen
171*59599516SKenneth E. Jansen	integer	function sparseloc( list, n, target )
172*59599516SKenneth E. Jansen
173*59599516SKenneth E. Jansenc-----------------------------------------------------------
174*59599516SKenneth E. Jansenc This function finds the location of the non-zero elements
175*59599516SKenneth E. Jansenc of the LHS matrix in the sparsely stored matrix
176*59599516SKenneth E. Jansenc lhsK(nflow*nflow,nnz*numnp)
177*59599516SKenneth E. Jansenc
178*59599516SKenneth E. Jansenc Nahid Razmara, Spring 2000. 	(Sparse Matrix)
179*59599516SKenneth E. Jansenc-----------------------------------------------------------
180*59599516SKenneth E. Jansen
181*59599516SKenneth E. Jansen	integer	list(n),	n,	target
182*59599516SKenneth E. Jansen	integer	rowvl,	rowvh,	rowv
183*59599516SKenneth E. Jansen
184*59599516SKenneth E. Jansenc
185*59599516SKenneth E. Jansenc.... Initialize
186*59599516SKenneth E. Jansenc
187*59599516SKenneth E. Jansen	rowvl = 1
188*59599516SKenneth E. Jansen	rowvh = n + 1
189*59599516SKenneth E. Jansenc
190*59599516SKenneth E. Jansenc.... do a binary search
191*59599516SKenneth E. Jansenc
192*59599516SKenneth E. Jansen100	if ( rowvh-rowvl .gt. 1 ) then
193*59599516SKenneth E. Jansen	    rowv = ( rowvh + rowvl ) / 2
194*59599516SKenneth E. Jansen	    if ( list(rowv) .gt. target ) then
195*59599516SKenneth E. Jansen		rowvh = rowv
196*59599516SKenneth E. Jansen	    else
197*59599516SKenneth E. Jansen		rowvl = rowv
198*59599516SKenneth E. Jansen	    endif
199*59599516SKenneth E. Jansen	    goto 100
200*59599516SKenneth E. Jansen	endif
201*59599516SKenneth E. Jansenc
202*59599516SKenneth E. Jansenc.... return
203*59599516SKenneth E. Jansenc
204*59599516SKenneth E. Jansen	sparseloc = rowvl
205*59599516SKenneth E. Jansenc
206*59599516SKenneth E. Jansen	return
207*59599516SKenneth E. Jansen	end
208*59599516SKenneth E. Jansen
209