FLA_UDdate_UT_unb_var1.c File Reference

(r)

Functions
FLA_Error	FLA_UDdate_UT_unb_var1 (FLA_Obj R, FLA_Obj C, FLA_Obj D, FLA_Obj T)

---

Function Documentation

FLA_Error FLA_UDdate_UT_unb_var1	(	FLA_Obj	R,
		FLA_Obj	C,
		FLA_Obj	D,
		FLA_Obj	T
	)

References FLA_Apply_HUD_UT(), FLA_Cont_with_1x3_to_1x2(), FLA_Cont_with_3x3_to_2x2(), FLA_Herk_external(), FLA_Househ3UD_UT(), FLA_MINUS_ONE, FLA_Obj_min_dim(), FLA_ONE, FLA_ONE_HALF, FLA_Part_1x2(), FLA_Part_2x2(), FLA_Repart_1x2_to_1x3(), FLA_Repart_2x2_to_3x3(), FLA_Scale_diag(), and FLA_Set_to_identity().

Referenced by FLA_UDdate_UT_internal().

{
  FLA_Obj RTL,   RTR,      R00,  r01,   R02, 
          RBL,   RBR,      r10t, rho11, r12t,
                           R20,  r21,   R22;

  FLA_Obj CL,    CR,       C0,  c1,  C2;

  FLA_Obj DL,    DR,       D0,  d1,  D2;

  FLA_Obj TTL,   TTR,      T00,  t01,   T02,
          TBL,   TBR,      t10t, tau11, w12t,
                           T20,  t21,   T22;

  FLA_Part_2x2( R,    &RTL, &RTR,
                      &RBL, &RBR,     0, 0, FLA_TL );

  FLA_Part_1x2( C,    &CL,  &CR,      0, FLA_LEFT );

  FLA_Part_1x2( D,    &DL,  &DR,      0, FLA_LEFT );

  FLA_Part_2x2( T,    &TTL, &TTR,
                      &TBL, &TBR,     0, 0, FLA_TL );

  while ( FLA_Obj_min_dim( RBR ) > 0 ){

    FLA_Repart_2x2_to_3x3( RTL, /**/ RTR,       &R00,  /**/ &r01,   &R02,
                        /* ************* */   /* ************************** */
                                                &r10t, /**/ &rho11, &r12t,
                           RBL, /**/ RBR,       &R20,  /**/ &r21,   &R22,
                           1, 1, FLA_BR );

    FLA_Repart_1x2_to_1x3( CL,  /**/ CR,        &C0, /**/ &c1, &C2,
                           1, FLA_RIGHT );

    FLA_Repart_1x2_to_1x3( DL,  /**/ DR,        &D0, /**/ &d1, &D2,
                           1, FLA_RIGHT );

    FLA_Repart_2x2_to_3x3( TTL, /**/ TTR,       &T00,  /**/ &t01,   &T02,
                        /* ************* */   /* ************************ */
                                                &t10t, /**/ &tau11, &w12t,
                           TBL, /**/ TBR,       &T20,  /**/ &t21,   &T22,
                           1, 1, FLA_BR );

    /*------------------------------------------------------------*/

    // Compute tau11, u1, and v2 from rho11, c1, and d1 such that tau11, u1,
    // and v1 determine an up/downdating UT Householder transform H such that
    // applying H from the left to the column vector consisting of rho11, c1,
    // and d1 annihilates the entries in c1 and d1 (and updates rho11).
    FLA_Househ3UD_UT( rho11,
                      c1,
                      d1, tau11 );

    // / r12t \       / r12t \ 
    // |  C2  | =  H' |  C2  | 
    // \  D2  /       \  D2  / 
    //
    // where H is formed from tau11, u1 (stored in c1) and v1 (stored in d1).
    FLA_Apply_HUD_UT( FLA_LEFT,
                      tau11, w12t,
                             r12t,
                      c1,    C2,
                      d1,    D2 );

    /*------------------------------------------------------------*/

    FLA_Cont_with_3x3_to_2x2( &RTL, /**/ &RTR,       R00,  r01,   /**/ R02,
                                                     r10t, rho11, /**/ r12t,
                            /* ************** */  /* ************************ */
                              &RBL, /**/ &RBR,       R20,  r21,   /**/ R22,
                              FLA_TL );

    FLA_Cont_with_1x3_to_1x2( &CL,  /**/ &CR,        C0, c1, /**/ C2,
                              FLA_LEFT );

    FLA_Cont_with_1x3_to_1x2( &DL,  /**/ &DR,        D0, d1, /**/ D2,
                              FLA_LEFT );

    FLA_Cont_with_3x3_to_2x2( &TTL, /**/ &TTR,       T00,  t01,   /**/ T02,
                                                     t10t, tau11, /**/ w12t,
                            /* ************** */  /* ********************** */
                              &TBL, /**/ &TBR,       T20,  t21,   /**/ T22,
                              FLA_TL );
  }

  // T = I + C' * C - D' * D;
  // T = striu( T ) + 0.5*diag( T );
  // NOTE: The only reason this 'herk' method of computing T works is because
  // up-and-downdating is used to up/downdate a system that is being solved
  // either by QR factorization, or the method of normal equations (Cholesky
  // factorization on A' * A), and in either case, R will have a real diagonal.

  FLA_Set_to_identity( T );

  FLA_Herk_external( FLA_UPPER_TRIANGULAR, FLA_CONJ_TRANSPOSE,
                     FLA_ONE, C, FLA_ONE, T );

  FLA_Herk_external( FLA_UPPER_TRIANGULAR, FLA_CONJ_TRANSPOSE,
                     FLA_MINUS_ONE, D, FLA_ONE, T );

  FLA_Scale_diag( FLA_NO_CONJUGATE, FLA_ONE_HALF, T );

  return FLA_SUCCESS;
}