libflame  revision_anchor
Functions
FLASH_Obj.c File Reference

(r)

Functions

FLA_Datatype FLASH_Obj_datatype (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)
 
dim_t FLASH_Obj_base_scalar_length (FLA_Obj H)
 
dim_t FLASH_Obj_base_scalar_width (FLA_Obj H)
 
FLA_Error FLASH_Obj_create (FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *b_mn, 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 *b_mn, 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, FLA_Obj *H_new)
 
FLA_Error FLASH_Obj_create_hier_conf_to_flat (FLA_Trans trans, FLA_Obj F, dim_t depth, dim_t *b_mn, 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_copy_of (FLA_Trans trans, FLA_Obj H_cur, FLA_Obj *H_new)
 
FLA_Error FLASH_Obj_create_hier_copy_of_flat (FLA_Obj F, dim_t depth, dim_t *b_mn, 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_without_buffer (FLA_Obj *H)
 
void FLASH_Obj_free_hierarchy (FLA_Obj *H)
 
void * FLASH_Obj_extract_buffer (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)
 
FLA_Error FLASH_Obj_attach_buffer (void *buffer, dim_t rs, dim_t cs, FLA_Obj *H)
 
FLA_Error FLASH_Obj_attach_buffer_hierarchy (FLA_Obj F, FLA_Obj *H)
 
void FLASH_print_struct (FLA_Obj H)
 
void FLASH_print_struct_helper (FLA_Obj H, int indent)
 

Function Documentation

◆ FLASH_Obj_attach_buffer()

FLA_Error FLASH_Obj_attach_buffer ( void *  buffer,
dim_t  rs,
dim_t  cs,
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_base_scalar_length(), FLASH_Obj_base_scalar_width(), and FLASH_Obj_datatype().

782 {
783  FLA_Obj flat_matrix;
784  dim_t m_base, n_base;
785  FLA_Datatype datatype;
786 
787  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
788  FLASH_Obj_attach_buffer_check( buffer, rs, cs, H );
789 
790  // Extract the scalar dimensions of the base object(s) and get its
791  // numerical datatype. (These fields will be set even if it has a NULL
792  // buffer, which it probably does since this function was just invoked.)
793  m_base = FLASH_Obj_base_scalar_length( *H );
794  n_base = FLASH_Obj_base_scalar_width( *H );
795  datatype = FLASH_Obj_datatype( *H );
796 
797  // Create a temporary conventional object and attach the given buffer.
798  // Segments of this buffer will be partitioned out to the various
799  // leaf-level matrices of the hierarchical matrix H.
800  FLA_Obj_create_without_buffer( datatype, m_base, n_base, &flat_matrix );
801  FLA_Obj_attach_buffer( buffer, rs, cs, &flat_matrix );
802 
803  // Recurse through the hierarchical matrix, assigning segments of
804  // flat_matrix to the various leaf-level matrices, similar to what
805  // we would do if we were creating the object outright.
806  FLASH_Obj_attach_buffer_hierarchy( flat_matrix, H );
807 
808  // Free the object (but don't free the buffer!).
809  FLA_Obj_free_without_buffer( &flat_matrix );
810 
811  return FLA_SUCCESS;
812 }
dim_t FLASH_Obj_base_scalar_length(FLA_Obj H)
Definition: FLASH_Obj.c:83
FLA_Error FLASH_Obj_attach_buffer_hierarchy(FLA_Obj F, FLA_Obj *H)
Definition: FLASH_Obj.c:815
unsigned long dim_t
Definition: FLA_type_defs.h:71
FLA_Error FLASH_Obj_attach_buffer_check(void *buffer, dim_t rs, dim_t cs, FLA_Obj *H)
Definition: FLASH_Obj_attach_buffer_check.c:13
dim_t FLASH_Obj_base_scalar_width(FLA_Obj H)
Definition: FLASH_Obj.c:113
FLA_Error FLA_Obj_attach_buffer(void *buffer, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition: FLA_Obj.c:522
FLA_Error FLA_Obj_create_without_buffer(FLA_Datatype datatype, dim_t m, dim_t n, FLA_Obj *obj)
Definition: FLA_Obj.c:362
Definition: FLA_type_defs.h:158
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int FLA_Datatype
Definition: FLA_type_defs.h:49
FLA_Error FLA_Obj_free_without_buffer(FLA_Obj *obj)
Definition: FLA_Obj.c:615
FLA_Datatype FLASH_Obj_datatype(FLA_Obj H)
Definition: FLASH_Obj.c:14

◆ FLASH_Obj_attach_buffer_hierarchy()

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_at_view(), FLA_Obj_col_stride(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_row_stride(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Part_2x1(), FLA_Repart_1x2_to_1x3(), FLA_Repart_2x1_to_3x1(), FLASH_Obj_attach_buffer_hierarchy_check(), FLASH_Obj_base_scalar_length(), and FLASH_Obj_base_scalar_width().

Referenced by FLASH_Obj_attach_buffer().

816 {
817  FLA_Obj FL, FR, F0, F1, F2;
818 
819  FLA_Obj HL, HR, H0, H1, H2;
820 
821  FLA_Obj F1T, F01,
822  F1B, F11,
823  F21;
824 
825  FLA_Obj H1T, H01,
826  H1B, H11,
827  H21;
828 
829  dim_t b_m, b_n;
830 
831  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
833 
834  if ( FLA_Obj_elemtype( *H ) == FLA_SCALAR )
835  {
836  // If we've recursed down to a leaf node, then we can simply attach
837  // the matrix buffer to the current leaf-level submatrix.
838  // Notice we use FLA_Obj_buffer_at_view() because we want to attach
839  // the buffer address referenced by the view of F.
841  FLA_Obj_row_stride( F ),
842  FLA_Obj_col_stride( F ), H );
843  }
844  else
845  {
846  FLA_Part_1x2( *H, &HL, &HR, 0, FLA_LEFT );
847 
848  FLA_Part_1x2( F, &FL, &FR, 0, FLA_LEFT );
849 
850  while ( FLA_Obj_width( HL ) < FLA_Obj_width( *H ) )
851  {
852 
853  FLA_Repart_1x2_to_1x3( HL, /**/ HR, &H0, /**/ &H1, &H2,
854  1, FLA_RIGHT );
855 
856  b_n = FLASH_Obj_base_scalar_width( H1 );
857 
858  FLA_Repart_1x2_to_1x3( FL, /**/ FR, &F0, /**/ &F1, &F2,
859  b_n, FLA_RIGHT );
860 
861  /*------------------------------------------------------------*/
862 
863  FLA_Part_2x1( H1, &H1T,
864  &H1B, 0, FLA_TOP );
865 
866  FLA_Part_2x1( F1, &F1T,
867  &F1B, 0, FLA_TOP );
868 
869  while ( FLA_Obj_length( H1T ) < FLA_Obj_length( H1 ) )
870  {
871 
872  FLA_Repart_2x1_to_3x1( H1T, &H01,
873  /* ** */ /* ** */
874  &H11,
875  H1B, &H21, 1, FLA_BOTTOM );
876 
877  b_m = FLASH_Obj_base_scalar_length( H11 );
878 
879  FLA_Repart_2x1_to_3x1( F1T, &F01,
880  /* ** */ /* ** */
881  &F11,
882  F1B, &F21, b_m, FLA_BOTTOM );
883 
884  /*------------------------------------------------------------*/
885 
887  FLASH_OBJ_PTR_AT( H11 ) );
888 
889  /*------------------------------------------------------------*/
890 
891  FLA_Cont_with_3x1_to_2x1( &H1T, H01,
892  H11,
893  /* ** */ /* ** */
894  &H1B, H21, FLA_TOP );
895 
896  FLA_Cont_with_3x1_to_2x1( &F1T, F01,
897  F11,
898  /* ** */ /* ** */
899  &F1B, F21, FLA_TOP );
900  }
901 
902  /*------------------------------------------------------------*/
903 
904  FLA_Cont_with_1x3_to_1x2( &HL, /**/ &HR, H0, H1, /**/ H2,
905  FLA_LEFT );
906 
907  FLA_Cont_with_1x3_to_1x2( &FL, /**/ &FR, F0, F1, /**/ F2,
908  FLA_LEFT );
909 
910  }
911  }
912 
913  return FLA_SUCCESS;
914 }
FLA_Error FLA_Repart_2x1_to_3x1(FLA_Obj AT, FLA_Obj *A0, FLA_Obj *A1, FLA_Obj AB, FLA_Obj *A2, dim_t mb, FLA_Side side)
Definition: FLA_View.c:226
dim_t FLASH_Obj_base_scalar_length(FLA_Obj H)
Definition: FLASH_Obj.c:83
FLA_Error FLA_Repart_1x2_to_1x3(FLA_Obj AL, FLA_Obj AR, FLA_Obj *A0, FLA_Obj *A1, FLA_Obj *A2, dim_t nb, FLA_Side side)
Definition: FLA_View.c:267
FLA_Error FLASH_Obj_attach_buffer_hierarchy(FLA_Obj F, FLA_Obj *H)
Definition: FLASH_Obj.c:815
unsigned long dim_t
Definition: FLA_type_defs.h:71
dim_t FLA_Obj_row_stride(FLA_Obj obj)
Definition: FLA_Query.c:167
dim_t FLASH_Obj_base_scalar_width(FLA_Obj H)
Definition: FLASH_Obj.c:113
FLA_Error FLA_Cont_with_3x1_to_2x1(FLA_Obj *AT, FLA_Obj A0, FLA_Obj A1, FLA_Obj *AB, FLA_Obj A2, FLA_Side side)
Definition: FLA_View.c:428
FLA_Error FLA_Obj_attach_buffer(void *buffer, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition: FLA_Obj.c:522
void * FLA_Obj_buffer_at_view(FLA_Obj obj)
Definition: FLA_Query.c:215
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_width(FLA_Obj obj)
Definition: FLA_Query.c:123
FLA_Error FLASH_Obj_attach_buffer_hierarchy_check(FLA_Obj F, FLA_Obj *H)
Definition: FLASH_Obj_attach_buffer_hierarchy_check.c:13
FLA_Error FLA_Cont_with_1x3_to_1x2(FLA_Obj *AL, FLA_Obj *AR, FLA_Obj A0, FLA_Obj A1, FLA_Obj A2, FLA_Side side)
Definition: FLA_View.c:475
FLA_Error FLA_Part_2x1(FLA_Obj A, FLA_Obj *A1, FLA_Obj *A2, dim_t mb, FLA_Side side)
Definition: FLA_View.c:76
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLA_Part_1x2(FLA_Obj A, FLA_Obj *A1, FLA_Obj *A2, dim_t nb, FLA_Side side)
Definition: FLA_View.c:110
dim_t FLA_Obj_col_stride(FLA_Obj obj)
Definition: FLA_Query.c:174
dim_t FLA_Obj_length(FLA_Obj obj)
Definition: FLA_Query.c:116
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_base_scalar_length()

dim_t FLASH_Obj_base_scalar_length ( FLA_Obj  H)

References FLA_Obj_view::base, FLA_Obj_base_buffer(), FLA_Obj_base_length(), FLA_Obj_col_stride(), FLA_Obj_elemtype(), FLA_Obj_row_stride(), and i.

Referenced by FLASH_Obj_attach_buffer(), FLASH_Obj_attach_buffer_check(), FLASH_Obj_attach_buffer_hierarchy(), FLASH_Obj_create_conf_to(), FLASH_Obj_scalar_length_tl(), FLASH_Part_create_1x2(), FLASH_Part_create_2x1(), and FLASH_Part_create_2x2().

84 {
85  FLA_Obj* buffer;
86  dim_t m;
87  dim_t rs, cs;
88  dim_t i;
89  dim_t m_base = 0;
90 
91  if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
92  return FLA_Obj_base_length( H );
93 
94  // Notice we use the base buffer since we are interested in the
95  // whole object, not just the part referened by the view.
96  buffer = FLA_Obj_base_buffer( H );
97  m = FLA_Obj_base_length( H );
98  rs = FLA_Obj_row_stride( H );
99  cs = FLA_Obj_col_stride( H );
100 
101  // Add up the row dimensions of all the base objects in the 0th
102  // column of objects.
103  for ( i = 0; i < m; ++i )
104  {
105  FLA_Obj hij = buffer[ i*rs + 0*cs ];
106 
107  m_base += (hij.base)->m_inner;
108  }
109 
110  return m_base;
111 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
dim_t FLA_Obj_row_stride(FLA_Obj obj)
Definition: FLA_Query.c:167
FLA_Base_obj * base
Definition: FLA_type_defs.h:168
void * FLA_Obj_base_buffer(FLA_Obj obj)
Definition: FLA_Query.c:210
dim_t FLA_Obj_base_length(FLA_Obj obj)
Definition: FLA_Query.c:192
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_col_stride(FLA_Obj obj)
Definition: FLA_Query.c:174
int i
Definition: bl1_axmyv2.c:145
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_base_scalar_width()

dim_t FLASH_Obj_base_scalar_width ( FLA_Obj  H)

References FLA_Obj_view::base, FLA_Obj_base_buffer(), FLA_Obj_base_width(), FLA_Obj_col_stride(), FLA_Obj_elemtype(), and FLA_Obj_row_stride().

Referenced by FLASH_Obj_attach_buffer(), FLASH_Obj_attach_buffer_check(), FLASH_Obj_attach_buffer_hierarchy(), FLASH_Obj_create_conf_to(), FLASH_Obj_scalar_width_tl(), FLASH_Part_create_1x2(), FLASH_Part_create_2x1(), and FLASH_Part_create_2x2().

114 {
115  FLA_Obj* buffer;
116  dim_t n;
117  dim_t rs, cs;
118  dim_t j;
119  dim_t n_base = 0;
120 
121  if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
122  return FLA_Obj_base_width( H );
123 
124  // Notice we use the base buffer since we are interested in the
125  // whole object, not just the part referened by the view.
126  buffer = FLA_Obj_base_buffer( H );
127  n = FLA_Obj_base_width( H );
128  rs = FLA_Obj_row_stride( H );
129  cs = FLA_Obj_col_stride( H );
130 
131  // Add up the column dimensions of all the base objects in the 0th
132  // row of objects.
133  for ( j = 0; j < n; ++j )
134  {
135  FLA_Obj hij = buffer[ 0*rs + j*cs ];
136 
137  n_base += (hij.base)->n_inner;
138  }
139 
140  return n_base;
141 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
dim_t FLA_Obj_row_stride(FLA_Obj obj)
Definition: FLA_Query.c:167
FLA_Base_obj * base
Definition: FLA_type_defs.h:168
void * FLA_Obj_base_buffer(FLA_Obj obj)
Definition: FLA_Query.c:210
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_base_width(FLA_Obj obj)
Definition: FLA_Query.c:198
dim_t FLA_Obj_col_stride(FLA_Obj obj)
Definition: FLA_Query.c:174
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_blocksizes()

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

References FLA_Check_error_level(), FLA_Obj_base_buffer(), FLA_Obj_base_length(), FLA_Obj_base_width(), FLA_Obj_elemtype(), and FLASH_Obj_blocksizes_check().

Referenced by FLASH_Obj_create_conf_to(), FLASH_Part_create_1x2(), FLASH_Part_create_2x1(), and FLASH_Part_create_2x2().

50 {
51  FLA_Elemtype elemtype;
52  FLA_Obj* buffer_H;
53  dim_t depth = 0;
54 
55  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
56  FLASH_Obj_blocksizes_check( H, b_m, b_n );
57 
58  // Recurse through the hierarchy to the first leaf node. We initialize
59  // the recursion here:
60  elemtype = FLA_Obj_elemtype( H );
61  buffer_H = ( FLA_Obj* ) FLA_Obj_base_buffer( H );
62 
63  while ( elemtype == FLA_MATRIX )
64  {
65  b_m[depth] = FLA_Obj_base_length( buffer_H[0] );
66  b_n[depth] = FLA_Obj_base_width( buffer_H[0] );
67  ++depth;
68 
69  // Get the element type of the top-leftmost underlying object. Also,
70  // get a pointer to the first element of the top-leftmost object and
71  // assume that it is of type FLA_Obj* in case elemtype is once again
72  // FLA_MATRIX.
73  elemtype = FLA_Obj_elemtype( buffer_H[0] );
74  buffer_H = ( FLA_Obj * ) FLA_Obj_base_buffer( buffer_H[0] );
75  }
76 
77  // At this point, the first depth elements of blocksizes have been filled
78  // with the blocksizes of H's various hierarchical levels. Return the
79  // matrix depth as a confirmation of how many blocksizes were found.
80  return depth;
81 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
int FLA_Elemtype
Definition: FLA_type_defs.h:50
void * FLA_Obj_base_buffer(FLA_Obj obj)
Definition: FLA_Query.c:210
dim_t FLA_Obj_base_length(FLA_Obj obj)
Definition: FLA_Query.c:192
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_base_width(FLA_Obj obj)
Definition: FLA_Query.c:198
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLASH_Obj_blocksizes_check(FLA_Obj H, dim_t *b_m, dim_t *b_n)
Definition: FLASH_Obj_blocksizes_check.c:13
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_create()

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

References FLASH_Obj_create_helper().

Referenced by FLASH_Obj_create_diag_panel(), and FLASH_Obj_create_hier_conf_to_flat().

144 {
145  FLASH_Obj_create_helper( FALSE, datatype, m, n, depth, b_mn, b_mn, H );
146 
147  return FLA_SUCCESS;
148 }
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)
Definition: FLASH_Obj.c:175

◆ FLASH_Obj_create_conf_to()

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

References FLA_Check_error_level(), FLA_free(), FLA_malloc(), FLASH_Obj_adjust_views(), FLASH_Obj_base_scalar_length(), FLASH_Obj_base_scalar_width(), FLASH_Obj_blocksizes(), FLASH_Obj_create_conf_to_check(), FLASH_Obj_create_ext(), FLASH_Obj_datatype(), FLASH_Obj_depth(), FLASH_Obj_scalar_col_offset(), FLASH_Obj_scalar_length(), FLASH_Obj_scalar_row_offset(), and FLASH_Obj_scalar_width().

Referenced by FLASH_CAQR_UT_inc_create_hier_matrices(), FLASH_Eig_gest(), and FLASH_Obj_create_copy_of().

407 {
408  FLA_Datatype datatype;
409  dim_t m_base, n_base;
410  dim_t m_view, n_view;
411  dim_t offm_scalar, offn_scalar;
412  dim_t depth;
413  dim_t* b_m;
414  dim_t* b_n;
415 
416  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
417  FLASH_Obj_create_conf_to_check( trans, H, H_new );
418 
419  // Acquire some properties of the hierarchical matrix object.
420  datatype = FLASH_Obj_datatype( H );
421  m_base = FLASH_Obj_base_scalar_length( H );
422  n_base = FLASH_Obj_base_scalar_width( H );
423  m_view = FLASH_Obj_scalar_length( H );
424  n_view = FLASH_Obj_scalar_width( H );
425  offm_scalar = FLASH_Obj_scalar_row_offset( H );
426  offn_scalar = FLASH_Obj_scalar_col_offset( H );
427  depth = FLASH_Obj_depth( H );
428 
429  // Allocate a pair of temporary arrays for the blocksizes, whose lengths
430  // are equal to the object's hierarchical depth.
431  b_m = ( dim_t* ) FLA_malloc( depth * sizeof( dim_t ) );
432  b_n = ( dim_t* ) FLA_malloc( depth * sizeof( dim_t ) );
433 
434  // Accumulate the blocksizes into the blocksize buffer.
435  FLASH_Obj_blocksizes( H, b_m, b_n );
436 
437  // Handle the transposition, if requested.
438  if ( trans == FLA_TRANSPOSE )
439  {
440  FLA_Check_error_code( FLA_NOT_YET_IMPLEMENTED );
441  }
442 
443  // Create the new hierarchical matrix object with the same base dimensions
444  // as the original object..
445  FLASH_Obj_create_ext( datatype, m_base, n_base, depth, b_m, b_n, H_new );
446 
447  // Adjust the hierarchical view of the new object to match that of the
448  // original object.
449  FLASH_Obj_adjust_views( FALSE, offm_scalar, offn_scalar, m_view, n_view, H, H_new );
450 
451  // Free the temporary blocksize buffers.
452  FLA_free( b_m );
453  FLA_free( b_n );
454 
455  return FLA_SUCCESS;
456 }
dim_t FLASH_Obj_base_scalar_length(FLA_Obj H)
Definition: FLASH_Obj.c:83
unsigned long dim_t
Definition: FLA_type_defs.h:71
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)
Definition: FLASH_Obj.c:151
FLA_Error FLASH_Obj_create_conf_to_check(FLA_Trans trans, FLA_Obj H_cur, FLA_Obj *H_new)
Definition: FLASH_Obj_create_conf_to_check.c:13
dim_t FLASH_Obj_base_scalar_width(FLA_Obj H)
Definition: FLASH_Obj.c:113
dim_t FLASH_Obj_scalar_length(FLA_Obj H)
Definition: FLASH_View.c:600
dim_t FLASH_Obj_scalar_row_offset(FLA_Obj H)
Definition: FLASH_View.c:693
void FLA_free(void *ptr)
Definition: FLA_Memory.c:247
void * FLA_malloc(size_t size)
Definition: FLA_Memory.c:111
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int FLA_Datatype
Definition: FLA_type_defs.h:49
dim_t FLASH_Obj_scalar_col_offset(FLA_Obj H)
Definition: FLASH_View.c:708
FLA_Error FLASH_Obj_adjust_views(FLA_Bool attach_buffer, dim_t offm, dim_t offn, dim_t m, dim_t n, FLA_Obj A, FLA_Obj *S)
Definition: FLASH_View.c:275
dim_t FLASH_Obj_scalar_width(FLA_Obj H)
Definition: FLASH_View.c:641
dim_t FLASH_Obj_blocksizes(FLA_Obj H, dim_t *b_m, dim_t *b_n)
Definition: FLASH_Obj.c:49
dim_t FLASH_Obj_depth(FLA_Obj H)
Definition: FLASH_Obj.c:20
FLA_Datatype FLASH_Obj_datatype(FLA_Obj H)
Definition: FLASH_Obj.c:14

◆ FLASH_Obj_create_copy_of()

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

References FLA_Obj_create_copy_of(), FLA_Obj_free(), FLASH_Copy(), FLASH_Obj_create_conf_to(), FLASH_Obj_create_flat_copy_of_hier(), and FLASH_Obj_hierarchify().

Referenced by FLASH_CAQR_UT_inc_solve(), FLASH_QR_UT_inc_solve(), FLASH_QR_UT_solve(), and FLASH_UDdate_UT_inc_update_rhs().

562 {
563  // Create a new object conformal to the current object.
564  FLASH_Obj_create_conf_to( trans, H_cur, H_new );
565 
566  // This is a workaround until we implement a FLASH_Copyt().
567  if ( trans == FLA_NO_TRANSPOSE || trans == FLA_CONJ_NO_TRANSPOSE )
568  {
569  // Copy the contents of the current object to the new object.
570  FLASH_Copy( H_cur, *H_new );
571 
572  // NOTE: we don't currently honor requests to conjugate!
573  // We could, if we had FLASH_Conj() implemented, but we don't
574  // currently.
575  }
576  else // if ( trans == FLA_TRANSPOSE || trans == FLA_CONJ_TRANSPOSE )
577  {
578  FLA_Obj F, F_trans;
579 
581  FLA_Obj_create_copy_of( trans, F, &F_trans );
582  FLASH_Obj_hierarchify( F_trans, *H_new );
583  FLA_Obj_free( &F );
584  FLA_Obj_free( &F_trans );
585  }
586 
587  return FLA_SUCCESS;
588 }
FLA_Error FLASH_Obj_create_conf_to(FLA_Trans trans, FLA_Obj H, FLA_Obj *H_new)
Definition: FLASH_Obj.c:406
FLA_Error FLA_Obj_free(FLA_Obj *obj)
Definition: FLA_Obj.c:588
FLA_Error FLASH_Obj_create_flat_copy_of_hier(FLA_Obj H, FLA_Obj *F)
Definition: FLASH_Obj.c:623
FLA_Error FLA_Obj_create_copy_of(FLA_Trans trans, FLA_Obj old, FLA_Obj *obj)
Definition: FLA_Obj.c:345
FLA_Error FLASH_Obj_hierarchify(FLA_Obj F, FLA_Obj H)
Definition: FLASH_Obj.c:773
Definition: FLA_type_defs.h:158
FLA_Error FLASH_Copy(FLA_Obj A, FLA_Obj B)
Definition: FLASH_Copy.c:15

◆ FLASH_Obj_create_ext()

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 
)

◆ FLASH_Obj_create_flat_conf_to_hier()

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().

528 {
529  FLA_Datatype datatype;
530  dim_t m_H, n_H;
531  dim_t m_F, n_F;
532 
533  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
535 
536  // Acquire the numerical datatype, length, and width of the
537  // hierarchical matrix object.
538  datatype = FLASH_Obj_datatype( H );
539  m_H = FLASH_Obj_scalar_length( H );
540  n_H = FLASH_Obj_scalar_width( H );
541 
542  // Factor in the transposition, if requested.
543  if ( trans == FLA_NO_TRANSPOSE )
544  {
545  m_F = m_H;
546  n_F = n_H;
547  }
548  else
549  {
550  m_F = n_H;
551  n_F = m_H;
552  }
553 
554  // Create the flat matrix object. Default to column-major storage.
555  FLA_Obj_create( datatype, m_F, n_F, 0, 0, F );
556 
557  return FLA_SUCCESS;
558 }
FLA_Error FLA_Obj_create(FLA_Datatype datatype, dim_t m, dim_t n, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition: FLA_Obj.c:55
unsigned long dim_t
Definition: FLA_type_defs.h:71
dim_t FLASH_Obj_scalar_length(FLA_Obj H)
Definition: FLASH_View.c:600
FLA_Error FLASH_Obj_create_flat_conf_to_hier_check(FLA_Trans trans, FLA_Obj H, FLA_Obj *F)
Definition: FLASH_Obj_create_flat_conf_to_hier_check.c:13
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int FLA_Datatype
Definition: FLA_type_defs.h:49
dim_t FLASH_Obj_scalar_width(FLA_Obj H)
Definition: FLASH_View.c:641
FLA_Datatype FLASH_Obj_datatype(FLA_Obj H)
Definition: FLASH_Obj.c:14

◆ FLASH_Obj_create_flat_copy_of_hier()

FLA_Error FLASH_Obj_create_flat_copy_of_hier ( FLA_Obj  H,
FLA_Obj F 
)

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

Referenced by FLA_LQ_UT_macro_task(), FLA_LU_piv_macro_task(), FLA_QR_UT_macro_task(), FLASH_Hermitianize(), FLASH_Max_elemwise_diff(), FLASH_Norm1(), FLASH_Obj_create_copy_of(), FLASH_Random_matrix(), FLASH_Random_spd_matrix(), FLASH_Set(), FLASH_Shift_diag(), and FLASH_Triangularize().

624 {
625  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
627 
628  // Create a flat object conformal to the hierarchical object.
629  FLASH_Obj_create_flat_conf_to_hier( FLA_NO_TRANSPOSE, H, F );
630 
631  // Flatten the hierarchical object's contents into the new flat object.
632  FLASH_Copy_hier_to_flat( 0, 0, H, *F );
633 
634  return FLA_SUCCESS;
635 }
FLA_Error FLASH_Obj_create_flat_conf_to_hier(FLA_Trans trans, FLA_Obj H, FLA_Obj *F)
Definition: FLASH_Obj.c:527
FLA_Error FLASH_Obj_create_flat_copy_of_hier_check(FLA_Obj H, FLA_Obj *F)
Definition: FLASH_Obj_create_flat_copy_of_hier_check.c:13
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLASH_Copy_hier_to_flat(dim_t i, dim_t j, FLA_Obj H, FLA_Obj F)
Definition: FLASH_Copy_other.c:110

◆ FLASH_Obj_create_helper()

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(), FLASH_Obj_create_hierarchy(), and i.

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

176 {
177  dim_t i;
178  FLA_Obj flat_matrix;
179 
180  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
181  FLASH_Obj_create_helper_check( without_buffer, datatype, m, n, depth, b_m, b_n, H );
182 
183  if ( depth == 0 )
184  {
185  // Base case: create a single contiguous matrix block. If we are
186  // creating an object with a buffer, then we use column-major order.
187  if ( without_buffer == FALSE )
188  FLA_Obj_create( datatype, m, n, 0, 0, H );
189  else
190  FLA_Obj_create_without_buffer( datatype, m, n, H );
191  }
192  else
193  {
194  // We need temporary arrays the same length as the blocksizes arrays.
195  dim_t* elem_sizes_m = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
196  dim_t* elem_sizes_n = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
197  dim_t* depth_sizes_m = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
198  dim_t* depth_sizes_n = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
199  dim_t* m_offsets = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
200  dim_t* n_offsets = ( dim_t * ) FLA_malloc( depth * sizeof( dim_t ) );
201 
202  // Fill two sets of arrays: elem_sizes_m/elem_sizes_n and depth_sizes_m/
203  // depth_sizes_n.
204  // - elem_sizes_m[i] will contain the number of numerical elements that span
205  // the row dimension of a block at the ith level of the hierarchy. This is
206  // just the product of all row blocksizes "internal" to and including the
207  // current blocking level. (The elem_sizes_n array tracks similar values
208  // in the column dimension.)
209  // - depth_sizes_m[i] is similar to elem_sizes_m[i]. The only difference is
210  // that instead of tracking the number of numerical elements in the row
211  // dimension, it tracks the number of "storage" blocks that span the m
212  // dimension of a block at the ith level, where the m dimension of a
213  // storage block is the block size given in b_m[depth-1], ie:
214  // the inner-most row dimension block size. (The depth_sizes_n array
215  // tracks similar values in the column dimension.)
216  elem_sizes_m[depth-1] = b_m[depth-1];
217  elem_sizes_n[depth-1] = b_n[depth-1];
218  depth_sizes_m[depth-1] = 1;
219  depth_sizes_n[depth-1] = 1;
220  for ( i = depth - 1; i > 0; --i )
221  {
222  elem_sizes_m[i-1] = elem_sizes_m[i] * b_m[i-1];
223  elem_sizes_n[i-1] = elem_sizes_n[i] * b_n[i-1];
224  depth_sizes_m[i-1] = depth_sizes_m[i] * b_m[i-1];
225  depth_sizes_n[i-1] = depth_sizes_n[i] * b_n[i-1];
226  }
227 
228  // Initialize the m_offsets and n_offsets arrays to zero.
229  for ( i = 0; i < depth; i++ )
230  {
231  m_offsets[i] = 0;
232  n_offsets[i] = 0;
233  }
234 
235  // Create a "flat" matrix object. All leaf-level child objects will refer
236  // to various offsets within this object's buffer. Whether we create the
237  // object with row- or column-major storage is moot, since either way it
238  // will be a 1-by-mn length matrix which we will partition through later
239  // on in FLASH_Obj_create_hierarchy(). Note that it is IMPORTANT that the
240  // matrix be 1-by-mn, and NOT m-by-n, since we want to use the 1x2
241  // partitioning routines to walk through it as we attach various parts of
242  // the buffer to the matrix hierarchy.
243  if ( without_buffer == FALSE )
244  FLA_Obj_create( datatype, 1, m*n, 0, 0, &flat_matrix );
245  else
246  FLA_Obj_create_without_buffer( datatype, m, n, &flat_matrix );
247 
248  // Recursively create the matrix hierarchy.
249  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 );
250 
251  // Free the flat_matrix object, but not its buffer. If we created a
252  // normal object with a buffer, we don't want to free the buffer because
253  // it is being used by the hierarchical objected we just created. If we
254  // created a bufferless object, we don't want to free the buffer because
255  // there was no buffer allocated in the first place.
256  FLA_Obj_free_without_buffer( &flat_matrix );
257 
258  // Free the local arrays.
259  FLA_free( elem_sizes_m );
260  FLA_free( elem_sizes_n );
261  FLA_free( depth_sizes_m );
262  FLA_free( depth_sizes_n );
263  FLA_free( m_offsets );
264  FLA_free( n_offsets );
265  }
266 
267  return FLA_SUCCESS;
268 }
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)
Definition: FLASH_Obj_create_helper_check.c:13
FLA_Error FLA_Obj_create(FLA_Datatype datatype, dim_t m, dim_t n, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition: FLA_Obj.c:55
unsigned long dim_t
Definition: FLA_type_defs.h:71
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)
Definition: FLASH_Obj.c:271
FLA_Error FLA_Obj_create_without_buffer(FLA_Datatype datatype, dim_t m, dim_t n, FLA_Obj *obj)
Definition: FLA_Obj.c:362
Definition: FLA_type_defs.h:158
void FLA_free(void *ptr)
Definition: FLA_Memory.c:247
void * FLA_malloc(size_t size)
Definition: FLA_Memory.c:111
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int i
Definition: bl1_axmyv2.c:145
FLA_Error FLA_Obj_free_without_buffer(FLA_Obj *obj)
Definition: FLA_Obj.c:615

◆ FLASH_Obj_create_hier_conf_to_flat()

FLA_Error FLASH_Obj_create_hier_conf_to_flat ( FLA_Trans  trans,
FLA_Obj  F,
dim_t  depth,
dim_t b_mn,
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().

460 {
461  FLA_Datatype datatype;
462  dim_t m_H, n_H;
463  dim_t m_F, n_F;
464 
465  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
466  FLASH_Obj_create_hier_conf_to_flat_check( trans, F, depth, b_mn, H );
467 
468  // Acquire the numerical datatype, length, and width of the flat matrix
469  // object.
470  datatype = FLA_Obj_datatype( F );
471  m_F = FLA_Obj_length( F );
472  n_F = FLA_Obj_width( F );
473 
474  // Factor in the transposition, if requested.
475  if ( trans == FLA_NO_TRANSPOSE )
476  {
477  m_H = m_F;
478  n_H = n_F;
479  }
480  else
481  {
482  m_H = n_F;
483  n_H = m_F;
484  }
485 
486  // Create the hierarchical matrix object.
487  FLASH_Obj_create( datatype, m_H, n_H, depth, b_mn, H );
488 
489  return FLA_SUCCESS;
490 }
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)
Definition: FLASH_Obj_create_hier_conf_to_flat_check.c:13
unsigned long dim_t
Definition: FLA_type_defs.h:71
FLA_Datatype FLA_Obj_datatype(FLA_Obj obj)
Definition: FLA_Query.c:13
dim_t FLA_Obj_width(FLA_Obj obj)
Definition: FLA_Query.c:123
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int FLA_Datatype
Definition: FLA_type_defs.h:49
FLA_Error FLASH_Obj_create(FLA_Datatype datatype, dim_t m, dim_t n, dim_t depth, dim_t *b_mn, FLA_Obj *H)
Definition: FLASH_Obj.c:143
dim_t FLA_Obj_length(FLA_Obj obj)
Definition: FLA_Query.c:116

◆ FLASH_Obj_create_hier_conf_to_flat_ext()

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().

494 {
495  FLA_Datatype datatype;
496  dim_t m_H, n_H;
497  dim_t m_F, n_F;
498 
499  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
500  FLASH_Obj_create_hier_conf_to_flat_ext_check( trans, F, depth, b_m, b_n, H );
501 
502  // Acquire the numerical datatype, length, and width of the flat matrix
503  // object.
504  datatype = FLA_Obj_datatype( F );
505  m_F = FLA_Obj_length( F );
506  n_F = FLA_Obj_width( F );
507 
508  // Factor in the transposition, if requested.
509  if ( trans == FLA_NO_TRANSPOSE )
510  {
511  m_H = m_F;
512  n_H = n_F;
513  }
514  else
515  {
516  m_H = n_F;
517  n_H = m_F;
518  }
519 
520  // Create the hierarchical matrix object.
521  FLASH_Obj_create_ext( datatype, m_H, n_H, depth, b_m, b_n, H );
522 
523  return FLA_SUCCESS;
524 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
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)
Definition: FLASH_Obj.c:151
FLA_Datatype FLA_Obj_datatype(FLA_Obj obj)
Definition: FLA_Query.c:13
dim_t FLA_Obj_width(FLA_Obj obj)
Definition: FLA_Query.c:123
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)
Definition: FLASH_Obj_create_hier_conf_to_flat_ext_check.c:13
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int FLA_Datatype
Definition: FLA_type_defs.h:49
dim_t FLA_Obj_length(FLA_Obj obj)
Definition: FLA_Query.c:116

◆ FLASH_Obj_create_hier_copy_of_flat()

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

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

Referenced by FLASH_CAQR_UT_inc_create_hier_matrices(), FLASH_LQ_UT_create_hier_matrices(), FLASH_LU_incpiv_create_hier_matrices(), FLASH_QR_UT_create_hier_matrices(), FLASH_QR_UT_inc_create_hier_matrices(), and FLASH_UDdate_UT_inc_create_hier_matrices().

592 {
593  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
594  FLASH_Obj_create_hier_copy_of_flat_check( F, depth, b_mn, H );
595 
596  // Create a hierarchical object conformal to the flat object.
597  FLASH_Obj_create_hier_conf_to_flat( FLA_NO_TRANSPOSE, F, depth, b_mn, H );
598 
599  // Initialize the contents of the hierarchical matrix object with the
600  // contents of the flat matrix object.
601  FLASH_Copy_flat_to_hier( F, 0, 0, *H );
602 
603  return FLA_SUCCESS;
604 }
FLA_Error FLASH_Obj_create_hier_conf_to_flat(FLA_Trans trans, FLA_Obj F, dim_t depth, dim_t *b_mn, FLA_Obj *H)
Definition: FLASH_Obj.c:459
FLA_Error FLASH_Obj_create_hier_copy_of_flat_check(FLA_Obj F, dim_t depth, dim_t *blocksizes, FLA_Obj *H)
Definition: FLASH_Obj_create_hier_copy_of_flat_check.c:13
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLASH_Copy_flat_to_hier(FLA_Obj F, dim_t i, dim_t j, FLA_Obj H)
Definition: FLASH_Copy_other.c:81

◆ FLASH_Obj_create_hier_copy_of_flat_ext()

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_flat_to_hier(), FLASH_Obj_create_hier_conf_to_flat_ext(), and FLASH_Obj_create_hier_copy_of_flat_ext_check().

608 {
609  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
610  FLASH_Obj_create_hier_copy_of_flat_ext_check( F, depth, b_m, b_n, H );
611 
612  // Create a hierarchical object conformal to the flat object.
613  FLASH_Obj_create_hier_conf_to_flat_ext( FLA_NO_TRANSPOSE, F, depth, b_m, b_n, H );
614 
615  // Initialize the contents of the hierarchical matrix object with the
616  // contents of the flat matrix object.
617  FLASH_Copy_flat_to_hier( F, 0, 0, *H );
618 
619  return FLA_SUCCESS;
620 }
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)
Definition: FLASH_Obj.c:493
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLASH_Copy_flat_to_hier(FLA_Obj F, dim_t i, dim_t j, FLA_Obj H)
Definition: FLASH_Copy_other.c:81
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)
Definition: FLASH_Obj_create_hier_copy_of_flat_ext_check.c:13

◆ FLASH_Obj_create_hierarchy()

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_at_view(), 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_check(), FLASH_Queue_set_block_size(), i, FLA_Obj_struct::id, FLA_Obj_struct::m_index, and FLA_Obj_struct::n_index.

Referenced by FLASH_Obj_create_helper().

272 {
273  dim_t i, j, b;
274  dim_t next_m, next_n;
275  dim_t num_m, num_n;
276  dim_t m_inner, n_inner;
277  dim_t elem_size_m_cur;
278  dim_t elem_size_n_cur;
279  FLA_Obj FL, FR, F0, F1, F2;
280  FLA_Obj* buffer_H;
281 
282  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
283  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 );
284 
285  if ( depth == 0 )
286  {
287  // If we're asked to create a zero-depth matrix, we interpret that as
288  // a request to create leaf-level objects using the remaining portion
289  // of the segment of the flat_matrix buffer that was passed in.
290  FLA_Obj_create_without_buffer( datatype, m, n, H );
291  FLA_Obj_attach_buffer( FLA_Obj_buffer_at_view( flat_matrix ), 1, m, H );
292 #ifdef FLA_ENABLE_SUPERMATRIX
293  FLASH_Queue_set_block_size( m * n * FLA_Obj_datatype_size( datatype ) );
294 #endif
295  H->base->id = id;
296 
297  // Fill in the m_index and n_index variables, which identify the
298  // location of the current leaf node, in units of storage blocks,
299  // within the overall matrix.
300  for ( i = 0; i < depth_overall; i++ )
301  {
302  H->base->m_index += m_offsets[i] * depth_sizes_m[i];
303  H->base->n_index += n_offsets[i] * depth_sizes_n[i];
304  }
305  }
306  else
307  {
308  // The "current" level's elem_size value. That is, the number of numerical
309  // scalar elements along one side of a full block on the current level,
310  // for the row and column dimensions.
311  elem_size_m_cur = elem_sizes_m[0];
312  elem_size_n_cur = elem_sizes_n[0];
313 
314  // Compute the number of rows and columns in the current hierarchical
315  // level of blocking.
316  num_m = m / elem_size_m_cur + ( (m % elem_size_m_cur) ? 1 : 0 );
317  num_n = n / elem_size_n_cur + ( (n % elem_size_n_cur) ? 1 : 0 );
318 
319  // The total number of scalar elements contained within/below this level
320  // of the hierarchy. (The edge cases are handled by the computation of
321  // next_m and next_n below, since they are passed in as the new m and n
322  // for the next recursive call.)
323  m_inner = m;
324  n_inner = n;
325 
326  // Create a matrix whose elements are FLA_Objs for the current level of
327  // blocking.
328  FLA_Obj_create_ext( datatype, FLA_MATRIX, num_m, num_n, m_inner, n_inner, 0, 0, H );
329 
330  if ( depth == depth_overall )
331  id = H->base->id;
332  else
333  H->base->id = id;
334 
335  // Grab the buffer from the new hierarchical object. This is an array of
336  // FLA_Objs.
337  buffer_H = ( FLA_Obj* ) FLA_Obj_buffer_at_view( *H );
338 
339  // Prepare to partition through the flat matrix so we can further allocate
340  // segments of it to the various hierarchical sub-matrices. (The second
341  // case occurs when the current function is called with a flat_matrix
342  // argument that was created without a buffer.)
343  if ( FLA_Obj_buffer_at_view( flat_matrix ) != NULL )
344  FLA_Part_1x2( flat_matrix, &FL, &FR, 0, FLA_LEFT );
345  else
346  FLA_Obj_create_without_buffer( datatype, 0, 0, &F1 );
347 
348  for ( j = 0; j < num_n; ++j )
349  {
350  // Determine the number of elements along the column dimension
351  // that will be contained within the submatrix referenced by
352  // the (i,j)th FLA_MATRIX element in the current matrix.
353  if ( j != num_n-1 || (n % elem_size_n_cur) == 0 )
354  next_n = elem_size_n_cur;
355  else
356  next_n = n % elem_size_n_cur;
357 
358  n_offsets[depth_overall-depth] = j;
359 
360  for ( i = 0; i < num_m; ++i )
361  {
362  // Determine the number of elements along the row dimension
363  // that will be contained within the submatrix referenced by
364  // the (i,j)th FLA_MATRIX element in the current matrix.
365  if ( i != num_m-1 || (m % elem_size_m_cur) == 0 )
366  next_m = elem_size_m_cur;
367  else
368  next_m = m % elem_size_m_cur;
369 
370  m_offsets[depth_overall-depth] = i;
371 
372  // Partition the next m*n elements from the flat matrix so we can
373  // "attach" them to the hierarchical matrices contained within the
374  // (i,j)th FLA_MATRIX object.
375  if ( FLA_Obj_buffer_at_view( flat_matrix ) != NULL )
376  {
377  b = min( FLA_Obj_width( FR ), next_m * next_n );
378  FLA_Repart_1x2_to_1x3( FL, /**/ FR, &F0, /**/ &F1, &F2,
379  b, FLA_RIGHT );
380  }
381 
382  // Recursively call ourselves, with the appropriate parameters for
383  // the next deeper level in the matrix hierarchy.
384  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 );
385 
386  // Continue with the repartitioning.
387  if ( FLA_Obj_buffer_at_view( flat_matrix ) != NULL )
388  {
389  FLA_Cont_with_1x3_to_1x2( &FL, /**/ &FR, F0, F1, /**/ F2,
390  FLA_LEFT );
391  }
392  }
393  }
394 
395  // Free the temporary flat matrix subpartition object, but only if it was
396  // created to begin with. Since it would have been created without a
397  // buffer, we must free it in a similar manner.
398  if ( FLA_Obj_buffer_at_view( flat_matrix ) == NULL )
400  }
401 
402  return FLA_SUCCESS;
403 }
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)
Definition: FLASH_Obj_create_hierarchy_check.c:13
FLA_Error FLA_Repart_1x2_to_1x3(FLA_Obj AL, FLA_Obj AR, FLA_Obj *A0, FLA_Obj *A1, FLA_Obj *A2, dim_t nb, FLA_Side side)
Definition: FLA_View.c:267
unsigned long dim_t
Definition: FLA_type_defs.h:71
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)
Definition: FLASH_Obj.c:271
FLA_Error FLA_Obj_create_ext(FLA_Datatype datatype, FLA_Elemtype elemtype, dim_t m, dim_t n, dim_t m_inner, dim_t n_inner, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition: FLA_Obj.c:64
FLA_Base_obj * base
Definition: FLA_type_defs.h:168
dim_t n_index
Definition: FLA_type_defs.h:135
FLA_Error FLA_Obj_attach_buffer(void *buffer, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition: FLA_Obj.c:522
unsigned long id
Definition: FLA_type_defs.h:133
FLA_Error FLA_Obj_create_without_buffer(FLA_Datatype datatype, dim_t m, dim_t n, FLA_Obj *obj)
Definition: FLA_Obj.c:362
void * FLA_Obj_buffer_at_view(FLA_Obj obj)
Definition: FLA_Query.c:215
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_width(FLA_Obj obj)
Definition: FLA_Query.c:123
FLA_Error FLA_Cont_with_1x3_to_1x2(FLA_Obj *AL, FLA_Obj *AR, FLA_Obj A0, FLA_Obj A1, FLA_Obj A2, FLA_Side side)
Definition: FLA_View.c:475
void FLASH_Queue_set_block_size(dim_t size)
Definition: FLASH_Queue.c:461
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLA_Part_1x2(FLA_Obj A, FLA_Obj *A1, FLA_Obj *A2, dim_t nb, FLA_Side side)
Definition: FLA_View.c:110
int i
Definition: bl1_axmyv2.c:145
FLA_Error FLA_Obj_free_without_buffer(FLA_Obj *obj)
Definition: FLA_Obj.c:615
dim_t FLA_Obj_datatype_size(FLA_Datatype datatype)
Definition: FLA_Query.c:61
dim_t m_index
Definition: FLA_type_defs.h:134

◆ FLASH_Obj_create_without_buffer()

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

References FLASH_Obj_create_helper().

160 {
161  FLASH_Obj_create_helper( TRUE, datatype, m, n, depth, b_mn, b_mn, H );
162 
163  return FLA_SUCCESS;
164 }
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)
Definition: FLASH_Obj.c:175

◆ FLASH_Obj_create_without_buffer_ext()

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().

Referenced by FLASH_Part_create_1x2(), FLASH_Part_create_2x1(), and FLASH_Part_create_2x2().

168 {
169  FLASH_Obj_create_helper( TRUE, datatype, m, n, depth, b_m, b_n, H );
170 
171  return FLA_SUCCESS;
172 }
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)
Definition: FLASH_Obj.c:175

◆ FLASH_Obj_datatype()

FLA_Datatype FLASH_Obj_datatype ( FLA_Obj  H)

◆ FLASH_Obj_depth()

dim_t FLASH_Obj_depth ( FLA_Obj  H)

References FLA_Obj_base_buffer(), and FLA_Obj_elemtype().

Referenced by FLASH_Apply_CAQ_UT_inc_create_workspace(), FLASH_Apply_pivots(), FLASH_Apply_Q_UT_create_workspace(), FLASH_Apply_Q_UT_inc_create_workspace(), FLASH_Apply_QUD_UT_inc_create_workspace(), FLASH_FS_incpiv(), FLASH_LQ_UT(), FLASH_LU_incpiv(), FLASH_LU_piv(), FLASH_Obj_create_conf_to(), FLASH_Part_create_1x2(), FLASH_Part_create_2x1(), FLASH_Part_create_2x2(), and FLASH_QR_UT().

21 {
22  FLA_Elemtype elemtype;
23  FLA_Obj* buffer_H;
24  dim_t depth = 0;
25 
26  // Recurse through the hierarchy to the first leaf node. We initialize
27  // the recursion here:
28  elemtype = FLA_Obj_elemtype( H );
29  buffer_H = ( FLA_Obj* ) FLA_Obj_base_buffer( H );
30 
31  while ( elemtype == FLA_MATRIX )
32  {
33  ++depth;
34 
35  // Get the element type of the top-leftmost underlying object. Also,
36  // get a pointer to the first element of the top-leftmost object and
37  // assume that it is of type FLA_Obj* in case elemtype is once again
38  // FLA_MATRIX.
39  elemtype = FLA_Obj_elemtype( buffer_H[0] );
40  buffer_H = ( FLA_Obj * ) FLA_Obj_base_buffer( buffer_H[0] );
41  }
42 
43  // At this point, the value of depth represents the depth of the matrix
44  // hierarchy.
45  return depth;
46 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
int FLA_Elemtype
Definition: FLA_type_defs.h:50
void * FLA_Obj_base_buffer(FLA_Obj obj)
Definition: FLA_Query.c:210
Definition: FLA_type_defs.h:158
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_extract_buffer()

void* FLASH_Obj_extract_buffer ( FLA_Obj  H)

References FLA_Obj_base_buffer(), and FLA_Obj_elemtype().

Referenced by FLASH_Obj_free().

742 {
743  FLA_Elemtype elemtype;
744  FLA_Obj* buffer_H;
745 
746  // Recurse through the hierarchy to the first leaf node to gain access
747  // to the address of the actual numerical data buffer (ie: the "flat"
748  // matrix object used in FLASH_Obj_create()). We initialize the search
749  // here:
750  elemtype = FLA_Obj_elemtype( H );
751  buffer_H = ( FLA_Obj* ) FLA_Obj_base_buffer( H );
752 
753  while ( elemtype == FLA_MATRIX )
754  {
755  elemtype = FLA_Obj_elemtype( buffer_H[0] );
756  buffer_H = ( FLA_Obj* ) FLA_Obj_base_buffer( buffer_H[0] );
757  }
758 
759  // At this point, the value in buffer_H is a pointer to the array that
760  // holds the numerical data.
761  return ( void* ) buffer_H;
762 }
int FLA_Elemtype
Definition: FLA_type_defs.h:50
void * FLA_Obj_base_buffer(FLA_Obj obj)
Definition: FLA_Query.c:210
Definition: FLA_type_defs.h:158
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_flatten()

FLA_Error FLASH_Obj_flatten ( FLA_Obj  H,
FLA_Obj  F 
)

References FLASH_Copy_hier_to_flat().

766 {
767  FLASH_Copy_hier_to_flat( 0, 0, H, F );
768 
769  return FLA_SUCCESS;
770 }
FLA_Error FLASH_Copy_hier_to_flat(dim_t i, dim_t j, FLA_Obj H, FLA_Obj F)
Definition: FLASH_Copy_other.c:110

◆ FLASH_Obj_free()

void FLASH_Obj_free ( FLA_Obj H)

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

Referenced by FLASH_CAQR_UT_inc_solve(), FLASH_Eig_gest(), FLASH_LQ_UT_solve(), FLASH_LU_incpiv_opt1(), FLASH_QR_UT_inc_opt1(), FLASH_QR_UT_inc_solve(), FLASH_QR_UT_solve(), FLASH_Random_spd_matrix(), and FLASH_UDdate_UT_inc_update_rhs().

639 {
640  FLA_Obj* buffer_H;
641 
642  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
644 
645  // Free the object according to whether it contains a hierarchy.
646  if ( FLA_Obj_elemtype( *H ) == FLA_MATRIX )
647  {
648  // Extract a pointer to the data buffer that was parititioned across all
649  // leaf-level submatrices.
650  buffer_H = ( FLA_Obj * ) FLASH_Obj_extract_buffer( *H );
651 
652  // Free the data buffer. This works because FLASH_Obj_extract_buffer()
653  // returns the starting address of the first element's buffer, which is
654  // also the starting address of the entire buffer.
655 #ifdef FLA_ENABLE_SCC
656  FLA_shfree( buffer_H );
657 #else
658  FLA_free( buffer_H );
659 #endif
660 
661  // All that remains now is to free the interior of the matrix hierarchy.
662  // This includes non-leaf buffers and their corresponding base objects
663  // as well as leaf-level base objects.
665  }
666  else
667  {
668  // If the matrix has no hierarchy, treat it like a flat object.
669  FLA_Obj_free( H );
670  }
671 }
FLA_Error FLA_Obj_free(FLA_Obj *obj)
Definition: FLA_Obj.c:588
void FLA_shfree(void *ptr)
Definition: FLA_Obj.c:27
Definition: FLA_type_defs.h:158
void FLA_free(void *ptr)
Definition: FLA_Memory.c:247
FLA_Error FLASH_Obj_free_check(FLA_Obj *H)
Definition: FLASH_Obj_free_check.c:13
void * FLASH_Obj_extract_buffer(FLA_Obj H)
Definition: FLASH_Obj.c:741
void FLASH_Obj_free_hierarchy(FLA_Obj *H)
Definition: FLASH_Obj.c:699
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_free_hierarchy()

void FLASH_Obj_free_hierarchy ( FLA_Obj H)

References FLA_Check_error_level(), FLA_Obj_base_buffer(), FLA_Obj_elemtype(), FLA_Obj_free(), FLA_Obj_free_without_buffer(), FLA_Obj_num_elem_alloc(), FLASH_Obj_free_hierarchy_check(), and i.

Referenced by FLASH_Obj_free(), and FLASH_Obj_free_without_buffer().

700 {
701  //dim_t m_H, n_H, rs, cs, i, j;
702  dim_t i;
703  dim_t n_elem_alloc;
704  FLA_Obj* buffer_H;
705 
706  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
708 
709  // If the element type of H is FLA_SCALAR, then it has no children to
710  // free, so free the base object. In order to avoid freeing the object's
711  // data buffer, which would have already been freed en masse by now if
712  // the calling function is FLASH_Obj_free(), we will call
713  // FLA_Obj_free_without_buffer().
714  if ( FLA_Obj_elemtype( *H ) == FLA_SCALAR )
715  {
717  return;
718  }
719  else
720  {
721  // Acquire the number of elements allocated when this node was
722  // created.
723  n_elem_alloc = FLA_Obj_num_elem_alloc( *H );
724 
725  // Acquire the array of objects contained inside of H.
726  buffer_H = ( FLA_Obj* ) FLA_Obj_base_buffer( *H );
727 
728  // For each allocated submatrix in H...
729  for ( i = 0; i < n_elem_alloc; ++i )
730  {
731  // Recurse with the ith element of the allocated buffer.
732  FLASH_Obj_free_hierarchy( &buffer_H[i] );
733  }
734 
735  // Finally, free the internal array of objects.
736  FLA_Obj_free( H );
737  }
738 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
FLA_Error FLA_Obj_free(FLA_Obj *obj)
Definition: FLA_Obj.c:588
dim_t FLA_Obj_num_elem_alloc(FLA_Obj obj)
Definition: FLA_Query.c:204
FLA_Error FLASH_Obj_free_hierarchy_check(FLA_Obj *H)
Definition: FLASH_Obj_free_hierachy_check.c:13
void * FLA_Obj_base_buffer(FLA_Obj obj)
Definition: FLA_Query.c:210
Definition: FLA_type_defs.h:158
void FLASH_Obj_free_hierarchy(FLA_Obj *H)
Definition: FLASH_Obj.c:699
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
int i
Definition: bl1_axmyv2.c:145
FLA_Error FLA_Obj_free_without_buffer(FLA_Obj *obj)
Definition: FLA_Obj.c:615
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_free_without_buffer()

void FLASH_Obj_free_without_buffer ( FLA_Obj H)

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

Referenced by FLASH_Part_free_1x2(), FLASH_Part_free_2x1(), and FLASH_Part_free_2x2().

675 {
676  if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
678 
679  // Free the object according to whether it contains a hierarchy.
680  if ( FLA_Obj_elemtype( *H ) == FLA_MATRIX )
681  {
682  // Skip freeing the numerical data buffer, since the object was
683  // presumably created with FLASH_Obj_create_without_buffer().
684 
685  // Free the interior of the matrix hierarchy. This includes non-leaf
686  // buffers and their corresponding base objects as well as leaf-level
687  // base objects.
689  }
690  else
691  {
692  // If the matrix has no hierarchy, treat it like a flat object with
693  // no internal data buffer.
695  }
696 }
void FLASH_Obj_free_hierarchy(FLA_Obj *H)
Definition: FLASH_Obj.c:699
unsigned int FLA_Check_error_level(void)
Definition: FLA_Check.c:18
FLA_Error FLA_Obj_free_without_buffer(FLA_Obj *obj)
Definition: FLA_Obj.c:615
FLA_Error FLASH_Obj_free_without_buffer_check(FLA_Obj *H)
Definition: FLASH_Obj_free_without_buffer_check.c:13
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_Obj_hierarchify()

FLA_Error FLASH_Obj_hierarchify ( FLA_Obj  F,
FLA_Obj  H 
)

References FLASH_Copy_flat_to_hier().

Referenced by FLASH_Hermitianize(), FLASH_Obj_create_copy_of(), FLASH_Random_matrix(), FLASH_Random_spd_matrix(), FLASH_Set(), FLASH_Shift_diag(), and FLASH_Triangularize().

774 {
775  FLASH_Copy_flat_to_hier( F, 0, 0, H );
776 
777  return FLA_SUCCESS;
778 }
FLA_Error FLASH_Copy_flat_to_hier(FLA_Obj F, dim_t i, dim_t j, FLA_Obj H)
Definition: FLASH_Copy_other.c:81

◆ FLASH_print_struct()

void FLASH_print_struct ( FLA_Obj  H)

References FLA_Obj_buffer_at_view(), FLA_Obj_col_stride(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_row_stride(), FLA_Obj_width(), FLASH_print_struct_helper(), and i.

918 {
919  dim_t m_H, n_H, rs, cs, i, j;
920  FLA_Obj* buffer_temp;
921 
922  m_H = FLA_Obj_length( H );
923  n_H = FLA_Obj_width( H );
924  rs = FLA_Obj_row_stride( H );
925  cs = FLA_Obj_col_stride( H );
926 
927  if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
929  else
930  {
931  for ( j = 0; j < n_H; ++j )
932  {
933  for ( i = 0; i < m_H; ++i )
934  {
935  buffer_temp = ( FLA_Obj* ) FLA_Obj_buffer_at_view( H );
936 
937  FLASH_print_struct_helper( buffer_temp[ j*cs + i*rs ], 0 );
938  }
939  }
940  }
941 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
dim_t FLA_Obj_row_stride(FLA_Obj obj)
Definition: FLA_Query.c:167
void * FLA_Obj_buffer_at_view(FLA_Obj obj)
Definition: FLA_Query.c:215
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_width(FLA_Obj obj)
Definition: FLA_Query.c:123
void FLASH_print_struct_helper(FLA_Obj H, int indent)
Definition: FLASH_Obj.c:944
dim_t FLA_Obj_col_stride(FLA_Obj obj)
Definition: FLA_Query.c:174
int i
Definition: bl1_axmyv2.c:145
dim_t FLA_Obj_length(FLA_Obj obj)
Definition: FLA_Query.c:116
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51

◆ FLASH_print_struct_helper()

void FLASH_print_struct_helper ( FLA_Obj  H,
int  indent 
)

References FLA_Obj_buffer_at_view(), FLA_Obj_col_stride(), FLA_Obj_datatype(), FLA_Obj_elemtype(), FLA_Obj_length(), FLA_Obj_row_stride(), FLA_Obj_width(), and i.

Referenced by FLASH_print_struct().

945 {
946  dim_t m_H, n_H, rs, cs, i, j, k;
947  FLA_Obj* buffer_temp;
948 
949  for ( i = 0; i < indent; ++i )
950  fprintf( stdout, " " );
951 
952  if ( FLA_Obj_elemtype( H ) == FLA_SCALAR )
953  {
954  fprintf( stdout, "LEAF (%3d | rs %3lu | cs %3lu | %3lu x %3lu | addr %p)\n",
955  FLA_Obj_datatype( H ),
957  FLA_Obj_length( H ), FLA_Obj_width( H ),
958  FLA_Obj_buffer_at_view( H ) );
959  fflush( stdout );
960  }
961  else
962  {
963  m_H = FLA_Obj_length( H );
964  n_H = FLA_Obj_width( H );
965  rs = FLA_Obj_row_stride( H );
966  cs = FLA_Obj_col_stride( H );
967 
968  fprintf( stdout, "MATRIX (%lux%lu):%d - %p\n",
969  m_H, n_H,
970  FLA_Obj_datatype( H ),
971  FLA_Obj_buffer_at_view( H ) );
972  fflush( stdout );
973 
974  for ( j = 0; j < n_H; ++j )
975  {
976  for ( i = 0; i < m_H; ++i )
977  {
978  for ( k = 0; k < indent; ++k )
979  fprintf( stdout, " " );
980 
981  buffer_temp = ( FLA_Obj* ) FLA_Obj_buffer_at_view( H );
982 
983  FLASH_print_struct_helper( buffer_temp[ j*cs + i*rs ],
984  indent + 1 );
985  }
986  }
987  }
988 }
unsigned long dim_t
Definition: FLA_type_defs.h:71
dim_t FLA_Obj_row_stride(FLA_Obj obj)
Definition: FLA_Query.c:167
FLA_Datatype FLA_Obj_datatype(FLA_Obj obj)
Definition: FLA_Query.c:13
void * FLA_Obj_buffer_at_view(FLA_Obj obj)
Definition: FLA_Query.c:215
Definition: FLA_type_defs.h:158
dim_t FLA_Obj_width(FLA_Obj obj)
Definition: FLA_Query.c:123
void FLASH_print_struct_helper(FLA_Obj H, int indent)
Definition: FLASH_Obj.c:944
dim_t FLA_Obj_col_stride(FLA_Obj obj)
Definition: FLA_Query.c:174
int i
Definition: bl1_axmyv2.c:145
dim_t FLA_Obj_length(FLA_Obj obj)
Definition: FLA_Query.c:116
FLA_Elemtype FLA_Obj_elemtype(FLA_Obj obj)
Definition: FLA_Query.c:51