7submodule(
linalg) linalg_sorting
12 module subroutine sort_dbl_array(x, ascend)
14 real(real64),
intent(inout),
dimension(:) :: x
15 logical,
intent(in),
optional :: ascend
19 integer(int32) :: n, info
22 if (
present(ascend))
then
34 call dlasrt(id, n, x, info)
38 module subroutine sort_dbl_array_ind(x, ind, ascend, err)
40 real(real64),
intent(inout),
dimension(:) :: x
41 integer(int32),
intent(inout),
dimension(:) :: ind
42 logical,
intent(in),
optional :: ascend
43 class(errors),
intent(inout),
optional,
target :: err
46 class(errors),
pointer :: errmgr
47 type(errors),
target :: deferr
48 character(len = :),
allocatable :: errmsg
54 if (
present(err))
then
59 if (
present(ascend))
then
66 if (
size(ind) /= n)
then
67 allocate(
character(len = 256) :: errmsg)
69 "Expected the tracking array to be of size ", n, &
70 ", but found an array of size ",
size(ind),
"."
71 call errmgr%report_error(
"sort_dbl_array_ind", trim(errmsg), &
78 call qsort_dbl_ind(dir, x, ind)
81100
format(a, i0, a, i0, a)
85 module subroutine sort_cmplx_array(x, ascend)
87 complex(real64),
intent(inout),
dimension(:) :: x
88 logical,
intent(in),
optional :: ascend
94 if (
present(ascend))
then
101 call qsort_cmplx(dir, x)
105 module subroutine sort_cmplx_array_ind(x, ind, ascend, err)
107 complex(real64),
intent(inout),
dimension(:) :: x
108 integer(int32),
intent(inout),
dimension(:) :: ind
109 logical,
intent(in),
optional :: ascend
110 class(errors),
intent(inout),
optional,
target :: err
113 class(errors),
pointer :: errmgr
114 type(errors),
target :: deferr
115 character(len = :),
allocatable :: errmsg
121 if (
present(err))
then
126 if (
present(ascend))
then
133 if (
size(ind) /= n)
then
134 allocate(
character(len = 256) :: errmsg)
136 "Expected the tracking array to be of size ", n, &
137 ", but found an array of size ",
size(ind),
"."
138 call errmgr%report_error(
"sort_cmplx_array_ind", trim(errmsg), &
145 call qsort_cmplx_ind(dir, x, ind)
148100
format(a, i0, a, i0, a)
152 module subroutine sort_eigen_cmplx(vals, vecs, ascend, err)
154 complex(real64),
intent(inout),
dimension(:) :: vals
155 complex(real64),
intent(inout),
dimension(:,:) :: vecs
156 logical,
intent(in),
optional :: ascend
157 class(errors),
intent(inout),
optional,
target :: err
160 class(errors),
pointer :: errmgr
161 type(errors),
target :: deferr
162 character(len = :),
allocatable :: errmsg
163 integer(int32) :: i, n, flag
165 integer(int32),
allocatable,
dimension(:) :: ind
168 if (
present(err))
then
173 if (
present(ascend))
then
181 if (
size(vecs, 1) /= n .or.
size(vecs, 2) /= n)
then
183 allocate(
character(len = 256) :: errmsg)
185 "Expected the eigenvector matrix to be of size ", n, &
186 "-by-", n,
", but found a matrix of size ",
size(vecs, 1), &
187 "-by-",
size(vecs, 2),
"."
188 call errmgr%report_error(
"sort_eigen_cmplx", trim(errmsg), &
193 allocate(ind(n), stat = flag)
195 call errmgr%report_error(
"sort_eigen_cmplx", &
196 "Insufficient memory available.", la_out_of_memory_error)
204 call qsort_cmplx_ind(dir, vals, ind)
211100
format(a, i0, a, i0, a, i0, a, i0, a)
215 module subroutine sort_eigen_dbl(vals, vecs, ascend, err)
217 real(real64),
intent(inout),
dimension(:) :: vals
218 real(real64),
intent(inout),
dimension(:,:) :: vecs
219 logical,
intent(in),
optional :: ascend
220 class(errors),
intent(inout),
optional,
target :: err
223 class(errors),
pointer :: errmgr
224 type(errors),
target :: deferr
225 character(len = :),
allocatable :: errmsg
226 integer(int32) :: i, n, flag
228 integer(int32),
allocatable,
dimension(:) :: ind
231 if (
present(err))
then
236 if (
present(ascend))
then
244 if (
size(vecs, 1) /= n .or.
size(vecs, 2) /= n)
then
246 allocate(
character(len = 256) :: errmsg)
248 "Expected the eigenvector matrix to be of size ", n, &
249 "-by-", n,
", but found a matrix of size ",
size(vecs, 1), &
250 "-by-",
size(vecs, 2),
"."
251 call errmgr%report_error(
"sort_eigen_dbl", trim(errmsg), &
256 allocate(ind(n), stat = flag)
258 call errmgr%report_error(
"sort_eigen_dbl", &
259 "Insufficient memory available.", la_out_of_memory_error)
267 call qsort_dbl_ind(dir, vals, ind)
274100
format(a, i0, a, i0, a, i0, a, i0, a)
293 recursive subroutine qsort_dbl_ind(ascend, x, ind)
295 logical,
intent(in) :: ascend
296 real(real64),
intent(inout),
dimension(:) :: x
297 integer(int32),
intent(inout),
dimension(:) :: ind
303 if (
size(x) > 1)
then
304 call dbl_partition_ind(ascend, x, ind, iq)
305 call qsort_dbl_ind(ascend, x(:iq-1), ind(:iq-1))
306 call qsort_dbl_ind(ascend, x(iq:), ind(iq:))
326 subroutine dbl_partition_ind(ascend, x, ind, marker)
328 logical,
intent(in) :: ascend
329 real(real64),
intent(inout),
dimension(:) :: x
330 integer(int32),
intent(inout),
dimension(:) :: ind
331 integer(int32),
intent(out) :: marker
334 integer(int32) :: i, j, itemp
335 real(real64) :: temp, pivot
346 if (x(j) <= pivot)
exit
351 if (x(i) >= pivot)
exit
363 else if (i == j)
then
376 if (x(j) >= pivot)
exit
381 if (x(i) <= pivot)
exit
393 else if (i == j)
then
419 recursive subroutine qsort_cmplx(ascend, x)
421 logical,
intent(in) :: ascend
422 complex(real64),
intent(inout),
dimension(:) :: x
428 if (
size(x) > 1)
then
429 call cmplx_partition(ascend, x, iq)
430 call qsort_cmplx(ascend, x(:iq-1))
431 call qsort_cmplx(ascend, x(iq:))
452 subroutine cmplx_partition(ascend, x, marker)
454 logical,
intent(in) :: ascend
455 complex(real64),
intent(inout),
dimension(:) :: x
456 integer(int32),
intent(out) :: marker
459 integer(int32) :: i, j
460 complex(real64) :: temp
461 real(real64) :: pivot
464 pivot = real(x(1), real64)
472 if (real(x(j), real64) <= pivot)
exit
477 if (real(x(i), real64) >= pivot)
exit
485 else if (i == j)
then
498 if (real(x(j), real64) >= pivot)
exit
503 if (real(x(i), real64) <= pivot)
exit
511 else if (i == j)
then
540 recursive subroutine qsort_cmplx_ind(ascend, x, ind)
542 logical,
intent(in) :: ascend
543 complex(real64),
intent(inout),
dimension(:) :: x
544 integer(int32),
intent(inout),
dimension(:) :: ind
550 if (
size(x) > 1)
then
551 call cmplx_partition_ind(ascend, x, ind, iq)
552 call qsort_cmplx_ind(ascend, x(:iq-1), ind(:iq-1))
553 call qsort_cmplx_ind(ascend, x(iq:), ind(iq:))
577 subroutine cmplx_partition_ind(ascend, x, ind, marker)
579 logical,
intent(in) :: ascend
580 complex(real64),
intent(inout),
dimension(:) :: x
581 integer(int32),
intent(inout),
dimension(:) :: ind
582 integer(int32),
intent(out) :: marker
585 integer(int32) :: i, j, itemp
586 complex(real64) :: temp
587 real(real64) :: pivot
590 pivot = real(x(1), real64)
598 if (real(x(j), real64) <= pivot)
exit
603 if (real(x(i), real64) >= pivot)
exit
615 else if (i == j)
then
628 if (real(x(j), real64) >= pivot)
exit
633 if (real(x(i), real64) <= pivot)
exit
645 else if (i == j)
then
Provides a set of common linear algebra routines.