1 /*
2
3 Copyright (C) 2008-2021 Michele Martone
4
5 This file is part of librsb.
6
7 librsb is free software; you can redistribute it and/or modify it
8 under the terms of the GNU Lesser General Public License as published
9 by the Free Software Foundation; either version 3 of the License, or
10 (at your option) any later version.
11
12 librsb is distributed in the hope that it will be useful, but WITHOUT
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15 License for more details.
16
17 You should have received a copy of the GNU Lesser General Public
18 License along with librsb; see the file COPYING.
19 If not, see <http://www.gnu.org/licenses/>.
20
21 */
22 /* @cond INNERDOC */
23 /**
24 * @file
25 * @author Michele Martone
26 * @brief
27 *
28 * Matrix cloning biotech.
29 * \internal
30 *
31 * */
32 #include "rsb_common.h"
33
34 //#define RSB_MTX_REASSIGN(OLD_MTXP,NEW_MTXP) {if(rsb_do_assign(NEW_MTXP,OLD_MTXP)) {RSB_PERR_GOTO(err,RSB_ERRM_ES);}
35 #define RSB_MTX_REASSIGN(OLD_MTXP,NEW_MTXP) { RSB_MTX_FREE(OLD_MTXP); (OLD_MTXP)=(NEW_MTXP); }
36
37 struct rsb_session_handle_t rsb_global_session_handle;
38
rsb__clone_area_with_extra(const void * src,size_t csize,size_t bsize,size_t esize)39 void * rsb__clone_area_with_extra(const void *src, size_t csize, size_t bsize, size_t esize)
40 {
41 /*!
42 * \ingroup gr_internals
43 * (m)allocates an area of bsize+csize+esize bytes and copies there csize bytes from src, at offset bsize
44 * */
45 rsb_byte_t * dst = NULL;
46
47 if(!src /* || esize < 0 */)
48 goto ret;
49 dst = rsb__malloc(csize+bsize+esize);
50 if(!dst)
51 goto ret;
52 rsb__memcpy(dst+bsize,src,csize);
53 ret:
54 return dst;
55 }
56
rsb__clone_area(const void * src,size_t size)57 void * rsb__clone_area(const void *src, size_t size)
58 {
59 /*!
60 * \ingroup gr_internals
61 * \param src the source data pointer
62 * \param size the amount of data to clone
63 * \return the cloned amount, or NULL in case of error
64 *
65 * allocates an area of size bytes and copies there data from src
66 * */
67 void * dst = NULL;
68
69 if(!src || size < 1)
70 goto ret;
71 dst = rsb__clone_area_with_extra(src,size,0,0);
72 ret:
73 return dst;
74 }
75
rsb__util_coo_alloc(void ** RSB_RESTRICT VAp,rsb_coo_idx_t ** RSB_RESTRICT IAp,rsb_coo_idx_t ** RSB_RESTRICT JAp,rsb_nnz_idx_t nnz,rsb_type_t typecode,rsb_bool_t do_calloc)76 rsb_err_t rsb__util_coo_alloc(void **RSB_RESTRICT VAp, rsb_coo_idx_t ** RSB_RESTRICT IAp, rsb_coo_idx_t ** RSB_RESTRICT JAp, rsb_nnz_idx_t nnz, rsb_type_t typecode, rsb_bool_t do_calloc)
77 {
78 rsb_err_t errval = RSB_ERR_NO_ERROR;
79 void *VA_ = NULL,*IA_ = NULL,*JA_ = NULL;
80
81 if( RSB_MATRIX_UNSUPPORTED_TYPE(typecode) )
82 {
83 errval = RSB_ERR_UNSUPPORTED_TYPE;
84 goto err;
85 }
86
87 if(do_calloc != RSB_BOOL_TRUE)
88 {
89 VA_ = rsb__malloc_vector((nnz),typecode),
90 IA_ = rsb__malloc(sizeof(rsb_coo_idx_t)*(nnz)),
91 JA_ = rsb__malloc(sizeof(rsb_coo_idx_t)*(nnz));
92 }
93 else
94 {
95 VA_ = rsb__calloc_vector((nnz),typecode),
96 IA_ = rsb__calloc(sizeof(rsb_coo_idx_t)*(nnz)),
97 JA_ = rsb__calloc(sizeof(rsb_coo_idx_t)*(nnz));
98 }
99
100 if(!VA_ || !IA_ || !JA_)
101 {
102 errval = RSB_ERR_ENOMEM;
103 RSB_PERR_GOTO(err,RSB_ERRM_ES);
104 }
105
106 *VAp = VA_;
107 *IAp = IA_;
108 *JAp = JA_;
109 goto done;
110 err:
111 RSB_CONDITIONAL_FREE(IA_);
112 RSB_CONDITIONAL_FREE(JA_);
113 RSB_CONDITIONAL_FREE(VA_);
114 done:
115 return errval;
116 }
117
rsb__util_coo_alloc_copy_and_stats(void ** RSB_RESTRICT VAp,rsb_coo_idx_t ** RSB_RESTRICT IAp,rsb_coo_idx_t ** RSB_RESTRICT JAp,const void * RSB_RESTRICT VA,const rsb_coo_idx_t * RSB_RESTRICT IA,const rsb_coo_idx_t * RSB_RESTRICT JA,rsb_coo_idx_t * RSB_RESTRICT mp,rsb_coo_idx_t * RSB_RESTRICT kp,rsb_nnz_idx_t nnz,rsb_nnz_idx_t ennz,rsb_type_t typecode,const rsb_coo_idx_t offi,const rsb_coo_idx_t offo,rsb_flags_t iflags,rsb_flags_t * RSB_RESTRICT flagsp)118 rsb_err_t rsb__util_coo_alloc_copy_and_stats(void **RSB_RESTRICT VAp, rsb_coo_idx_t ** RSB_RESTRICT IAp, rsb_coo_idx_t ** RSB_RESTRICT JAp, const void *RSB_RESTRICT VA, const rsb_coo_idx_t * RSB_RESTRICT IA, const rsb_coo_idx_t * RSB_RESTRICT JA, rsb_coo_idx_t*RSB_RESTRICT mp, rsb_coo_idx_t*RSB_RESTRICT kp, rsb_nnz_idx_t nnz, rsb_nnz_idx_t ennz, rsb_type_t typecode, const rsb_coo_idx_t offi, const rsb_coo_idx_t offo, rsb_flags_t iflags, rsb_flags_t*RSB_RESTRICT flagsp)
119 {
120 /*!
121 * Copies contents of a COO arrays triple to a freshly allocated COO arrays triple.
122 * Size is assumed to be nnz+ennz.
123 * Last ennz elements are not zeroed.
124 *
125 * Flags are determined: RSB_FLAG_UPPER_TRIANGULAR, RSB_FLAG_LOWER_TRIANGULAR.
126 *
127 * TODO: May implement input sanitization or zeroes detection.
128 * TODO: Check for nnz+ennz overflow.
129 */
130 rsb_err_t errval = RSB_ERR_NO_ERROR;
131 void *VA_ = NULL;
132 rsb_coo_idx_t *IA_ = NULL,*JA_ = NULL;
133
134 errval = rsb__util_coo_alloc((void**)(&VA_),&IA_,&JA_,nnz+ennz,typecode,RSB_BOOL_FALSE);
135 if(RSB_SOME_ERROR(errval))
136 goto err;
137
138 if(!VA && !IA && !JA)
139 goto nocopy; /* it's ok: alloc only semantics */
140 /* TODO: the following shall be made parallel */
141 if(mp || kp)
142 errval = rsb_util_coo_copy_and_stats(VA,IA,JA,VA_,IA_,JA_,mp,kp,nnz,typecode,offi,offo,iflags,flagsp);
143 else
144 {
145 errval = rsb_util_coo_copy(VA,IA,JA,VA_,IA_,JA_,nnz,typecode,offi,offo);
146 /* ... flags may not always be desired! */
147 /* if(flagsp)
148 (*flagsp)|=rsb__util_coo_determine_uplo_flags(IA_,JA_,nnz);*/
149 }
150 nocopy:
151 *VAp = VA_;
152 *IAp = IA_;
153 *JAp = JA_;
154 goto done;
155 err:
156 RSB_CONDITIONAL_FREE(IA_);
157 RSB_CONDITIONAL_FREE(JA_);
158 RSB_CONDITIONAL_FREE(VA_);
159 done:
160 return errval;
161 }
162
rsb__clone_area_parallel(const void * src,size_t size,size_t n)163 void * rsb__clone_area_parallel(const void *src, size_t size, size_t n)
164 {
165 /*!
166 * \ingroup gr_internals
167 * \param src the source data pointer
168 * \param size the amount of data to clone
169 * \return the cloned amount, or NULL in case of error
170 *
171 * allocates an area of size bytes and copies there data from src
172 * */
173 void * dst = NULL;
174
175 if(!src || size < 1)
176 goto ret;
177 dst = rsb__malloc(size*n);
178 if(!dst)
179 goto ret;
180 RSB_A_MEMCPY_parallel(dst,src,0,0,n,size);
181 ret:
182 return dst;
183 }
184
185
186 #if RSB_WANT_BITMAP
rsb__clone_options_t(const struct rsb_options_t * o,rsb_blk_idx_t M_b,rsb_blk_idx_t K_b)187 static void * rsb__clone_options_t(const struct rsb_options_t *o, rsb_blk_idx_t M_b, rsb_blk_idx_t K_b)
188 {
189 /*!
190 * \ingroup gr_internals
191 * \param o a valid option structure pointer
192 * \return the input pointer
193 *
194 * clones a rsb_options_t structure, deeply
195 *
196 * p.s.: the input rsb_options_t could be NULL. in that case it won't be cloned because there will be no need of doing so.
197 * */
198 struct rsb_options_t *no = NULL;
199
200 if(!o)
201 {RSB_PERR_GOTO(err,RSB_ERRM_ES);}
202
203 /* we allocate a new options structure */
204 if(! (no = rsb__clone_area(o,sizeof(*no))))
205 {RSB_PERR_GOTO(err,RSB_ERRM_ES);}
206
207 if( o->bitmap)
208 {
209 no->bitmap = rsb__clone_area(o->bitmap,RSB_BYTES_PER_BITMAP( M_b,K_b));
210 if(!no->bitmap)
211 {RSB_PERR_GOTO(err,RSB_ERRM_ES);}
212 }
213 return no;
214
215 err:
216 rsb__destroy_options_t(no);
217 return NULL;
218 }
219 #endif /* RSB_WANT_BITMAP */
220
221 #define RSB_ARE_MTX_COO_CONFORMANT(MTXAP,MTXBP) \
222 ( ( (MTXAP)->nnz == (MTXBP)->nnz ) && ( (MTXAP)->nr == (MTXBP)->nr ) && ( (MTXAP)->nc == (MTXBP)->nc ) && ( (MTXAP)->typecode == (MTXBP)->typecode ) )
223
rsb__mtx_shift_leaf_ptrs(struct rsb_mtx_t * RSB_RESTRICT mtxCp,const struct rsb_mtx_t * RSB_RESTRICT mtxAp,long smc)224 rsb_err_t rsb__mtx_shift_leaf_ptrs(struct rsb_mtx_t *RSB_RESTRICT mtxCp, const struct rsb_mtx_t *RSB_RESTRICT mtxAp, long smc)
225 {
226 /*
227 * Adjusts pointers displacements in the matrix tree.
228 * Please note that no pointer is being accessed.
229 */
230
231 rsb_err_t errval = RSB_ERR_NO_ERROR;
232 rsb_submatrix_idx_t n,si;
233
234 for(n=0;n<smc;++n)
235 if(mtxCp[n].nnz) /* If valid. FIXME: IF NOT (E.G. MERGED), SHALL BE COMPLETELY ZEROED. */
236 for(si=0;si<RSB_FOUR;++si)
237 if(mtxAp[n].sm[si])
238 {
239 RSB_PTR_SHIFT(mtxCp[n].sm[si],mtxAp,mtxCp,(struct rsb_mtx_t*));
240 //mtxCp[n].sm[si] = mtxCp+(mtxAp[n].sm[si]-mtxAp);
241 /* RSB_STDOUT("%03d/%03d: %p\n",n,si,mtxCp[n].sm[si]); */
242 }
243 return errval;
244 }
245
rsb__mtx_transplant_from_clone(struct rsb_mtx_t ** mtxDpp,struct rsb_mtx_t * mtxSp)246 rsb_err_t rsb__mtx_transplant_from_clone(struct rsb_mtx_t ** mtxDpp, struct rsb_mtx_t * mtxSp)
247 {
248 /*
249 Moves the inner content of mtxSp to mtxDp.
250 Shall free mtxSp at the end and not change the outer allocation status of mtxSp.
251 Shall work even if any of the two matrices is in-place.
252 Can only work if matrices match in size (nonzeroes, rows, columns, ...).
253
254 Untested.
255 * */
256 rsb_err_t errval = RSB_ERR_NO_ERROR;
257 struct rsb_mtx_t * mtxDp = NULL;
258 struct rsb_mtx_t *fsmS = NULL;
259 struct rsb_mtx_t *fsmD = NULL;
260 void * VS = NULL,* VD = NULL;
261 rsb_coo_idx_t * IS = NULL, *JS = NULL;
262 rsb_coo_idx_t * ID = NULL, *JD = NULL;
263 rsb_long_t smc = 0;
264
265 if( ( !mtxDpp) || ( !*mtxDpp) || (!mtxSp) )
266 {
267 errval = RSB_ERR_BADARGS;
268 RSB_PERR_GOTO(err, RSB_ERRM_E_MTXAP);
269 }
270
271 mtxDp = *mtxDpp;
272
273 if( ! ( RSB_ARE_MTX_COO_CONFORMANT( mtxSp, mtxDp ) ) )
274 {
275 errval = RSB_ERR_BADARGS;
276 RSB_PERR_GOTO(err, RSB_ERRM_ES);
277 }
278
279 #if 0
280 RSB_STDOUT(" ==== CLONING: ==== \n");
281 RSB_STDOUT("will transplant: \n");
282 RSB_STDOUT(RSB_PRINTF_MTX_SUMMARY_ARGS(mtxSp)),
283 RSB_STDOUT("to: \n");
284 RSB_STDOUT(RSB_PRINTF_MTX_SUMMARY_ARGS(mtxDp)),
285 RSB_STDOUT("\n");
286 RSB_STDOUT("S ip: : %x \n",(RSB_DO_FLAG_HAS(mtxSp->flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS)));
287 RSB_STDOUT("D ip: : %x \n",(RSB_DO_FLAG_HAS(mtxDp->flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS)));
288 #endif
289
290 fsmS = rsb__do_get_first_submatrix(mtxSp);
291 fsmD = rsb__do_get_first_submatrix(mtxDp);
292 VS = fsmS->VA, VD = fsmD->VA;
293 IS = fsmS->bpntr, JS = fsmS->bindx;
294 ID = fsmD->bpntr, JD = fsmD->bindx;
295
296 errval = rsb_util_coo_copy(VS, IS, JS, VD, ID, JD, mtxSp->nnz, mtxSp->typecode, 0, 0);
297 if(RSB_SOME_ERROR(errval))
298 goto err;
299
300 // correct mtxSp pointers recursively with the three offsets
301
302 /* get rid of the source COO arrays */
303 if(RSB_DO_FLAG_HAS(mtxSp->flags, RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS))
304 {
305 fsmS->VA = NULL;
306 fsmS->bindx = NULL;
307 fsmS->bpntr = NULL;
308 }
309 else
310 {
311 RSB_CONDITIONAL_FREE(fsmS->VA);
312 RSB_CONDITIONAL_FREE(fsmS->bindx);
313 RSB_CONDITIONAL_FREE(fsmS->bpntr);
314 }
315 smc = 1 + rsb__submatrices_max_ptr_diff(mtxSp);
316
317 if(RSB_DO_FLAG_HAS(mtxDp->flags, RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS))
318 RSB_DO_FLAG_ADD(mtxSp->flags, RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS);
319 else
320 ; /* mtxSp will retain its original flags with these regards */
321
322 /* set new matrix arrays by translating the submatrix pointers */
323 rsb__do_set_in_place_submatrices_offsets(mtxSp, smc, VD, ID, JD, mtxDp->el_size);
324
325 /* mtxDp->all_leaf_matrices shall be ok... */
326
327 /* free now unnecessary original destination matrix pointer */
328 fsmD->VA = NULL;
329 fsmD->bindx = NULL;
330 fsmD->bpntr = NULL;
331 RSB_MTX_FREE(mtxDp);
332
333 /* overwrite the output pointer */
334 mtxDp = mtxSp;
335 *mtxDpp = mtxDp;
336 mtxSp = NULL;
337
338 #if 0
339 RSB_STDOUT("obtained: \n");
340 RSB_STDOUT(RSB_PRINTF_MTX_SUMMARY_ARGS(mtxDp)),
341 RSB_STDOUT("\n");
342 #endif
343 err:
344 return errval;
345 }
346
347 #if 0
348 static rsb_err_t rsb_do_assign(struct rsb_mtx_t * mtxBp, const struct rsb_mtx_t * mtxAp)
349 {
350 rsb_submatrix_idx_t i,j;
351 rsb_err_t errval = RSB_ERR_NO_ERROR;
352 struct rsb_mtx_t * submatrix=NULL;
353
354 if(!mtxBp || !mtxAp)
355 goto err;
356
357 rsb__destroy_inner(mtxBp);
358
359 rsb__memcpy(mtxBp,mtxAp,sizeof(*mtxAp));
360 rsb__init_blank_pointers(mtxBp);
361
362 if(rsb__clone_inner(mtxAp,mtxBp)==NULL)
363 goto err;
364
365 RSB_SUBMATRIX_FOREACH(mtxAp,submatrix,i,j)
366 if(submatrix)
367 {
368 if((mtxBp->sm[i*2+j]=rsb__mtx_clone_simple(submatrix))==NULL)
369 goto err;
370 }
371
372 #if RSB_EXPERIMENTAL_SHOULD_TRAVERSE_RECURSIVE_MATRICES_AS_BLOCKS
373 if(mtxAp->all_leaf_matrices)
374 {
375 mtxBp->all_leaf_matrices=NULL;
376 errval = rsb__get_array_of_leaf_matrices(mtxBp,&mtxBp->all_leaf_matrices,&mtxBp->all_leaf_matrices_n);
377 if(RSB_SOME_ERROR(errval))
378 goto errr;
379 }
380 #endif /* RSB_EXPERIMENTAL_SHOULD_TRAVERSE_RECURSIVE_MATRICES_AS_BLOCKS */
381 goto errr;
382 err:
383 errval = RSB_ERR_GENERIC_ERROR;
384 errr:
385 return errval;
386 }
387 #endif
388
389
rsb__submatrices_max_ptr_diff_inner(const struct rsb_mtx_t * mtxRp,const struct rsb_mtx_t * mtxAp)390 static size_t rsb__submatrices_max_ptr_diff_inner(const struct rsb_mtx_t * mtxRp, const struct rsb_mtx_t * mtxAp)
391 {
392 /*!
393 * \ingroup gr_internals
394 * Note: this makes only sense if submatrices are allocated in one block.
395 */
396 size_t md = 0;
397 rsb_submatrix_idx_t i,j;
398 const struct rsb_mtx_t * submatrix;
399
400 RSB_SUBMATRIX_FOREACH(mtxAp,submatrix,i,j)
401 {
402 if(submatrix)
403 {
404 size_t sd = rsb__submatrices_max_ptr_diff_inner(mtxRp,submatrix);
405
406 md = RSB_MAX(md,sd);
407 md = RSB_MAX(md,(submatrix-mtxRp));
408 }
409 }
410 return md;
411 }
412
413 #if 0
414 size_t rsb__submatrices_max_ptr_diff_noholes(const struct rsb_mtx_t * mtxAp)
415 {
416 size_t md = 0;
417 rsb_submatrix_idx_t i,j;
418 const struct rsb_mtx_t * submatrix;
419
420 if(!mtxAp)
421 {
422 return 0;
423 }
424
425 RSB_SUBMATRIX_FOREACH(mtxAp,submatrix,i,j)
426 {
427 if(submatrix)
428 {
429 size_t sd = rsb__submatrices_max_ptr_diff(submatrix);
430 md = RSB_MAX(md,sd+(submatrix-mtxAp));
431 }
432 }
433 return md;
434 }
435 #endif /* 0 */
436
rsb__submatrices_max_ptr_diff(const struct rsb_mtx_t * mtxAp)437 size_t rsb__submatrices_max_ptr_diff(const struct rsb_mtx_t * mtxAp)
438 {
439 /*!
440 * \ingroup gr_internals
441 * Note: this makes only sense if submatrices are allocated in one block.
442 */
443 #if 0
444 return rsb__submatrices_max_ptr_diff_noholes(mtxAp);
445 #else
446 return rsb__submatrices_max_ptr_diff_inner(mtxAp, mtxAp);
447 #endif
448 }
449
rsb__clone_area_guided(void * RSB_RESTRICT dst,const void * RSB_RESTRICT src,size_t size,size_t nmemb,const struct rsb_mtx_t * RSB_RESTRICT mtxAp,const rsb_thread_t * RSB_RESTRICT cta,const rsb_thread_t nct,rsb_err_t * errvalp)450 static void * rsb__clone_area_guided(void * RSB_RESTRICT dst, const void *RSB_RESTRICT src, size_t size, size_t nmemb, const struct rsb_mtx_t *RSB_RESTRICT mtxAp, const rsb_thread_t * RSB_RESTRICT cta, const rsb_thread_t nct, rsb_err_t * errvalp)
451 {
452 /*
453 Initializes, eventually allocating and/or copying in parallel, using specified chunk sizes and array.
454 If dst supplied, will use it, otherwise will allocate one.
455 If src supplied, will use it, otherwise will only zero the arrays.
456 If mtxAp == NULL, then must also be cta == NULL && nct == 0.
457 Returns the either dst or the newly allocated area address.
458 */
459 rsb_err_t errval = RSB_ERR_BADARGS;
460
461 if(size*nmemb == 0)
462 {
463 RSB_PERR_GOTO(err,RSB_ERRM_ES);
464 }
465
466 if(mtxAp == NULL && (cta != NULL || nct != 0) )
467 {
468 RSB_PERR_GOTO(err,RSB_ERRM_ES);
469 }
470
471 if(dst == NULL && ( dst = rsb__malloc(size*nmemb) ) == NULL )
472 {
473 errval = RSB_ERR_ENOMEM;
474 RSB_PERR_GOTO(err,RSB_ERRM_ES);
475 }
476
477 #if RSB_WANT_OMP_RECURSIVE_KERNELS
478 if(mtxAp == NULL)
479 #endif /* RSB_WANT_OMP_RECURSIVE_KERNELS */
480 {
481 /* master thread */
482 RSB_A_MEMCPY(dst,src,0,0,nmemb,size);
483 errval = RSB_ERR_NO_ERROR;
484 goto done;
485 }
486
487 #if RSB_WANT_OMP_RECURSIVE_KERNELS
488 if(cta == NULL)
489 {
490 RSB_DEBUG_ASSERT((mtxAp)->all_leaf_matrices_n);
491 #pragma omp parallel shared(mtxAp)
492 {
493 rsb_submatrix_idx_t smi; /* submatrix index */
494 rsb_thread_t omt = omp_get_max_threads(), otn = omp_get_thread_num();
495 /* auto: each submatrix a round robin thread */
496 for(smi=0;smi<mtxAp->all_leaf_matrices_n;++smi) /* FIXME: make this an OpenMP-friendly macro */
497 if( ( smi % omt ) == otn )
498 {
499 const struct rsb_mtx_t * submatrix = mtxAp->all_leaf_matrices[smi].mtxlp;
500 size_t off = submatrix->nzoff;
501 rsb_nnz_idx_t nnz = submatrix->nnz;
502
503 if(src)
504 RSB_A_MEMCPY(dst,src,off,off,nnz,size);
505 else
506 RSB_A_BZERO(dst,off,nnz,size);
507 }
508 }
509 #pragma omp barrier
510 errval = RSB_ERR_NO_ERROR;
511 goto done;
512 }
513
514 #pragma omp parallel shared(mtxAp) RSB_NTC
515 {
516 /* guided: each submatrix a specified thread */
517 rsb_thread_t otn = omp_get_thread_num();
518 rsb_submatrix_idx_t cti; /* thread index */
519
520 for(cti=0;cti<nct;++cti)
521 if( otn == cta[cti] )
522 {
523 const struct rsb_mtx_t * submatrix = mtxAp->all_leaf_matrices[cti].mtxlp;
524 size_t off = submatrix->nzoff;
525 rsb_nnz_idx_t nnz = submatrix->nnz;
526
527 if(src)
528 RSB_A_MEMCPY(dst,src,off,off,nnz,size);
529 else
530 RSB_A_BZERO(dst,off,nnz,size);
531 }
532 }
533 #pragma omp barrier
534 #endif /* RSB_WANT_OMP_RECURSIVE_KERNELS */
535 errval = RSB_ERR_NO_ERROR;
536 goto done;
537 err:
538 dst = NULL;
539 done:
540 /* FIXME: errval unused so far */
541 RSB_CONDITIONAL_ERRPSET(errvalp,errval);
542 return dst;
543 }
544
rsb__clone_simple_extra(const struct rsb_mtx_t * mtxAp,rsb_submatrix_idx_t esmc)545 struct rsb_mtx_t *rsb__clone_simple_extra(const struct rsb_mtx_t *mtxAp, rsb_submatrix_idx_t esmc)
546 {
547 /*!
548 * \ingroup gr_internals
549 *
550 * Clones a whole matrix, retaining the same submatrices structure.
551 * TODO: need better error handling.
552 * FIXME : Unfinished: only the RSB_FLAG_ASSEMBLED_IN_COO_ARRAYS case is handled.
553 * TODO: rename from rsb__clone_simple_extra to rsb__mtx_clone_simple_extra.
554 */
555 struct rsb_mtx_t *mtxCp = NULL;
556 rsb_flags_t flags = RSB_FLAG_NOFLAGS;
557 #if RSB_ALLOW_INTERNAL_GETENVS
558 rsb_time_t ct = RSB_TIME_ZERO;
559 rsb_time_t mact = RSB_TIME_ZERO;
560 #endif /* RSB_ALLOW_INTERNAL_GETENVS */
561
562 if(!mtxAp)
563 {
564 RSB_PERR_GOTO(nerr,RSB_ERRM_ES);
565 }
566
567 flags = mtxAp->flags;
568 #if RSB_ALLOW_INTERNAL_GETENVS
569 ct = -rsb_time();
570 #endif /* RSB_ALLOW_INTERNAL_GETENVS */
571 RSB_DO_FLAG_DEL(flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS);
572 if(RSB_DO_FLAG_HAS(flags,RSB_FLAG_ASSEMBLED_IN_COO_ARRAYS))
573 {
574 if(rsb__is_root_matrix(mtxAp))
575 {
576 /* this is a trick: fitting the whole recursive matrix in three arrays */
577 /* void *IA = NULL,*JA = NULL,*VA = NULL; */
578 const rsb_nnz_idx_t nnz = mtxAp->nnz;
579 /* rsb_bool_t is_bio = rsb__do_is_matrix_binary_loaded(mtxAp);*/ /* binary I/O matrix (20120930 FIXME why is this unused ?) */
580 /* rsb_long_t smc = rsb__submatrices(mtxAp); */
581 rsb_long_t smc = 1 + rsb__submatrices_max_ptr_diff(mtxAp);
582 const struct rsb_mtx_t *fsm = rsb__do_get_first_submatrix(mtxAp);
583 rsb_err_t errval = RSB_ERR_NO_ERROR;
584 void * VA = NULL;
585 rsb_coo_idx_t * bpntr = NULL, * bindx = NULL;
586
587 RSB_DEBUG_ASSERT(1+rsb__submatrices_max_ptr_diff(mtxAp)>=rsb__submatrices(mtxAp));
588
589 /* RSB_STDOUT("MAX PTR DIFF: %d SUBM COUNT:%d\n",1+rsb__submatrices_max_ptr_diff(mtxAp), rsb__submatrices(mtxAp)); */
590
591 /* mtxCp = rsb__clone_area(mtxAp,sizeof(struct rsb_mtx_t)*(smc+esmc)); */
592 mtxCp = rsb__clone_area_with_extra(mtxAp,sizeof(struct rsb_mtx_t)*(smc),0,sizeof(struct rsb_mtx_t)*(esmc));
593
594 if(!mtxCp)
595 {
596 RSB_PERR_GOTO(nerr,RSB_ERRM_PAL);
597 }
598
599 errval = rsb__mtx_shift_leaf_ptrs(mtxCp, mtxAp, smc);
600
601 mtxCp->all_leaf_matrices = NULL;
602 #if 0
603 errval = rsb__get_array_of_leaf_matrices(mtxCp,&(mtxCp->all_leaf_matrices),&(mtxCp->all_leaf_matrices_n));
604 #else
605 mtxCp->all_leaf_matrices = rsb__clone_area_with_extra(mtxAp->all_leaf_matrices,sizeof(mtxAp->all_leaf_matrices[0])*(mtxCp->all_leaf_matrices_n),0,0);
606 if( mtxCp->all_leaf_matrices == NULL )
607 {
608 errval = RSB_ERR_ENOMEM;
609 }
610 else
611 {
612 rsb_submatrix_idx_t si;
613 for(si=0;si<mtxCp->all_leaf_matrices_n;++si)
614 RSB_PTR_SHIFT(mtxCp->all_leaf_matrices[si].mtxlp,mtxAp,mtxCp,(struct rsb_mtx_t*));
615 for(si=0;si<mtxCp->all_leaf_matrices_n;++si)
616 RSB_DO_FLAG_DEL(mtxCp->all_leaf_matrices[si].mtxlp->flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS);
617 RSB_DO_FLAG_DEL(mtxCp->flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS);
618 }
619 #endif
620 if(RSB_SOME_ERROR(errval))
621 {
622 RSB_PERR_GOTO(nerr,RSB_ERRM_NAOL);
623 }
624 /* FIXME: and what when nnz==0 and VA!=NULL ? */
625 mtxCp->bindx = NULL; mtxCp->bpntr = NULL; mtxCp->VA = NULL;
626 if(nnz)
627 {
628 #if RSB_ALLOW_INTERNAL_GETENVS
629 mact = -rsb_time();
630 #endif /* RSB_ALLOW_INTERNAL_GETENVS */
631 RSB_ASSERT( fsm->el_size );
632 RSB_ASSERT( fsm->bindx );
633 RSB_ASSERT( fsm->bpntr );
634 RSB_ASSERT( fsm->VA );
635
636 #if RSB_WANT_SM_TO_THREAD_MOD_MAPPING
637 if( 1 /*rsb__util_atoi(getenv("RSB_CLONE_SERIAL")) */ )
638 {
639 /* rsb_thread_t nct = (rsb__util_atoi(getenv("RSB_CLONE_SERIAL"))) - 1; */
640 rsb_thread_t nct = 0;
641 bindx = rsb__clone_area_guided(NULL,fsm->bindx,sizeof(rsb_coo_idx_t),nnz,mtxAp,NULL,nct,&errval);
642 bpntr = rsb__clone_area_guided(NULL,fsm->bpntr,sizeof(rsb_coo_idx_t),nnz,mtxAp,NULL,nct,&errval);
643 VA = rsb__clone_area_guided(NULL,fsm->VA ,mtxAp->el_size, nnz,mtxAp,NULL,nct,&errval);
644 }
645 else
646 #endif
647 {
648 bindx = rsb__clone_area_parallel(fsm->bindx,sizeof(rsb_coo_idx_t),nnz);
649 bpntr = rsb__clone_area_parallel(fsm->bpntr,sizeof(rsb_coo_idx_t),nnz);
650 VA = rsb__clone_area_parallel(fsm->VA,mtxAp->el_size,nnz);
651 }
652
653 #if RSB_ALLOW_INTERNAL_GETENVS
654 mact += rsb_time();
655 #endif /* RSB_ALLOW_INTERNAL_GETENVS */
656 if(!bindx || !bpntr || !VA || !mtxAp->el_size)
657 {
658 RSB_ASSERT( mtxCp->el_size );
659 RSB_ASSERT( bindx );
660 RSB_ASSERT( bpntr );
661 RSB_ASSERT( VA );
662 RSB_PERR_GOTO(ierr,RSB_ERRM_NNTC);
663 }
664 }
665 #if !RSB_ALLOW_EMPTY_MATRICES
666 else
667 if(RSB_DO_FLAG_HAS(flags,RSB_FLAG_UNIT_DIAG_IMPLICIT))
668 {
669 // ok
670 }
671 else
672 {
673 RSB_PERR_GOTO(ierr,RSB_ERRM_NDIANN);
674 }
675 #endif /* RSB_ALLOW_EMPTY_MATRICES */
676 rsb__do_set_in_place_submatrices_offsets(mtxCp,smc,VA,bpntr,bindx,mtxCp->el_size);
677
678 if(!smc)
679 mtxCp->bindx = bindx,
680 mtxCp->bpntr = bpntr,
681 mtxCp->VA = VA;
682
683 /* note: the cloned matrix won't have the is_bio property */
684 goto ret;
685 ierr:
686 RSB_ERROR(RSB_ERRM_ES);
687 if(mtxCp)
688 {
689 RSB_CONDITIONAL_FREE(mtxCp->bpntr);
690 RSB_CONDITIONAL_FREE(mtxCp->bindx);
691 RSB_CONDITIONAL_FREE(mtxCp->VA);
692 RSB_CONDITIONAL_FREE(mtxCp->all_leaf_matrices);
693 }
694 RSB_CONDITIONAL_FREE(mtxCp);
695 goto ret;
696 }
697 else
698 {
699 RSB_PERR_GOTO(ret,"no cloning possible for a non root matrix!\n");
700 /* no cloning for a non root */
701 }
702 }
703 else
704 {
705 /* we allocate a new matrix structure */
706 mtxCp = rsb__clone_area(mtxAp,sizeof(*mtxCp));
707
708 if(!mtxCp)
709 {RSB_PERR_GOTO(err,RSB_ERRM_ES);}
710
711 rsb__init_blank_pointers(mtxCp);
712
713 RSB_MTX_REASSIGN(mtxCp,(struct rsb_mtx_t*)mtxAp);
714 goto ret;
715 }
716 err:
717 RSB_MTX_FREE(mtxCp);
718 nerr:
719 RSB_CONDITIONAL_FREE(mtxCp);
720 ret:
721 #if RSB_ALLOW_INTERNAL_GETENVS
722 ct += rsb_time();
723 if( rsb__util_atoi(getenv("RSB_MTX_CLONE_STATS") ) != 0)
724 if(mtxCp)
725 {
726 size_t szv = rsb__get_sizeof(mtxCp);
727 RSB_STDOUT("Cloned a %zd nnz, %zd bytes matrix in %0.2lgs (%0.3lg MiB/s x 2 = r+w); of which %0.2lgs for the main arrays.\n",
728 (size_t)(mtxCp->nnz),szv,ct,(((rsb_time_t)szv)/ct)/RSB_MEGABYTE,mact);
729 }
730 #endif /* RSB_ALLOW_INTERNAL_GETENV */
731 if(mtxCp)
732 RSB_DO_FLAG_DEL(mtxCp->flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS);
733 return mtxCp;
734 }
735
rsb__mtx_clone_simple(const struct rsb_mtx_t * mtxAp)736 struct rsb_mtx_t *rsb__mtx_clone_simple(const struct rsb_mtx_t *mtxAp)
737 {
738 /* TODO: rename from rsb__mtx_clone_simple to rsb__mtx_clone_simple */
739 return rsb__clone_simple_extra(mtxAp, 0);
740 }
741
rsb__clone_coo(const struct rsb_mtx_t * mtxAp,rsb_trans_t transA,const void * alphap,rsb_type_t typecode,struct rsb_coo_matrix_t * dcoop,rsb_flags_t flags)742 rsb_err_t rsb__clone_coo(const struct rsb_mtx_t * mtxAp, rsb_trans_t transA, const void *alphap, rsb_type_t typecode, struct rsb_coo_matrix_t*dcoop, rsb_flags_t flags/*, rsb_extff_t eflags*/)
743 {
744 /*
745 TODO: may integrate here fortran indices handling and so on
746 TODO: missing checks for indices overflow
747 TODO: shall juggle appropriately with 'sorted' flags
748 */
749 rsb_flags_t cflags = RSB_FLAG_NOFLAGS;
750 rsb_err_t errval = RSB_ERR_NO_ERROR;
751 rsb_nnz_idx_t dels = 0;
752 rsb_coo_idx_t ioff,joff;
753 struct rsb_coo_matrix_t dcoo,scoo;
754 rsb_bool_t expsymm = RSB_BOOL_FALSE, expherm = RSB_BOOL_FALSE;
755
756 RSB_BZERO_P(&scoo);
757 RSB_BZERO_P(&dcoo);
758 ioff = joff = ( flags & RSB_FLAG_FORTRAN_INDICES_INTERFACE )?1:0;
759 scoo.nr = dcoo.nr = mtxAp->nr;
760 scoo.nc = dcoo.nc = mtxAp->nc;
761 dcoo.nnz = scoo.nnz = mtxAp->nnz;
762 expsymm = (RSB_DO_FLAG_HAVE_XOR(flags,mtxAp->flags,RSB_FLAG_SYMMETRIC) && !RSB_DO_FLAG_HAS(mtxAp->flags,RSB_FLAG_DIAGONAL));
763 expherm = (RSB_DO_FLAG_HAVE_XOR(flags,mtxAp->flags,RSB_FLAG_HERMITIAN) && !RSB_DO_FLAG_HAS(mtxAp->flags,RSB_FLAG_DIAGONAL));
764 if(expsymm || expherm)
765 dcoo.nnz *= 2;/* of course, this is overkill in the case of a diagonal matrix */
766 if(RSB_DO_FLAG_HAVE_XOR(flags,mtxAp->flags,RSB_FLAG_UNIT_DIAG_IMPLICIT))
767 dels = RSB_MIN(dcoo.nr,dcoo.nc);
768 dcoo.nnz += dels;
769 //if(dels)
770 // RSB_STDOUT("on diag: %d\n",dels);
771
772 scoo.typecode = mtxAp->typecode;
773 dcoo.typecode = typecode;
774 if(dcoo.nnz>0)
775 {
776 if(rsb__allocate_coo_matrix_t(&dcoo)!=&dcoo)
777 {
778 errval = RSB_ERR_INTERNAL_ERROR;
779 RSB_PERR_GOTO(ierr,RSB_ERRM_PAL);
780 }
781 if(rsb__allocate_coo_matrix_t(&scoo)!=&scoo)
782 {
783 errval = RSB_ERR_INTERNAL_ERROR;
784 RSB_PERR_GOTO(ierr,RSB_ERRM_PAL);
785 }
786 errval = rsb__do_get_coo_noalloc(mtxAp,scoo.VA,dcoo.IA,dcoo.JA,NULL,/*mtxAp->*/flags);
787 if(RSB_SOME_ERROR(errval))
788 {
789 RSB_PERR_GOTO(ierr,RSB_ERRM_NL);
790 }
791 errval = rsb__do_copy_converted_scaled(scoo.VA,dcoo.VA,alphap,mtxAp->typecode,typecode,mtxAp->nnz,transA);
792 if(RSB_SOME_ERROR(errval))
793 {
794 RSB_PERR_GOTO(ierr,RSB_ERRM_NL);
795 }
796 if(expsymm || expherm)
797 RSB_COO_MEMCPY(dcoo.VA,dcoo.IA,dcoo.JA,dcoo.VA,dcoo.JA,dcoo.IA,scoo.nnz,0,scoo.nnz,RSB_SIZEOF(typecode));
798 if(expherm)
799 rsb__util_do_conjugate(((rsb_byte_t*)(dcoo.VA))+(RSB_SIZEOF(typecode)*scoo.nnz),typecode,scoo.nnz);
800 if(RSB_DO_FLAG_HAVE_XOR(flags,mtxAp->flags,RSB_FLAG_UNIT_DIAG_IMPLICIT))
801 rsb__do_fill_with_diag(dcoo.VA,dcoo.IA,dcoo.JA,ioff,joff,dcoo.nnz-dels,typecode,dels);
802 rsb__destroy_coo_matrix_t(&scoo);
803 RSB_BZERO_P(&scoo);
804 }
805
806 if(RSB_DO_FLAG_HAS(flags,RSB_FLAG_UPPER_TRIANGULAR))
807 RSB_DO_FLAG_ADD(cflags,RSB_FLAG_UPPER_TRIANGULAR);
808 if(RSB_DO_FLAG_HAS(flags,RSB_FLAG_LOWER_TRIANGULAR))
809 RSB_DO_FLAG_ADD(cflags,RSB_FLAG_LOWER_TRIANGULAR);
810 if(RSB_DO_FLAG_HAVE_XOR(flags,mtxAp->flags,RSB_FLAG_UNIT_DIAG_IMPLICIT))
811 RSB_DO_FLAG_ADD(cflags,RSB_FLAG_UNIT_DIAG_IMPLICIT);
812 if((cflags != RSB_FLAG_NOFLAGS) || expsymm || expherm)
813 {
814 rsb__util_sort_row_major_inner(dcoo.VA,dcoo.IA,dcoo.JA,dcoo.nnz,dcoo.nr,dcoo.nc,typecode,flags);
815 RSB_DO_FLAG_ADD(flags,RSB_FLAG_SORTED_INPUT);
816 dcoo.nnz = rsb_weed_out_duplicates(dcoo.IA,dcoo.JA,dcoo.VA,dcoo.nnz,typecode,flags);
817 errval = rsb__do_cleanup_nnz(dcoo.VA,dcoo.IA,dcoo.JA,dcoo.nnz,0,0,dcoo.nr,dcoo.nc,&dcoo.nnz,dcoo.typecode,cflags); /* FIXME: are we using roff,coff well here ? */
818 }
819 if(RSB_SOME_ERROR(errval))
820 {
821 errval = RSB_ERR_INTERNAL_ERROR;
822 RSB_PERR_GOTO(ierr,RSB_ERRM_CP);
823 }
824
825 if(RSB_DOES_TRANSPOSE(transA))
826 rsb__transpose_coo_matrix_t(&dcoo);
827 *dcoop = dcoo;
828 ierr:
829 return errval;
830 }
831
rsb__mtx_clone(struct rsb_mtx_t ** mtxBpp,rsb_type_t typecode,rsb_trans_t transA,const void * alphap,const struct rsb_mtx_t * mtxAp,rsb_flags_t flags)832 rsb_err_t rsb__mtx_clone(struct rsb_mtx_t ** mtxBpp, rsb_type_t typecode, rsb_trans_t transA, const void *alphap, const struct rsb_mtx_t * mtxAp, rsb_flags_t flags)
833 {
834 /*!
835 * \ingroup gr_internals
836 * clones a rsb_mtx_t structure, deeply
837 * TODO: rename rsb__mtx_clone -> rsb_cln ?
838 * This routine may/shall be optimized in plenty of ways, in the future.
839 * */
840 rsb_err_t errval = RSB_ERR_NO_ERROR;
841 struct rsb_mtx_t * mtxCp = NULL;
842
843 if( typecode != RSB_NUMERICAL_TYPE_SAME_TYPE && RSB_MATRIX_UNSUPPORTED_TYPE(typecode) )
844 {
845 errval = RSB_ERR_UNSUPPORTED_TYPE;
846 RSB_PERR_GOTO(err,RSB_ERRM_ES);
847 }
848
849 if( (!mtxAp) || (!mtxBpp) )
850 {
851 errval = RSB_ERR_BADARGS;
852 RSB_PERR_GOTO(err,"user did not supply a valid matrix pointer\n");
853 }
854
855 if(RSB_DO_FLAG_HAS(flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS))
856 {
857 errval = RSB_ERR_BADARGS;
858 RSB_PERR_GOTO(err,"user supplied wrong flags (RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS is illegal here)\n");
859 }
860
861 if( flags == RSB_FLAG_IDENTICAL_FLAGS )
862 flags = mtxAp->flags;
863 if(typecode == RSB_NUMERICAL_TYPE_SAME_TYPE)
864 typecode = mtxAp->typecode;
865 RSB_DO_FLAG_DEL(flags,RSB_FLAG_FORTRAN_INDICES_INTERFACE);/* unnecessary here */
866 RSB_DO_FLAG_DEL(flags,RSB_FLAG_NON_ROOT_MATRIX);
867
868 /* what about RSB_FLAG_DISCARD_ZEROS ? */
869 /* what about many other structural flags ? */
870
871 if(((!alphap) || RSB_IS_ELEMENT_ONE(alphap,mtxAp->typecode)) &&
872 (flags == mtxAp->flags ) && /* FIXME: this condition is unnecessarily strict and cause inefficiencies */
873 /* RSB_DOES_NOT_TRANSPOSE(transA) && */ (typecode==mtxAp->typecode) )
874 {
875 if( (*mtxBpp) != mtxAp)
876 mtxCp = rsb__mtx_clone_simple(mtxAp);
877 else
878 mtxCp = *mtxBpp;
879 if( transA == RSB_TRANSPOSITION_C )
880 errval = rsb__do_transpose(&mtxCp,RSB_BOOL_TRUE);
881 else
882 if( transA == RSB_TRANSPOSITION_T )
883 errval = rsb__do_transpose(&mtxCp,RSB_BOOL_FALSE);
884 if( (*mtxBpp) == mtxAp)
885 {
886 *mtxBpp = mtxCp;
887 goto ok;
888 }
889 }
890 else
891 {
892 struct rsb_coo_matrix_t dcoo;
893 RSB_BZERO_P(&dcoo);
894 #if 0
895 struct rsb_coo_matrix_t scoo;
896 RSB_BZERO_P(&scoo);
897 scoo.nr = dcoo.nr = mtxAp->nr;
898 scoo.nc = dcoo.nc = mtxAp->nc;
899 dcoo.nnz = scoo.nnz = mtxAp->nnz;
900 scoo.typecode = mtxAp->typecode;
901 dcoo.typecode = typecode;
902 if(mtxAp->nnz>0)
903 {
904 if(rsb__allocate_coo_matrix_t(&dcoo)!=&dcoo)
905 {
906 errval = RSB_ERR_INTERNAL_ERROR;
907 RSB_PERR_GOTO(ierr,RSB_ERRM_PAL);
908 }
909 if(rsb__allocate_coo_matrix_t(&scoo)!=&scoo)
910 {
911 errval = RSB_ERR_INTERNAL_ERROR;
912 RSB_PERR_GOTO(ierr,RSB_ERRM_PAL);
913 }
914 errval = rsb__do_get_coo_noalloc(mtxAp,scoo.VA,dcoo.IA,dcoo.JA,NULL,mtxAp->flags);
915 rsb__do_copy_converted_scaled(scoo.VA,dcoo.VA,alphap,mtxAp->typecode,typecode,mtxAp->nnz,transA);
916 rsb__destroy_coo_matrix_t(&scoo);
917 RSB_BZERO_P(&scoo);
918 }
919 if(RSB_DOES_TRANSPOSE(transA))
920 rsb__transpose_coo_matrix_t(&dcoo);
921 #else
922 errval = rsb__clone_coo(mtxAp,transA,alphap,typecode,&dcoo,flags);
923 if(RSB_SOME_ERROR(errval))
924 {
925 RSB_PERR_GOTO(ierr,RSB_ERRM_NL);
926 }
927 #endif
928 mtxCp = rsb__do_mtx_alloc_from_coo_inplace(dcoo.VA,dcoo.IA,dcoo.JA,dcoo.nnz,dcoo.typecode,dcoo.nr,dcoo.nc,RSB_DEFAULT_ROW_BLOCKING,RSB_DEFAULT_COL_BLOCKING,flags,NULL);
929 if(mtxCp)
930 RSB_DO_FLAG_DEL(mtxCp->flags,RSB_FLAG_EXTERNALLY_ALLOCATED_ARRAYS);
931 ierr:
932 if(!mtxCp)
933 rsb__destroy_coo_matrix_t(&dcoo);
934 }
935
936 if( (*mtxBpp) == NULL )
937 {
938 *mtxBpp = mtxCp;
939 }
940 else
941 {
942 RSB_MTX_REASSIGN(*mtxBpp,mtxCp);
943 }
944 ok:
945 err:
946 return errval;
947 }
948
949 #if 0
950 void * rsb__clone_inner(const struct rsb_mtx_t *mtxAp, struct rsb_mtx_t *mtxCp)
951 {
952 /*!
953 * \ingroup gr_internals
954 * clones a rsb_mtx_t structure, deeply
955 *
956 * \param matrix valid matrix pointer (to an empty mtxAp)
957 * \param mtxCp valid matrix pointer
958 * \return a pointer to the cloned structure (mtxCp) in case of success, NULL otherwise
959 *
960 * \note matrix flags are largely ignored in this function.
961 **/
962
963 if(!mtxAp || !mtxCp)
964 {RSB_PERR_GOTO(err,RSB_ERRM_ES);}
965
966 #if RSB_WANT_BITMAP
967 /* we allocate a new options structure */
968 mtxCp->options = rsb__clone_options_t(mtxAp->options,mtxAp->M_b,mtxAp->K_b);
969
970 if(! mtxCp->options && mtxAp->options )
971 {RSB_PERR_GOTO(err_opt,RSB_ERRM_ES);}
972 #endif /* RSB_WANT_BITMAP */
973 if( mtxAp->rpntr && (mtxAp->flags & RSB_FLAG_OWN_PARTITIONING_ARRAYS))
974 {
975 mtxCp->rpntr = rsb__clone_area(mtxAp->rpntr,sizeof(rsb_coo_idx_t)*(mtxAp->M_b+1));
976 if(!mtxCp->rpntr)
977 {RSB_PERR_GOTO(err_rpntr,RSB_ERRM_ES);}
978 }
979 else
980 mtxCp->rpntr = mtxAp->rpntr;
981
982 if( mtxAp->cpntr && (mtxAp->flags & RSB_FLAG_OWN_PARTITIONING_ARRAYS))
983 {
984 mtxCp->cpntr = rsb__clone_area(mtxAp->cpntr,sizeof(rsb_coo_idx_t)*(mtxAp->K_b+1));
985 if(!mtxCp->cpntr)
986 {RSB_PERR_GOTO(err_cpntr,RSB_ERRM_ES);}
987 }
988 else
989 mtxCp->cpntr = mtxAp->cpntr;
990
991 if( mtxAp->bindx)
992 {
993 mtxCp->bindx = rsb__clone_area(mtxAp->bindx,sizeof(rsb_nnz_idx_t)*(mtxAp->block_count+1));
994 if(!mtxCp->bindx)
995 {RSB_PERR_GOTO(err_bindx,RSB_ERRM_ES);}
996 }
997
998 if( mtxAp->indptr)
999 {
1000 mtxCp->indptr = rsb__clone_area(mtxAp->indptr,sizeof(rsb_nnz_idx_t)*(mtxAp->block_count+1));
1001 if(!mtxCp->indptr)
1002 {RSB_PERR_GOTO(err_indptr,RSB_ERRM_ES);}
1003 }
1004
1005 if( mtxAp->bpntr)
1006 {
1007 mtxCp->bpntr = rsb__clone_area(mtxAp->bpntr,sizeof(rsb_nnz_idx_t)*(mtxAp->Mdim+1));
1008 if(!mtxCp->bpntr)
1009 {RSB_PERR_GOTO(err_bpntr,RSB_ERRM_ES);}
1010 }
1011
1012 #if RSB_WANT_BITMAP
1013 if( mtxAp->VA)
1014 {
1015 mtxCp->VA = rsb__clone_area(mtxAp->VA,(RSB_TOTAL_BLOCK_BYTES(mtxAp,mtxAp->options)));
1016 if(!mtxCp->VA)
1017 {RSB_PERR_GOTO(err_va,RSB_ERRM_ES);}
1018 }
1019 #endif
1020 goto ret;
1021
1022 #if RSB_WANT_BITMAP
1023 err_va:
1024 if( mtxAp->VA)
1025 RSB_CONDITIONAL_FREE(mtxCp->VA);
1026 #endif
1027 err_bpntr:
1028 if( mtxAp->bpntr )
1029 RSB_CONDITIONAL_FREE(mtxCp->bpntr);
1030 err_indptr:
1031 if( mtxAp->indptr )
1032 RSB_CONDITIONAL_FREE(mtxCp->indptr);
1033 err_bindx:
1034 if( mtxAp->bindx )
1035 RSB_CONDITIONAL_FREE(mtxCp->bindx);
1036 err_cpntr:
1037 if( mtxAp->cpntr && (mtxAp->flags & RSB_FLAG_OWN_PARTITIONING_ARRAYS))
1038 RSB_CONDITIONAL_FREE(mtxCp->cpntr);
1039 err_rpntr:
1040 if( mtxAp->rpntr && (mtxAp->flags & RSB_FLAG_OWN_PARTITIONING_ARRAYS))
1041 RSB_CONDITIONAL_FREE(mtxCp->rpntr);
1042 #if RSB_WANT_BITMAP
1043 err_opt:
1044 RSB_CONDITIONAL_FREE(mtxCp->options);
1045 #endif /* RSB_WANT_BITMAP */
1046 err:
1047 mtxCp = NULL;
1048 ret:
1049 return mtxCp;
1050 }
1051 #endif
1052
1053 /* @endcond */
1054