|
| 1 | + |
| 2 | +function compress_ids(cell_to_bgcell::AbstractVector{<:Integer}) |
| 3 | + nccells = 0 |
| 4 | + c = -1 |
| 5 | + for bgcell in cell_to_bgcell |
| 6 | + if c != bgcell |
| 7 | + nccells += 1 |
| 8 | + c = bgcell |
| 9 | + end |
| 10 | + end |
| 11 | + ccell_to_first_cell = zeros(Int32,nccells+1) |
| 12 | + ccell = 0 |
| 13 | + c = -1 |
| 14 | + for (cell,bgcell) in enumerate(cell_to_bgcell) |
| 15 | + if c != bgcell |
| 16 | + ccell += 1 |
| 17 | + c = bgcell |
| 18 | + ccell_to_first_cell[ccell] = cell |
| 19 | + end |
| 20 | + end |
| 21 | + ccell_to_first_cell[end] = length(cell_to_bgcell)+1 |
| 22 | + ccell_to_first_cell |
| 23 | +end |
| 24 | + |
| 25 | +function compress_ids(cell_ids,cell_to_bgcell::AbstractVector{<:Integer}) |
| 26 | + ccell_to_first_cell = compress_ids(cell_to_bgcell) |
| 27 | + nccells = length(ccell_to_first_cell)-1 |
| 28 | + ccell_to_cell = view(ccell_to_first_cell,1:nccells) |
| 29 | + ccell_ids = view(cell_ids,ccell_to_cell) |
| 30 | + ccell_ids |
| 31 | +end |
| 32 | + |
| 33 | +function compress_contributions( |
| 34 | + cell_to_mat, |
| 35 | + cell_to_bgcell::AbstractVector{<:Integer}, |
| 36 | + ccell_to_first_cell::AbstractVector{<:Integer}=compress_ids(cell_to_bgcell)) |
| 37 | + |
| 38 | + nccells = length(ccell_to_first_cell)-1 |
| 39 | + ccell_to_cell = view(ccell_to_first_cell,1:nccells) |
| 40 | + ccell_to_mat = CompressedCellArray(cell_to_mat,ccell_to_first_cell) |
| 41 | + ccell_to_mat |
| 42 | +end |
| 43 | + |
| 44 | +struct CompressedCellArray{T,A,B} <: AbstractVector{T} |
| 45 | + cell_to_array::A |
| 46 | + ccell_to_first_cell::B |
| 47 | + function CompressedCellArray(cell_to_array::A, ccell_to_first_cell::B) where {A,B} |
| 48 | + T = eltype(cell_to_array) |
| 49 | + new{T,A,B}(cell_to_array,ccell_to_first_cell) |
| 50 | + end |
| 51 | +end |
| 52 | + |
| 53 | +Base.size(a::CompressedCellArray) = (length(a.ccell_to_first_cell)-1,) |
| 54 | +Base.IndexStyle(::Type{<:CompressedCellArray}) = IndexLinear() |
| 55 | + |
| 56 | +function Base.getindex(a::CompressedCellArray,ccell::Integer) |
| 57 | + cell_first = a.ccell_to_first_cell[ccell] |
| 58 | + cell_last = a.ccell_to_first_cell[ccell+1]-1 |
| 59 | + array = a.cell_to_array[cell_first] |
| 60 | + cell_range = (cell_first+1):cell_last |
| 61 | + _compress!(array,a.cell_to_array,cell_range) |
| 62 | +end |
| 63 | + |
| 64 | +@inline function _compress!(array,cell_to_array,cell_range) |
| 65 | + r = copy(array) |
| 66 | + for cell in cell_range |
| 67 | + carray = cell_to_array[cell] |
| 68 | + for k in eachindex(r) |
| 69 | + r[k] = r[k] + carray[k] |
| 70 | + end |
| 71 | + end |
| 72 | + r |
| 73 | +end |
| 74 | + |
| 75 | +@inline function _compress!(matvec::Tuple,cell_to_matvec,cell_range) |
| 76 | + mat, vec = matvec |
| 77 | + rm = copy(mat) |
| 78 | + rv = copy(vec) |
| 79 | + for cell in cell_range |
| 80 | + cmat, cvec = cell_to_matvec[cell] |
| 81 | + for k in eachindex(mat) |
| 82 | + rm[k] = rm[k] + cmat[k] |
| 83 | + end |
| 84 | + for i in eachindex(vec) |
| 85 | + rv[i] = rv[i] + cvec[i] |
| 86 | + end |
| 87 | + end |
| 88 | + (rm,rv) |
| 89 | +end |
| 90 | + |
| 91 | +function array_cache(a::CompressedCellArray) |
| 92 | + array = testitem(a.cell_to_array) |
| 93 | + _compress_cache(array,a) |
| 94 | +end |
| 95 | + |
| 96 | +function _compress_cache(array,a) |
| 97 | + CachedArray(copy(array)), array_cache(a.cell_to_array) |
| 98 | +end |
| 99 | + |
| 100 | +function _compress_cache(matvec::Tuple,a) |
| 101 | + mat, vec = matvec |
| 102 | + cmatvec = (CachedArray(copy(mat)), CachedArray(copy(vec))) |
| 103 | + cmatvec, array_cache(a.cell_to_array) |
| 104 | +end |
| 105 | + |
| 106 | +function getindex!(cache,a::CompressedCellArray,ccell::Integer) |
| 107 | + c, ccache = cache |
| 108 | + cell_first = a.ccell_to_first_cell[ccell] |
| 109 | + cell_last = a.ccell_to_first_cell[ccell+1]-1 |
| 110 | + array = getindex!(ccache,a.cell_to_array,cell_first) |
| 111 | + cell_range = (cell_first+1):cell_last |
| 112 | + _compress!(c,array,ccache,a.cell_to_array,cell_range) |
| 113 | +end |
| 114 | + |
| 115 | +@inline function _compress!(c,array,ccache,cell_to_array,cell_range) |
| 116 | + setsize!(c,size(array)) |
| 117 | + r = c.array |
| 118 | + copyto!(r,array) |
| 119 | + for cell in cell_range |
| 120 | + carray = getindex!(ccache,cell_to_array,cell) |
| 121 | + for k in eachindex(array) |
| 122 | + r[k] = r[k] + carray[k] |
| 123 | + end |
| 124 | + end |
| 125 | + r |
| 126 | +end |
| 127 | + |
| 128 | +@inline function _compress!(c,matvec::Tuple,ccache,cell_to_matvec,cell_range) |
| 129 | + cm, cv = c |
| 130 | + mat, vec = matvec |
| 131 | + setsize!(cm,size(mat)) |
| 132 | + setsize!(cv,size(vec)) |
| 133 | + rm = cm.array |
| 134 | + rv = cv.array |
| 135 | + copyto!(rm,mat) |
| 136 | + copyto!(rv,vec) |
| 137 | + for cell in cell_range |
| 138 | + cmat, cvec = getindex!(ccache,cell_to_matvec,cell) |
| 139 | + for k in eachindex(mat) |
| 140 | + rm[k] = rm[k] + cmat[k] |
| 141 | + end |
| 142 | + for i in eachindex(vec) |
| 143 | + rv[i] = rv[i] + cvec[i] |
| 144 | + end |
| 145 | + end |
| 146 | + (rm, rv) |
| 147 | +end |
0 commit comments