1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
|
subroutine spreshape(A_m, A_n, A_mnel, A_icol, A_R, A_I,
$ B_m, B_n, B_mnel, B_icol, B_R, B_I,
$ nel, it, ij, ptr, p)
*
* PURPOSE
* reshape a sparse matrix, ie, computes the
* sparse matrix B such that : B = matrix(A, m_B, n_B)
* with m_B * n_B = m_A * n_A (this condition being
* verified at the calling level)
*
* depending upon the argument it, A may be (and so B) :
* a/ a boolean sparse (it = -1). In this case the arrays
* R_A, I_A, R_B, I_B are not used.
* b/ a real sparse (it = 0). The arrays I_A and I_B are not used.
*
* c/ a complex sparse (it = 1)
*
* ptr : work array of size max(A_n,B_m)+1
* p : work arrays of size (nel)
* ij: work array of size (2,nel)
*
* AUTHOR
* Bruno Pincon
*
implicit none
integer A_m, A_n, A_mnel(*), A_icol(*)
integer B_m, B_n, B_mnel(*), B_icol(*)
integer nel, it, ij(2,*), ptr(*), p(*)
double precision A_R(*), A_I(*), B_R(*), B_I(*)
* internal vars
integer i,j, k, l, ka, kb, num
* 1/ stack the matrix indices by column order
call iset(A_n+1, 0, ptr, 1)
do k = 1, nel
j = A_icol(k)
ptr(j+1) = ptr(j+1) + 1
enddo
ptr(1) = 1
do j = 2, A_n
ptr(j) = ptr(j) + ptr(j-1)
enddo
ka = 0
do i = 1, A_m
do l = 1, A_mnel(i)
ka = ka + 1
j = A_icol(ka)
k = ptr(j)
ij(1,k) = i
ij(2,k) = j
p(k) = ka
ptr(j) = ptr(j) + 1
enddo
enddo
* 2/ applies the reshaping onto the indices and computes B_mnel
call iset(B_m, 0, B_mnel, 1) ! init B_mnel to 0
do k = 1, nel
num = ij(1,k) + (ij(2,k)-1)*A_m - 1
ij(2,k) = num/B_m + 1
ij(1,k) = num+1 - B_m*(ij(2,k)-1)
B_mnel(ij(1,k)) = B_mnel(ij(1,k)) + 1
enddo
* 3/ computes the others parts of B
call sz2ptr(B_mnel, B_m, ptr)
do k = 1, nel
i = ij(1,k)
kb = ptr(i)
B_icol(kb) = ij(2,k)
if ( it .ge. 0 ) B_R(kb) = A_R(p(k))
if ( it .eq. 1 ) B_I(kb) = A_I(p(k))
ptr(i) = ptr(i) + 1
enddo
end
|
