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