7submodule(
linalg) linalg_sorting
14 module subroutine sort_dbl_array(x, ascend)
16 real(real64),
intent(inout),
dimension(:) :: x
17 logical,
intent(in),
optional :: ascend
21 integer(int32) :: n, info
24 if (
present(ascend))
then
36 call dlasrt(id, n, x, info)
40 module subroutine sort_dbl_array_ind(x, ind, ascend, err)
42 real(real64),
intent(inout),
dimension(:) :: x
43 integer(int32),
intent(inout),
dimension(:) :: ind
44 logical,
intent(in),
optional :: ascend
45 class(errors),
intent(inout),
optional,
target :: err
48 class(errors),
pointer :: errmgr
49 type(errors),
target :: deferr
50 character(len = :),
allocatable :: errmsg
56 if (
present(err))
then
61 if (
present(ascend))
then
68 if (
size(ind) /= n)
then
69 allocate(
character(len = 256) :: errmsg)
71 "Expected the tracking array to be of size ", n, &
72 ", but found an array of size ",
size(ind),
"."
73 call errmgr%report_error(
"sort_dbl_array_ind", trim(errmsg), &
80 call qsort_dbl_ind(dir, x, ind)
83100
format(a, i0, a, i0, a)
87 module subroutine sort_cmplx_array(x, ascend)
89 complex(real64),
intent(inout),
dimension(:) :: x
90 logical,
intent(in),
optional :: ascend
96 if (
present(ascend))
then
103 call qsort_cmplx(dir, x)
107 module subroutine sort_cmplx_array_ind(x, ind, ascend, err)
109 complex(real64),
intent(inout),
dimension(:) :: x
110 integer(int32),
intent(inout),
dimension(:) :: ind
111 logical,
intent(in),
optional :: ascend
112 class(errors),
intent(inout),
optional,
target :: err
115 class(errors),
pointer :: errmgr
116 type(errors),
target :: deferr
117 character(len = :),
allocatable :: errmsg
123 if (
present(err))
then
128 if (
present(ascend))
then
135 if (
size(ind) /= n)
then
136 allocate(
character(len = 256) :: errmsg)
138 "Expected the tracking array to be of size ", n, &
139 ", but found an array of size ",
size(ind),
"."
140 call errmgr%report_error(
"sort_cmplx_array_ind", trim(errmsg), &
147 call qsort_cmplx_ind(dir, x, ind)
150100
format(a, i0, a, i0, a)
154 module subroutine sort_eigen_cmplx(vals, vecs, ascend, err)
156 complex(real64),
intent(inout),
dimension(:) :: vals
157 complex(real64),
intent(inout),
dimension(:,:) :: vecs
158 logical,
intent(in),
optional :: ascend
159 class(errors),
intent(inout),
optional,
target :: err
162 class(errors),
pointer :: errmgr
163 type(errors),
target :: deferr
164 character(len = :),
allocatable :: errmsg
165 integer(int32) :: i, n, flag
167 integer(int32),
allocatable,
dimension(:) :: ind
170 if (
present(err))
then
175 if (
present(ascend))
then
183 if (
size(vecs, 1) /= n .or.
size(vecs, 2) /= n)
then
185 allocate(
character(len = 256) :: errmsg)
187 "Expected the eigenvector matrix to be of size ", n, &
188 "-by-", n,
", but found a matrix of size ",
size(vecs, 1), &
189 "-by-",
size(vecs, 2),
"."
190 call errmgr%report_error(
"sort_eigen_cmplx", trim(errmsg), &
195 allocate(ind(n), stat = flag)
197 call errmgr%report_error(
"sort_eigen_cmplx", &
198 "Insufficient memory available.", la_out_of_memory_error)
206 call qsort_cmplx_ind(dir, vals, ind)
213100
format(a, i0, a, i0, a, i0, a, i0, a)
217 module subroutine sort_eigen_dbl(vals, vecs, ascend, err)
219 real(real64),
intent(inout),
dimension(:) :: vals
220 real(real64),
intent(inout),
dimension(:,:) :: vecs
221 logical,
intent(in),
optional :: ascend
222 class(errors),
intent(inout),
optional,
target :: err
225 class(errors),
pointer :: errmgr
226 type(errors),
target :: deferr
227 character(len = :),
allocatable :: errmsg
228 integer(int32) :: i, n, flag
230 integer(int32),
allocatable,
dimension(:) :: ind
233 if (
present(err))
then
238 if (
present(ascend))
then
246 if (
size(vecs, 1) /= n .or.
size(vecs, 2) /= n)
then
248 allocate(
character(len = 256) :: errmsg)
250 "Expected the eigenvector matrix to be of size ", n, &
251 "-by-", n,
", but found a matrix of size ",
size(vecs, 1), &
252 "-by-",
size(vecs, 2),
"."
253 call errmgr%report_error(
"sort_eigen_dbl", trim(errmsg), &
258 allocate(ind(n), stat = flag)
260 call errmgr%report_error(
"sort_eigen_dbl", &
261 "Insufficient memory available.", la_out_of_memory_error)
269 call qsort_dbl_ind(dir, vals, ind)
276100
format(a, i0, a, i0, a, i0, a, i0, a)
295 recursive subroutine qsort_dbl_ind(ascend, x, ind)
297 logical,
intent(in) :: ascend
298 real(real64),
intent(inout),
dimension(:) :: x
299 integer(int32),
intent(inout),
dimension(:) :: ind
305 if (
size(x) > 1)
then
306 call dbl_partition_ind(ascend, x, ind, iq)
307 call qsort_dbl_ind(ascend, x(:iq-1), ind(:iq-1))
308 call qsort_dbl_ind(ascend, x(iq:), ind(iq:))
328 subroutine dbl_partition_ind(ascend, x, ind, marker)
330 logical,
intent(in) :: ascend
331 real(real64),
intent(inout),
dimension(:) :: x
332 integer(int32),
intent(inout),
dimension(:) :: ind
333 integer(int32),
intent(out) :: marker
336 integer(int32) :: i, j, itemp
337 real(real64) :: temp, pivot
348 if (x(j) <= pivot)
exit
353 if (x(i) >= pivot)
exit
365 else if (i == j)
then
378 if (x(j) >= pivot)
exit
383 if (x(i) <= pivot)
exit
395 else if (i == j)
then
421 recursive subroutine qsort_cmplx(ascend, x)
423 logical,
intent(in) :: ascend
424 complex(real64),
intent(inout),
dimension(:) :: x
430 if (
size(x) > 1)
then
431 call cmplx_partition(ascend, x, iq)
432 call qsort_cmplx(ascend, x(:iq-1))
433 call qsort_cmplx(ascend, x(iq:))
454 subroutine cmplx_partition(ascend, x, marker)
456 logical,
intent(in) :: ascend
457 complex(real64),
intent(inout),
dimension(:) :: x
458 integer(int32),
intent(out) :: marker
461 integer(int32) :: i, j
462 complex(real64) :: temp
463 real(real64) :: pivot
466 pivot = real(x(1), real64)
474 if (real(x(j), real64) <= pivot)
exit
479 if (real(x(i), real64) >= pivot)
exit
487 else if (i == j)
then
500 if (real(x(j), real64) >= pivot)
exit
505 if (real(x(i), real64) <= pivot)
exit
513 else if (i == j)
then
542 recursive subroutine qsort_cmplx_ind(ascend, x, ind)
544 logical,
intent(in) :: ascend
545 complex(real64),
intent(inout),
dimension(:) :: x
546 integer(int32),
intent(inout),
dimension(:) :: ind
552 if (
size(x) > 1)
then
553 call cmplx_partition_ind(ascend, x, ind, iq)
554 call qsort_cmplx_ind(ascend, x(:iq-1), ind(:iq-1))
555 call qsort_cmplx_ind(ascend, x(iq:), ind(iq:))
579 subroutine cmplx_partition_ind(ascend, x, ind, marker)
581 logical,
intent(in) :: ascend
582 complex(real64),
intent(inout),
dimension(:) :: x
583 integer(int32),
intent(inout),
dimension(:) :: ind
584 integer(int32),
intent(out) :: marker
587 integer(int32) :: i, j, itemp
588 complex(real64) :: temp
589 real(real64) :: pivot
592 pivot = real(x(1), real64)
600 if (real(x(j), real64) <= pivot)
exit
605 if (real(x(i), real64) >= pivot)
exit
617 else if (i == j)
then
630 if (real(x(j), real64) >= pivot)
exit
635 if (real(x(i), real64) <= pivot)
exit
647 else if (i == j)
then
Provides a set of common linear algebra routines.