FLASH_main_prototypes.h File Reference

(r)

Go to the source code of this file.

Functions

FLA_Error FLASH_Obj_blocksizes_check (FLA_Obj H, dim_t *b_m, dim_t *b_n)
FLA_Error FLASH_Obj_create_helper_check (FLA_Bool without_buffer, FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_hierarchy_check (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *elem_sizes_m, dim_t *elem_sizes_n, FLA_Obj flat_matrix, FLA_Obj *H, unsigned long id, dim_t depth_overall, dim_t *depth_sizes_m, dim_t *depth_sizes_n, dim_t *m_offsets, dim_t *n_offsets)
FLA_Error FLASH_Obj_create_conf_to_check (FLA_Trans trans, FLA_Obj H_cur, FLA_Obj *H_new)
FLA_Error FLASH_Obj_create_hier_conf_to_flat_check (FLA_Trans trans, FLA_Obj F, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
FLA_Error FLASH_Obj_create_hier_conf_to_flat_ext_check (FLA_Trans trans, FLA_Obj F, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_flat_conf_to_hier_check (FLA_Trans trans, FLA_Obj H, FLA_Obj *F)
FLA_Error FLASH_Obj_create_hier_copy_of_flat_check (FLA_Obj F, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
FLA_Error FLASH_Obj_create_hier_copy_of_flat_ext_check (FLA_Obj F, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_flat_copy_of_hier_check (FLA_Obj H, FLA_Obj *F)
FLA_Error FLASH_Obj_free_check (FLA_Obj *H)
FLA_Error FLASH_Obj_free_without_buffer_check (FLA_Obj *H)
FLA_Error FLASH_Obj_free_hierarchy_check (FLA_Obj *H)
FLA_Error FLASH_Obj_attach_buffer_check (void *buffer, dim_t ldim, FLA_Obj *H)
FLA_Error FLASH_Obj_attach_buffer_hierarchy_check (FLA_Obj F, FLA_Obj *H)
FLA_Error FLASH_Axpy_submatrix_to_global (FLA_Obj alpha, dim_t m, dim_t n, void *buffer, dim_t ldim, dim_t i, dim_t j, FLA_Obj H)
FLA_Error FLASH_Axpy_global_to_submatrix (FLA_Obj alpha, dim_t i, dim_t j, FLA_Obj H, dim_t m, dim_t n, void *buffer, dim_t ldim)
FLA_Error FLASH_Axpy_subobject_to_global (FLA_Obj alpha, FLA_Obj F, dim_t i, dim_t j, FLA_Obj H)
FLA_Error FLASH_Axpy_global_to_subobject (FLA_Obj alpha, dim_t i, dim_t j, FLA_Obj H, FLA_Obj F)
FLA_Error FLASH_Axpy_hierarchy (int direction, FLA_Obj alpha, FLA_Obj F, dim_t i, dim_t j, FLA_Obj *H)
FLA_Error FLASH_Axpy_hierarchy_r (int direction, FLA_Obj alpha, FLA_Obj F, dim_t i, dim_t j, FLA_Obj *H, FLASH_Obj_queue *queue)
FLA_Error FLASH_Copy_submatrix_to_global (dim_t m, dim_t n, void *buffer, dim_t ldim, dim_t i, dim_t j, FLA_Obj H)
FLA_Error FLASH_Copy_global_to_submatrix (dim_t i, dim_t j, FLA_Obj H, dim_t m, dim_t n, void *buffer, dim_t ldim)
FLA_Error FLASH_Copy_subobject_to_global (FLA_Obj F, dim_t i, dim_t j, FLA_Obj H)
FLA_Error FLASH_Copy_global_to_subobject (dim_t i, dim_t j, FLA_Obj H, FLA_Obj F)
FLA_Error FLASH_Copy_hierarchy (int direction, FLA_Obj F, dim_t i, dim_t j, FLA_Obj *H)
FLA_Error FLASH_Copy_hierarchy_r (int direction, FLA_Obj F, dim_t i, dim_t j, FLA_Obj *H, FLASH_Obj_queue *queue)
FLA_Datatype FLASH_Obj_datatype (FLA_Obj H)
dim_t FLASH_Obj_scalar_length (FLA_Obj H)
dim_t FLASH_Obj_scalar_width (FLA_Obj H)
dim_t FLASH_Obj_depth (FLA_Obj H)
dim_t FLASH_Obj_blocksizes (FLA_Obj H, dim_t *b_m, dim_t *b_n)
FLA_Error FLASH_Obj_create (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
FLA_Error FLASH_Obj_create_ext (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_without_buffer (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
FLA_Error FLASH_Obj_create_without_buffer_ext (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_helper (FLA_Bool without_buffer, FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_hierarchy (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *elem_sizes_m, dim_t *elem_sizes_n, FLA_Obj flat_matrix, FLA_Obj *H, unsigned long id, dim_t depth_overall, dim_t *depth_sizes_m, dim_t *depth_sizes_n, dim_t *m_offsets, dim_t *n_offsets)
FLA_Error FLASH_Obj_create_conf_to (FLA_Trans trans, FLA_Obj H_cur, FLA_Obj *H_new)
FLA_Error FLASH_Obj_create_hier_conf_to_flat (FLA_Trans trans, FLA_Obj F, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
FLA_Error FLASH_Obj_create_hier_conf_to_flat_ext (FLA_Trans trans, FLA_Obj F, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_flat_conf_to_hier (FLA_Trans trans, FLA_Obj H, FLA_Obj *F)
FLA_Error FLASH_Obj_create_hier_copy_of_flat (FLA_Obj F, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
FLA_Error FLASH_Obj_create_hier_copy_of_flat_ext (FLA_Obj F, dim_t depth, dim_t *b_m, dim_t *b_n, FLA_Obj *H)
FLA_Error FLASH_Obj_create_flat_copy_of_hier (FLA_Obj H, FLA_Obj *F)
void FLASH_Obj_free (FLA_Obj *H)
void FLASH_Obj_free_hierarchy (FLA_Obj *H)
void FLASH_Obj_free_without_buffer (FLA_Obj *H)
FLA_Error FLASH_Obj_attach_buffer (void *buffer, dim_t ldim, FLA_Obj *H)
FLA_Error FLASH_Obj_attach_buffer_hierarchy (FLA_Obj F, FLA_Obj *H)
FLA_Error FLASH_Obj_flatten (FLA_Obj H, FLA_Obj F)
FLA_Error FLASH_Obj_hierarchify (FLA_Obj F, FLA_Obj H)
void * FLASH_Obj_extract_buffer (FLA_Obj H)
FLA_Error FLASH_Obj_show (char *header, FLA_Obj H, char *elem_format, char *footer)
void FLASH_print_struct (FLA_Obj H)
void FLASH_print_struct_helper (FLA_Obj H, int indent)
void FLASH_Obj_exec (void *func, FLASH_Obj_queue *queue)
void * FLASH_Obj_exec_parallel (void *arg)
void FLASH_Obj_push (int direction, FLA_Obj alpha, FLA_Obj F, FLA_Obj H, FLASH_Obj_queue *queue)


Function Documentation

FLA_Error FLASH_Axpy_global_to_submatrix ( FLA_Obj  alpha,
dim_t  i,
dim_t  j,
FLA_Obj  H,
dim_t  m,
dim_t  n,
void *  buffer,
dim_t  ldim 
)

References FLA_Check_consistent_object_datatype(), FLA_Check_error_level(), FLA_Check_if_scalar(), FLA_Check_submatrix_dims_and_offset(), FLA_Obj_attach_buffer(), FLA_Obj_create_without_buffer(), FLA_Obj_free_without_buffer(), FLASH_Axpy_hierarchy(), and FLASH_Obj_datatype().

00073 {
00074     FLA_Obj      flat_matrix;
00075     FLA_Datatype datatype;
00076     FLA_Error    e_val;
00077 
00078     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00079     {
00080         e_val = FLA_Check_if_scalar( alpha );
00081         FLA_Check_error_code( e_val );
00082 
00083         e_val = FLA_Check_consistent_object_datatype( alpha, H );
00084         FLA_Check_error_code( e_val );
00085 
00086         e_val = FLA_Check_submatrix_dims_and_offset( m, n, i, j, H );
00087         FLA_Check_error_code( e_val );
00088     }
00089 
00090     // Acquire the datatype from the hierarchical matrix object.
00091     datatype = FLASH_Obj_datatype( H );
00092 
00093     // Create a temporary conventional matrix object of the requested datatype
00094     // and dimensions and attach the given buffer containing the incoming data.
00095     FLA_Obj_create_without_buffer( datatype, m, n, &flat_matrix );
00096     FLA_Obj_attach_buffer( buffer, ldim, &flat_matrix );
00097 
00098     // Recurse through H, adding in the corresponding elements of flat_matrix,
00099     // starting at the (i,j) element offset.
00100     FLASH_Axpy_hierarchy( FLA_HIER_TO_FLAT, alpha, flat_matrix, i, j, &H );
00101 
00102     // Free the object (but don't free the buffer!).
00103     FLA_Obj_free_without_buffer( &flat_matrix );
00104 
00105     return FLA_SUCCESS;
00106 }

FLA_Error FLASH_Axpy_global_to_subobject ( FLA_Obj  alpha,
dim_t  i,
dim_t  j,
FLA_Obj  H,
FLA_Obj  F 
)

References FLASH_Axpy_hierarchy().

00118 {
00119     FLASH_Axpy_hierarchy( FLA_HIER_TO_FLAT, alpha, F, i, j, &H );
00120 
00121     return FLA_SUCCESS;
00122 }

FLA_Error FLASH_Axpy_hierarchy ( int  direction,
FLA_Obj  alpha,
FLA_Obj  F,
dim_t  i,
dim_t  j,
FLA_Obj H 
)

References FLASH_Axpy_hierarchy(), FLASH_Axpy_hierarchy_r(), FLASH_Obj_exec(), FLASH_Obj_queue_s::head, FLASH_Obj_queue_s::n_tasks, and FLASH_Obj_queue_s::tail.

Referenced by FLASH_Axpy_global_to_submatrix(), FLASH_Axpy_global_to_subobject(), FLASH_Axpy_hierarchy(), FLASH_Axpy_submatrix_to_global(), FLASH_Axpy_subobject_to_global(), and FLASH_Obj_exec_parallel().

00126 {
00127     FLA_Error e_val;
00128     FLASH_Obj_queue queue;
00129 
00130     queue.n_tasks = 0;
00131     queue.head = NULL;
00132     queue.tail = NULL;
00133 
00134     e_val = FLASH_Axpy_hierarchy_r( direction, alpha, F, i, j, H, &queue );
00135 
00136 #ifdef FLA_ENABLE_MULTITHREADING
00137     FLASH_Obj_exec( ( void* ) FLASH_Axpy_hierarchy, &queue );
00138 #endif
00139     return e_val;
00140 }

FLA_Error FLASH_Axpy_hierarchy_r ( int  direction,
FLA_Obj  alpha,
FLA_Obj  F,
dim_t  i,
dim_t  j,
FLA_Obj H,
FLASH_Obj_queue queue 
)

References FLA_Axpy_external(), FLA_Cont_with_1x3_to_1x2(), FLA_Cont_with_3x1_to_2x1(), FLA_Obj_buffer(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Part_2x1(), FLA_Part_2x2(), FLA_Repart_1x2_to_1x3(), FLA_Repart_2x1_to_3x1(), FLA_Repart_2x2_to_3x3(), FLASH_Axpy_hierarchy_r(), FLASH_Obj_push(), FLASH_Obj_scalar_length(), FLASH_Obj_scalar_width(), and FLASH_Queue_get_num_threads().

Referenced by FLASH_Axpy_hierarchy(), and FLASH_Axpy_hierarchy_r().

00144 {
00145     dim_t    m_scal, n_scal;
00146     dim_t    row_offset, col_offset;
00147     dim_t    b_m, b_n;
00148     FLA_Obj* buffer;
00149     FLA_Obj* h11p;
00150 
00151     FLA_Obj HTL,   HTR,    H00, H01, H02,
00152             HBL,   HBR,    H10, H11, H12,
00153                            H20, H21, H22;
00154 
00155     FLA_Obj F1T,             F01,
00156             F1B,             F11,
00157                              F21;
00158 
00159     FLA_Obj hT,              h01,
00160             hB,              h11,
00161                              h21;
00162 
00163     FLA_Obj FL,    FR,       F0,  F1,  F2;
00164 
00165     FLA_Obj HBR_L,  HBR_R,     h0,  h1,  h2;
00166 
00167     // Once we get down to a submatrix whose elements are scalars, we are down
00168     // to our base case.
00169     if ( FLA_Obj_elemtype( *H ) == FLA_SCALAR )
00170     {
00171         dim_t m_flat, n_flat;
00172 
00173         // Grab the length and width of the conventional matrix argument.
00174         m_flat = FLA_Obj_length( F );
00175         n_flat = FLA_Obj_width( F );
00176 
00177         // Advance to the (i,j) offset of the hierarchical matrix.
00178         FLA_Part_2x2( *H,   &HTL, &HTR,
00179                             &HBL, &HBR,     i, j, FLA_TL );
00180 
00181         // Repartition the destination submatrix, H, such that H11 is conformal
00182         // to dimensions of the flat conventional matrix F.
00183         FLA_Repart_2x2_to_3x3( HTL, /**/ HTR,       &H00, /**/ &H01, &H02,
00184                             /* ************* */   /* ******************** */
00185                                                     &H10, /**/ &H11, &H12,
00186                                HBL, /**/ HBR,       &H20, /**/ &H21, &H22,
00187                                m_flat, n_flat, FLA_BR );
00188 
00189 #ifdef FLA_ENABLE_MULTITHREADING
00190         if ( FLASH_Queue_get_num_threads() > 1 )
00191         {
00192             FLASH_Obj_push( direction, alpha, F, H11, queue );
00193         }   
00194         else
00195 #endif
00196         {
00197             // Depending on which top-level function invoked us, we either add the
00198             // the source data in the flat matrix to the leaf-level submatrix of the
00199             // hierarchical matrix, or add the data in the hierarchical submatrix to
00200             // the appropriate place in the flat matrix.
00201             if      ( direction == FLA_FLAT_TO_HIER )
00202             {
00203                 FLA_Axpy_external( alpha, F, H11 );
00204             }
00205             else if ( direction == FLA_HIER_TO_FLAT )
00206             {
00207                 FLA_Axpy_external( alpha, H11, F );
00208             }
00209         }
00210 
00211         // Return, since we're done for this leaf.
00212         return FLA_SUCCESS;
00213     }
00214 
00215     // If the above conditional did not hold, then we proceed.
00216 
00217     // Grab the buffer of H, whose elemtype we know must be FLA_MATRIX. The
00218     // elements of this buffer are objects.
00219     buffer = ( FLA_Obj* ) FLA_Obj_buffer( *H );
00220 
00221     // Extract the overall scalar dimensions of the first (top-left) matrix
00222     // referenced in the hierarchical matrix's buffer.
00223     m_scal = FLASH_Obj_scalar_length( buffer[0] );
00224     n_scal = FLASH_Obj_scalar_width( buffer[0] );
00225 
00226     // The i and j offsets given to the current function are relative to the
00227     // entire matrix, which is represented by H. Since the elemtype of H is
00228     // FLA_MATRIX, we now compute the number of submatrix elements we will skip
00229     // in both the m and n dimensions to advance ourselves to the proper
00230     // hierarchical submatrix within H, the one that contains the destination
00231     // for the data (or at least the top-left most portion of the data) in the
00232     // raw buffer matrix. Notice that integer division is performed here.
00233     b_m = i / m_scal;
00234     b_n = j / n_scal;
00235 
00236     // We use the indices computed above to advance through H such that the
00237     // (i,j) element represented by the (i,j) values passed to the current
00238     // function is situated at the top-left of HBR.
00239     FLA_Part_2x2( *H,   &HTL, &HTR,
00240                         &HBL, &HBR,     b_m, b_n, FLA_TL );
00241 
00242     // Prepare to march laterally through HBR.
00243     FLA_Part_1x2( HBR,     &HBR_L,  &HBR_R,    0, FLA_LEFT);
00244 
00245     // Prepare to march laterally through the conventional flat matrix.
00246     FLA_Part_1x2( F,       &FL,  &FR,          0, FLA_LEFT );
00247 
00248     // Iterate laterally through F.
00249     while ( FLA_Obj_width( FL ) < FLA_Obj_width( F ) )
00250     {
00251         // Handle the initial edge case and "interior" cases differently.
00252         if ( FLA_Obj_width( FL ) == 0 )
00253         {
00254             // Compute the width of the matrix that we will recurse into.
00255             // Since this is the first iteration of the while loop, this
00256             // blocksize may be smaller than subsequent full blocks if the
00257             // j offset passed to this routine is non-zero.
00258             b_n = n_scal - (j % n_scal);
00259 
00260             // Reduce the blocksize if the width of F is smaller than the
00261             // value computed above.
00262             b_n = min( FLA_Obj_width( FR ), b_n );
00263 
00264             // If the j offset is non-zero, we must modulo it to the
00265             // hierarchical blocksize for the element submatrices of H
00266             // in preparation for our recursion.
00267             col_offset = j % n_scal;
00268         }
00269         else
00270         {
00271             // The width of the matrix into which we will recursve may be
00272             // assumed to be the same as the width of the element submatrices
00273             // of H. 
00274             b_n = n_scal;
00275 
00276             // Reduce the blocksize if the width of F is smaller than the
00277             // element submatrices.
00278             b_n = min( FLA_Obj_width( FR ), b_n );
00279 
00280             // The column offset will be zero because we will be starting from
00281             // the left edge of the submatrix.
00282             col_offset = 0;
00283         }
00284 
00285         // Expose one column of hierarchical submatrices within HBR.
00286         FLA_Repart_1x2_to_1x3( HBR_L,  /**/ HBR_R,        &h0, /**/ &h1, &h2,
00287                                1, FLA_RIGHT );
00288 
00289         // Expose one panel of F of width b_n, computed above.
00290         FLA_Repart_1x2_to_1x3( FL,  /**/ FR,        &F0, /**/ &F1, &F2,
00291                                b_n, FLA_RIGHT );
00292 
00293         /*------------------------------------------------------------*/
00294 
00295         // Prepare to march vertically through h1 (the first submatrix column
00296         // of HBR).
00297         FLA_Part_2x1( h1,    &hT,
00298                              &hB,            0, FLA_TOP );
00299 
00300         // Prepare to march vertically through F.
00301         FLA_Part_2x1( F1,    &F1T,
00302                              &F1B,           0, FLA_TOP );
00303 
00304         // Iterate vertically through the current column-panel of F, F1.
00305         while ( FLA_Obj_length( F1T ) < FLA_Obj_length( F1 ) )
00306         {
00307 
00308             if ( FLA_Obj_length( F1T ) == 0 )
00309             {
00310                 // Compute the length of the matrix that we will recurse into.
00311                 // Since this is the first iteration of the while loop, this
00312                 // blocksize may be smaller than subsequent full blocks if the
00313                 // i offset passed to this routine is non-zero.
00314                 b_m = m_scal - (i % m_scal);
00315 
00316                 // Reduce the blocksize if the length of F is smaller than the
00317                 // value computed above.
00318                 b_m = min( FLA_Obj_length( F1B ), b_m );
00319 
00320                 // If the i offset is non-zero, we must modulo it to the
00321                 // hierarchical blocksize for the element submatrices of H
00322                 // in preparation for our recursion.
00323                 row_offset = i % m_scal;
00324             }
00325             else
00326             {
00327                 // The length of the matrix into which we will recursve may be
00328                 // assumed to be the same as the length of the element submatrices
00329                 // of H. 
00330                 b_m = m_scal;
00331 
00332                 // Reduce the blocksize if the length of the first panel of
00333                 // F is smaller than the the element submatrices.
00334                 b_m = min( FLA_Obj_length( F1B ), b_m );
00335 
00336                 // The column offset will be zero because we will be starting from
00337                 // the top edge of the submatrix.
00338                 row_offset = 0;
00339             }
00340 
00341             // Expose one hierarchical submatrix within the current column of h1.
00342             FLA_Repart_2x1_to_3x1( hT,                &h01,
00343                                 /* ** */            /* ** */
00344                                                       &h11,
00345                                    hB,                &h21,        1, FLA_BOTTOM );
00346 
00347             // Expose one block of F1 (the current panel of F) of length b_m,
00348             // computed above.
00349             FLA_Repart_2x1_to_3x1( F1T,               &F01,
00350                                 /* ** */            /* ** */
00351                                                       &F11,
00352                                    F1B,               &F21,      b_m, FLA_BOTTOM );
00353 
00354             /*------------------------------------------------------------*/
00355 
00356             // This is a pointer to the FLA_Obj element referenced by h11.
00357             h11p = ( FLA_Obj* ) FLA_Obj_buffer( h11 );
00358 
00359             // Recurse with new (i,j) offset parameters, a narrowed block of
00360             // F, and a corresponding hierarchical block of H.
00361             FLASH_Axpy_hierarchy_r( direction, alpha, F11, row_offset, col_offset, h11p, queue );
00362 
00363             /*------------------------------------------------------------*/
00364 
00365             FLA_Cont_with_3x1_to_2x1( &hT,                h01,
00366                                                           h11,
00367                                    /* ** */            /* ** */
00368                                       &hB,                h21,     FLA_TOP );
00369 
00370             FLA_Cont_with_3x1_to_2x1( &F1T,               F01,
00371                                                           F11,
00372                                     /* ** */           /* ** */
00373                                       &F1B,               F21,     FLA_TOP );
00374         }
00375 
00376         /*------------------------------------------------------------*/
00377 
00378         FLA_Cont_with_1x3_to_1x2( &FL,  /**/ &FR,        F0, F1, /**/ F2,
00379                                   FLA_LEFT );
00380 
00381         FLA_Cont_with_1x3_to_1x2( &HBR_L,  /**/ &HBR_R,        h0, h1, /**/ h2,
00382                                   FLA_LEFT );
00383     }
00384 
00385     return FLA_SUCCESS;
00386 
00387 }

FLA_Error FLASH_Axpy_submatrix_to_global ( FLA_Obj  alpha,
dim_t  m,
dim_t  n,
void *  buffer,
dim_t  ldim,
dim_t  i,
dim_t  j,
FLA_Obj  H 
)

References FLA_Check_consistent_object_datatype(), FLA_Check_error_level(), FLA_Check_if_scalar(), FLA_Check_submatrix_dims_and_offset(), FLA_Obj_attach_buffer(), FLA_Obj_create_without_buffer(), FLA_Obj_free_without_buffer(), FLASH_Axpy_hierarchy(), and FLASH_Obj_datatype().

00036 {
00037     FLA_Obj      flat_matrix;
00038     FLA_Datatype datatype;
00039     FLA_Error    e_val;
00040 
00041     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00042     {
00043         e_val = FLA_Check_if_scalar( alpha );
00044         FLA_Check_error_code( e_val );
00045 
00046         e_val = FLA_Check_consistent_object_datatype( alpha, H );
00047         FLA_Check_error_code( e_val );
00048 
00049         e_val = FLA_Check_submatrix_dims_and_offset( m, n, i, j, H );
00050         FLA_Check_error_code( e_val );
00051     }
00052 
00053     // Acquire the datatype from the hierarchical matrix object.
00054     datatype = FLASH_Obj_datatype( H );
00055 
00056     // Create a temporary conventional matrix object of the requested datatype
00057     // and dimensions and attach the given buffer containing the incoming data.
00058     FLA_Obj_create_without_buffer( datatype, m, n, &flat_matrix );
00059     FLA_Obj_attach_buffer( buffer, ldim, &flat_matrix );
00060 
00061     // Recurse through H, adding in the corresponding elements of flat_matrix,
00062     // starting at the (i,j) element offset.
00063     FLASH_Axpy_hierarchy( FLA_FLAT_TO_HIER, alpha, flat_matrix, i, j, &H );
00064 
00065     // Free the object (but don't free the buffer!).
00066     FLA_Obj_free_without_buffer( &flat_matrix );
00067 
00068     return FLA_SUCCESS;
00069 }

FLA_Error FLASH_Axpy_subobject_to_global ( FLA_Obj  alpha,
FLA_Obj  F,
dim_t  i,
dim_t  j,
FLA_Obj  H 
)

References FLASH_Axpy_hierarchy().

00110 {
00111     FLASH_Axpy_hierarchy( FLA_FLAT_TO_HIER, alpha, F, i, j, &H );
00112 
00113     return FLA_SUCCESS;
00114 }

FLA_Error FLASH_Copy_global_to_submatrix ( dim_t  i,
dim_t  j,
FLA_Obj  H,
dim_t  m,
dim_t  n,
void *  buffer,
dim_t  ldim 
)

References FLA_Check_error_level(), FLA_Check_submatrix_dims_and_offset(), FLA_Obj_attach_buffer(), FLA_Obj_create_without_buffer(), FLA_Obj_free_without_buffer(), FLASH_Copy_hierarchy(), and FLASH_Obj_datatype().

00067 {
00068     FLA_Obj      flat_matrix;
00069     FLA_Datatype datatype;
00070     FLA_Error    e_val;
00071 
00072     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00073     {
00074         e_val = FLA_Check_submatrix_dims_and_offset( m, n, i, j, H );
00075         FLA_Check_error_code( e_val );
00076     }
00077 
00078     // Acquire the datatype from the hierarchical matrix object.
00079     datatype = FLASH_Obj_datatype( H );
00080 
00081     // Create a temporary conventional matrix object of the requested datatype
00082     // and dimensions and attach the given buffer containing the incoming data.
00083     FLA_Obj_create_without_buffer( datatype, m, n, &flat_matrix );
00084     FLA_Obj_attach_buffer( buffer, ldim, &flat_matrix );
00085 
00086     // Recurse through H, adding in the corresponding elements of flat_matrix,
00087     // starting at the (i,j) element offset.
00088     FLASH_Copy_hierarchy( FLA_HIER_TO_FLAT, flat_matrix, i, j, &H );
00089 
00090     // Free the object (but don't free the buffer!).
00091     FLA_Obj_free_without_buffer( &flat_matrix );
00092 
00093     return FLA_SUCCESS;
00094 }

FLA_Error FLASH_Copy_global_to_subobject ( dim_t  i,
dim_t  j,
FLA_Obj  H,
FLA_Obj  F 
)

References FLASH_Copy_hierarchy().

Referenced by FLASH_Obj_create_flat_copy_of_hier(), and FLASH_Obj_flatten().

00106 {
00107     FLASH_Copy_hierarchy( FLA_HIER_TO_FLAT, F, i, j, &H );
00108 
00109     return FLA_SUCCESS;
00110 }

FLA_Error FLASH_Copy_hierarchy ( int  direction,
FLA_Obj  F,
dim_t  i,
dim_t  j,
FLA_Obj H 
)

References FLASH_Copy_hierarchy(), FLASH_Copy_hierarchy_r(), FLASH_Obj_exec(), FLASH_Obj_queue_s::head, FLASH_Obj_queue_s::n_tasks, and FLASH_Obj_queue_s::tail.

Referenced by FLASH_Copy_global_to_submatrix(), FLASH_Copy_global_to_subobject(), FLASH_Copy_hierarchy(), FLASH_Copy_submatrix_to_global(), FLASH_Copy_subobject_to_global(), and FLASH_Obj_exec_parallel().

00114 {
00115     FLA_Error e_val;
00116     FLASH_Obj_queue queue;
00117 
00118     queue.n_tasks = 0;
00119     queue.head = NULL;
00120     queue.tail = NULL;
00121 
00122     e_val = FLASH_Copy_hierarchy_r( direction, F, i, j, H, &queue );
00123 
00124 #ifdef FLA_ENABLE_MULTITHREADING
00125     FLASH_Obj_exec( ( void* ) FLASH_Copy_hierarchy, &queue );
00126 #endif
00127     return e_val;
00128 }

FLA_Error FLASH_Copy_hierarchy_r ( int  direction,
FLA_Obj  F,
dim_t  i,
dim_t  j,
FLA_Obj H,
FLASH_Obj_queue queue 
)

References FLA_Cont_with_1x3_to_1x2(), FLA_Cont_with_3x1_to_2x1(), FLA_Copy_external(), FLA_Obj_buffer(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_width(), FLA_ONE, FLA_Part_1x2(), FLA_Part_2x1(), FLA_Part_2x2(), FLA_Repart_1x2_to_1x3(), FLA_Repart_2x1_to_3x1(), FLA_Repart_2x2_to_3x3(), FLASH_Copy_hierarchy_r(), FLASH_Obj_push(), FLASH_Obj_scalar_length(), FLASH_Obj_scalar_width(), and FLASH_Queue_get_num_threads().

Referenced by FLASH_Copy_hierarchy(), and FLASH_Copy_hierarchy_r().

00132 {
00133     dim_t    m_scal, n_scal;
00134     dim_t    row_offset, col_offset;
00135     dim_t    b_m, b_n;
00136     FLA_Obj* buffer;
00137     FLA_Obj* h11p;
00138 
00139     FLA_Obj HTL,   HTR,    H00, H01, H02,
00140             HBL,   HBR,    H10, H11, H12,
00141                            H20, H21, H22;
00142 
00143     FLA_Obj F1T,             F01,
00144             F1B,             F11,
00145                              F21;
00146 
00147     FLA_Obj hT,              h01,
00148             hB,              h11,
00149                              h21;
00150 
00151     FLA_Obj FL,    FR,       F0,  F1,  F2;
00152 
00153     FLA_Obj HBR_L,  HBR_R,     h0,  h1,  h2;
00154 
00155     // Once we get down to a submatrix whose elements are scalars, we are down
00156     // to our base case.
00157     if ( FLA_Obj_elemtype( *H ) == FLA_SCALAR )
00158     {
00159         dim_t m_flat, n_flat;
00160 
00161         // Grab the length and width of the conventional matrix argument.
00162         m_flat = FLA_Obj_length( F );
00163         n_flat = FLA_Obj_width( F );
00164 
00165         // Advance to the (i,j) offset of the hierarchical matrix.
00166         FLA_Part_2x2( *H,   &HTL, &HTR,
00167                             &HBL, &HBR,     i, j, FLA_TL );
00168 
00169         // Repartition the destination submatrix, H, such that H11 is conformal
00170         // to dimensions of the flat conventional matrix F.
00171         FLA_Repart_2x2_to_3x3( HTL, /**/ HTR,       &H00, /**/ &H01, &H02,
00172                             /* ************* */   /* ******************** */
00173                                                     &H10, /**/ &H11, &H12,
00174                                HBL, /**/ HBR,       &H20, /**/ &H21, &H22,
00175                                m_flat, n_flat, FLA_BR );
00176 
00177 #ifdef FLA_ENABLE_MULTITHREADING
00178         if ( FLASH_Queue_get_num_threads() > 1 )
00179         {
00180             FLASH_Obj_push( direction, FLA_ONE, F, H11, queue );
00181         }   
00182         else
00183 #endif  
00184         {
00185             // Depending on which top-level function invoked us, we either add the
00186             // the source data in the flat matrix to the leaf-level submatrix of the
00187             // hierarchical matrix, or add the data in the hierarchical submatrix to
00188             // the appropriate place in the flat matrix.
00189             if      ( direction == FLA_FLAT_TO_HIER )
00190             {
00191                 FLA_Copy_external( F, H11 );
00192             }
00193             else if ( direction == FLA_HIER_TO_FLAT )
00194             {
00195                 FLA_Copy_external( H11, F );
00196             }
00197         }
00198 
00199         // Return, since we're done for this leaf.
00200         return FLA_SUCCESS;
00201     }
00202 
00203     // If the above conditional did not hold, then we proceed.
00204 
00205     // Grab the buffer of H, whose elemtype we know must be FLA_MATRIX. The
00206     // elements of this buffer are objects.
00207     buffer = ( FLA_Obj* ) FLA_Obj_buffer( *H );
00208 
00209     // Extract the overall scalar dimensions of the first (top-left) matrix
00210     // referenced in the hierarchical matrix's buffer.
00211     m_scal = FLASH_Obj_scalar_length( buffer[0] );
00212     n_scal = FLASH_Obj_scalar_width( buffer[0] );
00213 
00214     // The i and j offsets given to the current function are relative to the
00215     // entire matrix, which is represented by H. Since the elemtype of H is
00216     // FLA_MATRIX, we now compute the number of submatrix elements we will skip
00217     // in both the m and n dimensions to advance ourselves to the proper
00218     // hierarchical submatrix within H, the one that contains the destination
00219     // for the data (or at least the top-left most portion of the data) in the
00220     // raw buffer matrix. Notice that integer division is performed here.
00221     b_m = i / m_scal;
00222     b_n = j / n_scal;
00223 
00224     // We use the indices computed above to advance through H such that the
00225     // (i,j) element represented by the (i,j) values passed to the current
00226     // function is situated at the top-left of HBR.
00227     FLA_Part_2x2( *H,   &HTL, &HTR,
00228                         &HBL, &HBR,     b_m, b_n, FLA_TL );
00229 
00230     // Prepare to march laterally through HBR.
00231     FLA_Part_1x2( HBR,     &HBR_L,  &HBR_R,    0, FLA_LEFT);
00232 
00233     // Prepare to march laterally through the conventional flat matrix.
00234     FLA_Part_1x2( F,       &FL,  &FR,          0, FLA_LEFT );
00235 
00236     // Iterate laterally through F.
00237     while ( FLA_Obj_width( FL ) < FLA_Obj_width( F ) )
00238     {
00239         // Handle the initial edge case and "interior" cases differently.
00240         if ( FLA_Obj_width( FL ) == 0 )
00241         {
00242             // Compute the width of the matrix that we will recurse into.
00243             // Since this is the first iteration of the while loop, this
00244             // blocksize may be smaller than subsequent full blocks if the
00245             // j offset passed to this routine is non-zero.
00246             b_n = n_scal - (j % n_scal);
00247 
00248             // Reduce the blocksize if the width of F is smaller than the
00249             // value computed above.
00250             b_n = min( FLA_Obj_width( FR ), b_n );
00251 
00252             // If the j offset is non-zero, we must modulo it to the
00253             // hierarchical blocksize for the element submatrices of H
00254             // in preparation for our recursion.
00255             col_offset = j % n_scal;
00256         }
00257         else
00258         {
00259             // The width of the matrix into which we will recursve may be
00260             // assumed to be the same as the width of the element submatrices
00261             // of H. 
00262             b_n = n_scal;
00263 
00264             // Reduce the blocksize if the width of F is smaller than the
00265             // element submatrices.
00266             b_n = min( FLA_Obj_width( FR ), b_n );
00267 
00268             // The column offset will be zero because we will be starting from
00269             // the left edge of the submatrix.
00270             col_offset = 0;
00271         }
00272 
00273         // Expose one column of hierarchical submatrices within HBR.
00274         FLA_Repart_1x2_to_1x3( HBR_L,  /**/ HBR_R,        &h0, /**/ &h1, &h2,
00275                                1, FLA_RIGHT );
00276 
00277         // Expose one panel of F of width b_n, computed above.
00278         FLA_Repart_1x2_to_1x3( FL,  /**/ FR,        &F0, /**/ &F1, &F2,
00279                                b_n, FLA_RIGHT );
00280 
00281         /*------------------------------------------------------------*/
00282 
00283         // Prepare to march vertically through h1 (the first submatrix column
00284         // of HBR).
00285         FLA_Part_2x1( h1,    &hT,
00286                              &hB,            0, FLA_TOP );
00287 
00288         // Prepare to march vertically through F.
00289         FLA_Part_2x1( F1,    &F1T,
00290                              &F1B,           0, FLA_TOP );
00291 
00292         // Iterate vertically through the current column-panel of F, F1.
00293         while ( FLA_Obj_length( F1T ) < FLA_Obj_length( F1 ) )
00294         {
00295 
00296             if ( FLA_Obj_length( F1T ) == 0 )
00297             {
00298                 // Compute the length of the matrix that we will recurse into.
00299                 // Since this is the first iteration of the while loop, this
00300                 // blocksize may be smaller than subsequent full blocks if the
00301                 // i offset passed to this routine is non-zero.
00302                 b_m = m_scal - (i % m_scal);
00303 
00304                 // Reduce the blocksize if the length of F is smaller than the
00305                 // value computed above.
00306                 b_m = min( FLA_Obj_length( F1B ), b_m );
00307 
00308                 // If the i offset is non-zero, we must modulo it to the
00309                 // hierarchical blocksize for the element submatrices of H
00310                 // in preparation for our recursion.
00311                 row_offset = i % m_scal;
00312             }
00313             else
00314             {
00315                 // The length of the matrix into which we will recursve may be
00316                 // assumed to be the same as the length of the element submatrices
00317                 // of H. 
00318                 b_m = m_scal;
00319 
00320                 // Reduce the blocksize if the length of the first panel of
00321                 // F is smaller than the the element submatrices.
00322                 b_m = min( FLA_Obj_length( F1B ), b_m );
00323 
00324                 // The column offset will be zero because we will be starting from
00325                 // the top edge of the submatrix.
00326                 row_offset = 0;
00327             }
00328 
00329             // Expose one hierarchical submatrix within the current column of h1.
00330             FLA_Repart_2x1_to_3x1( hT,                &h01,
00331                                 /* ** */            /* ** */
00332                                                       &h11,
00333                                    hB,                &h21,        1, FLA_BOTTOM );
00334 
00335             // Expose one block of F1 (the current panel of F) of length b_m,
00336             // computed above.
00337             FLA_Repart_2x1_to_3x1( F1T,               &F01,
00338                                 /* ** */            /* ** */
00339                                                       &F11,
00340                                    F1B,               &F21,      b_m, FLA_BOTTOM );
00341 
00342             /*------------------------------------------------------------*/
00343 
00344             // This is a pointer to the FLA_Obj element referenced by h11.
00345             h11p = ( FLA_Obj* ) FLA_Obj_buffer( h11 );
00346 
00347             // Recurse with new (i,j) offset parameters, a narrowed block of
00348             // F, and a corresponding hierarchical block of H.
00349             FLASH_Copy_hierarchy_r( direction, F11, row_offset, col_offset, h11p, queue );
00350 
00351             /*------------------------------------------------------------*/
00352 
00353             FLA_Cont_with_3x1_to_2x1( &hT,                h01,
00354                                                           h11,
00355                                    /* ** */            /* ** */
00356                                       &hB,                h21,     FLA_TOP );
00357 
00358             FLA_Cont_with_3x1_to_2x1( &F1T,               F01,
00359                                                           F11,
00360                                     /* ** */           /* ** */
00361                                       &F1B,               F21,     FLA_TOP );
00362         }
00363 
00364         /*------------------------------------------------------------*/
00365 
00366         FLA_Cont_with_1x3_to_1x2( &FL,  /**/ &FR,        F0, F1, /**/ F2,
00367                                   FLA_LEFT );
00368 
00369         FLA_Cont_with_1x3_to_1x2( &HBR_L,  /**/ &HBR_R,        h0, h1, /**/ h2,
00370                                   FLA_LEFT );
00371     }
00372 
00373     return FLA_SUCCESS;
00374 
00375 }

FLA_Error FLASH_Copy_submatrix_to_global ( dim_t  m,
dim_t  n,
void *  buffer,
dim_t  ldim,
dim_t  i,
dim_t  j,
FLA_Obj  H 
)

References FLA_Check_error_level(), FLA_Check_submatrix_dims_and_offset(), FLA_Obj_attach_buffer(), FLA_Obj_create_without_buffer(), FLA_Obj_free_without_buffer(), FLASH_Copy_hierarchy(), and FLASH_Obj_datatype().

00036 {
00037     FLA_Obj      flat_matrix;
00038     FLA_Datatype datatype;
00039     FLA_Error    e_val;
00040 
00041     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00042     {
00043         e_val = FLA_Check_submatrix_dims_and_offset( m, n, i, j, H );
00044         FLA_Check_error_code( e_val );
00045     }
00046 
00047     // Acquire the datatype from the hierarchical matrix object.
00048     datatype = FLASH_Obj_datatype( H );
00049 
00050     // Create a temporary conventional matrix object of the requested datatype
00051     // and dimensions and attach the given buffer containing the incoming data.
00052     FLA_Obj_create_without_buffer( datatype, m, n, &flat_matrix );
00053     FLA_Obj_attach_buffer( buffer, ldim, &flat_matrix );
00054 
00055     // Recurse through H, adding in the corresponding elements of flat_matrix,
00056     // starting at the (i,j) element offset.
00057     FLASH_Copy_hierarchy( FLA_FLAT_TO_HIER, flat_matrix, i, j, &H );
00058 
00059     // Free the object (but don't free the buffer!).
00060     FLA_Obj_free_without_buffer( &flat_matrix );
00061 
00062     return FLA_SUCCESS;
00063 }

FLA_Error FLASH_Copy_subobject_to_global ( FLA_Obj  F,
dim_t  i,
dim_t  j,
FLA_Obj  H 
)

References FLASH_Copy_hierarchy().

Referenced by FLASH_Obj_create_hier_copy_of_flat(), FLASH_Obj_create_hier_copy_of_flat_ext(), and FLASH_Obj_hierarchify().

00098 {
00099     FLASH_Copy_hierarchy( FLA_FLAT_TO_HIER, F, i, j, &H );
00100 
00101     return FLA_SUCCESS;
00102 }

FLA_Error FLASH_Obj_attach_buffer ( void *  buffer,
dim_t  ldim,
FLA_Obj H 
)

References FLA_Check_error_level(), FLA_Obj_attach_buffer(), FLA_Obj_create_without_buffer(), FLA_Obj_free_without_buffer(), FLASH_Obj_attach_buffer_check(), FLASH_Obj_attach_buffer_hierarchy(), FLASH_Obj_datatype(), FLASH_Obj_scalar_length(), and FLASH_Obj_scalar_width().

00812 {
00813     FLA_Obj      flat_matrix;
00814     dim_t        m_H, n_H;
00815     FLA_Datatype datatype;
00816 
00817     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00818         FLASH_Obj_attach_buffer_check( buffer, ldim, H );
00819 
00820     // Extract the overall scalar dimensions of the hierarchical matrix and
00821     // its base numerical datatype. (These fields will be set even if it has
00822     // a NULL buffer, which presumably it does given that this function was
00823     // just invoked.)
00824     m_H      = FLASH_Obj_scalar_length( *H );
00825     n_H      = FLASH_Obj_scalar_width( *H );
00826     datatype = FLASH_Obj_datatype( *H );
00827 
00828     // Create a temporary conventional object and attach the given buffer.
00829     // Segments of this buffer will be partitioned out to the various
00830     // leaf-level matrices of the hierarchical matrix H.
00831     FLA_Obj_create_without_buffer( datatype, m_H, n_H, &flat_matrix );
00832     FLA_Obj_attach_buffer( buffer, ldim, &flat_matrix );
00833 
00834     // Recurse through the hierarchical matrix, assigning segments of
00835     // flat_matrix to the various leaf-level matrices, similar to what
00836     // we would do if we were creating the object outright.
00837     FLASH_Obj_attach_buffer_hierarchy( flat_matrix, H );
00838 
00839     // Free the object (but don't free the buffer!).
00840     FLA_Obj_free_without_buffer( &flat_matrix );
00841 
00842     return FLA_SUCCESS;
00843 }

FLA_Error FLASH_Obj_attach_buffer_check ( void *  buffer,
dim_t  ldim,
FLA_Obj H 
)

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_attach_buffer().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( H );
00040   FLA_Check_error_code( e_val );
00041 
00042   return FLA_SUCCESS;
00043 }

FLA_Error FLASH_Obj_attach_buffer_hierarchy ( FLA_Obj  F,
FLA_Obj H 
)

References FLA_Check_error_level(), FLA_Cont_with_1x3_to_1x2(), FLA_Cont_with_3x1_to_2x1(), FLA_Obj_attach_buffer(), FLA_Obj_buffer(), FLA_Obj_elemtype(), FLA_Obj_ldim(), FLA_Obj_length(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Part_2x1(), FLA_Repart_1x2_to_1x3(), FLA_Repart_2x1_to_3x1(), FLASH_Obj_attach_buffer_hierarchy(), FLASH_Obj_attach_buffer_hierarchy_check(), FLASH_Obj_scalar_length(), and FLASH_Obj_scalar_width().

Referenced by FLASH_Obj_attach_buffer(), and FLASH_Obj_attach_buffer_hierarchy().

00847 {
00848     FLA_Obj FL,    FR,       F0,  F1,  F2;
00849 
00850     FLA_Obj HL,    HR,       H0,  H1,  H2;
00851 
00852     FLA_Obj F1T,              F01,
00853             F1B,              F11,
00854                               F21;
00855 
00856     FLA_Obj H1T,              H01,
00857             H1B,              H11,
00858                               H21;
00859 
00860     FLA_Obj* H11p;
00861 
00862     dim_t b_m, b_n;
00863 
00864     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00865         FLASH_Obj_attach_buffer_hierarchy_check( F, H );
00866 
00867     if ( FLA_Obj_elemtype( *H ) == FLA_SCALAR )
00868     {
00869         // If we've recursed down to a leaf node, then we can simply attach the
00870         // matrix buffer to the current leaf-level submatrix.
00871         FLA_Obj_attach_buffer( FLA_Obj_buffer( F ), FLA_Obj_ldim( F ), H );
00872     }
00873     else
00874     {
00875         FLA_Part_1x2( *H,    &HL,  &HR,      0, FLA_LEFT );
00876 
00877         FLA_Part_1x2(  F,    &FL,  &FR,      0, FLA_LEFT );
00878 
00879         while ( FLA_Obj_width( HL ) < FLA_Obj_width( *H ) )
00880         {
00881 
00882             FLA_Repart_1x2_to_1x3( HL,  /**/ HR,        &H0, /**/ &H1, &H2,
00883                                    1, FLA_RIGHT );
00884 
00885             b_n = FLASH_Obj_scalar_width( H1 );
00886 
00887             FLA_Repart_1x2_to_1x3( FL,  /**/ FR,        &F0, /**/ &F1, &F2,
00888                                    b_n, FLA_RIGHT );
00889 
00890             /*------------------------------------------------------------*/
00891 
00892             FLA_Part_2x1( H1,    &H1T,
00893                                  &H1B,            0, FLA_TOP );
00894 
00895             FLA_Part_2x1( F1,    &F1T,
00896                                  &F1B,            0, FLA_TOP );
00897 
00898             while ( FLA_Obj_length( H1T ) < FLA_Obj_length( H1 ) )
00899             {
00900 
00901                 FLA_Repart_2x1_to_3x1( H1T,               &H01,
00902                                     /* ** */            /* ** */
00903                                                           &H11,
00904                                        H1B,               &H21,        1, FLA_BOTTOM );
00905 
00906                 b_m = FLASH_Obj_scalar_length( H11 );
00907 
00908                 FLA_Repart_2x1_to_3x1( F1T,               &F01,
00909                                     /* ** */            /* ** */
00910                                                           &F11,
00911                                        F1B,               &F21,      b_m, FLA_BOTTOM );
00912 
00913                 /*------------------------------------------------------------*/
00914 
00915                 H11p = ( FLA_Obj* ) FLA_Obj_buffer( H11 );
00916 
00917                 FLASH_Obj_attach_buffer_hierarchy( F11, H11p );
00918 
00919                 /*------------------------------------------------------------*/
00920 
00921                 FLA_Cont_with_3x1_to_2x1( &H1T,               H01,
00922                                                               H11,
00923                                         /* ** */           /* ** */
00924                                           &H1B,               H21,     FLA_TOP );
00925 
00926                 FLA_Cont_with_3x1_to_2x1( &F1T,               F01,
00927                                                               F11,
00928                                         /* ** */           /* ** */
00929                                           &F1B,               F21,     FLA_TOP );
00930             }
00931 
00932             /*------------------------------------------------------------*/
00933 
00934             FLA_Cont_with_1x3_to_1x2( &HL,  /**/ &HR,        H0, H1, /**/ H2,
00935                                       FLA_LEFT );
00936 
00937             FLA_Cont_with_1x3_to_1x2( &FL,  /**/ &FR,        F0, F1, /**/ F2,
00938                                       FLA_LEFT );
00939 
00940         }
00941     }
00942 
00943     return FLA_SUCCESS;
00944 }

FLA_Error FLASH_Obj_attach_buffer_hierarchy_check ( FLA_Obj  F,
FLA_Obj H 
)

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_attach_buffer_hierarchy().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( H );
00040   FLA_Check_error_code( e_val );
00041 
00042   return FLA_SUCCESS;
00043 }

dim_t FLASH_Obj_blocksizes ( FLA_Obj  H,
dim_t b_m,
dim_t b_n 
)

References FLA_Check_error_level(), FLA_Obj_buffer(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_width(), and FLASH_Obj_blocksizes_check().

Referenced by FLASH_Obj_create_conf_to().

00072 {
00073     FLA_Elemtype elemtype;
00074     FLA_Obj*     buffer_H;
00075     dim_t        depth = 0;
00076 
00077     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00078         FLASH_Obj_blocksizes_check( H, b_m, b_n );
00079 
00080     // Recurse through the hierarchy to the first leaf node. We initialize
00081     // the recursion here:
00082     elemtype = FLA_Obj_elemtype( H );
00083     buffer_H = ( FLA_Obj* ) FLA_Obj_buffer( H );
00084 
00085     while ( elemtype == FLA_MATRIX )
00086     {
00087         b_m[depth] = FLA_Obj_length( buffer_H[0] );
00088         b_n[depth] = FLA_Obj_width( buffer_H[0] );
00089         ++depth;
00090 
00091         // Get the element type of the top-leftmost underlying object. Also,
00092         // get a pointer to the first element of the top-leftmost object and
00093         // assume that it is of type FLA_Obj* in case elemtype is once again
00094         // FLA_MATRIX.
00095         elemtype = FLA_Obj_elemtype( buffer_H[0] );
00096         buffer_H = ( FLA_Obj * ) FLA_Obj_buffer( buffer_H[0] );
00097     }
00098 
00099     // At this point, the first depth elements of blocksizes have been filled
00100     // with the blocksizes of H's various hierarchical levels. Return the
00101     // matrix depth as a confirmation of how many blocksizes were found.
00102     return depth;
00103 }

FLA_Error FLASH_Obj_blocksizes_check ( FLA_Obj  H,
dim_t b_m,
dim_t b_n 
)

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_blocksizes().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( b_m );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( b_n );
00043   FLA_Check_error_code( e_val );
00044 
00045   return FLA_SUCCESS;
00046 }

FLA_Error FLASH_Obj_create ( FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t blocksizes,
FLA_Obj H 
)

References FLASH_Obj_create_helper().

Referenced by FLA_QR_UT_inc_create_U_panel(), and FLASH_Obj_create_hier_conf_to_flat().

00178 {
00179     FLASH_Obj_create_helper( FALSE, datatype, m, n, depth, blocksizes, blocksizes, H );
00180 
00181     return FLA_SUCCESS;
00182 }

FLA_Error FLASH_Obj_create_conf_to ( FLA_Trans  trans,
FLA_Obj  H_cur,
FLA_Obj H_new 
)

References FLA_Check_error_level(), FLA_free(), FLA_malloc(), FLASH_Obj_blocksizes(), FLASH_Obj_create_conf_to_check(), FLASH_Obj_create_ext(), FLASH_Obj_datatype(), FLASH_Obj_depth(), FLASH_Obj_scalar_length(), and FLASH_Obj_scalar_width().

Referenced by FLA_QR_UT_inc_create_U().

00439 {
00440     FLA_Datatype datatype;
00441     dim_t        m_H_cur, n_H_cur;
00442     dim_t        m_H_new, n_H_new;
00443     dim_t        d_H;
00444     dim_t*       b_m;
00445     dim_t*       b_n;
00446 
00447     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00448         FLASH_Obj_create_conf_to_check( trans, H_cur, H_new );
00449 
00450     // Acquire the numerical datatype, length, width, and depth of the
00451     // hierarchical matrix object.
00452     datatype = FLASH_Obj_datatype( H_cur );
00453     m_H_cur  = FLASH_Obj_scalar_length( H_cur );
00454     n_H_cur  = FLASH_Obj_scalar_width( H_cur );
00455     d_H      = FLASH_Obj_depth( H_cur );
00456 
00457     // Allocate a pair of temporary arrays for the blocksizes, whose length is
00458     // equal to the object's hierarchical depth.
00459     b_m = ( dim_t * ) FLA_malloc( d_H * sizeof( dim_t ) );
00460     b_n = ( dim_t * ) FLA_malloc( d_H * sizeof( dim_t ) );
00461 
00462     // Accumulate the blocksizes into a buffer.
00463     FLASH_Obj_blocksizes( H_cur, b_m, b_n );
00464 
00465     // Factor in the transposition, if requested.
00466     if ( trans == FLA_NO_TRANSPOSE )
00467     {
00468         m_H_new = m_H_cur;
00469         n_H_new = n_H_cur;
00470     }
00471     else
00472     {
00473         m_H_new = n_H_cur;
00474         n_H_new = m_H_cur;
00475     }
00476 
00477     // Create the new hierarchical matrix object.
00478     FLASH_Obj_create_ext( datatype, m_H_new, n_H_new, d_H, b_m, b_n, H_new ); 
00479 
00480     // Free the temporary blocksize arrays.
00481     FLA_free( b_m );
00482     FLA_free( b_n );
00483 
00484     return FLA_SUCCESS;
00485 }

FLA_Error FLASH_Obj_create_conf_to_check ( FLA_Trans  trans,
FLA_Obj  H_cur,
FLA_Obj H_new 
)

References FLA_Check_null_pointer(), and FLA_Check_valid_real_trans().

Referenced by FLASH_Obj_create_conf_to().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_valid_real_trans( trans );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( H_new );
00043   FLA_Check_error_code( e_val );
00044 
00045   return FLA_SUCCESS;
00046 }

FLA_Error FLASH_Obj_create_ext ( FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

FLA_Error FLASH_Obj_create_flat_conf_to_hier ( FLA_Trans  trans,
FLA_Obj  H,
FLA_Obj F 
)

References FLA_Check_error_level(), FLA_Obj_create(), FLASH_Obj_create_flat_conf_to_hier_check(), FLASH_Obj_datatype(), FLASH_Obj_scalar_length(), and FLASH_Obj_scalar_width().

Referenced by FLASH_Obj_create_flat_copy_of_hier().

00557 {
00558     FLA_Datatype datatype;
00559     dim_t        m_H, n_H;
00560     dim_t        m_F, n_F;
00561 
00562     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00563         FLASH_Obj_create_flat_conf_to_hier_check( trans, H, F );
00564 
00565     // Acquire the numerical datatype, length, and width of the
00566     // hierarchical matrix object.
00567     datatype = FLASH_Obj_datatype( H );
00568     m_H      = FLASH_Obj_scalar_length( H );
00569     n_H      = FLASH_Obj_scalar_width( H );
00570 
00571     // Factor in the transposition, if requested.
00572     if ( trans == FLA_NO_TRANSPOSE )
00573     {
00574         m_F = m_H;
00575         n_F = n_H;
00576     }
00577     else
00578     {
00579         m_F = n_H;
00580         n_F = m_H;
00581     }
00582 
00583     // Create the flat matrix object.
00584     FLA_Obj_create( datatype, m_F, n_F, F );
00585 
00586     return FLA_SUCCESS;
00587 }

FLA_Error FLASH_Obj_create_flat_conf_to_hier_check ( FLA_Trans  trans,
FLA_Obj  H,
FLA_Obj F 
)

References FLA_Check_null_pointer(), and FLA_Check_valid_real_trans().

Referenced by FLASH_Obj_create_flat_conf_to_hier().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_valid_real_trans( trans );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( F );
00043   FLA_Check_error_code( e_val );
00044 
00045   return FLA_SUCCESS;
00046 }

FLA_Error FLASH_Obj_create_flat_copy_of_hier ( FLA_Obj  H,
FLA_Obj F 
)

References FLA_Check_error_level(), FLASH_Copy_global_to_subobject(), FLASH_Obj_create_flat_conf_to_hier(), and FLASH_Obj_create_flat_copy_of_hier_check().

Referenced by FLASH_Max_elemwise_diff(), FLASH_Norm1(), FLASH_Obj_set_to_scalar(), FLASH_Obj_shift_diagonal(), FLASH_Obj_show(), FLASH_Random_matrix(), FLASH_Random_spd_matrix(), and FLASH_Triangularize().

00623 {
00624     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00625         FLASH_Obj_create_flat_copy_of_hier_check( H, F );
00626 
00627     // Create a flat object conformal to the hierarchical object.
00628     FLASH_Obj_create_flat_conf_to_hier( FLA_NO_TRANSPOSE, H, F );
00629 
00630     // Flatten the hierarchical object's contents into the new flat object.
00631     FLASH_Copy_global_to_subobject( 0, 0, H, *F );
00632 
00633     return FLA_SUCCESS;
00634 }

FLA_Error FLASH_Obj_create_flat_copy_of_hier_check ( FLA_Obj  H,
FLA_Obj F 
)

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_create_flat_copy_of_hier().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( F );
00040   FLA_Check_error_code( e_val );
00041 
00042   return FLA_SUCCESS;
00043 }

FLA_Error FLASH_Obj_create_helper ( FLA_Bool  without_buffer,
FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLA_Check_error_level(), FLA_free(), FLA_malloc(), FLA_Obj_create(), FLA_Obj_create_without_buffer(), FLA_Obj_free_without_buffer(), FLASH_Obj_create_helper_check(), and FLASH_Obj_create_hierarchy().

Referenced by FLASH_Obj_create(), FLASH_Obj_create_ext(), FLASH_Obj_create_without_buffer(), and FLASH_Obj_create_without_buffer_ext().

00210 {
00211     dim_t     i;
00212     FLA_Obj   flat_matrix;
00213 
00214     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00215         FLASH_Obj_create_helper_check( without_buffer, datatype, m, n, depth, b_m, b_n, H );
00216 
00217     if ( depth == 0 )
00218     {
00219         // Base case: create a single contiguous matrix block.
00220         if ( without_buffer == FALSE )
00221             FLA_Obj_create( datatype, m, n, H );
00222         else
00223             FLA_Obj_create_without_buffer( datatype, m, n, H );
00224     }
00225     else
00226     {
00227         // We need temporary arrays the same length as the blocksizes arrays.
00228         dim_t* elem_sizes_m  = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
00229         dim_t* elem_sizes_n  = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
00230         dim_t* depth_sizes_m = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
00231         dim_t* depth_sizes_n = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
00232         dim_t* m_offsets     = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
00233         dim_t* n_offsets     = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
00234         
00235         // Fill two sets of arrays: elem_sizes_m/elem_sizes_n and depth_sizes_m/
00236         // depth_sizes_n.
00237         //  - elem_sizes_m[i] will contain the number of numerical elements that span
00238         //    the row dimension of a block at the ith level of the hierarchy. This is
00239         //    just the product of all row blocksizes "internal" to and including the
00240         //    current blocking level. (The elem_sizes_n array tracks similar values
00241         //    in the column dimension.)
00242         //  - depth_sizes_m[i] is similar to elem_sizes_m[i]. The only difference is
00243         //    that instead of tracking the number of numerical elements in the row
00244         //    dimension, it tracks the number of "storage" blocks that span the m
00245         //    dimension of a block at the ith level, where the m dimension of a
00246         //    storage block is the block size given in b_m[depth-1], ie:
00247         //    the inner-most row dimension block size. (The depth_sizes_n array
00248         //    tracks similar values in the column dimension.)
00249         elem_sizes_m[depth-1]  = b_m[depth-1];
00250         elem_sizes_n[depth-1]  = b_n[depth-1];
00251         depth_sizes_m[depth-1] = 1;
00252         depth_sizes_n[depth-1] = 1;
00253         for ( i = depth - 1; i > 0; --i )
00254         {
00255             elem_sizes_m[i-1]  = elem_sizes_m[i]  * b_m[i-1];
00256             elem_sizes_n[i-1]  = elem_sizes_n[i]  * b_n[i-1];
00257             depth_sizes_m[i-1] = depth_sizes_m[i] * b_m[i-1];
00258             depth_sizes_n[i-1] = depth_sizes_n[i] * b_n[i-1];
00259         }
00260     
00261         // Initialize the m_offsets and n_offsets arrays to zero.
00262         for ( i = 0; i < depth; i++ )
00263         {
00264             m_offsets[i] = 0;
00265             n_offsets[i] = 0;
00266         }
00267 
00268         // Create a "flat" matrix object. All leaf-level child objects will refer
00269         // to various offsets within this object's buffer.
00270         // Notice that the normal invocation requests a 1 x m*n flat matrix,
00271         // whereas the without_buffer line requests an m x n matrix. This is
00272         // simply how Tze originally implemented the corresponding cases.
00273         // There may be a subtle difference between these two methods--one
00274         // which I can't see at the time of this writing.
00275         if ( without_buffer == FALSE )
00276             FLA_Obj_create( datatype, 1, m*n, &flat_matrix );
00277         else
00278             FLA_Obj_create_without_buffer( datatype, m, n, &flat_matrix );
00279         
00280         // Recursively create the matrix hierarchy.
00281         FLASH_Obj_create_hierarchy( datatype, m, n, depth, elem_sizes_m, elem_sizes_n, flat_matrix, H, 0, depth, depth_sizes_m, depth_sizes_n, m_offsets, n_offsets );
00282         
00283         // Free the flat_matrix object, but not its buffer. If we created a
00284         // normal object with a buffer, we don't want to free the buffer because
00285         // it is being used by the hierarchical objected we just created. If we
00286         // created a bufferless object, we don't want to free the buffer because
00287         // there was no buffer allocated in the first place.
00288         FLA_Obj_free_without_buffer( &flat_matrix );
00289         
00290         // Free the local arrays.
00291         FLA_free( elem_sizes_m );
00292         FLA_free( elem_sizes_n );
00293         FLA_free( depth_sizes_m );
00294         FLA_free( depth_sizes_n );
00295         FLA_free( m_offsets );
00296         FLA_free( n_offsets );
00297     }
00298 
00299     return FLA_SUCCESS;
00300 }

FLA_Error FLASH_Obj_create_helper_check ( FLA_Bool  without_buffer,
FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLA_Check_null_pointer(), and FLA_Check_valid_datatype().

Referenced by FLASH_Obj_create_helper().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_valid_datatype( datatype );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( b_m );
00043   FLA_Check_error_code( e_val );
00044 
00045   e_val = FLA_Check_null_pointer( b_n );
00046   FLA_Check_error_code( e_val );
00047 
00048   e_val = FLA_Check_null_pointer( H );
00049   FLA_Check_error_code( e_val );
00050 
00051   // A value of depth < 0 should cause an error.
00052 
00053   // Values of m < 1, n < 1 should cause an error. (or < 0?)
00054 
00055   // First depth entries in blocksize_m, _n should be checked; values < 1 should cause error.
00056 
00057   return FLA_SUCCESS;
00058 }

FLA_Error FLASH_Obj_create_hier_conf_to_flat ( FLA_Trans  trans,
FLA_Obj  F,
dim_t  depth,
dim_t blocksizes,
FLA_Obj H 
)

References FLA_Check_error_level(), FLA_Obj_datatype(), FLA_Obj_length(), FLA_Obj_width(), FLASH_Obj_create(), and FLASH_Obj_create_hier_conf_to_flat_check().

Referenced by FLASH_Obj_create_hier_copy_of_flat().

00489 {
00490     FLA_Datatype datatype;
00491     dim_t        m_H, n_H;
00492     dim_t        m_F, n_F;
00493 
00494     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00495         FLASH_Obj_create_hier_conf_to_flat_check( trans, F, depth, blocksizes, H );
00496 
00497     // Acquire the numerical datatype, length, and width of the flat matrix
00498     // object.
00499     datatype = FLA_Obj_datatype( F );
00500     m_F      = FLA_Obj_length( F );
00501     n_F      = FLA_Obj_width( F );
00502 
00503     // Factor in the transposition, if requested.
00504     if ( trans == FLA_NO_TRANSPOSE )
00505     {
00506         m_H = m_F;
00507         n_H = n_F;
00508     }
00509     else
00510     {
00511         m_H = n_F;
00512         n_H = m_F;
00513     }
00514 
00515     // Create the hierarchical matrix object.
00516     FLASH_Obj_create( datatype, m_H, n_H, depth, blocksizes, H );
00517 
00518     return FLA_SUCCESS;
00519 }

FLA_Error FLASH_Obj_create_hier_conf_to_flat_check ( FLA_Trans  trans,
FLA_Obj  F,
dim_t  depth,
dim_t blocksizes,
FLA_Obj H 
)

References FLA_Check_null_pointer(), and FLA_Check_valid_real_trans().

Referenced by FLASH_Obj_create_hier_conf_to_flat().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_valid_real_trans( trans );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( blocksizes );
00043   FLA_Check_error_code( e_val );
00044 
00045   e_val = FLA_Check_null_pointer( H );
00046   FLA_Check_error_code( e_val );
00047 
00048   // A value of depth < 0 should cause an error.
00049 
00050   // First depth entries in blocksize should be checked; values < 1 should cause error.
00051 
00052   return FLA_SUCCESS;
00053 }

FLA_Error FLASH_Obj_create_hier_conf_to_flat_ext ( FLA_Trans  trans,
FLA_Obj  F,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLA_Check_error_level(), FLA_Obj_datatype(), FLA_Obj_length(), FLA_Obj_width(), FLASH_Obj_create_ext(), and FLASH_Obj_create_hier_conf_to_flat_ext_check().

Referenced by FLASH_Obj_create_hier_copy_of_flat_ext().

00523 {
00524     FLA_Datatype datatype;
00525     dim_t        m_H, n_H;
00526     dim_t        m_F, n_F;
00527 
00528     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00529         FLASH_Obj_create_hier_conf_to_flat_ext_check( trans, F, depth, b_m, b_n, H );
00530 
00531     // Acquire the numerical datatype, length, and width of the flat matrix
00532     // object.
00533     datatype = FLA_Obj_datatype( F );
00534     m_F      = FLA_Obj_length( F );
00535     n_F      = FLA_Obj_width( F );
00536 
00537     // Factor in the transposition, if requested.
00538     if ( trans == FLA_NO_TRANSPOSE )
00539     {
00540         m_H = m_F;
00541         n_H = n_F;
00542     }
00543     else
00544     {
00545         m_H = n_F;
00546         n_H = m_F;
00547     }
00548 
00549     // Create the hierarchical matrix object.
00550     FLASH_Obj_create_ext( datatype, m_H, n_H, depth, b_m, b_n, H );
00551 
00552     return FLA_SUCCESS;
00553 }

FLA_Error FLASH_Obj_create_hier_conf_to_flat_ext_check ( FLA_Trans  trans,
FLA_Obj  F,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLA_Check_null_pointer(), and FLA_Check_valid_real_trans().

Referenced by FLASH_Obj_create_hier_conf_to_flat_ext().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_valid_real_trans( trans );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( b_m );
00043   FLA_Check_error_code( e_val );
00044 
00045   e_val = FLA_Check_null_pointer( b_n );
00046   FLA_Check_error_code( e_val );
00047 
00048   e_val = FLA_Check_null_pointer( H );
00049   FLA_Check_error_code( e_val );
00050 
00051   // A value of depth < 0 should cause an error.
00052 
00053   // First depth entries in blocksize should be checked; values < 1 should cause error.
00054 
00055   return FLA_SUCCESS;
00056 }

FLA_Error FLASH_Obj_create_hier_copy_of_flat ( FLA_Obj  F,
dim_t  depth,
dim_t blocksizes,
FLA_Obj H 
)

References FLA_Check_error_level(), FLASH_Copy_subobject_to_global(), FLASH_Obj_create_hier_conf_to_flat(), and FLASH_Obj_create_hier_copy_of_flat_check().

Referenced by FLASH_LU_incpiv_create_hier_matrices(), and FLASH_QR_UT_inc_create_hier_matrices().

00591 {
00592     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00593         FLASH_Obj_create_hier_copy_of_flat_check( F, depth, blocksizes, H );
00594 
00595     // Create a hierarchical object conformal to the flat object.
00596     FLASH_Obj_create_hier_conf_to_flat( FLA_NO_TRANSPOSE, F, depth, blocksizes, H );
00597 
00598     // Initialize the contents of the hierarchical matrix object with the
00599     // contents of the flat matrix object.
00600     FLASH_Copy_subobject_to_global( F, 0, 0, *H );
00601     
00602     return FLA_SUCCESS;
00603 }

FLA_Error FLASH_Obj_create_hier_copy_of_flat_check ( FLA_Obj  F,
dim_t  depth,
dim_t blocksizes,
FLA_Obj H 
)

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_create_hier_copy_of_flat().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( blocksizes );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( H );
00043   FLA_Check_error_code( e_val );
00044 
00045   // A value of depth < 0 should cause an error.
00046 
00047   // First depth entries in blocksize should be checked; values < 1 should cause error.
00048 
00049   return FLA_SUCCESS;
00050 }

FLA_Error FLASH_Obj_create_hier_copy_of_flat_ext ( FLA_Obj  F,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLA_Check_error_level(), FLASH_Copy_subobject_to_global(), FLASH_Obj_create_hier_conf_to_flat_ext(), and FLASH_Obj_create_hier_copy_of_flat_ext_check().

00607 {
00608     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00609         FLASH_Obj_create_hier_copy_of_flat_ext_check( F, depth, b_m, b_n, H );
00610 
00611     // Create a hierarchical object conformal to the flat object.
00612     FLASH_Obj_create_hier_conf_to_flat_ext( FLA_NO_TRANSPOSE, F, depth, b_m, b_n, H );
00613 
00614     // Initialize the contents of the hierarchical matrix object with the
00615     // contents of the flat matrix object.
00616     FLASH_Copy_subobject_to_global( F, 0, 0, *H );
00617     
00618     return FLA_SUCCESS;
00619 }

FLA_Error FLASH_Obj_create_hier_copy_of_flat_ext_check ( FLA_Obj  F,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_create_hier_copy_of_flat_ext().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( b_m );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( b_n );
00043   FLA_Check_error_code( e_val );
00044 
00045   e_val = FLA_Check_null_pointer( H );
00046   FLA_Check_error_code( e_val );
00047 
00048   // A value of depth < 0 should cause an error.
00049 
00050   // First depth entries in blocksize should be checked; values < 1 should cause error.
00051 
00052   return FLA_SUCCESS;
00053 }

FLA_Error FLASH_Obj_create_hierarchy ( FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t elem_sizes_m,
dim_t elem_sizes_n,
FLA_Obj  flat_matrix,
FLA_Obj H,
unsigned long  id,
dim_t  depth_overall,
dim_t depth_sizes_m,
dim_t depth_sizes_n,
dim_t m_offsets,
dim_t n_offsets 
)

References FLA_Obj_view::base, FLA_Check_error_level(), FLA_Cont_with_1x3_to_1x2(), FLA_Obj_attach_buffer(), FLA_Obj_buffer(), FLA_Obj_create_ext(), FLA_Obj_create_without_buffer(), FLA_Obj_datatype_size(), FLA_Obj_free_without_buffer(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Repart_1x2_to_1x3(), FLASH_Obj_create_hierarchy(), FLASH_Obj_create_hierarchy_check(), FLASH_Queue_set_block_size(), FLA_Obj_struct::id, FLA_Obj_struct::m_index, and FLA_Obj_struct::n_index.

Referenced by FLASH_Obj_create_helper(), and FLASH_Obj_create_hierarchy().

00304 {
00305     dim_t    i, j, b;
00306     dim_t    next_m, next_n;
00307     dim_t    num_m, num_n;
00308     dim_t    m_inner, n_inner;
00309     dim_t    elem_size_m_cur;
00310     dim_t    elem_size_n_cur;
00311     FLA_Obj  FL, FR, F0, F1, F2;
00312     FLA_Obj* buffer_H;
00313 
00314     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00315         FLASH_Obj_create_hierarchy_check( datatype, m, n, depth, elem_sizes_m, elem_sizes_n, flat_matrix, H, id, depth_overall, depth_sizes_m, depth_sizes_n, m_offsets, n_offsets );
00316 
00317     if ( depth == 0 )
00318     {
00319         // If we're asked to create a zero-depth matrix, we interpret that as
00320         // a request to create leaf-level objects using the remaining portion
00321         // of the segment of the flat_matrix buffer that was passed in.
00322         FLA_Obj_create_without_buffer( datatype, m, n, H );
00323         FLA_Obj_attach_buffer( FLA_Obj_buffer( flat_matrix ), m, H );
00324 #ifdef FLA_ENABLE_SUPERMATRIX
00325         FLASH_Queue_set_block_size( m * n * FLA_Obj_datatype_size( datatype ) );
00326 #endif
00327         H->base->id = id;
00328 
00329         // Fill in the m_index and n_index variables, which identify the
00330         // location of the current leaf node, in units of storage blocks,
00331         // within the overall matrix.
00332         for ( i = 0; i < depth_overall; i++ )
00333         {
00334             H->base->m_index += m_offsets[i] * depth_sizes_m[i];
00335             H->base->n_index += n_offsets[i] * depth_sizes_n[i];
00336         }
00337     }
00338     else
00339     {
00340         // The "current" level's elem_size value. That is, the number of numerical
00341         // scalar elements along one side of a full block on the current level,
00342         // for the row and column dimensions.
00343         elem_size_m_cur = elem_sizes_m[0];
00344         elem_size_n_cur = elem_sizes_n[0];
00345 
00346         // Compute the number of rows and columns in the current hierarchical
00347         // level of blocking.
00348         num_m = m / elem_size_m_cur + ( (m % elem_size_m_cur) ? 1 : 0 );
00349         num_n = n / elem_size_n_cur + ( (n % elem_size_n_cur) ? 1 : 0 );
00350 
00351         // The total number of scalar elements contained within/below this level
00352         // of the hierarchy. (The edge cases are handled by the computation of
00353         // next_m and next_n below, since they are passed in as the new m and n
00354         // for the next recursive call.)
00355         m_inner = m;
00356         n_inner = n;
00357         
00358         // Create a matrix whose elements are FLA_Objs for the current level of
00359         // blocking.
00360         FLA_Obj_create_ext( datatype, FLA_MATRIX, num_m, num_n, m_inner, n_inner, H );
00361 
00362         if ( depth == depth_overall )
00363             id = H->base->id;
00364         else
00365             H->base->id = id;
00366 
00367         // Grab the buffer from the new hierarchical object. This is an array of
00368         // FLA_Objs.
00369         buffer_H = ( FLA_Obj* ) FLA_Obj_buffer( *H );
00370         
00371         // Prepare to partition through the flat matrix so we can further allocate
00372         // segments of it to the various hierarchical sub-matrices. (The second
00373         // case occurs when the current function is called with a flat_matrix
00374         // argument that was created without a buffer.)
00375         if ( FLA_Obj_buffer( flat_matrix ) != NULL )
00376             FLA_Part_1x2( flat_matrix, &FL, &FR, 0, FLA_LEFT );
00377         else
00378             FLA_Obj_create_without_buffer( datatype, 0, 0, &F1 );
00379         
00380         for ( j = 0; j < num_n; ++j )
00381         {
00382             // Determine the number of elements along the column dimension
00383             // that will be contained within the submatrix referenced by
00384             // the (i,j)th FLA_MATRIX element in the current matrix.
00385             if ( j != num_n-1 || (n % elem_size_n_cur) == 0 )
00386                 next_n = elem_size_n_cur;
00387             else
00388                 next_n = n % elem_size_n_cur;
00389             
00390             n_offsets[depth_overall-depth] = j;
00391                         
00392             for ( i = 0; i < num_m; ++i )
00393             {           
00394                 // Determine the number of elements along the row dimension
00395                 // that will be contained within the submatrix referenced by
00396                 // the (i,j)th FLA_MATRIX element in the current matrix.
00397                 if ( i != num_m-1 || (m % elem_size_m_cur) == 0 )
00398                     next_m = elem_size_m_cur;
00399                 else
00400                     next_m = m % elem_size_m_cur;
00401 
00402                 m_offsets[depth_overall-depth] = i;
00403                 
00404                 // Partition the next m*n elements from the flat matrix so we can
00405                 // "attach" them to the hierarchical matrices contained within the
00406                 // (i,j)th FLA_MATRIX object.
00407                 if ( FLA_Obj_buffer( flat_matrix ) != NULL )
00408                 {
00409                     b = min( FLA_Obj_width( FR ), next_m * next_n );
00410                     FLA_Repart_1x2_to_1x3( FL,  /**/ FR,        &F0, /**/ &F1, &F2,
00411                                            b, FLA_RIGHT );
00412                 }
00413                 
00414                 // Recursively call ourselves, with the appropriate parameters for
00415                 // the next deeper level in the matrix hierarchy.
00416                 FLASH_Obj_create_hierarchy( datatype, next_m, next_n, depth-1, &elem_sizes_m[1], &elem_sizes_n[1], F1, &buffer_H[j*num_m + i], id, depth_overall, depth_sizes_m, depth_sizes_n, m_offsets, n_offsets );
00417 
00418                 // Continue with the repartitioning.
00419                 if ( FLA_Obj_buffer( flat_matrix ) != NULL )
00420                 {
00421                     FLA_Cont_with_1x3_to_1x2( &FL,  /**/ &FR,        F0, F1, /**/ F2,
00422                                               FLA_LEFT );
00423                 }
00424             }
00425         }
00426 
00427         // Free the temporary flat matrix subpartition object, but only if it was
00428         // created to begin with. Since it would have been created without a
00429         // buffer, we must free it in a similar manner.
00430         if ( FLA_Obj_buffer( flat_matrix ) == NULL )
00431             FLA_Obj_free_without_buffer( &F1 );
00432     }
00433     
00434     return FLA_SUCCESS;
00435 }

FLA_Error FLASH_Obj_create_hierarchy_check ( FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t elem_sizes_m,
dim_t elem_sizes_n,
FLA_Obj  flat_matrix,
FLA_Obj H,
unsigned long  id,
dim_t  depth_overall,
dim_t depth_sizes_m,
dim_t depth_sizes_n,
dim_t m_offsets,
dim_t n_offsets 
)

References FLA_Check_null_pointer(), and FLA_Check_valid_datatype().

Referenced by FLASH_Obj_create_hierarchy().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_valid_datatype( datatype );
00040   FLA_Check_error_code( e_val );
00041 
00042   e_val = FLA_Check_null_pointer( elem_sizes_m );
00043   FLA_Check_error_code( e_val );
00044 
00045   e_val = FLA_Check_null_pointer( elem_sizes_n );
00046   FLA_Check_error_code( e_val );
00047 
00048   e_val = FLA_Check_null_pointer( H );
00049   FLA_Check_error_code( e_val );
00050 
00051   e_val = FLA_Check_null_pointer( depth_sizes_m );
00052   FLA_Check_error_code( e_val );
00053 
00054   e_val = FLA_Check_null_pointer( depth_sizes_n );
00055   FLA_Check_error_code( e_val );
00056 
00057   e_val = FLA_Check_null_pointer( m_offsets );
00058   FLA_Check_error_code( e_val );
00059 
00060   e_val = FLA_Check_null_pointer( n_offsets );
00061   FLA_Check_error_code( e_val );
00062 
00063   // A value of depth < 0 should cause an error.
00064 
00065   // Values of m < 1, n < 1 should cause an error. (or < 0?)
00066 
00067   // First depth entries in depth_sizes_m,_n elem_sizes_m,_n m_,n_offsets should be checked; values < 1 should cause error.
00068 
00069   return FLA_SUCCESS;
00070 }

FLA_Error FLASH_Obj_create_without_buffer ( FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t blocksizes,
FLA_Obj H 
)

References FLASH_Obj_create_helper().

00194 {
00195     FLASH_Obj_create_helper( TRUE, datatype, m, n, depth, blocksizes, blocksizes, H );
00196 
00197     return FLA_SUCCESS;
00198 }

FLA_Error FLASH_Obj_create_without_buffer_ext ( FLA_Datatype  datatype,
dim_t  m,
dim_t  n,
dim_t  depth,
dim_t b_m,
dim_t b_n,
FLA_Obj H 
)

References FLASH_Obj_create_helper().

00202 {
00203     FLASH_Obj_create_helper( TRUE, datatype, m, n, depth, b_m, b_n, H );
00204 
00205     return FLA_SUCCESS;
00206 }

FLA_Datatype FLASH_Obj_datatype ( FLA_Obj  H  ) 

dim_t FLASH_Obj_depth ( FLA_Obj  H  ) 

References FLA_Obj_buffer(), and FLA_Obj_elemtype().

Referenced by FLASH_Apply_Q_UT_inc_create_workspace(), FLASH_FS_incpiv(), FLASH_LU_incpiv(), and FLASH_Obj_create_conf_to().

00043 {
00044     FLA_Elemtype elemtype;
00045     FLA_Obj*     buffer_H;
00046     dim_t        depth = 0;
00047 
00048     // Recurse through the hierarchy to the first leaf node. We initialize
00049     // the recursion here:
00050     elemtype = FLA_Obj_elemtype( H );
00051     buffer_H  = ( FLA_Obj* ) FLA_Obj_buffer( H );
00052 
00053     while ( elemtype == FLA_MATRIX )
00054     {
00055         ++depth;
00056 
00057         // Get the element type of the top-leftmost underlying object. Also,
00058         // get a pointer to the first element of the top-leftmost object and
00059         // assume that it is of type FLA_Obj* in case elemtype is once again
00060         // FLA_MATRIX.
00061         elemtype = FLA_Obj_elemtype( buffer_H[0] );
00062         buffer_H = ( FLA_Obj * ) FLA_Obj_buffer( buffer_H[0] );
00063     }
00064 
00065     // At this point, the value of depth represents the depth of the matrix
00066     // hierarchy.
00067     return depth;
00068 }

void FLASH_Obj_exec ( void *  func,
FLASH_Obj_queue queue 
)

References FLASH_Thread_s::args, FLA_Check_error_level(), FLA_Check_pthread_create_result(), FLA_Check_pthread_join_result(), FLA_Lock_destroy(), FLA_Lock_init(), FLASH_Obj_exec_parallel(), FLASH_Queue_get_num_threads(), FLASH_Obj_variables::func, FLASH_Thread_s::id, FLASH_Obj_variables::lock, FLASH_Obj_queue_s::n_tasks, and FLASH_Obj_variables::queue.

Referenced by FLASH_Axpy_hierarchy(), and FLASH_Copy_hierarchy().

00054 {
00055    int i;
00056    int n_threads = FLASH_Queue_get_num_threads();
00057    void* (*thread_entry_point)( void* );
00058 #ifdef FLA_ENABLE_WINDOWS_BUILD
00059    FLASH_Thread* thread;
00060 #endif
00061 
00062    // The structre to save all the execution variables.
00063    FLASH_Obj_vars args;
00064    
00065    // Return if there are no tasks.
00066    if ( queue->n_tasks == 0 )
00067       return;
00068 
00069    // Determine which function to send threads to.
00070    thread_entry_point = FLASH_Obj_exec_parallel;
00071 
00072    // Save the queue and function pointer.
00073    args.func  = func;
00074    args.queue = queue;
00075 
00076    // Initialize the lock.
00077    FLA_Lock_init( &(args.lock) );
00078 
00079    // Allocate the thread structures array. Here, an array of FLASH_Thread
00080    // structures of length n_threads is allocated and the fields of each
00081    // structure set to appropriate values.
00082 #ifdef FLA_ENABLE_WINDOWS_BUILD
00083    thread = ( FLASH_Thread* ) _alloca( n_threads * sizeof( FLASH_Thread ) );
00084 #else
00085    FLASH_Thread thread[n_threads];
00086 #endif
00087 
00088    // Initialize the thread structures array.
00089    for ( i = 0; i < n_threads; i++ )
00090    {
00091       // Save the thread's identifier.
00092       thread[i].id = i;
00093 
00094       // Save the pointer to the necessary variables with the thread.
00095       thread[i].args = ( void* ) &args;
00096 
00097       // The pthread object, if it was even compiled into the FLASH_Thread
00098       // structure, will be initialized by the pthread implementation when we
00099       // call pthread_create() and does not need to be touched at this time.
00100    }
00101 
00102 #if FLA_MULTITHREADING_MODEL == FLA_OPENMP
00103 
00104    // An OpenMP parallel for region spawns n_threads threads. Each thread
00105    // executes the work function with a different FLASH_Thread argument.
00106    // An implicit synchronization point exists at the end of the curly
00107    // brace scope.
00108    #pragma omp parallel for \
00109            private( i ) \
00110            shared( thread, n_threads, thread_entry_point ) \
00111            schedule( static, 1 ) \
00112            num_threads( n_threads )
00113    for ( i = 0; i < n_threads; ++i )
00114    {
00115       thread_entry_point( ( void* ) &thread[i] );
00116    }
00117 
00118 #elif FLA_MULTITHREADING_MODEL == FLA_PTHREADS
00119 
00120    // Create each POSIX thread needed in addition to the main thread.
00121    for ( i = 1; i < n_threads; i++ )
00122    {
00123       int pthread_e_val;
00124 
00125       // Create thread i with default attributes.
00126       pthread_e_val = pthread_create( &(thread[i].pthread_obj),
00127                                       NULL,
00128                                       thread_entry_point,
00129                                       ( void* ) &thread[i] );
00130       if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00131       {
00132          FLA_Error e_val = FLA_Check_pthread_create_result( pthread_e_val );
00133          FLA_Check_error_code( e_val );
00134       }
00135    }
00136 
00137    // The main thread is assigned the role of thread 0. Here we manually
00138    // execute it as a worker thread.
00139    thread_entry_point( ( void* ) &thread[0] );
00140 
00141    // Wait for non-main threads to finish.
00142    for ( i = 1; i < n_threads; i++ )
00143    {
00144       // These two variables are declared local to this for loop since this
00145       // is the only place they are needed, and since they would show up as
00146       // unused variables if FLA_MULTITHREADING_MODEL == FLA_PTHREADS.
00147       // Strangely, the Intel compiler produces code that results in an
00148       // "unaligned access" runtime message if thread_status is declared as
00149       // an int. Declaring it as a long or void* appears to force the
00150       // compiler (not surprisingly) into aligning it to an 8-byte boundary.
00151       int   pthread_e_val;
00152       void* thread_status;
00153 
00154       // Wait for thread i to invoke its respective pthread_exit().
00155       // The return value passed to pthread_exit() is provided to us
00156       // via status, if one was given.
00157       pthread_e_val = pthread_join( thread[i].pthread_obj,
00158                                     ( void** ) &thread_status );
00159 
00160       if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00161       {
00162          FLA_Error e_val = FLA_Check_pthread_join_result( pthread_e_val );
00163          FLA_Check_error_code( e_val );
00164       }
00165    }
00166 
00167 #endif
00168 
00169    // Destroy the lock.
00170    FLA_Lock_destroy( &(args.lock) );
00171 
00172    return;
00173 }

void* FLASH_Obj_exec_parallel ( void *  arg  ) 

References FLASH_Obj_task_s::alpha, FLASH_Thread_s::args, FLASH_Obj_task_s::direction, FLASH_Obj_task_s::F, FLA_Axpy_external(), FLA_Copy_external(), FLA_free(), FLA_Lock_acquire(), FLA_Lock_release(), FLASH_Axpy_hierarchy(), FLASH_Copy_hierarchy(), FLASH_Obj_variables::func, FLASH_Obj_task_s::H, FLASH_Obj_queue_s::head, FLASH_Thread_s::id, FLASH_Obj_variables::lock, FLASH_Obj_queue_s::n_tasks, FLASH_Obj_task_s::next_task, FLASH_Obj_task_s::prev_task, FLASH_Obj_variables::queue, and FLASH_Obj_queue_s::tail.

Referenced by FLASH_Obj_exec().

00177 {
00178    FLASH_Thread*   me;
00179    FLASH_Obj_vars* args;
00180    FLASH_Obj_task* t;
00181    FLA_Bool        condition = TRUE;
00182 
00183    // Interpret the thread argument as what it really is--a pointer to an
00184    // FLASH_Thread structure.
00185    me = ( FLASH_Thread* ) arg;
00186 
00187    // Extract the variables from the current thread.
00188    args = ( FLASH_Obj_vars* ) me->args;
00189    
00190    // Loop until all the tasks have committed.
00191    while ( condition )
00192    {
00193       t = NULL;
00194 
00195       FLA_Lock_acquire( &(args->lock) ); // L ***
00196 
00197       // If there are no more tasks, stop execution.
00198       if ( args->queue->n_tasks == 0 )
00199       {
00200          condition = FALSE;
00201       }
00202       else
00203       {
00204          // Dequeue the first task.
00205          t = args->queue->head;
00206 
00207          if ( args->queue->n_tasks == 1 )
00208          {
00209             // Clear the queue of its only task.
00210             args->queue->head = NULL;
00211             args->queue->tail = NULL;
00212          }
00213          else
00214          {
00215             // Adjust pointers in queue.
00216             args->queue->head = t->next_task;
00217             args->queue->head->prev_task = NULL;
00218          }
00219 
00220          // Decrement number of tasks on queue.
00221          args->queue->n_tasks--;
00222       }
00223 
00224       FLA_Lock_release( &(args->lock) ); // L ***
00225 
00226       // Execute the task.
00227       if ( t != NULL )
00228       {
00229          // FLASH_Axpy
00230          if      ( args->func == (void *) FLASH_Axpy_hierarchy )
00231          {
00232             if      ( t->direction == FLA_FLAT_TO_HIER )
00233             {
00234                FLA_Axpy_external( t->alpha, t->F, t->H );
00235             }
00236             else if ( t->direction == FLA_HIER_TO_FLAT )
00237             {
00238                FLA_Axpy_external( t->alpha, t->H, t->F );
00239             }
00240          }
00241          // FLASH_Copy
00242          else if ( args->func == (void *) FLASH_Copy_hierarchy )
00243          {
00244             if      ( t->direction == FLA_FLAT_TO_HIER )
00245             {
00246                FLA_Copy_external( t->F, t->H );
00247             }
00248             else if ( t->direction == FLA_HIER_TO_FLAT )
00249             {
00250                FLA_Copy_external( t->H, t->F );
00251             }
00252          }
00253          else
00254          {
00255             FLA_Check_error_code( FLA_NOT_YET_IMPLEMENTED );
00256          }
00257 
00258          // Free the task structure in parallel.
00259          FLA_free( t );
00260       }
00261    }
00262 
00263 #if FLA_MULTITHREADING_MODEL == FLA_PTHREADS
00264    // If this is a non-main thread, then exit with a zero (normal) error code.
00265    // The main thread cannot call pthread_exit() because this routine never
00266    // returns. The main thread must proceed so it can oversee the joining of
00267    // the exited non-main pthreads.
00268    if ( me->id != 0 )
00269       pthread_exit( ( void* ) NULL );
00270 #endif
00271 
00272    return ( void* ) NULL;
00273 }

void* FLASH_Obj_extract_buffer ( FLA_Obj  H  ) 

References FLA_Obj_buffer(), and FLA_Obj_elemtype().

Referenced by FLASH_Obj_free().

00745 {
00746     FLA_Elemtype elemtype;
00747     FLA_Obj*     buffer_H;
00748 
00749     // Recurse through the hierarchy to the first leaf node to gain access
00750     // to the address of the actual numerical data buffer (ie: the "flat"
00751     // matrix object used in FLASH_Obj_create()). We initialize the search
00752     // here:
00753     elemtype = FLA_Obj_elemtype( H );
00754     buffer_H  = ( FLA_Obj* ) FLA_Obj_buffer( H );
00755 
00756     while ( elemtype == FLA_MATRIX )
00757     {
00758         elemtype = FLA_Obj_elemtype( buffer_H[0] );
00759         buffer_H  = ( FLA_Obj * ) FLA_Obj_buffer( buffer_H[0] );
00760     }
00761 
00762     // At this point, the value in buffer_H is a pointer to the array that
00763     // holds the numerical data.
00764     return ( void* ) buffer_H;
00765 }

FLA_Error FLASH_Obj_flatten ( FLA_Obj  H,
FLA_Obj  F 
)

References FLASH_Copy_global_to_subobject().

00769 {
00770     FLASH_Copy_global_to_subobject( 0, 0, H, F );
00771 
00772     return FLA_SUCCESS;
00773 }

void FLASH_Obj_free ( FLA_Obj H  ) 

References FLA_Check_error_level(), FLA_free(), FLA_Obj_elemtype(), FLA_Obj_free(), FLASH_Obj_extract_buffer(), FLASH_Obj_free_check(), and FLASH_Obj_free_hierarchy().

Referenced by FLA_QR_UT_inc_free_U_panel(), and FLASH_Random_spd_matrix().

00638 {
00639     FLA_Obj* buffer_H;
00640 
00641     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00642         FLASH_Obj_free_check( H );
00643 
00644     // Free the object accordin got whether it contains a hierarchy.
00645     if ( FLA_Obj_elemtype( *H ) == FLA_MATRIX )
00646     {
00647         // Extract a pointer to the data buffer that was parititioned across all
00648         // leaf-level submatrices.
00649         buffer_H = ( FLA_Obj * ) FLASH_Obj_extract_buffer( *H );
00650     
00651         // Free the data buffer. This works because FLASH_Obj_extract_buffer()
00652         // returns the starting address of the first element's buffer, which is
00653         // also the starting address of the entire buffer.
00654         FLA_free( buffer_H );
00655     
00656         // All that remains now is to free the interior of the matrix hierarchy.
00657         // This includes non-leaf buffers and their corresponding base objects
00658         // as well as leaf-level base objects.
00659         FLASH_Obj_free_hierarchy( H );
00660 
00661         // Finally, free the top-level FLASH object itself.
00662         FLA_Obj_free( H );
00663     }
00664     else
00665     {
00666         // If the matrix has no hierarchy, treat it like a flat object.
00667         FLA_Obj_free( H );
00668     }
00669 }

FLA_Error FLASH_Obj_free_check ( FLA_Obj H  ) 

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_free().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( H );
00040   FLA_Check_error_code( e_val );
00041 
00042   return FLA_SUCCESS;
00043 }

void FLASH_Obj_free_hierarchy ( FLA_Obj H  ) 

References FLA_Check_error_level(), FLA_Obj_buffer(), FLA_Obj_elemtype(), FLA_Obj_free_without_buffer(), FLA_Obj_length(), FLA_Obj_width(), FLASH_Obj_free_hierarchy(), and FLASH_Obj_free_hierarchy_check().

Referenced by FLASH_Obj_free(), FLASH_Obj_free_hierarchy(), and FLASH_Obj_free_without_buffer().

00701 {
00702     dim_t    i,j, m_H, n_H;
00703     FLA_Obj* buffer_H;
00704     FLA_Obj* obj_ij_ptr;
00705 
00706     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00707         FLASH_Obj_free_hierarchy_check( H );
00708 
00709     m_H = FLA_Obj_length( *H );
00710     n_H = FLA_Obj_width( *H );
00711 
00712     // If the element type of H is FLA_SCALAR, then it has no children to free,
00713     // so we can return. (The base object of H will be freed by the caller,
00714     // whether that is a recursive instance of the current function,
00715     // FLASH_Obj_free_hierarchy(), or FLASH_Obj_free() or
00716     // FLASH_Obj_free_without_buffer()).
00717     if ( FLA_Obj_elemtype( *H ) == FLA_SCALAR ) return;
00718 
00719     // Since we didn't return, we know that H is of type FLA_MATRIX. Here, we
00720     // acquire the array of objects contained inside of H.
00721     buffer_H = ( FLA_Obj* ) FLA_Obj_buffer( *H );
00722 
00723     // For each submatrix in H...
00724     for ( j = 0; j < n_H; ++j )
00725     {
00726         for ( i = 0; i < m_H; ++i )
00727         {
00728             obj_ij_ptr = &buffer_H[ j*m_H + i ];
00729 
00730             // Recurse. The object at element (i,j) may also be of element type
00731             // FLA_MATRIX. If so, we will need to free its children.
00732             FLASH_Obj_free_hierarchy( obj_ij_ptr );
00733 
00734             // Free the base object at element (i,j). In order to avoid freeing
00735             // the object's data buffer, which would have already been freed en
00736             // masse by now if the calling function is FLASH_Obj_free(), we
00737             // will call FLA_Obj_free_without_buffer().
00738             FLA_Obj_free_without_buffer( obj_ij_ptr );
00739         }
00740     }
00741 }

FLA_Error FLASH_Obj_free_hierarchy_check ( FLA_Obj H  ) 

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_free_hierarchy().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( H );
00040   FLA_Check_error_code( e_val );
00041 
00042   return FLA_SUCCESS;
00043 }

void FLASH_Obj_free_without_buffer ( FLA_Obj H  ) 

References FLA_Check_error_level(), FLA_Obj_elemtype(), FLA_Obj_free(), FLA_Obj_free_without_buffer(), FLASH_Obj_free_hierarchy(), and FLASH_Obj_free_without_buffer_check().

00673 {
00674     if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
00675         FLASH_Obj_free_without_buffer_check( H );
00676 
00677     // Free the object according to whether it contains a hierarchy.
00678     if ( FLA_Obj_elemtype( *H ) == FLA_MATRIX )
00679     {
00680         // Skip freeing the numerical data buffer, since the object was
00681         // presumably created with FLASH_Obj_create_without_buffer().
00682 
00683         // Free the interior of the matrix hierarchy. This includes non-leaf
00684         // buffers and their corresponding base objects as well as leaf-level
00685         // base objects.
00686         FLASH_Obj_free_hierarchy( H );
00687 
00688         // Finally, free the top-level FLASH object itself.
00689         FLA_Obj_free( H );
00690     }
00691     else
00692     {
00693         // If the matrix has no hierarchy, treat it like a flat object with
00694         // no internal data buffer.
00695         FLA_Obj_free_without_buffer( H );
00696     }
00697 }

FLA_Error FLASH_Obj_free_without_buffer_check ( FLA_Obj H  ) 

References FLA_Check_null_pointer().

Referenced by FLASH_Obj_free_without_buffer().

00036 {
00037   FLA_Error e_val;
00038 
00039   e_val = FLA_Check_null_pointer( H );
00040   FLA_Check_error_code( e_val );
00041 
00042   return FLA_SUCCESS;
00043 }

FLA_Error FLASH_Obj_hierarchify ( FLA_Obj  F,
FLA_Obj  H 
)

void FLASH_Obj_push ( int  direction,
FLA_Obj  alpha,
FLA_Obj  F,
FLA_Obj  H,
FLASH_Obj_queue queue 
)

References FLASH_Obj_task_s::alpha, FLASH_Obj_task_s::direction, FLASH_Obj_task_s::F, FLA_malloc(), FLASH_Obj_task_s::H, FLASH_Obj_queue_s::head, FLASH_Obj_queue_s::n_tasks, FLASH_Obj_task_s::next_task, FLASH_Obj_task_s::prev_task, and FLASH_Obj_queue_s::tail.

Referenced by FLASH_Axpy_hierarchy_r(), and FLASH_Copy_hierarchy_r().

00278 {
00279    FLASH_Obj_task* t = (FLASH_Obj_task *) FLA_malloc( sizeof(FLASH_Obj_task) );
00280 
00281    // Initialize the task elements.
00282    t->direction = direction;
00283    t->alpha     = alpha;
00284    t->F         = F;
00285    t->H         = H;
00286    t->prev_task = NULL;
00287    t->next_task = NULL;
00288   
00289    // Add the task to the tail of the queue (and the head if queue is empty).
00290    if ( queue->n_tasks == 0 )
00291    {
00292       queue->head = t;
00293       queue->tail = t;
00294    }
00295    else
00296    {
00297       t->prev_task = queue->tail;
00298       queue->tail->next_task = t;
00299       queue->tail            = t;
00300    }
00301 
00302    // Increment the number of tasks.
00303    queue->n_tasks++;
00304 
00305    return;
00306 }

dim_t FLASH_Obj_scalar_length ( FLA_Obj  H  ) 

References FLA_Obj_view::base, FLA_Cont_with_3x1_to_2x1(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Part_2x1(), FLA_Repart_2x1_to_3x1(), and FLA_Obj_struct::m_inner.

Referenced by FLA_Check_submatrix_dims_and_offset(), FLA_QR_UT_inc_create_U_panel(), FLASH_Apply_Q_UT(), FLASH_Apply_Q_UT_inc(), FLASH_Apply_Q_UT_inc_create_workspace(), FLASH_Apply_Q_UT_UD(), FLASH_Axpy_hierarchy_r(), FLASH_Copy_hierarchy_r(), FLASH_Obj_attach_buffer(), FLASH_Obj_attach_buffer_hierarchy(), FLASH_Obj_create_conf_to(), FLASH_Obj_create_flat_conf_to_hier(), FLASH_QR_UT_inc_noopt(), FLASH_QR_UT_inc_opt1(), and FLASH_QR_UT_UD().

00107 {
00108     FLA_Obj  HT,              H0,
00109              HB,              H1,
00110                               H2;
00111     FLA_Obj* H1p;
00112 
00113     dim_t b = 0;
00114 
00115     if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
00116         return FLA_Obj_length( H );
00117 
00118     FLA_Part_2x1( H,    &HT,
00119                         &HB,            0, FLA_TOP );
00120 
00121     while ( FLA_Obj_length( HT ) < FLA_Obj_length( H ) )
00122     {
00123         FLA_Repart_2x1_to_3x1( HT,                &H0,
00124                             /* ** */            /* ** */
00125                                                   &H1,
00126                                HB,                &H2,        1, FLA_BOTTOM );
00127 
00128         /*------------------------------------------------------------*/
00129 
00130         H1p = FLASH_OBJ_PTR_AT( H1 );
00131         b += H1p->base->m_inner;
00132 
00133         /*------------------------------------------------------------*/
00134 
00135         FLA_Cont_with_3x1_to_2x1( &HT,                H0,
00136                                                       H1,
00137                                 /* ** */           /* ** */
00138                                   &HB,                H2,     FLA_TOP );
00139   }
00140 
00141   return b;
00142 }

dim_t FLASH_Obj_scalar_width ( FLA_Obj  H  ) 

References FLA_Obj_view::base, FLA_Cont_with_1x3_to_1x2(), FLA_Obj_elemtype(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Repart_1x2_to_1x3(), and FLA_Obj_struct::n_inner.

Referenced by FLA_Check_submatrix_dims_and_offset(), FLASH_Apply_Q_UT_inc_create_workspace(), FLASH_Axpy_hierarchy_r(), FLASH_Copy_hierarchy_r(), FLASH_FS_incpiv(), FLASH_LU_incpiv(), FLASH_Obj_attach_buffer(), FLASH_Obj_attach_buffer_hierarchy(), FLASH_Obj_create_conf_to(), and FLASH_Obj_create_flat_conf_to_hier().

00146 {
00147     FLA_Obj  HL,    HR,       H0,  H1,  H2;
00148     FLA_Obj* H1p;
00149 
00150     dim_t b = 0;
00151 
00152     if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
00153         return FLA_Obj_width( H );
00154 
00155     FLA_Part_1x2( H,    &HL,  &HR,      0, FLA_LEFT );
00156 
00157     while ( FLA_Obj_width( HL ) < FLA_Obj_width( H ) )
00158     {
00159         FLA_Repart_1x2_to_1x3( HL,  /**/ HR,        &H0, /**/ &H1, &H2,
00160                                1, FLA_RIGHT );
00161 
00162         /*------------------------------------------------------------*/
00163 
00164         H1p = FLASH_OBJ_PTR_AT( H1 );
00165         b += H1p->base->n_inner;
00166 
00167         /*------------------------------------------------------------*/
00168 
00169         FLA_Cont_with_1x3_to_1x2( &HL,  /**/ &HR,        H0, H1, /**/ H2,
00170                                   FLA_LEFT );
00171     }
00172 
00173     return b;
00174 }

FLA_Error FLASH_Obj_show ( char *  header,
FLA_Obj  H,
char *  elem_format,
char *  footer 
)

References FLA_Obj_elemtype(), FLA_Obj_free(), FLA_Obj_show(), and FLASH_Obj_create_flat_copy_of_hier().

00785 {
00786     FLA_Obj F;
00787     
00788     // Show the object according to whether it contains a hierarchy.
00789     if ( FLA_Obj_elemtype( H ) == FLA_MATRIX )
00790     {
00791         // Create a temporary flat object conformal to and copied from the
00792         // hierarchical object.
00793         FLASH_Obj_create_flat_copy_of_hier( H, &F );
00794     
00795         // We can now display the temporary matrix in a straightforward manner.
00796         FLA_Obj_show( header, F, elem_format, footer );
00797     
00798         // Free the temporary flat matrix.
00799         FLA_Obj_free( &F );
00800     }
00801     else
00802     {
00803         // Display the flat matrix.
00804         FLA_Obj_show( header, H, elem_format, footer );
00805     }
00806 
00807     return FLA_SUCCESS;
00808 }

void FLASH_print_struct ( FLA_Obj  H  ) 

References FLA_Obj_buffer(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_width(), and FLASH_print_struct_helper().

00948 {
00949     dim_t    i, j, m_H, n_H;
00950     FLA_Obj* buffer_temp;
00951 
00952     m_H = FLA_Obj_length( H );
00953     n_H = FLA_Obj_width( H );
00954 
00955     if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
00956         FLASH_print_struct_helper( H, 0 );
00957     else
00958     {
00959         for ( j = 0; j < n_H; ++j )
00960         {
00961             for ( i = 0; i < m_H; ++i )
00962             {
00963                 buffer_temp = ( FLA_Obj* ) FLA_Obj_buffer( H );
00964 
00965                 FLASH_print_struct_helper( buffer_temp[ j*m_H + i ], 0 );
00966             }
00967         }
00968     }
00969 }

void FLASH_print_struct_helper ( FLA_Obj  H,
int  indent 
)

References FLA_Obj_buffer(), FLA_Obj_datatype(), FLA_Obj_elemtype(), FLA_Obj_ldim(), FLA_Obj_length(), FLA_Obj_width(), and FLASH_print_struct_helper().

Referenced by FLASH_print_struct(), and FLASH_print_struct_helper().

00973 {
00974     dim_t    i, j, k, m_H, n_H;
00975     FLA_Obj* buffer_temp;
00976 
00977     for ( i = 0; i < indent; ++i )
00978         fprintf( stdout, "  " );
00979 
00980     if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
00981     {
00982         fprintf( stdout, "LEAF (%3d | ldim %3d | %3d x %3d | addr %p)\n",
00983                          FLA_Obj_datatype( H ), FLA_Obj_ldim( H ),
00984                          FLA_Obj_length( H ), FLA_Obj_width( H ),
00985                          FLA_Obj_buffer( H ) );
00986         fflush( stdout );
00987     }
00988     else
00989     {
00990         m_H = FLA_Obj_length( H );
00991         n_H = FLA_Obj_width( H );
00992         
00993         fprintf( stdout, "MATRIX (%dx%d):%d - %p\n",
00994                          m_H, n_H,
00995                          FLA_Obj_datatype( H ),
00996                          FLA_Obj_buffer( H ) );
00997         fflush( stdout );
00998         
00999         for ( j = 0; j < n_H; ++j )
01000         {
01001             for ( i = 0; i < m_H; ++i )
01002             {
01003                 for ( k = 0; k < indent; ++k )
01004                     fprintf( stdout, "  " );
01005                 
01006                 buffer_temp = ( FLA_Obj* ) FLA_Obj_buffer( H );
01007 
01008                 FLASH_print_struct_helper( buffer_temp[ j*m_H + i ],
01009                                            indent + 1);
01010             }
01011         }
01012     }
01013 }


Generated on Mon Jul 6 05:45:51 2009 for libflame by  doxygen 1.5.9