1! This software is Copyright (c) 2012-2015 The Regents of the 2! University of California. All Rights Reserved. 3! 4! Permission to copy and modify this software and its documentation 5! for educational, research and non-profit purposes, without fee, 6! and without a written agreement is hereby granted, provided that 7! the above copyright notice, this paragraph and the following three 8! paragraphs appear in all copies. 9! 10! Permission to make commercial use of this software may be obtained 11! by contacting: 12! 13! Technology Transfer Office 14! 9500 Gilman Drive, Mail Code 0910 15! University of California 16! La Jolla, CA 92093-0910 17! (858) 534-5815 18! invent@ucsd.edu 19! 20! This software program and documentation are copyrighted by The 21! Regents of the University of California. The software program and 22! documentation are supplied "as is", without any accompanying 23! services from The Regents. The Regents does not warrant that the 24! operation of the program will be uninterrupted or error-free. The 25! end-user understands that the program was developed for research 26! purposes and is advised not to rely exclusively on the program for 27! any reason. 28! 29! IN NO EVENT SHALL THE UNIVERSITY OF CALIFORNIA BE LIABLE TO ANY 30! PARTY FOR DIRECT, INDIRECT, SPECIAL, INCIDENTAL, OR CONSEQUENTIAL 31! DAMAGES, INCLUDING LOST PROFITS, ARISING OUT OF THE USE OF THIS 32! SOFTWARE AND ITS DOCUMENTATION, EVEN IF THE UNIVERSITY OF 33! CALIFORNIA HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 34! THE UNIVERSITY OF CALIFORNIA SPECIFICALLY DISCLAIMS ANY 35! WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES 36! OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE 37! SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS" BASIS, AND THE 38! UNIVERSITY OF CALIFORNIA HAS NO OBLIGATIONS TO PROVIDE 39! MAINTENANCE, SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS. 40 41 SUBROUTINE svLS_LHS_CREATE(lhs, commu, gnNo, nNo, nnz, gNodes, 42 2 rowPtr, colPtr, nFaces) 43 44 INCLUDE "svLS_STD.h" 45 46 TYPE(svLS_lhsType), INTENT(INOUT) :: lhs 47 TYPE(svLS_commuType), INTENT(IN) :: commu 48 INTEGER, INTENT(IN) :: gnNo, nNo, nnz 49 INTEGER, INTENT(IN) :: gNodes(nNo), rowPtr(nNo+1), colPtr(nnz) 50 INTEGER, INTENT(IN) :: nFaces 51 52 INTEGER i, j, k, a, Ac, ai, s, e, nTasks, tF, maxnNo, ierr, 53 2 stat(MPI_STATUS_SIZE) 54 55 INTEGER comm 56 INTEGER, ALLOCATABLE :: aNodes(:,:), gtlPtr(:), ltg(:), 57 2 part(:), sCount(:), disp(:) 58 59 IF (lhs%foC) THEN 60 PRINT *, "LHS is not free" 61 PRINT *, "You may use svLS_LHS_FREE to free this structure" 62 END IF 63 64 lhs%foC = .TRUE. 65 lhs%gnNo = gnNo 66 lhs%nNo = nNo 67 lhs%nnz = nnz 68 lhs%commu = commu 69 lhs%nFaces = nFaces 70 71 nTasks = commu%nTasks 72 comm = commu%comm 73 tF = commu%tF 74 75 ALLOCATE (lhs%colPtr(nnz), lhs%rowPtr(2,nNo), lhs%diagPtr(nNo), 76 2 lhs%map(nNo), lhs%cS(nTasks), lhs%face(nFaces)) 77 78 IF (nTasks .EQ. 1) THEN 79 DO i=1, nnz 80 lhs%colPtr(i) = colPtr(i) 81 END DO 82 DO Ac=1, nNo 83 s = rowPtr(Ac) 84 e = rowPtr(Ac+1) - 1 85 DO i=s, e 86 a = colPtr(i) 87 IF (Ac .EQ. a) THEN 88 lhs%diagPtr(Ac) = i 89 EXIT 90 END IF 91 END DO 92 93 lhs%rowPtr(1,Ac) = s 94 lhs%rowPtr(2,Ac) = e 95 96 lhs%map(Ac) = Ac 97 END DO 98 99 lhs%mynNo = nNo 100 RETURN 101 END IF 102 103 CALL MPI_ALLREDUCE (nNo, maxnNo, 1, mpint, MPI_MAX, comm, ierr) 104 105 ALLOCATE(aNodes(maxnNo,nTasks), part(maxnNo), sCount(nTasks), 106 2 disp(nTasks), gtlPtr(gnNo), ltg(nNo)) 107 108 part = 0 109 part(1:nNo) = gNodes 110 111 DO i=1, nTasks 112 disp(i) = (i-1)*maxnNo 113 sCount(i) = maxnNo 114 END DO 115 116 CALL MPI_ALLGATHERV(part, maxnNo, mpint, aNodes, sCount, disp, 117 2 mpint, comm, ierr) 118 119 gtlPtr = 0 120 DO a=1, nNo 121 Ac = gNodes(a) 122 gtlPtr(Ac) = a 123 END DO 124 125 DO i=nTasks, 1, -1 126 IF (i .EQ. tF) CYCLE 127 128 DO a=1, maxnNo 129 Ac = aNodes(a,i) 130 IF (Ac .EQ. 0) EXIT 131 ai = gtlPtr(Ac) 132 IF (ai .NE. 0) THEN 133 IF (aNodes(ai,tF) .NE. 0) THEN 134 IF (i .GT. tF) aNodes(ai,tF) = 0 135 ELSE 136 aNodes(a,i) = 0 137 END IF 138 ELSE 139 aNodes(a,i) = 0 140 END IF 141 END DO 142 END DO 143 144 j = 1 145 lhs%cS(1)%ptr = 1 146 DO i=1, nTasks 147 lhs%cS(i)%n = 0 148 lhs%cS(i)%ptr = j 149 IF (i.NE.tF .AND. i.NE.1) THEN 150 lhs%cS(i)%ptr = lhs%cS(i-1)%ptr + lhs%cS(i-1)%n 151 END IF 152 153 DO a=1, maxnNo 154 Ac = aNodes(a,i) 155 IF (Ac .NE. 0) THEN 156 lhs%cS(i)%n = lhs%cS(i)%n + 1 157 ai = gtlPtr(Ac) 158 IF (i.GT.tF .OR. aNodes(ai,tF).NE.0) THEN 159 ltg(j) = Ac 160 j = j + 1 161 aNodes(ai,tF) = 0 162 END IF 163 END IF 164 END DO 165 166 IF (i .LT. tF) THEN 167 lhs%cS(i)%tag = nTasks*i + tF 168 ELSE 169 lhs%cS(i)%tag = nTasks*tF + i 170 END IF 171 IF (lhs%cS(i)%n .EQ. 0) lhs%cS(i)%tag = 0 172 END DO 173 174 lhs%cS(tF)%tag = 0 175 lhs%mynNo = lhs%cS(tF)%ptr + lhs%cS(tF)%n - 1 176 177 gtlPtr = 0 178 DO a=1, nNo 179 Ac = ltg(a) 180 gtlPtr(Ac) = a 181 END DO 182 DO a=1, nNo 183 Ac = gNodes(a) 184 lhs%map(a) = gtlPtr(Ac) 185 END DO 186 187 DEALLOCATE(aNodes, part, gtlPtr, sCount, disp) 188 189 DO a=1, nNo 190 Ac = lhs%map(a) 191 lhs%rowPtr(1,Ac) = rowPtr(a) 192 lhs%rowPtr(2,Ac) = rowPtr(a+1) - 1 193 END DO 194 195 DO i=1, nnz 196 lhs%colPtr(i) = lhs%map(colPtr(i)) 197 END DO 198 199 DO Ac=1, nNo 200 DO i=lhs%rowPtr(1,Ac), lhs%rowPtr(2,Ac) 201 a = lhs%colPtr(i) 202 IF (Ac .EQ. a) THEN 203 lhs%diagPtr(Ac) = i 204 EXIT 205 END IF 206 END DO 207 END DO 208 209 IF (tF .NE. 1) THEN 210 i = tF - 1 211 i = lhs%cS(i)%ptr + lhs%cS(i)%n - 1 212 ALLOCATE(part(i)) 213 END IF 214 215 DO i=1, nTasks 216 lhs%cS(i)%nBl = 0 217 IF (lhs%cS(i)%tag .NE. 0) THEN 218 s = lhs%cS(i)%ptr 219 e = s + lhs%cS(i)%n - 1 220 IF (i .LT. tF) THEN 221 CALL MPI_RECV(part(s:e), lhs%cS(i)%n, mpint, i-1, 222 2 lhs%cS(i)%tag, comm, stat, ierr) 223 224 k = 0 225 DO j=s, e 226 k = k + 1 227 IF (part(j).NE.ltg(k) .OR. j.EQ.s) THEN 228 lhs%cS(i)%nBl = lhs%cS(i)%nBl + 1 229 DO k=1, lhs%cS(tF)%ptr 230 IF (part(j) .EQ. ltg(k)) EXIT 231 END DO 232 END IF 233 END DO 234 a = lhs%cS(i)%nBl 235 ALLOCATE(lhs%cS(i)%blPtr(a), lhs%cS(i)%blN(a)) 236 237 k = 0 238 a = 0 239 DO j=s, e 240 k = k + 1 241 IF (part(j).NE.ltg(k) .OR. j.EQ.s) THEN 242 a = a + 1 243 lhs%cS(i)%blN(a) = 1 244 DO k=1, lhs%cS(tF)%ptr 245 IF (part(j) .EQ. ltg(k)) THEN 246 lhs%cS(i)%blPtr(a) = k 247 EXIT 248 END IF 249 END DO 250 ELSE 251 lhs%cS(i)%blN(a) = lhs%cS(i)%blN(a) + 1 252 END IF 253 END DO 254 ELSE 255 CALL MPI_SEND(ltg(s:e), lhs%cS(i)%n, mpint, i-1, 256 2 lhs%cS(i)%tag, comm, stat, ierr) 257 END IF 258 END IF 259 END DO 260 261 IF (ALLOCATED(part)) DEALLOCATE(part) 262 DEALLOCATE(ltg) 263 264 RETURN 265 END SUBROUTINE svLS_LHS_CREATE 266 267!==================================================================== 268 269 SUBROUTINE svLS_LHS_FREE(lhs) 270 271 INCLUDE "svLS_STD.h" 272 273 TYPE(svLS_lhsType), INTENT(INOUT) :: lhs 274 275 INTEGER faIn, i 276 277 IF (.NOT.lhs%foC) THEN 278 PRINT *, 'Cannot free LHS' 279 PRINT *, 'It is not created yet' 280 STOP 281 END IF 282 283 DO faIn = 1, lhs%nFaces 284 IF (lhs%face(faIn)%foC) CALL svLS_BC_FREE(lhs, faIn) 285 END DO 286 287 DO i=1, lhs%commu%nTasks 288 IF (ALLOCATED(lhs%cS(i)%blPtr)) THEN 289 DEALLOCATE(lhs%cS(i)%blPtr, lhs%cS(i)%blN) 290 END IF 291 END DO 292 293 lhs%foC = .FALSE. 294 lhs%gnNo = 0 295 lhs%nNo = 0 296 lhs%nnz = 0 297 lhs%nFaces = 0 298 299 DEALLOCATE (lhs%colPtr, lhs%rowPtr, lhs%diagPtr, lhs%map, lhs%cS, 300 2 lhs%face) 301 302 RETURN 303 END SUBROUTINE svLS_LHS_FREE 304 305!==================================================================== 306 307 SUBROUTINE svLS_LHS_CREATE_C(pLHS, commu, gnNo, nNo, nnz, gNodes, 308 2 rowPtr, colPtr, nFaces) 309 310 INCLUDE "svLS_STD.h" 311 312 TYPE(svLS_lhsType), POINTER, INTENT(OUT) :: pLHS 313 TYPE(svLS_commuType), INTENT(IN) :: commu 314 INTEGER, INTENT(IN) :: gnNo, nNo, nnz 315 INTEGER, INTENT(IN) :: gNodes(nNo), rowPtr(nNo+1), colPtr(nnz) 316 INTEGER, INTENT(IN) :: nFaces 317 318 TYPE(svLS_lhsType), TARGET, SAVE :: lhs 319 320 CALL svLS_LHS_CREATE(lhs, commu, gnNo, nNo, nnz, gNodes, 321 2 rowPtr, colPtr, nFaces) 322 323 pLHS => lhs 324 325 RETURN 326 END SUBROUTINE svLS_LHS_CREATE_C 327