Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

cposvxx (3p)

Name

cposvxx - compute the solution to a complex system of linear equations A*X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices

Synopsis

SUBROUTINE  CPOSVXX(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, EQUED, S, B,
LDB, X, LDX, RCOND, RPVGRW, BERR, N_ERR_BNDS,  ERR_BNDS_NORM,
ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, RWORK, INFO)


CHARACTER*1 EQUED, FACT, UPLO

INTEGER INFO, LDA, LDAF, LDB, LDX, N, NRHS, NPARAMS, N_ERR_BNDS

REAL RCOND, RPVGRW

COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), WORK(*), X(LDX,*)

REAL    S(*),   PARAMS(*),   BERR(*),RWORK(*),   ERR_BNDS_NORM(NRHS,*),
ERR_BNDS_COMP(NRHS,*)


SUBROUTINE CPOSVXX_64(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, EQUED,  S,
B,   LDB,   X,   LDX,   RCOND,   RPVGRW,   BERR,  N_ERR_BNDS,
ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS,  WORK,  RWORK,
INFO)


CHARACTER*1 EQUED, FACT, UPLO

INTEGER*8 INFO, LDA, LDAF, LDB, LDX, N, NRHS, NPARAMS, N_ERR_BNDS

REAL RCOND, RPVGRW

COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), WORK(*), X(LDX,*)

REAL    S(*),   PARAMS(*),   BERR(*),RWORK(*),   ERR_BNDS_NORM(NRHS,*),
ERR_BNDS_COMP(NRHS,*)


F95 INTERFACE
SUBROUTINE POSVXX(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, EQUED,  S,  B,
LDB,  X, LDX, RCOND, RPVGRW, BERR, N_ERR_BNDS, ERR_BNDS_NORM,
ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, RWORK, INFO)


REAL, DIMENSION(:,:) :: ERR_BNDS_NORM, ERR_BNDS_COMP

INTEGER :: N, NRHS, LDA, LDAF, LDB, LDX, N_ERR_BNDS, NPARAMS, INFO

CHARACTER(LEN=1) :: FACT, UPLO, EQUED

REAL, DIMENSION(:) :: S, BERR, PARAMS, RWORK

COMPLEX, DIMENSION(:,:) :: A, AF, B, X

COMPLEX, DIMENSION(:) :: WORK

REAL :: RCOND, RPVGRW


SUBROUTINE POSVXX_64(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF,  EQUED,  S,
B,   LDB,   X,   LDX,   RCOND,   RPVGRW,   BERR,  N_ERR_BNDS,
ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS,  WORK,  RWORK,
INFO)


REAL, DIMENSION(:,:) :: ERR_BNDS_NORM, ERR_BNDS_COMP

INTEGER(8) :: N, NRHS, LDA, LDAF, LDB, LDX, N_ERR_BNDS, NPARAMS, INFO

CHARACTER(LEN=1) :: FACT, UPLO, EQUED

REAL, DIMENSION(:) :: S, BERR, PARAMS, RWORK

COMPLEX, DIMENSION(:,:) :: A, AF, B, X

COMPLEX, DIMENSION(:) :: WORK

REAL :: RCOND, RPVGRW


C INTERFACE
#include <sunperf.h>

void  cposvxx  (char fact, char uplo, int n, int nrhs, floatcomplex *a,
int lda, floatcomplex *af, int ldaf, char *equed,  float  *s,
floatcomplex  *b,  int  ldb,  floatcomplex *x, int ldx, float
*rcond, float *rpvgrw, float  *berr,  int  n_err_bnds,  float
*err_bnds_norm,  float  *err_bnds_comp,  int  nparams,  float
*params, int *info);


void cposvxx_64 (char fact, char uplo, long n, long nrhs,  floatcomplex
*a, long lda, floatcomplex *af, long ldaf, char *equed, float
*s, floatcomplex *b, long ldb,  floatcomplex  *x,  long  ldx,
float  *rcond,  float  *rpvgrw, float *berr, long n_err_bnds,
float *err_bnds_norm,  float  *err_bnds_comp,  long  nparams,
float *params, long *info);

Description

Oracle Solaris Studio Performance Library                          cposvxx(3P)



NAME
       cposvxx  - compute the solution to a complex system of linear equations
       A*X = B, where A is an N-by-N symmetric positive definite matrix and  X
       and B are N-by-NRHS matrices


SYNOPSIS
       SUBROUTINE  CPOSVXX(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, EQUED, S, B,
                 LDB, X, LDX, RCOND, RPVGRW, BERR, N_ERR_BNDS,  ERR_BNDS_NORM,
                 ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, RWORK, INFO)


       CHARACTER*1 EQUED, FACT, UPLO

       INTEGER INFO, LDA, LDAF, LDB, LDX, N, NRHS, NPARAMS, N_ERR_BNDS

       REAL RCOND, RPVGRW

       COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), WORK(*), X(LDX,*)

       REAL    S(*),   PARAMS(*),   BERR(*),RWORK(*),   ERR_BNDS_NORM(NRHS,*),
                 ERR_BNDS_COMP(NRHS,*)


       SUBROUTINE CPOSVXX_64(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, EQUED,  S,
                 B,   LDB,   X,   LDX,   RCOND,   RPVGRW,   BERR,  N_ERR_BNDS,
                 ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS,  WORK,  RWORK,
                 INFO)


       CHARACTER*1 EQUED, FACT, UPLO

       INTEGER*8 INFO, LDA, LDAF, LDB, LDX, N, NRHS, NPARAMS, N_ERR_BNDS

       REAL RCOND, RPVGRW

       COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), WORK(*), X(LDX,*)

       REAL    S(*),   PARAMS(*),   BERR(*),RWORK(*),   ERR_BNDS_NORM(NRHS,*),
                 ERR_BNDS_COMP(NRHS,*)


   F95 INTERFACE
       SUBROUTINE POSVXX(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, EQUED,  S,  B,
                 LDB,  X, LDX, RCOND, RPVGRW, BERR, N_ERR_BNDS, ERR_BNDS_NORM,
                 ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, RWORK, INFO)


       REAL, DIMENSION(:,:) :: ERR_BNDS_NORM, ERR_BNDS_COMP

       INTEGER :: N, NRHS, LDA, LDAF, LDB, LDX, N_ERR_BNDS, NPARAMS, INFO

       CHARACTER(LEN=1) :: FACT, UPLO, EQUED

       REAL, DIMENSION(:) :: S, BERR, PARAMS, RWORK

       COMPLEX, DIMENSION(:,:) :: A, AF, B, X

       COMPLEX, DIMENSION(:) :: WORK

       REAL :: RCOND, RPVGRW


       SUBROUTINE POSVXX_64(FACT, UPLO, N, NRHS, A, LDA, AF, LDAF,  EQUED,  S,
                 B,   LDB,   X,   LDX,   RCOND,   RPVGRW,   BERR,  N_ERR_BNDS,
                 ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS,  WORK,  RWORK,
                 INFO)


       REAL, DIMENSION(:,:) :: ERR_BNDS_NORM, ERR_BNDS_COMP

       INTEGER(8) :: N, NRHS, LDA, LDAF, LDB, LDX, N_ERR_BNDS, NPARAMS, INFO

       CHARACTER(LEN=1) :: FACT, UPLO, EQUED

       REAL, DIMENSION(:) :: S, BERR, PARAMS, RWORK

       COMPLEX, DIMENSION(:,:) :: A, AF, B, X

       COMPLEX, DIMENSION(:) :: WORK

       REAL :: RCOND, RPVGRW


   C INTERFACE
       #include <sunperf.h>

       void  cposvxx  (char fact, char uplo, int n, int nrhs, floatcomplex *a,
                 int lda, floatcomplex *af, int ldaf, char *equed,  float  *s,
                 floatcomplex  *b,  int  ldb,  floatcomplex *x, int ldx, float
                 *rcond, float *rpvgrw, float  *berr,  int  n_err_bnds,  float
                 *err_bnds_norm,  float  *err_bnds_comp,  int  nparams,  float
                 *params, int *info);


       void cposvxx_64 (char fact, char uplo, long n, long nrhs,  floatcomplex
                 *a, long lda, floatcomplex *af, long ldaf, char *equed, float
                 *s, floatcomplex *b, long ldb,  floatcomplex  *x,  long  ldx,
                 float  *rcond,  float  *rpvgrw, float *berr, long n_err_bnds,
                 float *err_bnds_norm,  float  *err_bnds_comp,  long  nparams,
                 float *params, long *info);


PURPOSE
       cposvxx  uses  the  Cholesky  factorization A = U**T*U or A = L*L**T to
       compute the solution to a complex system of linear equations A * X = B,
       where A is an N-by-N symmetric positive definite matrix and X and B are
       N-by-NRHS matrices.

       If requested, both normwise and maximum componentwise error bounds  are
       returned.  CPOSVXX  will return a solution with a tiny guaranteed error
       (O(eps) where eps is the working machine precision) unless  the  matrix
       is  very ill-conditioned, in which case a warning is returned. Relevant
       condition numbers also are calculated and returned.

       CPOSVXX accepts user-provided factorizations and equilibration factors;
       see  the  definitions  of  the  FACT  and  EQUED options.  Solving with
       refinement and using a factorization from a previous CPOSVXX call  will
       also  produce  a solution with either O(eps) errors or warnings, but we
       cannot make that claim for  general  user-provided  factorizations  and
       equilibration  factors  if  they  differ from what CPOSVXX would itself
       produce.


ARGUMENTS
       FACT (input)
                 FACT is CHARACTER*1
                 Specifies whether or not the factored form of the matrix A is
                 supplied on entry, and if not, whether the matrix A should be
                 equilibrated before it is factored.
                 = 'F':  On entry, AF contains the factored form of A.
                 If EQUED is not 'N', the matrix A has been equilibrated  with
                 scaling factors given by S.  A and AF are not modified.
                 = 'N':  The matrix A will be copied to AF and factored.
                 =  'E':  The matrix A will be equilibrated if necessary, then
                 copied to AF and factored.


       UPLO (input)
                 UPLO is CHARACTER*1
                 = 'U':  Upper triangle of A is stored;
                 = 'L':  Lower triangle of A is stored.


       N (input)
                 N is INTEGER
                 The number of linear equations, i.e., the order of the matrix
                 A.  N >= 0.


       NRHS (input)
                 NRHS is INTEGER
                 The  number  of right hand sides, i.e., the number of columns
                 of the matrices B and X.  NRHS >= 0.


       A (input/output)
                 A is COMPLEX array, dimension (LDA,N)
                 On entry, the symmetric matrix A, except if FACT  =  'F'  and
                 EQUED = diag(S)*A*diag(S).  If UPLO = 'U', the leading N-by-N
                 upper triangular part of A contains the upper triangular part
                 of  the matrix A, and the strictly lower triangular part of A
                 is not referenced.  If UPLO = 'L', the leading  N-by-N  lower
                 triangular  part  of  A contains the lower triangular part of
                 the matrix A, and the strictly upper triangular part of A  is
                 not  referenced.   A is not modified if FACT = 'F' or 'N', or
                 if FACT = 'E' and EQUED =
                 On exit, if FACT = 'E' and EQUED = 'Y', A is  overwritten  by
                 diag(S)*A*diag(S).


       LDA (input)
                 LDA is INTEGER
                 The leading dimension of the array A.
                 LDA >= max(1,N).


       AF (input/output)
                 AF is COMPLEX array, dimension (LDAF,N)
                 If FACT = 'F', then AF is an input argument and on entry con-
                 tains the triangular factor U or L from the Cholesky  factor-
                 ization  A=U**T*U  or A=L*L**T, in the same storage format as
                 A. If EQUED .ne. 'N', then AF is the  factored  form  of  the
                 equilibrated matrix diag(S)*A*diag(S).
                 If  FACT  =  'N',  then  AF is an output argument and on exit
                 returns the triangular factor U or L from the  Cholesky  fac-
                 torization A=U**T*U or A=L*L**T of the original matrix A.
                 If  FACT  =  'E',  then  AF is an output argument and on exit
                 returns the triangular factor U or L from the  Cholesky  fac-
                 torization  A=U**T*U or A=L*L**T of the equilibrated matrix A
                 (see the description of A for the form  of  the  equilibrated
                 matrix).


       LDAF (input)
                 LDAF is INTEGER
                 The leading dimension of the array AF.
                 LDAF >= max(1,N).


       EQUED (input/output)
                 EQUED is CHARACTER*1
                 Specifies the form of equilibration that was done.
                 = 'N':  No equilibration (always true if FACT = 'N').
                 =  'Y':   Both row and column equilibration, i.e., A has been
                 replaced by diag(S) * A * diag(S).  EQUED is an  input  argu-
                 ment if FACT = 'F'; otherwise, it is an output argument.


       S (input/output)
                 S  is  REAL array, dimension (N) The row scale factors for A.
                 If EQUED = 'Y', A is multiplied on  the  left  and  right  by
                 diag(S). S is an input argument if FACT = = 'Y', each element
                 of S must be positive. If S is output, each element of S is a
                 power  of  the radix. If S is input, each element of S should
                 be a power of the radix to ensure  a  reliable  solution  and
                 error  estimates.  Scaling  by  powers  of the radix does not
                 cause rounding errors unless the result underflows  or  over-
                 flows.   Rounding errors during scaling lead to refining with
                 a matrix that is not equivalent to the input matrix,  produc-
                 ing error estimates that may not be reliable.


       B (input/output)
                 B is COMPLEX array, dimension (LDB,NRHS)
                 On entry, the N-by-NRHS right hand side matrix B.
                 On exit,
                 if EQUED = 'N', B is not modified;
                 if EQUED = 'Y', B is overwritten by diag(S)*B;


       LDB (input)
                 LDB is INTEGER
                 The leading dimension of the array B.
                 LDB >= max(1,N).


       X (output)
                 X is COMPLEX array, dimension (LDX,NRHS)
                 If  INFO = 0, the N-by-NRHS solution matrix X to the original
                 system of equations.  Note that A and B are modified on  exit
                 if  EQUED .ne. 'N', and the solution to the equilibrated sys-
                 tem is inv(diag(S))*X.


       LDX (input)
                 LDX is INTEGER
                 The leading dimension of the array X.
                 LDX >= max(1,N).


       RCOND (output)
                 RCOND is REAL
                 Reciprocal scaled condition number. This is  an  estimate  of
                 the  reciprocal  Skeel condition number of the matrix A after
                 equilibration (if done). If this is  less  than  the  machine
                 precision  (in particular, if it is zero), the matrix is sin-
                 gular to working precision. Note that the error may still  be
                 small  even  if  this  number  is  very  small and the matrix
                 appears ill- conditioned.


       RPVGRW (output)
                 RPVGRW is REAL
                 Reciprocal pivot growth. On exit, this contains the  recipro-
                 cal  pivot  growth  factor norm(A)/norm(U). The "max absolute
                 element" norm is used. If this is much less than 1, then  the
                 stability  of  the  LU  factorization  of  the (equilibrated)
                 matrix A could be poor.  This also means that the solution X,
                 estimated  condition numbers, and error bounds could be unre-
                 liable. If factorization fails with 0<INFO<=N, then this con-
                 tains the reciprocal pivot growth factor for the leading INFO
                 columns of A.


       BERR (output)
                 BERR is REAL array, dimension (NRHS)
                 Componentwise relative backward error. This is the component-
                 wise  relative  backward  error  of each solution vector X(j)
                 (i.e., the smallest relative change in any element of A or  B
                 that makes X(j) an exact solution).


       N_ERR_BNDS (input)
                 N_ERR_BNDS is INTEGER
                 Number of error bounds to return for each right hand side and
                 each type (normwise or componentwise). See ERR_BNDS_NORM  and
                 ERR_BNDS_COMP below.


       ERR_BNDS_NORM (output)
                 ERR_BNDS_NORM is REAL array, dimension (NRHS, N_ERR_BNDS)
                 For  each  right-hand  side,  this array contains information
                 about various error bounds and condition numbers  correspond-
                 ing  to the normwise relative error, which is defined as fol-
                 lows:
                 Normwise relative error in the ith solution vector:

                           max_j (abs(XTRUE(j,i) - X(j,i)))
                          ------------------------------
                                max_j abs(X(j,i))

                 The array is indexed by the  type  of  error  information  as
                 described  below.  There  currently are up to three pieces of
                 information returned.
                 The first index in ERR_BNDS_NORM(i,:) corresponds to the  ith
                 right-hand side.
                 The second index in ERR_BNDS_NORM(:,err) contains the follow-
                 ing three fields:
                 err = 1 "Trust/don't trust" boolean. Trust the answer if  the
                 reciprocal  condition  number  is  less  than  the  threshold
                 sqrt(n) * slamch('Epsilon').
                 err = 2  "Guaranteed"  error  bound:  The  estimated  forward
                 error,  almost  certainly  within  a factor of 10 of the true
                 error so long as the next entry is greater than the threshold
                 sqrt(n)  * slamch('Epsilon'). This error bound should only be
                 trusted if the previous boolean is true.
                 err = 3   Reciprocal  condition  number:  Estimated  normwise
                 reciprocal  condition  number.   Compared  with the threshold
                 sqrt(n) * slamch('Epsilon') to determine if the  error  esti-
                 mate is "guaranteed". These reciprocal condition numbers are
                 1  /  (norm(Z^{-1},inf) * norm(Z,inf)) for some appropriately
                 scaled matrix Z.
                 Let Z = S*A, where S scales each row by a power of the  radix
                 so all absolute row sums of Z are approximately 1.
                 See  Lapack  Working  Note  165 for further details and extra
                 cautions.


       ERR_BNDS_COMP (output)
                 ERR_BNDS_COMP is REAL array, dimension (NRHS, N_ERR_BNDS)
                 For each right-hand side,  this  array  contains  information
                 about  various error bounds and condition numbers correspond-
                 ing to the componentwise relative error, which is defined  as
                 follows:
                 Componentwise relative error in the ith solution vector:

                              abs(XTRUE(j,i) - X(j,i))
                       max_j ----------------------
                                   abs(X(j,i))

                 The  array  is indexed by the right-hand side i (on which the
                 componentwise relative error depends), and the type of  error
                 information  as  described  below.  There currently are up to
                 three pieces of  information  returned  for  each  right-hand
                 side. If componentwise accuracy is not requested (PARAMS(3) =
                 0.0), then ERR_BNDS_COMP is not accessed. If N_ERR_BNDS  .LT.
                 3,   then  at  most  the  first  (:,N_ERR_BNDS)  entries  are
                 returned.
                 The first index in ERR_BNDS_COMP(i,:) corresponds to the  ith
                 right-hand side.
                 The second index in ERR_BNDS_COMP(:,err) contains the follow-
                 ing three fields:
                 err = 1 "Trust/don't trust" boolean. Trust the answer if  the
                 reciprocal  condition  number  is  less  than  the  threshold
                 sqrt(n) * slamch('Epsilon').
                 err = 2  "Guaranteed"  error  bound:  The  estimated  forward
                 error,  almost  certainly  within  a factor of 10 of the true
                 error so long as the next entry is greater than the threshold
                 sqrt(n)  * slamch('Epsilon'). This error bound should only be
                 trusted if the previous boolean is true.
                 err = 3  Reciprocal condition number: Estimated componentwise
                 reciprocal  condition  number.   Compared  with the threshold
                 sqrt(n) * slamch('Epsilon') to determine if the  error  esti-
                 mate is "guaranteed". These reciprocal condition numbers are
                 1  /  (norm(Z^{-1},inf) * norm(Z,inf)) for some appropriately
                 scaled matrix Z.
                 Let Z = S*(A*diag(x)), where x is the solution for  the  cur-
                 rent  right-hand side and S scales each row of A*diag(x) by a
                 power of the radix so all absolute row sums of Z are approxi-
                 mately 1.
                 See  Lapack  Working  Note  165 for further details and extra
                 cautions.


       NPARAMS (input)
                 NPARAMS is INTEGER
                 Specifies the number of parameters set in PARAMS. If .LE.  0,
                 the  PARAMS  array is never referenced and default values are
                 used.


       PARAMS (input/output)
                 PARAMS is REAL array, dimension NPARAMS  Specifies  algorithm
                 parameters.  If an entry is .LT. 0.0, then that entry will be
                 filled with default value  used  for  that  parameter.   Only
                 positions  up  to NPARAMS are accessed; defaults are used for
                 higher-numbered parameters.
                 PARAMS(LA_LINRX_ITREF_I = 1) : Whether to  perform  iterative
                 refinement or not.
                 Default: 1.0
                 =  0.0  : No refinement is performed, and no error bounds are
                 computed.
                 = 1.0 : Use the double-precision refinement algorithm, possi-
                 bly with doubled-single computations if the compilation envi-
                 ronment does not support DOUBLE PRECISION.
                 (other values are reserved for future use)
                 PARAMS(LA_LINRX_ITHRESH_I = 2) : Maximum number  of  residual
                 computations allowed for refinement.
                 Default: 10
                 Aggressive:  Set  to 100 to permit convergence using approxi-
                 mate
                 factorizations or factorizations other than LU. If  the  fac-
                 torization  uses a technique other than Gaussian elimination,
                 the guarantees in  err_bnds_norm  and  err_bnds_comp  may  no
                 longer be trustworthy.
                 PARAMS(LA_LINRX_CWISE_I  =  3) : Flag determining if the code
                 will attempt to find a solution with small componentwise rel-
                 ative  error  in the double-precision algorithm.  Positive is
                 true, 0.0 is false.
                 Default: 1.0 (attempt componentwise convergence)


       WORK (output)
                 WORK is COMPLEX array, dimension (2*N)


       RWORK (output)
                 RWORK is REAL array, dimension (2*N)


       INFO (output)
                 INFO is INTEGER
                 = 0:  Successful exit. The solution to every right-hand  side
                 is guaranteed.
                 < 0:  If INFO = -i, the i-th argument had an illegal value
                 >  0 and <= N:  U(INFO,INFO) is exactly zero.  The factoriza-
                 tion has been completed, but the factor U is  exactly  singu-
                 lar,  so the solution and error bounds could not be computed.
                 RCOND = 0 is returned.
                 = N+J: The solution corresponding to the Jth right-hand  side
                 is  not  guaranteed.  The  solutions  corresponding  to other
                 right- hand sides K with K > J may not be guaranteed as well,
                 but  only  the  first  such right-hand side is reported. If a
                 small componentwise error is not requested (PARAMS(3) =  0.0)
                 then  the  Jth  right-hand  side is the first with a normwise
                 error bound that is not guaranteed (the smallest J such  that
                 ERR_BNDS_NORM(J,1)  =  0.0). By default (PARAMS(3) = 1.0) the
                 Jth right-hand side is the first with either  a  normwise  or
                 componentwise  error bound that is not guaranteed (the small-
                 est  J  such  that  either  ERR_BNDS_NORM(J,1)   =   0.0   or
                 ERR_BNDS_COMP(J,1)    =   0.0).   See   the   definition   of
                 ERR_BNDS_NORM(:,1) and ERR_BNDS_COMP(:,1).
                 To get information about all of the  right-hand  sides  check
                 ERR_BNDS_NORM or ERR_BNDS_COMP.



                                  7 Nov 2015                       cposvxx(3P)