GSP
Quick Navigator

Search Site

Unix VPS
A - Starter
B - Basic
C - Preferred
D - Commercial
MPS - Dedicated
Previous VPSs
* Sign Up! *

Support
Contact Us
Online Help
Handbooks
Domain Status
Man Pages

FAQ
Virtual Servers
Pricing
Billing
Technical

Network
Facilities
Connectivity
Topology Map

Miscellaneous
Server Agreement
Year 2038
Credits
 

USA Flag

 

 

Man Pages
slarfb_gett.f(3) LAPACK slarfb_gett.f(3)

slarfb_gett.f


subroutine slarfb_gett (IDENT, M, N, K, T, LDT, A, LDA, B, LDB, WORK, LDWORK)
SLARFB_GETT

SLARFB_GETT

Purpose:

 SLARFB_GETT applies a real Householder block reflector H from the
 left to a real (K+M)-by-N  "triangular-pentagonal" matrix
 composed of two block matrices: an upper trapezoidal K-by-N matrix A
 stored in the array A, and a rectangular M-by-(N-K) matrix B, stored
 in the array B. The block reflector H is stored in a compact
 WY-representation, where the elementary reflectors are in the
 arrays A, B and T. See Further Details section.

Parameters

IDENT

          IDENT is CHARACTER*1
          If IDENT = not 'I', or not 'i', then V1 is unit
             lower-triangular and stored in the left K-by-K block of
             the input matrix A,
          If IDENT = 'I' or 'i', then  V1 is an identity matrix and
             not stored.
          See Further Details section.

M

          M is INTEGER
          The number of rows of the matrix B.
          M >= 0.

N

          N is INTEGER
          The number of columns of the matrices A and B.
          N >= 0.

K

          K is INTEGER
          The number or rows of the matrix A.
          K is also order of the matrix T, i.e. the number of
          elementary reflectors whose product defines the block
          reflector. 0 <= K <= N.

T

          T is REAL array, dimension (LDT,K)
          The upper-triangular K-by-K matrix T in the representation
          of the block reflector.

LDT

          LDT is INTEGER
          The leading dimension of the array T. LDT >= K.

A

          A is REAL array, dimension (LDA,N)
          On entry:
           a) In the K-by-N upper-trapezoidal part A: input matrix A.
           b) In the columns below the diagonal: columns of V1
              (ones are not stored on the diagonal).
          On exit:
            A is overwritten by rectangular K-by-N product H*A.
          See Further Details section.

LDA

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

B

          B is REAL array, dimension (LDB,N)
          On entry:
            a) In the M-by-(N-K) right block: input matrix B.
            b) In the M-by-N left block: columns of V2.
          On exit:
            B is overwritten by rectangular M-by-N product H*B.
          See Further Details section.

LDB

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

WORK

          WORK is REAL array,
          dimension (LDWORK,max(K,N-K))

LDWORK

          LDWORK is INTEGER
          The leading dimension of the array WORK. LDWORK>=max(1,K).

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors:

 November 2020, Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

Further Details:

    (1) Description of the Algebraic Operation.
    The matrix A is a K-by-N matrix composed of two column block
    matrices, A1, which is K-by-K, and A2, which is K-by-(N-K):
    A = ( A1, A2 ).
    The matrix B is an M-by-N matrix composed of two column block
    matrices, B1, which is M-by-K, and B2, which is M-by-(N-K):
    B = ( B1, B2 ).
    Perform the operation:
       ( A_out ) := H * ( A_in ) = ( I - V * T * V**T ) * ( A_in ) =
       ( B_out )        ( B_in )                          ( B_in )
                  = ( I - ( V1 ) * T * ( V1**T, V2**T ) ) * ( A_in )
                          ( V2 )                            ( B_in )
     On input:
    a) ( A_in )  consists of two block columns:
       ( B_in )
       ( A_in ) = (( A1_in ) ( A2_in )) = (( A1_in ) ( A2_in ))
       ( B_in )   (( B1_in ) ( B2_in ))   ((     0 ) ( B2_in )),
       where the column blocks are:
       (  A1_in )  is a K-by-K upper-triangular matrix stored in the
                   upper triangular part of the array A(1:K,1:K).
       (  B1_in )  is an M-by-K rectangular ZERO matrix and not stored.
       ( A2_in )  is a K-by-(N-K) rectangular matrix stored
                  in the array A(1:K,K+1:N).
       ( B2_in )  is an M-by-(N-K) rectangular matrix stored
                  in the array B(1:M,K+1:N).
    b) V = ( V1 )
           ( V2 )
       where:
       1) if IDENT == 'I',V1 is a K-by-K identity matrix, not stored;
       2) if IDENT != 'I',V1 is a K-by-K unit lower-triangular matrix,
          stored in the lower-triangular part of the array
          A(1:K,1:K) (ones are not stored),
       and V2 is an M-by-K rectangular stored the array B(1:M,1:K),
                 (because on input B1_in is a rectangular zero
                  matrix that is not stored and the space is
                  used to store V2).
    c) T is a K-by-K upper-triangular matrix stored
       in the array T(1:K,1:K).
    On output:
    a) ( A_out ) consists of two  block columns:
       ( B_out )
       ( A_out ) = (( A1_out ) ( A2_out ))
       ( B_out )   (( B1_out ) ( B2_out )),
       where the column blocks are:
       ( A1_out )  is a K-by-K square matrix, or a K-by-K
                   upper-triangular matrix, if V1 is an
                   identity matrix. AiOut is stored in
                   the array A(1:K,1:K).
       ( B1_out )  is an M-by-K rectangular matrix stored
                   in the array B(1:M,K:N).
       ( A2_out )  is a K-by-(N-K) rectangular matrix stored
                   in the array A(1:K,K+1:N).
       ( B2_out )  is an M-by-(N-K) rectangular matrix stored
                   in the array B(1:M,K+1:N).
    The operation above can be represented as the same operation
    on each block column:
       ( A1_out ) := H * ( A1_in ) = ( I - V * T * V**T ) * ( A1_in )
       ( B1_out )        (     0 )                          (     0 )
       ( A2_out ) := H * ( A2_in ) = ( I - V * T * V**T ) * ( A2_in )
       ( B2_out )        ( B2_in )                          ( B2_in )
    If IDENT != 'I':
       The computation for column block 1:
       A1_out: = A1_in - V1*T*(V1**T)*A1_in
       B1_out: = - V2*T*(V1**T)*A1_in
       The computation for column block 2, which exists if N > K:
       A2_out: = A2_in - V1*T*( (V1**T)*A2_in + (V2**T)*B2_in )
       B2_out: = B2_in - V2*T*( (V1**T)*A2_in + (V2**T)*B2_in )
    If IDENT == 'I':
       The operation for column block 1:
       A1_out: = A1_in - V1*T**A1_in
       B1_out: = - V2*T**A1_in
       The computation for column block 2, which exists if N > K:
       A2_out: = A2_in - T*( A2_in + (V2**T)*B2_in )
       B2_out: = B2_in - V2*T*( A2_in + (V2**T)*B2_in )
    (2) Description of the Algorithmic Computation.
    In the first step, we compute column block 2, i.e. A2 and B2.
    Here, we need to use the K-by-(N-K) rectangular workspace
    matrix W2 that is of the same size as the matrix A2.
    W2 is stored in the array WORK(1:K,1:(N-K)).
    In the second step, we compute column block 1, i.e. A1 and B1.
    Here, we need to use the K-by-K square workspace matrix W1
    that is of the same size as the as the matrix A1.
    W1 is stored in the array WORK(1:K,1:K).
    NOTE: Hence, in this routine, we need the workspace array WORK
    only of size WORK(1:K,1:max(K,N-K)) so it can hold both W2 from
    the first step and W1 from the second step.
    Case (A), when V1 is unit lower-triangular, i.e. IDENT != 'I',
    more computations than in the Case (B).
    if( IDENT != 'I' ) then
     if ( N > K ) then
       (First Step - column block 2)
       col2_(1) W2: = A2
       col2_(2) W2: = (V1**T) * W2 = (unit_lower_tr_of_(A1)**T) * W2
       col2_(3) W2: = W2 + (V2**T) * B2 = W2 + (B1**T) * B2
       col2_(4) W2: = T * W2
       col2_(5) B2: = B2 - V2 * W2 = B2 - B1 * W2
       col2_(6) W2: = V1 * W2 = unit_lower_tr_of_(A1) * W2
       col2_(7) A2: = A2 - W2
     else
       (Second Step - column block 1)
       col1_(1) W1: = A1
       col1_(2) W1: = (V1**T) * W1 = (unit_lower_tr_of_(A1)**T) * W1
       col1_(3) W1: = T * W1
       col1_(4) B1: = - V2 * W1 = - B1 * W1
       col1_(5) square W1: = V1 * W1 = unit_lower_tr_of_(A1) * W1
       col1_(6) square A1: = A1 - W1
     end if
    end if
    Case (B), when V1 is an identity matrix, i.e. IDENT == 'I',
    less computations than in the Case (A)
    if( IDENT == 'I' ) then
     if ( N > K ) then
       (First Step - column block 2)
       col2_(1) W2: = A2
       col2_(3) W2: = W2 + (V2**T) * B2 = W2 + (B1**T) * B2
       col2_(4) W2: = T * W2
       col2_(5) B2: = B2 - V2 * W2 = B2 - B1 * W2
       col2_(7) A2: = A2 - W2
     else
       (Second Step - column block 1)
       col1_(1) W1: = A1
       col1_(3) W1: = T * W1
       col1_(4) B1: = - V2 * W1 = - B1 * W1
       col1_(6) upper-triangular_of_(A1): = A1 - W1
     end if
    end if
    Combine these cases (A) and (B) together, this is the resulting
    algorithm:
    if ( N > K ) then
      (First Step - column block 2)
      col2_(1)  W2: = A2
      if( IDENT != 'I' ) then
        col2_(2)  W2: = (V1**T) * W2
                      = (unit_lower_tr_of_(A1)**T) * W2
      end if
      col2_(3)  W2: = W2 + (V2**T) * B2 = W2 + (B1**T) * B2]
      col2_(4)  W2: = T * W2
      col2_(5)  B2: = B2 - V2 * W2 = B2 - B1 * W2
      if( IDENT != 'I' ) then
        col2_(6)    W2: = V1 * W2 = unit_lower_tr_of_(A1) * W2
      end if
      col2_(7) A2: = A2 - W2
    else
    (Second Step - column block 1)
      col1_(1) W1: = A1
      if( IDENT != 'I' ) then
        col1_(2) W1: = (V1**T) * W1
                    = (unit_lower_tr_of_(A1)**T) * W1
      end if
      col1_(3) W1: = T * W1
      col1_(4) B1: = - V2 * W1 = - B1 * W1
      if( IDENT != 'I' ) then
        col1_(5) square W1: = V1 * W1 = unit_lower_tr_of_(A1) * W1
        col1_(6_a) below_diag_of_(A1): =  - below_diag_of_(W1)
      end if
      col1_(6_b) up_tr_of_(A1): = up_tr_of_(A1) - up_tr_of_(W1)
    end if

Definition at line 390 of file slarfb_gett.f.

Generated automatically by Doxygen for LAPACK from the source code.
Mon Jun 28 2021 Version 3.10.0

Search for    or go to Top of page |  Section 3 |  Main Index

Powered by GSP Visit the GSP FreeBSD Man Page Interface.
Output converted with ManDoc.