Developer Reference

  • 2021.1
  • 12/04/2020
  • Public Content
Contents

p?pttrsv

Solves a single triangular linear system via frontsolve or backsolve where the triangular matrix is a factor of a tridiagonal matrix computed by
p?pttrf
.

Syntax

void
pspttrsv
(
char
*uplo
,
MKL_INT
*n
,
MKL_INT
*nrhs
,
float
*d
,
float
*e
,
MKL_INT
*ja
,
MKL_INT
*desca
,
float
*b
,
MKL_INT
*ib
,
MKL_INT
*descb
,
float
*af
,
MKL_INT
*laf
,
float
*work
,
MKL_INT
*lwork
,
MKL_INT
*info
);
void
pdpttrsv
(
char
*uplo
,
MKL_INT
*n
,
MKL_INT
*nrhs
,
double
*d
,
double
*e
,
MKL_INT
*ja
,
MKL_INT
*desca
,
double
*b
,
MKL_INT
*ib
,
MKL_INT
*descb
,
double
*af
,
MKL_INT
*laf
,
double
*work
,
MKL_INT
*lwork
,
MKL_INT
*info
);
void
pcpttrsv
(
char
*uplo
,
char
*trans
,
MKL_INT
*n
,
MKL_INT
*nrhs
,
float
*d
,
MKL_Complex8
*e
,
MKL_INT
*ja
,
MKL_INT
*desca
,
MKL_Complex8
*b
,
MKL_INT
*ib
,
MKL_INT
*descb
,
MKL_Complex8
*af
,
MKL_INT
*laf
,
MKL_Complex8
*work
,
MKL_INT
*lwork
,
MKL_INT
*info
);
void
pzpttrsv
(
char
*uplo
,
char
*trans
,
MKL_INT
*n
,
MKL_INT
*nrhs
,
double
*d
,
MKL_Complex16
*e
,
MKL_INT
*ja
,
MKL_INT
*desca
,
MKL_Complex16
*b
,
MKL_INT
*ib
,
MKL_INT
*descb
,
MKL_Complex16
*af
,
MKL_INT
*laf
,
MKL_Complex16
*work
,
MKL_INT
*lwork
,
MKL_INT
*info
);
Include Files
  • mkl_scalapack.h
Description
The
p?pttrsv
function
solves a tridiagonal triangular system of linear equations
A
(1:
n
,
ja
:
ja
+
n
-1)*
X
=
B
(
jb
:
jb+n
-1, 1:
nrhs
)
or
A
(1:
n
,
ja
:
ja
+
n
-1)
T
*
X
=
B
(
jb
:
jb+n
-1, 1:
nrhs
)
for real flavors,
A
(1:
n
,
ja
:
ja
+
n
-1)
H
*
X
=
B
(
jb
:
jb+n
-1, 1:
nrhs
)
for complex flavors,
where
A
(1:
n
,
ja
:
ja
+
n
-1)
is a tridiagonal triangular matrix factor produced by the Cholesky factorization code
p?pttrf
and is stored in
A
(1:
n
,
ja
:
ja
+
n
-1)
and
af
. The matrix stored in
A
(1:
n
,
ja
:
ja
+
n
-1)
is either upper or lower triangular according to
uplo
.
The function
p?pttrf
must be called first.
Input Parameters
uplo
(global) Must be
'U'
or
'L'
.
If
uplo
=
'U'
, upper triangle of
A
(1:
n
,
ja
:
ja
+
n
-1)
is stored;
If
uplo
=
'L'
, lower triangle of
A
(1:
n
,
ja
:
ja
+
n
-1)
is stored.
trans
(global) Must be
'N'
or
'C'
.
If
trans
=
'N'
, solve with
A
(1:
n
,
ja
:
ja
+
n
-1)
;
If
trans
=
'C'
(for complex flavors), solve with conjugate transpose (
A
(1:
n
,
ja
:
ja
+
n
-1)
)
H
.
n
(global)
The number of rows and columns to be operated on, that is, the order of the distributed submatrix
A
(1:
n
,
ja
:
ja
+
n
-1)
.
n
0
.
nrhs
(global)
The number of right hand sides; the number of columns of the distributed submatrix
B
(
jb
:
jb+n
-1, 1:
nrhs
);
nrhs
0
.
d
(local)
Pointer to the local part of the global vector storing the main diagonal of the matrix; must be of size
nb_a
.
e
(local)
Pointer to the local part of the global vector
du
storing the upper diagonal of the matrix; must be of size
nb_a
. Globally,
du
(
n
) is not referenced, and
du
must be aligned with
d
.
ja
(global) The index in the global matrix
A
that points to the start of the matrix to be operated on (which may be either all of
A
or a submatrix of
A
).
desca
(global and local) array of size
dlen_
. The array descriptor for the distributed matrix
A
.
If
1D type (
dtype_a
= 501 or 502)
, then
dlen
7
;
If
2D type (
dtype_a
= 1)
, then
dlen
9
.
Contains information on mapping of
A
to memory. See ScaLAPACK manual for full description and options.
b
(local)
Pointer into the local memory to an array of local lead size
lld_b
nb
.
On entry, this array contains the local pieces of the right hand sides
B
(
jb
:
jb+n
-1, 1:
nrhs
)
.
ib
(global) The row index in the global matrix
B
that points to the first row of the matrix to be operated on (which may be either all of
B
or a submatrix of
B
).
descb
(global and local) array of size
dlen_
. The array descriptor for the distributed matrix
B
.
If
1D type (
dtype_b
= 502)
, then
dlen
7
;
If
2D type (
dtype_b
= 1)
, then
dlen
9
.
Contains information on mapping of
B
to memory. See ScaLAPACK manual for full description and options.
laf
(local)
The size of user-input auxiliary fill-in space
af
. Must be
laf
(
nb
+2
*bw
)
*bw
.
If
laf
is not large enough, an error code will be returned and the minimum acceptable size will be returned in
af
[0]
.
work
(local)
The array
work
is a temporary workspace array of size
lwork
. This space may be overwritten in between
function calls
.
lwork
(local or global) The size of the user-input workspace
work
, must be at least
lwork
(10+2*min(100,
nrhs
))*
npcol
+4*
nrhs
. If
lwork
is too small, the minimal acceptable size will be returned in
work
[0]
and an error code is returned.
Output Parameters
d
,
e
(local).
On exit, these arrays contain information on the factors of the matrix.
af
(local)
The array
af
is of size
laf
. It contains auxiliary fill-in space. The fill-in space is created in a call to the factorization
function
p?pbtrf
and is stored in
af
. If a linear system is to be solved using
p?pttrs
after the factorization
function
,
af
must not be altered after the factorization.
b
On exit, this array contains the local piece of the solutions distributed matrix
X
.
work
[0]
On exit,
work
[0]
contains the minimum value of
lwork
.
info
(local)
= 0
: successful exit
< 0
: if the
i
-th argument is an array and the
j
-th entry
, indexed
j
-1,
had an illegal value,
then
info
= - (
i
*100 +
j
),
if the
i
-th argument is a scalar and had an illegal value,
then
info
= -
i
.

Product and Performance Information

1

Performance varies by use, configuration and other factors. Learn more at www.Intel.com/PerformanceIndex.