[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: Help with [] for defined types
From: |
David Bateman |
Subject: |
Re: Help with [] for defined types |
Date: |
Tue, 3 Aug 2004 15:29:20 +0200 |
User-agent: |
Mutt/1.4.1i |
Hi Andy,
According to Andy Adler <address@hidden> (on 08/02/04):
> I've done this, and it seems to be being called twice, the first time
> with dv=[0,0] and the second time with the right output size, but
>
> echo 's=sparse(diag([1,3],1)); p=[s,s]' | octave -qfH
> DEBUG:sparse( void)
> installing sparse type at type-id = 40
> install sm_sm
> DEBUG:sparse - matrix_to_sparse
> DEBUG:sparse( SuperMatrix A)
> DEBUG:sparse_resize <--- calls sparse([1;2],[2;3],[1;3],0,0)
> DEBUG:sparse - assemble_sparse
> error: sparse row index out of range <-- error because rows,cols= 0,0
> DEBUG:sparse( SuperMatrix A)
> DEBUG:sparse_resize <--- calls sparse([],[],[],3,6)
> DEBUG:sparse( SuperMatrix A)
> DEBUG:sparse destructor
> DEBUG:sparse destructor
> error: evaluating assignment expression near line 1, column 28
>
> I don't understand why octave_sparse::resize is being called twice here
> Am I missing something?
The explanation is in the comment in tree_matrix::rvalue, which states
// The line below might seem crazy, since we take a copy
// of the first argument, resize it to be empty and then resize
// it to be full. This is done since it means that there is no
// recopying of data, as would happen if we used a single resize.
// It should be noted that resize operation is also significantly
// slower than the do_cat_op function, so it makes sense to have an
// empty matrix and copy all data.
//
// We might also start with a empty octave_value using
// ctmp = octave_value_typeinfo::lookup_type
// (tmp.begin() -> begin() -> type_name());
// and then directly resize. However, for some types there might be
// some additional setup needed, and so this should be avoided.
The error comes about because, resize is supposed to be destructive. So
you should modifiy your ASSEMBLE_SPARSE macro to be something like
#define ASSEMBLE_SPARSE2( TYPX, RESIZE) \
int nnz= MAX( ridxA.length(), cidxA.length() ); \
TYPX* coefX= (TYPX*)oct_sparse_malloc( nnz * sizeof(TYPX)); \
int * ridxX= (int *)oct_sparse_malloc( nnz * sizeof(int) ); \
int * cidxX= (int *)oct_sparse_malloc( (n+1)* sizeof(int)); cidxX[0]= 0; \
\
bool ri_scalar = (ridxA.length() == 1); \
bool ci_scalar = (cidxA.length() == 1); \
bool cf_scalar = (coefA.length() == 1); \
\
sort_idxl* sidx = (sort_idxl *)oct_sparse_malloc( nnz* sizeof(sort_idxl) );\
int actual_nnz=0; \
OCTAVE_QUIT; \
for (int i=0; i<nnz; i++) { \
double rowidx= ( ri_scalar ? ridxA(0) : ridxA(i) ) - 1; \
double colidx= ( ci_scalar ? cidxA(0) : cidxA(i) ) - 1; \
if (rowidx >= m || rowidx < 0) \
if (! RESIZE ) \
error("sparse row index out of range"); \
else if (colidx >= n || colidx < 0) \
if (! RESIZE ) \
error("sparse column index out of range"); \
if (error_state) { actual_nnz=0; break; } \
if ( coefA(cf_scalar?0:i) != 0. ) { \
sidx[actual_nnz].val = (long) ( rowidx + m*colidx ); \
sidx[actual_nnz].idx = i; \
actual_nnz++; \
} \
} \
<snip the rest of the old ASSEMBLE_SPARSE>
#define ASSEMBLE_SPARSE( TYPX ) ASSEMBLE_SPARSE2 (TYPX, false)
You can then call ASSEMBLE_SPARSE2 in the resize operations
> I've implemented these, but they don't seem to be being called.
> Definitely, octave_sparse concat (const octave_sparse& ra,
> const octave_sparse& rb, const Array<int>& ra_idx);
> is not called.
>
> Also, I'm not sure I understand what ra_idx is and how I'm supposed to
> interpret it.
As the failing resize operations in tree_matrix::rvalue is followed be
if (error_state)
goto done;
This is not surprising. Fix up resize so that it is destructive and you'll
probably get to the code that calls concat. Then another set of bugs will
probably appear :-)
Cheers
David
--
David Bateman address@hidden
Motorola CRM +33 1 69 35 48 04 (Ph)
Parc Les Algorithmes, Commune de St Aubin +33 1 69 35 77 01 (Fax)
91193 Gif-Sur-Yvette FRANCE
The information contained in this communication has been classified as:
[x] General Business Information
[ ] Motorola Internal Use Only
[ ] Motorola Confidential Proprietary