Skip site navigation (1)Skip section navigation (2)

FreeBSD Manual Pages

  
 
  

home | help
PSSTEIN(l)			       )			    PSSTEIN(l)

NAME
       PSSTEIN - compute the eigenvectors of a symmetric tridiagonal matrix in
       parallel, using inverse iteration

SYNOPSIS
       SUBROUTINE PSSTEIN( N, D, E, M, W, IBLOCK, ISPLIT, ORFAC,  Z,  IZ,  JZ,
			   DESCZ,  WORK, LWORK,	IWORK, LIWORK, IFAIL, ICLUSTR,
			   GAP,	INFO )

	   INTEGER	   INFO, IZ, JZ, LIWORK, LWORK,	M, N

	   REAL		   ORFAC

	   INTEGER	   DESCZ( * ), IBLOCK( * ), ICLUSTR( * ), IFAIL( *  ),
			   ISPLIT( * ),	IWORK( * )

	   REAL		   D( *	), E( *	), GAP(	* ), W(	* ), WORK( * ),	Z( * )

PURPOSE
       PSSTEIN	computes the eigenvectors of a symmetric tridiagonal matrix in
       parallel, using inverse iteration. The eigenvectors found correspond to
       user specified eigenvalues. PSSTEIN does	not orthogonalize vectors that
       are on different	processes. The extent  of  orthogonalization  is  con-
       trolled	by the input parameter LWORK.  Eigenvectors that are to	be or-
       thogonalized are	computed by the	same process. PSSTEIN decides  on  the
       allocation of work among	the processes and then calls SSTEIN2 (modified
       LAPACK routine) on each individual process. If  insufficient  workspace
       is allocated, the expected orthogonalization may	not be done.

       Note : If the eigenvectors obtained are not orthogonal, increase
	      LWORK and	run the	code again.

       Notes
       =====

       Each  global data object	is described by	an associated description vec-
       tor.  This vector stores	the information	required to establish the map-
       ping between an object element and its corresponding process and	memory
       location.

       Let A be	a generic term for any 2D block	 cyclicly  distributed	array.
       Such a global array has an associated description vector	DESCA.	In the
       following comments, the character _ should be read as  "of  the	global
       array".

       NOTATION	       STORED IN      EXPLANATION
       ---------------	--------------	--------------------------------------
       DTYPE_A(global) DESCA( DTYPE_ )The descriptor type.  In this case,
				      DTYPE_A =	1.
       CTXT_A (global) DESCA( CTXT_ ) The BLACS	context	handle,	indicating
				      the BLACS	process	grid A is distribu-
				      ted over.	The context itself is glo-
				      bal, but the handle (the integer
				      value) may vary.
       M_A    (global) DESCA( M_ )    The number of rows in the	global
				      array A.
       N_A    (global) DESCA( N_ )    The number of columns in the global
				      array A.
       MB_A   (global) DESCA( MB_ )   The blocking factor used to distribute
				      the rows of the array.
       NB_A   (global) DESCA( NB_ )   The blocking factor used to distribute
				      the columns of the array.
       RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
				      row  of  the  array  A  is  distributed.
       CSRC_A (global) DESCA( CSRC_ ) The process column over which the
				      first column of the array	A is
				      distributed.
       LLD_A  (local)  DESCA( LLD_ )  The leading dimension of the local
				      array.  LLD_A >= MAX(1,LOCr(M_A)).

       Let K be	the number of rows or columns of a distributed matrix, and as-
       sume that its process grid has dimension	r x c.
       LOCr( K ) denotes the number of elements	of K that a process would  re-
       ceive if	K were distributed over	the r processes	of its process column.
       Similarly, LOCc(	K ) denotes the	number of elements of K	that a process
       would receive if	K were distributed over	the c processes	of its process
       row.
       The  values  of	LOCr()	and LOCc() may be determined via a call	to the
       ScaLAPACK tool function,	NUMROC:
	       LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
	       LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).  An	 upper
       bound for these quantities may be computed by:
	       LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
	       LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A

ARGUMENTS
       P = NPROW * NPCOL is the	total number of	processes

       N       (global input) INTEGER
	       The order of the	tridiagonal matrix T.  N >= 0.

       D       (global input) REAL array, dimension (N)
	       The n diagonal elements of the tridiagonal matrix T.

       E       (global input) REAL array, dimension (N-1)
	       The (n-1) off-diagonal elements of the tridiagonal matrix T.

       M       (global input) INTEGER
	       The total number	of eigenvectors	to be found. 0 <= M <= N.

       W       (global input/global output) REAL array,	dim (M)
	       On input, the first M elements of W contain all the eigenvalues
	       for which eigenvectors are  to  be  computed.  The  eigenvalues
	       should  be grouped by split-off block and ordered from smallest
	       to largest within the block (The	output array  W	 from  PSSTEBZ
	       with  ORDER='b'	is expected here). This	array should be	repli-
	       cated on	all processes.	On output, the first M	elements  con-
	       tain the	input eigenvalues in ascending order.

	       Note  : To obtain orthogonal vectors, it	is best	if eigenvalues
	       are computed to highest accuracy	( this can be done by  setting
	       ABSTOL  to  the underflow threshold = SLAMCH('U') --- ABSTOL is
	       an input	parameter to PSSTEBZ )

       IBLOCK  (global input) INTEGER array, dimension (N)
	       The submatrix indices associated	with the corresponding	eigen-
	       values  in W -- 1 for eigenvalues belonging to the first	subma-
	       trix from the top, 2 for	those belonging	to the	second	subma-
	       trix,  etc.  (The  output array IBLOCK from PSSTEBZ is expected
	       here).

       ISPLIT  (global input) INTEGER array, dimension (N)
	       The splitting points, at	which T	breaks	up  into  submatrices.
	       The  first  submatrix  consists of rows/columns 1 to ISPLIT(1),
	       the second of rows/columns ISPLIT(1)+1 through ISPLIT(2), etc.,
	       and  the	 NSPLIT-th consists of rows/columns ISPLIT(NSPLIT-1)+1
	       through ISPLIT(NSPLIT)=N	(The output array ISPLIT from  PSSTEBZ
	       is expected here.)

       ORFAC   (global input) REAL
	       ORFAC  specifies	 which	eigenvectors should be orthogonalized.
	       Eigenvectors that correspond to eigenvalues  which  are	within
	       ORFAC*||T||  of	each other are to be orthogonalized.  However,
	       if the workspace	is insufficient	(see  LWORK),  this  tolerance
	       may  be	decreased  until all eigenvectors to be	orthogonalized
	       can be stored in	one process.   No  orthogonalization  will  be
	       done if ORFAC equals zero.  A default value of 10^-3 is used if
	       ORFAC is	negative.  ORFAC should	be identical on	all processes.

       Z       (local output) REAL array,
	       dimension (DESCZ(DLEN_),	N/npcol	+ NB) Z	contains the  computed
	       eigenvectors  associated	 with  the  specified eigenvalues. Any
	       vector which fails to converge is set to	 its  current  iterate
	       after  MAXITS iterations	( See SSTEIN2 ).  On output, Z is dis-
	       tributed	across the P processes in block	cyclic format.

       IZ      (global input) INTEGER
	       Z's global row index, which points to the beginning of the sub-
	       matrix which is to be operated on.

       JZ      (global input) INTEGER
	       Z's  global  column index, which	points to the beginning	of the
	       submatrix which is to be	operated on.

       DESCZ   (global and local input)	INTEGER	array of dimension DLEN_.
	       The array descriptor for	the distributed	matrix Z.

       WORK    (local workspace/global output) REAL array,
	       dimension ( LWORK ) On output, WORK(1) gives a lower  bound  on
	       the  workspace  (  LWORK	) that guarantees the user desired or-
	       thogonalization (see ORFAC).  Note that this  may  overestimate
	       the minimum workspace needed.

       LWORK   (local input) integer
	       LWORK   controls	 the  extent of	orthogonalization which	can be
	       done. The number	of eigenvectors	for which storage is allocated
	       on  each	 process  is  NVEC = floor(( LWORK- max(5*N,NP00*MQ00)
	       )/N).  Eigenvectors corresponding  to  eigenvalue  clusters  of
	       size NVEC - ceil(M/P) + 1 are guaranteed	to be orthogonal ( the
	       orthogonality is	similar	to that	obtained from SSTEIN2).	  Note
	       :   LWORK   must	 be  no	 smaller  than:	 max(5*N,NP00*MQ00)  +
	       ceil(M/P)*N, and	should have the	same input value on  all  pro-
	       cesses.	 It  is	 the minimum value of LWORK input on different
	       processes that is significant.

	       If LWORK	= -1, then LWORK is global input and a workspace query
	       is assumed; the routine only calculates the minimum and optimal
	       size for	all work arrays. Each of these values is  returned  in
	       the  first  entry of the	corresponding work array, and no error
	       message is issued by PXERBLA.

       IWORK   (local workspace/global output) INTEGER array,
	       dimension ( 3*N+P+1 ) On	return,	IWORK(1) contains  the	amount
	       of integer workspace required.  On return, the IWORK(2) through
	       IWORK(P+2) indicate the eigenvectors computed by	each  process.
	       Process	I  computes  eigenvectors  indexed  IWORK(I+2)+1 thru'
	       IWORK(I+3).

       LIWORK  (local input) INTEGER
	       Size of array IWORK.  Must be >=	3*N + P	+ 1

	       If LIWORK = -1, then LIWORK is global  input  and  a  workspace
	       query  is  assumed; the routine only calculates the minimum and
	       optimal size for	all work arrays. Each of these values  is  re-
	       turned  in the first entry of the corresponding work array, and
	       no error	message	is issued by PXERBLA.

       IFAIL   (global output) integer array, dimension	(M)
	       On normal exit, all elements of IFAIL are zero.	If one or more
	       eigenvectors  fail  to  converge	after MAXITS iterations	(as in
	       SSTEIN),	then INFO > 0 is returned.  If	mod(INFO,M+1)>0,  then
	       for  I=1	to mod(INFO,M+1), the eigenvector corresponding	to the
	       eigenvalue W(IFAIL(I)) failed to	converge ( W refers to the ar-
	       ray of eigenvalues on output ).

	       ICLUSTR	(global	 output)  integer  array, dimension (2*P) This
	       output array contains indices of	eigenvectors corresponding  to
	       a  cluster  of eigenvalues that could not be orthogonalized due
	       to insufficient workspace (see LWORK, ORFAC and	INFO).	Eigen-
	       vectors	 corresponding	to  clusters  of  eigenvalues  indexed
	       ICLUSTR(2*I-1) to ICLUSTR(2*I), I = 1 to	INFO/(M+1), could  not
	       be orthogonalized due to	lack of	workspace. Hence the eigenvec-
	       tors corresponding to these clusters  may  not  be  orthogonal.
	       ICLUSTR	is  a  zero  terminated	 array --- ( ICLUSTR(2*K).NE.0
	       .AND. ICLUSTR(2*K+1).EQ.0 ) if and only if K is the  number  of
	       clusters.

       GAP     (global output) REAL array, dimension (P)
	       This  output  array  contains the gap between eigenvalues whose
	       eigenvectors could not be  orthogonalized.  The	INFO/M	output
	       values  in this array correspond	to the INFO/(M+1) clusters in-
	       dicated by the array ICLUSTR. As	a result, the dot product  be-
	       tween  eigenvectors corresponding to the	I^th cluster may be as
	       high as ( O(n)*macheps )	/ GAP(I).

       INFO    (global output) INTEGER
	       = 0:  successful	exit
	       < 0:  If	the i-th argument is an	array and the j-entry  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.  < 0	:   if
	       INFO = -I, the I-th argument had	an illegal value
	       > 0 :  if mod(INFO,M+1) = I, then I eigenvectors	failed to con-
	       verge in	MAXITS iterations. Their indices are stored in the ar-
	       ray  IFAIL.  if INFO/(M+1) = I, then eigenvectors corresponding
	       to I clusters of	eigenvalues could not be orthogonalized	due to
	       insufficient  workspace.	The indices of the clusters are	stored
	       in the array ICLUSTR.

ScaLAPACK version 1.7		13 August 2001			    PSSTEIN(l)

NAME | SYNOPSIS | PURPOSE | ARGUMENTS

Want to link to this manual page? Use this URL:
<https://www.freebsd.org/cgi/man.cgi?query=psstein&manpath=FreeBSD+12.0-RELEASE+and+Ports>

home | help