LAPACK  3.8.0 LAPACK: Linear Algebra PACKage

## ◆ cpftrf()

 subroutine cpftrf ( character TRANSR, character UPLO, integer N, complex, dimension( 0: * ) A, integer INFO )

CPFTRF

Purpose:
``` CPFTRF computes the Cholesky factorization of a complex Hermitian
positive definite matrix A.

The factorization has the form
A = U**H * U,  if UPLO = 'U', or
A = L  * L**H,  if UPLO = 'L',
where U is an upper triangular matrix and L is lower triangular.

This is the block version of the algorithm, calling Level 3 BLAS.```
Parameters
 [in] TRANSR ``` TRANSR is CHARACTER*1 = 'N': The Normal TRANSR of RFP A is stored; = 'C': The Conjugate-transpose TRANSR of RFP A is stored.``` [in] UPLO ``` UPLO is CHARACTER*1 = 'U': Upper triangle of RFP A is stored; = 'L': Lower triangle of RFP A is stored.``` [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension ( N*(N+1)/2 ); On entry, the Hermitian matrix A in RFP format. RFP format is described by TRANSR, UPLO, and N as follows: If TRANSR = 'N' then RFP A is (0:N,0:k-1) when N is even; k=N/2. RFP A is (0:N-1,0:k) when N is odd; k=N/2. IF TRANSR = 'C' then RFP is the Conjugate-transpose of RFP A as defined when TRANSR = 'N'. The contents of RFP A are defined by UPLO as follows: If UPLO = 'U' the RFP A contains the nt elements of upper packed A. If UPLO = 'L' the RFP A contains the elements of lower packed A. The LDA of RFP A is (N+1)/2 when TRANSR = 'C'. When TRANSR is 'N' the LDA is N+1 when N is even and N is odd. See the Note below for more details. On exit, if INFO = 0, the factor U or L from the Cholesky factorization RFP A = U**H*U or RFP A = L*L**H.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, the leading minor of order i is not positive definite, and the factorization could not be completed. Further Notes on RFP Format: ============================ We first consider Standard Packed Format when N is even. We give an example where N = 6. AP is Upper AP is Lower 00 01 02 03 04 05 00 11 12 13 14 15 10 11 22 23 24 25 20 21 22 33 34 35 30 31 32 33 44 45 40 41 42 43 44 55 50 51 52 53 54 55 Let TRANSR = 'N'. RFP holds AP as follows: For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last three columns of AP upper. The lower triangle A(4:6,0:2) consists of conjugate-transpose of the first three columns of AP upper. For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first three columns of AP lower. The upper triangle A(0:2,0:2) consists of conjugate-transpose of the last three columns of AP lower. To denote conjugate we place -- above the element. This covers the case N even and TRANSR = 'N'. RFP A RFP A -- -- -- 03 04 05 33 43 53 -- -- 13 14 15 00 44 54 -- 23 24 25 10 11 55 33 34 35 20 21 22 -- 00 44 45 30 31 32 -- -- 01 11 55 40 41 42 -- -- -- 02 12 22 50 51 52 Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate- transpose of RFP A above. One therefore gets: RFP A RFP A -- -- -- -- -- -- -- -- -- -- 03 13 23 33 00 01 02 33 00 10 20 30 40 50 -- -- -- -- -- -- -- -- -- -- 04 14 24 34 44 11 12 43 44 11 21 31 41 51 -- -- -- -- -- -- -- -- -- -- 05 15 25 35 45 55 22 53 54 55 22 32 42 52 We next consider Standard Packed Format when N is odd. We give an example where N = 5. AP is Upper AP is Lower 00 01 02 03 04 00 11 12 13 14 10 11 22 23 24 20 21 22 33 34 30 31 32 33 44 40 41 42 43 44 Let TRANSR = 'N'. RFP holds AP as follows: For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last three columns of AP upper. The lower triangle A(3:4,0:1) consists of conjugate-transpose of the first two columns of AP upper. For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first three columns of AP lower. The upper triangle A(0:1,1:2) consists of conjugate-transpose of the last two columns of AP lower. To denote conjugate we place -- above the element. This covers the case N odd and TRANSR = 'N'. RFP A RFP A -- -- 02 03 04 00 33 43 -- 12 13 14 10 11 44 22 23 24 20 21 22 -- 00 33 34 30 31 32 -- -- 01 11 44 40 41 42 Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate- transpose of RFP A above. One therefore gets: RFP A RFP A -- -- -- -- -- -- -- -- -- 02 12 22 00 01 00 10 20 30 40 50 -- -- -- -- -- -- -- -- -- 03 13 23 33 11 33 11 21 31 41 51 -- -- -- -- -- -- -- -- -- 04 14 24 34 44 43 44 22 32 42 52```
Date
December 2016

Definition at line 213 of file cpftrf.f.

213 *
214 * -- LAPACK computational routine (version 3.7.0) --
215 * -- LAPACK is a software package provided by Univ. of Tennessee, --
216 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
217 * December 2016
218 *
219 * .. Scalar Arguments ..
220  CHARACTER transr, uplo
221  INTEGER n, info
222 * ..
223 * .. Array Arguments ..
224  COMPLEX a( 0: * )
225 *
226 * =====================================================================
227 *
228 * .. Parameters ..
229  REAL one
230  COMPLEX cone
231  parameter( one = 1.0e+0, cone = ( 1.0e+0, 0.0e+0 ) )
232 * ..
233 * .. Local Scalars ..
234  LOGICAL lower, nisodd, normaltransr
235  INTEGER n1, n2, k
236 * ..
237 * .. External Functions ..
238  LOGICAL lsame
239  EXTERNAL lsame
240 * ..
241 * .. External Subroutines ..
242  EXTERNAL xerbla, cherk, cpotrf, ctrsm
243 * ..
244 * .. Intrinsic Functions ..
245  INTRINSIC mod
246 * ..
247 * .. Executable Statements ..
248 *
249 * Test the input parameters.
250 *
251  info = 0
252  normaltransr = lsame( transr, 'N' )
253  lower = lsame( uplo, 'L' )
254  IF( .NOT.normaltransr .AND. .NOT.lsame( transr, 'C' ) ) THEN
255  info = -1
256  ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo, 'U' ) ) THEN
257  info = -2
258  ELSE IF( n.LT.0 ) THEN
259  info = -3
260  END IF
261  IF( info.NE.0 ) THEN
262  CALL xerbla( 'CPFTRF', -info )
263  RETURN
264  END IF
265 *
266 * Quick return if possible
267 *
268  IF( n.EQ.0 )
269  \$ RETURN
270 *
271 * If N is odd, set NISODD = .TRUE.
272 * If N is even, set K = N/2 and NISODD = .FALSE.
273 *
274  IF( mod( n, 2 ).EQ.0 ) THEN
275  k = n / 2
276  nisodd = .false.
277  ELSE
278  nisodd = .true.
279  END IF
280 *
281 * Set N1 and N2 depending on LOWER
282 *
283  IF( lower ) THEN
284  n2 = n / 2
285  n1 = n - n2
286  ELSE
287  n1 = n / 2
288  n2 = n - n1
289  END IF
290 *
291 * start execution: there are eight cases
292 *
293  IF( nisodd ) THEN
294 *
295 * N is odd
296 *
297  IF( normaltransr ) THEN
298 *
299 * N is odd and TRANSR = 'N'
300 *
301  IF( lower ) THEN
302 *
303 * SRPA for LOWER, NORMAL and N is odd ( a(0:n-1,0:n1-1) )
304 * T1 -> a(0,0), T2 -> a(0,1), S -> a(n1,0)
305 * T1 -> a(0), T2 -> a(n), S -> a(n1)
306 *
307  CALL cpotrf( 'L', n1, a( 0 ), n, info )
308  IF( info.GT.0 )
309  \$ RETURN
310  CALL ctrsm( 'R', 'L', 'C', 'N', n2, n1, cone, a( 0 ), n,
311  \$ a( n1 ), n )
312  CALL cherk( 'U', 'N', n2, n1, -one, a( n1 ), n, one,
313  \$ a( n ), n )
314  CALL cpotrf( 'U', n2, a( n ), n, info )
315  IF( info.GT.0 )
316  \$ info = info + n1
317 *
318  ELSE
319 *
320 * SRPA for UPPER, NORMAL and N is odd ( a(0:n-1,0:n2-1)
321 * T1 -> a(n1+1,0), T2 -> a(n1,0), S -> a(0,0)
322 * T1 -> a(n2), T2 -> a(n1), S -> a(0)
323 *
324  CALL cpotrf( 'L', n1, a( n2 ), n, info )
325  IF( info.GT.0 )
326  \$ RETURN
327  CALL ctrsm( 'L', 'L', 'N', 'N', n1, n2, cone, a( n2 ), n,
328  \$ a( 0 ), n )
329  CALL cherk( 'U', 'C', n2, n1, -one, a( 0 ), n, one,
330  \$ a( n1 ), n )
331  CALL cpotrf( 'U', n2, a( n1 ), n, info )
332  IF( info.GT.0 )
333  \$ info = info + n1
334 *
335  END IF
336 *
337  ELSE
338 *
339 * N is odd and TRANSR = 'C'
340 *
341  IF( lower ) THEN
342 *
343 * SRPA for LOWER, TRANSPOSE and N is odd
344 * T1 -> A(0,0) , T2 -> A(1,0) , S -> A(0,n1)
345 * T1 -> a(0+0) , T2 -> a(1+0) , S -> a(0+n1*n1); lda=n1
346 *
347  CALL cpotrf( 'U', n1, a( 0 ), n1, info )
348  IF( info.GT.0 )
349  \$ RETURN
350  CALL ctrsm( 'L', 'U', 'C', 'N', n1, n2, cone, a( 0 ), n1,
351  \$ a( n1*n1 ), n1 )
352  CALL cherk( 'L', 'C', n2, n1, -one, a( n1*n1 ), n1, one,
353  \$ a( 1 ), n1 )
354  CALL cpotrf( 'L', n2, a( 1 ), n1, info )
355  IF( info.GT.0 )
356  \$ info = info + n1
357 *
358  ELSE
359 *
360 * SRPA for UPPER, TRANSPOSE and N is odd
361 * T1 -> A(0,n1+1), T2 -> A(0,n1), S -> A(0,0)
362 * T1 -> a(n2*n2), T2 -> a(n1*n2), S -> a(0); lda = n2
363 *
364  CALL cpotrf( 'U', n1, a( n2*n2 ), n2, info )
365  IF( info.GT.0 )
366  \$ RETURN
367  CALL ctrsm( 'R', 'U', 'N', 'N', n2, n1, cone, a( n2*n2 ),
368  \$ n2, a( 0 ), n2 )
369  CALL cherk( 'L', 'N', n2, n1, -one, a( 0 ), n2, one,
370  \$ a( n1*n2 ), n2 )
371  CALL cpotrf( 'L', n2, a( n1*n2 ), n2, info )
372  IF( info.GT.0 )
373  \$ info = info + n1
374 *
375  END IF
376 *
377  END IF
378 *
379  ELSE
380 *
381 * N is even
382 *
383  IF( normaltransr ) THEN
384 *
385 * N is even and TRANSR = 'N'
386 *
387  IF( lower ) THEN
388 *
389 * SRPA for LOWER, NORMAL, and N is even ( a(0:n,0:k-1) )
390 * T1 -> a(1,0), T2 -> a(0,0), S -> a(k+1,0)
391 * T1 -> a(1), T2 -> a(0), S -> a(k+1)
392 *
393  CALL cpotrf( 'L', k, a( 1 ), n+1, info )
394  IF( info.GT.0 )
395  \$ RETURN
396  CALL ctrsm( 'R', 'L', 'C', 'N', k, k, cone, a( 1 ), n+1,
397  \$ a( k+1 ), n+1 )
398  CALL cherk( 'U', 'N', k, k, -one, a( k+1 ), n+1, one,
399  \$ a( 0 ), n+1 )
400  CALL cpotrf( 'U', k, a( 0 ), n+1, info )
401  IF( info.GT.0 )
402  \$ info = info + k
403 *
404  ELSE
405 *
406 * SRPA for UPPER, NORMAL, and N is even ( a(0:n,0:k-1) )
407 * T1 -> a(k+1,0) , T2 -> a(k,0), S -> a(0,0)
408 * T1 -> a(k+1), T2 -> a(k), S -> a(0)
409 *
410  CALL cpotrf( 'L', k, a( k+1 ), n+1, info )
411  IF( info.GT.0 )
412  \$ RETURN
413  CALL ctrsm( 'L', 'L', 'N', 'N', k, k, cone, a( k+1 ),
414  \$ n+1, a( 0 ), n+1 )
415  CALL cherk( 'U', 'C', k, k, -one, a( 0 ), n+1, one,
416  \$ a( k ), n+1 )
417  CALL cpotrf( 'U', k, a( k ), n+1, info )
418  IF( info.GT.0 )
419  \$ info = info + k
420 *
421  END IF
422 *
423  ELSE
424 *
425 * N is even and TRANSR = 'C'
426 *
427  IF( lower ) THEN
428 *
429 * SRPA for LOWER, TRANSPOSE and N is even (see paper)
430 * T1 -> B(0,1), T2 -> B(0,0), S -> B(0,k+1)
431 * T1 -> a(0+k), T2 -> a(0+0), S -> a(0+k*(k+1)); lda=k
432 *
433  CALL cpotrf( 'U', k, a( 0+k ), k, info )
434  IF( info.GT.0 )
435  \$ RETURN
436  CALL ctrsm( 'L', 'U', 'C', 'N', k, k, cone, a( k ), n1,
437  \$ a( k*( k+1 ) ), k )
438  CALL cherk( 'L', 'C', k, k, -one, a( k*( k+1 ) ), k, one,
439  \$ a( 0 ), k )
440  CALL cpotrf( 'L', k, a( 0 ), k, info )
441  IF( info.GT.0 )
442  \$ info = info + k
443 *
444  ELSE
445 *
446 * SRPA for UPPER, TRANSPOSE and N is even (see paper)
447 * T1 -> B(0,k+1), T2 -> B(0,k), S -> B(0,0)
448 * T1 -> a(0+k*(k+1)), T2 -> a(0+k*k), S -> a(0+0)); lda=k
449 *
450  CALL cpotrf( 'U', k, a( k*( k+1 ) ), k, info )
451  IF( info.GT.0 )
452  \$ RETURN
453  CALL ctrsm( 'R', 'U', 'N', 'N', k, k, cone,
454  \$ a( k*( k+1 ) ), k, a( 0 ), k )
455  CALL cherk( 'L', 'N', k, k, -one, a( 0 ), k, one,
456  \$ a( k*k ), k )
457  CALL cpotrf( 'L', k, a( k*k ), k, info )
458  IF( info.GT.0 )
459  \$ info = info + k
460 *
461  END IF
462 *
463  END IF
464 *
465  END IF
466 *
467  RETURN
468 *
469 * End of CPFTRF
470 *
subroutine ctrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
CTRSM
Definition: ctrsm.f:182
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
CHERK
Definition: cherk.f:175
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cpotrf(UPLO, N, A, LDA, INFO)
CPOTRF
Definition: cpotrf.f:109
Here is the call graph for this function:
Here is the caller graph for this function: