1 /* Copyright (C) 2005 The Scalable Software Infrastructure Project. All rights reserved.
2
3 Redistribution and use in source and binary forms, with or without
4 modification, are permitted provided that the following conditions are met:
5 1. Redistributions of source code must retain the above copyright
6 notice, this list of conditions and the following disclaimer.
7 2. Redistributions in binary form must reproduce the above copyright
8 notice, this list of conditions and the following disclaimer in the
9 documentation and/or other materials provided with the distribution.
10 3. Neither the name of the project nor the names of its contributors
11 may be used to endorse or promote products derived from this software
12 without specific prior written permission.
13
14 THIS SOFTWARE IS PROVIDED BY THE SCALABLE SOFTWARE INFRASTRUCTURE PROJECT
15 ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
16 TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
17 PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE SCALABLE SOFTWARE INFRASTRUCTURE
18 PROJECT BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,
19 OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
20 SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
21 INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
22 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
23 ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
24 POSSIBILITY OF SUCH DAMAGE.
25 */
26
27 #ifdef HAVE_CONFIG_H
28 #include "lis_config.h"
29 #else
30 #ifdef HAVE_CONFIG_WIN_H
31 #include "lis_config_win.h"
32 #endif
33 #endif
34
35 #include <stdio.h>
36 #include <stdlib.h>
37 #ifdef HAVE_MALLOC_H
38 #include <malloc.h>
39 #endif
40 #include <string.h>
41 #include <stdarg.h>
42 #include <math.h>
43 #ifdef _OPENMP
44 #include <omp.h>
45 #endif
46 #ifdef USE_MPI
47 #include <mpi.h>
48 #endif
49 #include "lislib.h"
50
51 /************************************************
52 * lis_matrix_malloc
53 * lis_matrix_realloc
54 ************************************************/
55
56 #undef __FUNC__
57 #define __FUNC__ "lis_matrix_create_rco"
lis_matrix_create_rco(LIS_INT local_n,LIS_INT global_n,LIS_Comm comm,LIS_INT annz,LIS_INT * nnz,LIS_MATRIX * Amat)58 LIS_INT lis_matrix_create_rco(LIS_INT local_n, LIS_INT global_n, LIS_Comm comm, LIS_INT annz, LIS_INT *nnz, LIS_MATRIX *Amat)
59 {
60 LIS_INT nprocs,my_rank;
61 LIS_INT is,ie,err;
62 LIS_INT i,k;
63 LIS_INT *ranges;
64
65 LIS_DEBUG_FUNC_IN;
66
67 *Amat = NULL;
68
69 if( global_n>0 && local_n>global_n )
70 {
71 LIS_SETERR2(LIS_ERR_ILL_ARG,"local n(=%D) is larger than global n(=%D)\n",local_n,global_n);
72 return LIS_ERR_ILL_ARG;
73 }
74 if( local_n<0 || global_n<0 )
75 {
76 LIS_SETERR2(LIS_ERR_ILL_ARG,"local n(=%D) or global n(=%D) are less than 0\n",local_n,global_n);
77 return LIS_ERR_ILL_ARG;
78 }
79 if( local_n==0 && global_n==0 )
80 {
81 LIS_SETERR2(LIS_ERR_ILL_ARG,"local n(=%D) and global n(=%D) are 0\n",local_n,global_n);
82 return LIS_ERR_ILL_ARG;
83 }
84
85 *Amat = (LIS_MATRIX)lis_malloc( sizeof(struct LIS_MATRIX_STRUCT),"lis_matrix_create_rco::Amat" );
86 if( NULL==*Amat )
87 {
88 LIS_SETERR_MEM(sizeof(struct LIS_MATRIX_STRUCT));
89 return LIS_OUT_OF_MEMORY;
90 }
91 lis_matrix_init(Amat);
92
93 err = lis_ranges_create(comm,&local_n,&global_n,&ranges,&is,&ie,&nprocs,&my_rank);
94 if( err )
95 {
96 lis_matrix_destroy(*Amat);
97 *Amat = NULL;
98 return err;
99 }
100 (*Amat)->ranges = ranges;
101
102 (*Amat)->w_nnz = (LIS_INT *)lis_malloc(local_n*sizeof(LIS_INT),"lis_matrix_create_rco::Amat->w_nnz");
103 if( (*Amat)->w_nnz==NULL )
104 {
105 LIS_SETERR_MEM(local_n*sizeof(LIS_INT));
106 return LIS_OUT_OF_MEMORY;
107 }
108 if( nnz==NULL )
109 {
110 (*Amat)->w_annz = annz;
111 for(k=0;k<local_n;k++) (*Amat)->w_nnz[k] = (*Amat)->w_annz;
112 }
113 else
114 {
115 i = 0;
116 for(k=0;k<local_n;k++)
117 {
118 (*Amat)->w_nnz[k] = nnz[k];
119 i += nnz[k];
120 }
121 (*Amat)->w_annz = i/local_n;
122 }
123 err = lis_matrix_malloc_rco(local_n,(*Amat)->w_nnz,&(*Amat)->w_row,&(*Amat)->w_index,&(*Amat)->w_value);
124 if( err )
125 {
126 lis_free((*Amat)->w_nnz);
127 return err;
128 }
129 (*Amat)->status = LIS_MATRIX_ASSEMBLING;
130
131 (*Amat)->n = local_n;
132 (*Amat)->gn = global_n;
133 (*Amat)->np = local_n;
134 (*Amat)->comm = comm;
135 (*Amat)->my_rank = my_rank;
136 (*Amat)->nprocs = nprocs;
137 (*Amat)->is = is;
138 (*Amat)->ie = ie;
139
140 LIS_DEBUG_FUNC_OUT;
141 return LIS_SUCCESS;
142 }
143
144 #undef __FUNC__
145 #define __FUNC__ "lis_matrix_malloc_rco"
lis_matrix_malloc_rco(LIS_INT n,LIS_INT nnz[],LIS_INT ** row,LIS_INT *** index,LIS_SCALAR *** value)146 LIS_INT lis_matrix_malloc_rco(LIS_INT n, LIS_INT nnz[], LIS_INT **row, LIS_INT ***index, LIS_SCALAR ***value)
147 {
148 LIS_INT i,j;
149 LIS_INT *w_row,**w_index;
150 LIS_SCALAR **w_value;
151
152 LIS_DEBUG_FUNC_IN;
153
154 w_row = NULL;
155 w_index = NULL;
156 w_value = NULL;
157
158 w_row = (LIS_INT *)lis_malloc( n*sizeof(LIS_INT),"lis_matrix_malloc_rco::w_row" );
159 if( w_row==NULL )
160 {
161 LIS_SETERR_MEM(n*sizeof(LIS_INT));
162 return LIS_OUT_OF_MEMORY;
163 }
164 w_index = (LIS_INT **)lis_malloc( n*sizeof(LIS_INT *),"lis_matrix_malloc_rco::w_index" );
165 if( w_index==NULL )
166 {
167 LIS_SETERR_MEM(n*sizeof(LIS_INT *));
168 lis_free2(3,w_row,w_index,w_value);
169 return LIS_OUT_OF_MEMORY;
170 }
171 w_value = (LIS_SCALAR **)lis_malloc( n*sizeof(LIS_SCALAR *),"lis_matrix_malloc_rco::w_value" );
172 if( w_value==NULL )
173 {
174 LIS_SETERR_MEM(n*sizeof(LIS_SCALAR *));
175 lis_free2(3,w_row,w_index,w_value);
176 return LIS_OUT_OF_MEMORY;
177 }
178 if( nnz!=NULL )
179 {
180 for(i=0;i<n;i++)
181 {
182 w_index[i] = NULL;
183 w_value[i] = NULL;
184 if( nnz[i]==0 ) continue;
185 w_index[i] = (LIS_INT *)lis_malloc( nnz[i]*sizeof(LIS_INT),"lis_matrix_malloc_rco::w_index[i]" );
186 if( w_index[i]==NULL )
187 {
188 LIS_SETERR_MEM(nnz[i]*sizeof(LIS_INT));
189 break;
190 }
191 w_value[i] = (LIS_SCALAR *)lis_malloc( nnz[i]*sizeof(LIS_SCALAR),"lis_matrix_malloc_rco::w_value[i]" );
192 if( w_value[i]==NULL )
193 {
194 LIS_SETERR_MEM(nnz[i]*sizeof(LIS_SCALAR));
195 break;
196 }
197 }
198 if(i<n)
199 {
200 for(j=0;j<i;j++)
201 {
202 if( w_index[i] ) lis_free(w_index[i]);
203 if( w_value[i] ) lis_free(w_value[i]);
204 }
205 lis_free2(3,w_row,w_index,w_value);
206 return LIS_OUT_OF_MEMORY;
207 }
208 }
209 #ifdef _OPENMP
210 #pragma omp parallel for private(i)
211 #endif
212 for(i=0;i<n;i++) w_row[i] = 0;
213 *row = w_row;
214 *index = w_index;
215 *value = w_value;
216
217 LIS_DEBUG_FUNC_OUT;
218 return LIS_SUCCESS;
219 }
220
221 #undef __FUNC__
222 #define __FUNC__ "lis_matrix_realloc_rco"
lis_matrix_realloc_rco(LIS_INT row,LIS_INT nnz,LIS_INT *** index,LIS_SCALAR *** value)223 LIS_INT lis_matrix_realloc_rco(LIS_INT row, LIS_INT nnz, LIS_INT ***index, LIS_SCALAR ***value)
224 {
225 LIS_INT **w_index;
226 LIS_SCALAR **w_value;
227
228 LIS_DEBUG_FUNC_IN;
229
230 w_index = *index;
231 w_value = *value;
232
233 w_index[row] = (LIS_INT *)lis_realloc(w_index[row],nnz*sizeof(LIS_INT));
234 if( w_index[row]==NULL )
235 {
236 LIS_SETERR_MEM(nnz*sizeof(LIS_INT));
237 return LIS_OUT_OF_MEMORY;
238 }
239 w_value[row] = (LIS_SCALAR *)lis_realloc(w_value[row],nnz*sizeof(LIS_SCALAR));
240 if( w_value[row]==NULL )
241 {
242 LIS_SETERR_MEM(nnz*sizeof(LIS_SCALAR));
243 return LIS_OUT_OF_MEMORY;
244 }
245 *index = w_index;
246 *value = w_value;
247
248 LIS_DEBUG_FUNC_OUT;
249 return LIS_SUCCESS;
250 }
251
252 #undef __FUNC__
253 #define __FUNC__ "lis_matrix_convert_rco2csr"
lis_matrix_convert_rco2csr(LIS_MATRIX Ain,LIS_MATRIX Aout)254 LIS_INT lis_matrix_convert_rco2csr(LIS_MATRIX Ain, LIS_MATRIX Aout)
255 {
256 LIS_INT i,j,k,n,nnz,err;
257 LIS_INT *ptr,*index;
258 LIS_SCALAR *value;
259
260 LIS_DEBUG_FUNC_IN;
261
262 ptr = NULL;
263 index = NULL;
264 value = NULL;
265
266 n = Ain->n;
267 nnz = 0;
268 #ifdef _OPENMP
269 #pragma omp parallel for reduction(+:nnz) private(i)
270 #endif
271 for(i=0;i<n;i++)
272 {
273 nnz += Ain->w_row[i];
274 }
275
276 err = lis_matrix_malloc_csr(n,nnz,&ptr,&index,&value);
277 if( err )
278 {
279 return err;
280 }
281
282 #ifdef _NUMA
283 #pragma omp parallel for private(i)
284 for(i=0;i<n+1;i++) ptr[i] = 0;
285 #else
286 ptr[0] = 0;
287 #endif
288 for(i=0;i<n;i++)
289 {
290 ptr[i+1] = ptr[i] + Ain->w_row[i];
291 }
292 #ifdef _OPENMP
293 #pragma omp parallel for private(i,j,k)
294 #endif
295 for(i=0;i<n;i++)
296 {
297 k = ptr[i];
298 for(j=0;j<Ain->w_row[i];j++)
299 {
300 index[k] = Ain->w_index[i][j];
301 value[k] = Ain->w_value[i][j];
302 k++;
303 }
304 }
305
306 err = lis_matrix_set_csr(nnz,ptr,index,value,Aout);
307 if( err )
308 {
309 lis_free2(3,ptr,index,value);
310 return err;
311 }
312 err = lis_matrix_assemble(Aout);
313 if( err )
314 {
315 lis_matrix_storage_destroy(Aout);
316 return err;
317 }
318
319 LIS_DEBUG_FUNC_OUT;
320 return LIS_SUCCESS;
321 }
322
323 #undef __FUNC__
324 #define __FUNC__ "lis_matrix_convert_rco2bsr"
lis_matrix_convert_rco2bsr(LIS_MATRIX Ain,LIS_MATRIX Aout)325 LIS_INT lis_matrix_convert_rco2bsr(LIS_MATRIX Ain, LIS_MATRIX Aout)
326 {
327 LIS_INT i,j,k,n,gn,nnz,bnnz,nr,nc,bnr,bnc,err;
328 LIS_INT ii,jj,kk,bj,jpos,ij,kv,bi;
329 LIS_INT *iw,*iw2;
330 LIS_INT *bptr,*bindex;
331 LIS_SCALAR *value;
332
333 LIS_DEBUG_FUNC_IN;
334
335 bnr = Ain->conv_bnr;
336 bnc = Ain->conv_bnc;
337 n = Ain->n;
338 gn = Ain->gn;
339 nr = 1 + (n-1)/bnr;
340 nc = 1 + (gn-1)/bnc;
341 bptr = NULL;
342 bindex = NULL;
343 value = NULL;
344 iw = NULL;
345 iw2 = NULL;
346
347
348 bptr = (LIS_INT *)lis_malloc( (nr+1)*sizeof(LIS_INT),"lis_matrix_convert_rco2bsr::bptr" );
349 if( bptr==NULL )
350 {
351 LIS_SETERR_MEM((nr+1)*sizeof(LIS_INT));
352 lis_free2(5,bptr,bindex,value,iw,iw2);
353 return LIS_OUT_OF_MEMORY;
354 }
355
356 #ifdef _OPENMP
357 #pragma omp parallel private(i,k,ii,j,bj,kk,ij,jj,iw,iw2,kv,jpos)
358 #endif
359 {
360 iw = (LIS_INT *)lis_malloc( nc*sizeof(LIS_INT),"lis_matrix_convert_rco2bsr::iw" );
361 iw2 = (LIS_INT *)lis_malloc( nc*sizeof(LIS_INT),"lis_matrix_convert_rco2bsr::iw2" );
362 memset(iw,0,nc*sizeof(LIS_INT));
363
364 #ifdef _OPENMP
365 #pragma omp for
366 #endif
367 for(i=0;i<nr;i++)
368 {
369 k = 0;
370 kk = bnr*i;
371 jj = 0;
372 for(ii=0;ii+kk<n&&ii<bnr;ii++)
373 {
374 for(j=0;j<Ain->w_row[kk+ii];j++)
375 {
376 bj = Ain->w_index[kk+ii][j]/bnc;
377 jpos = iw[bj];
378 if( jpos==0 )
379 {
380 iw[bj] = 1;
381 iw2[jj] = bj;
382 jj++;
383 }
384 }
385 }
386 for(bj=0;bj<jj;bj++)
387 {
388 k++;
389 ii = iw2[bj];
390 iw[ii]=0;
391 }
392 bptr[i+1] = k;
393 }
394 lis_free(iw);
395 lis_free(iw2);
396 }
397
398 bptr[0] = 0;
399 for(i=0;i<nr;i++)
400 {
401 bptr[i+1] += bptr[i];
402 }
403 bnnz = bptr[nr];
404 nnz = bnnz*bnr*bnc;
405
406 bindex = (LIS_INT *)lis_malloc( bnnz*sizeof(LIS_INT),"lis_matrix_convert_rco2bsr::bindex" );
407 if( bindex==NULL )
408 {
409 LIS_SETERR_MEM((nr+1)*sizeof(LIS_INT));
410 lis_free2(3,bptr,bindex,value);
411 return LIS_OUT_OF_MEMORY;
412 }
413 value = (LIS_SCALAR *)lis_malloc( nnz*sizeof(LIS_SCALAR),"lis_matrix_convert_rco2bsr::value" );
414 if( value==NULL )
415 {
416 LIS_SETERR_MEM(nnz*sizeof(LIS_SCALAR));
417 lis_free2(3,bptr,bindex,value);
418 return LIS_OUT_OF_MEMORY;
419 }
420
421 /* convert bsr */
422 #ifdef _OPENMP
423 #pragma omp parallel private(bi,i,ii,k,j,bj,jpos,kv,kk,ij,jj,iw)
424 #endif
425 {
426 iw = (LIS_INT *)lis_malloc( nc*sizeof(LIS_INT),"lis_matrix_convert_rco2bsr::iw" );
427 memset(iw,0,nc*sizeof(LIS_INT));
428
429 #ifdef _OPENMP
430 #pragma omp for
431 #endif
432 for(bi=0;bi<nr;bi++)
433 {
434 i = bi*bnr;
435 ii = 0;
436 kk = bptr[bi];
437 while( i+ii<n && ii<=bnr-1 )
438 {
439 for( k=0;k<Ain->w_row[i+ii];k++)
440 {
441 j = Ain->w_index[i+ii][k];
442 bj = j/bnc;
443 j = j%bnc;
444 jpos = iw[bj];
445 if( jpos==0 )
446 {
447 kv = kk * bnr * bnc;
448 iw[bj] = kv+1;
449 bindex[kk] = bj;
450 for(jj=0;jj<bnr*bnc;jj++) value[kv+jj] = 0.0;
451 ij = j*bnr + ii;
452 value[kv+ij] = Ain->w_value[i+ii][k];
453 kk = kk+1;
454 }
455 else
456 {
457 ij = j*bnr + ii;
458 value[jpos+ij-1] = Ain->w_value[i+ii][k];
459 }
460 }
461 ii = ii+1;
462 }
463 for(j=bptr[bi];j<bptr[bi+1];j++)
464 {
465 iw[bindex[j]] = 0;
466 }
467 }
468 lis_free(iw);
469 }
470
471 err = lis_matrix_set_bsr(bnr,bnc,bnnz,bptr,bindex,value,Aout);
472 if( err )
473 {
474 lis_free2(3,bptr,bindex,value);
475 return err;
476 }
477 err = lis_matrix_assemble(Aout);
478 if( err )
479 {
480 lis_matrix_storage_destroy(Aout);
481 return err;
482 }
483 LIS_DEBUG_FUNC_OUT;
484 return LIS_SUCCESS;
485 }
486
487 #undef __FUNC__
488 #define __FUNC__ "lis_matrix_convert_rco2csc"
lis_matrix_convert_rco2csc(LIS_MATRIX Ain,LIS_MATRIX Aout)489 LIS_INT lis_matrix_convert_rco2csc(LIS_MATRIX Ain, LIS_MATRIX Aout)
490 {
491 LIS_INT i,j,k,l,n,nnz,err;
492 LIS_INT *ptr,*index,*iw;
493 LIS_SCALAR *value;
494
495 LIS_DEBUG_FUNC_IN;
496
497 ptr = NULL;
498 index = NULL;
499 value = NULL;
500 iw = NULL;
501 n = Ain->n;
502
503
504 iw = (LIS_INT *)lis_malloc(n*sizeof(LIS_INT),"lis_matrix_convert_rco2csc::iw");
505 if( iw==NULL )
506 {
507 LIS_SETERR_MEM(n*sizeof(LIS_INT));
508 lis_free2(4,ptr,index,value,iw);
509 return LIS_OUT_OF_MEMORY;
510 }
511 ptr = (LIS_INT *)lis_malloc((n+1)*sizeof(LIS_INT),"lis_matrix_convert_rco2csc::ptr");
512 if( ptr==NULL )
513 {
514 LIS_SETERR_MEM((n+1)*sizeof(LIS_INT));
515 lis_free2(4,ptr,index,value,iw);
516 return LIS_OUT_OF_MEMORY;
517 }
518
519 for(i=0;i<n;i++) iw[i] = 0;
520 for(i=0;i<n;i++)
521 {
522 for(j=0;j<Ain->w_row[i];j++)
523 {
524 iw[Ain->w_index[i][j]]++;
525 }
526 }
527 ptr[0] = 0;
528 for(i=0;i<n;i++)
529 {
530 ptr[i+1] = ptr[i] + iw[i];
531 iw[i] = ptr[i];
532 }
533 nnz = ptr[n];
534
535 index = (LIS_INT *)lis_malloc( nnz*sizeof(LIS_INT),"lis_matrix_convert_rco2csc::index" );
536 if( index==NULL )
537 {
538 LIS_SETERR_MEM(nnz*sizeof(LIS_INT));
539 lis_free2(4,ptr,index,value,iw);
540 return LIS_OUT_OF_MEMORY;
541 }
542 value = (LIS_SCALAR *)lis_malloc( nnz*sizeof(LIS_SCALAR),"lis_matrix_convert_rco2csc::value" );
543 if( value==NULL )
544 {
545 LIS_SETERR_MEM(nnz*sizeof(LIS_SCALAR));
546 lis_free2(4,ptr,index,value,iw);
547 return LIS_OUT_OF_MEMORY;
548 }
549
550 for(i=0;i<n;i++)
551 {
552 for(j=0;j<Ain->w_row[i];j++)
553 {
554 k = Ain->w_index[i][j];
555 l = iw[k];
556 value[l] = Ain->w_value[i][j];
557 index[l] = i;
558 iw[k]++;
559 }
560 }
561
562 err = lis_matrix_set_csc(nnz,ptr,index,value,Aout);
563 if( err )
564 {
565 lis_free2(4,ptr,index,value,iw);
566 return err;
567 }
568 err = lis_matrix_assemble(Aout);
569 if( err )
570 {
571 lis_matrix_storage_destroy(Aout);
572 return err;
573 }
574
575 lis_free(iw);
576
577 LIS_DEBUG_FUNC_OUT;
578 return LIS_SUCCESS;
579 }
580
581 /*
582 #undef __FUNC__
583 #define __FUNC__ "lis_matrix_malloc_brco"
584 LIS_INT lis_matrix_malloc_brco(LIS_INT n, LIS_INT **row, LIS_INT ***index, LIS_SCALAR ***value)
585 {
586 LIS_INT i,j;
587 LIS_INT *w_row,**w_index;
588 LIS_SCALAR **w_value;
589
590 LIS_DEBUG_FUNC_IN;
591
592 w_row = NULL;
593 w_index = NULL;
594 w_value = NULL;
595
596 w_row = (LIS_INT *)lis_malloc( n*sizeof(LIS_INT),"lis_matrix_malloc_brco::w_row" );
597 if( w_row==NULL )
598 {
599 LIS_SETERR_MEM(n*sizeof(LIS_INT));
600 return LIS_OUT_OF_MEMORY;
601 }
602 w_index = (LIS_INT **)lis_malloc( n*sizeof(LIS_INT *),"lis_matrix_malloc_brco::w_index" );
603 if( w_index==NULL )
604 {
605 LIS_SETERR_MEM(n*sizeof(LIS_INT *));
606 lis_free2(3,w_row,w_index,w_value);
607 return LIS_OUT_OF_MEMORY;
608 }
609 w_value = (LIS_SCALAR **)lis_malloc( n*sizeof(LIS_SCALAR *),"lis_matrix_malloc_brco::w_value" );
610 if( w_value==NULL )
611 {
612 LIS_SETERR_MEM(n*sizeof(LIS_SCALAR *));
613 lis_free2(3,w_row,w_index,w_value);
614 return LIS_OUT_OF_MEMORY;
615 }
616 #pragma omp parallel for private(i)
617 for(i=0;i<n;i++)
618 {
619 w_row[i] = 0;
620 w_index[i] = NULL;
621 }
622 *row = w_row;
623 *index = w_index;
624 *value = w_value;
625
626 LIS_DEBUG_FUNC_OUT;
627 return LIS_SUCCESS;
628 }
629
630 #undef __FUNC__
631 #define __FUNC__ "lis_matrix_malloc_vrco"
632 LIS_INT lis_matrix_malloc_vrco(LIS_INT n, LIS_INT **nnz, LIS_INT **row, LIS_INT ***index, LIS_SCALAR ****value)
633 {
634 LIS_INT i,j;
635 LIS_INT *w_nnz, *w_row, **w_index;
636 LIS_SCALAR ***w_value;
637
638 LIS_DEBUG_FUNC_IN;
639
640 w_nnz = NULL;
641 w_row = NULL;
642 w_index = NULL;
643 w_value = NULL;
644
645 w_nnz = (LIS_INT *)lis_malloc( n*sizeof(LIS_INT),"lis_matrix_malloc_vrco::w_ptr" );
646 if( w_nnz==NULL )
647 {
648 LIS_SETERR_MEM(n*sizeof(LIS_INT));
649 return LIS_OUT_OF_MEMORY;
650 }
651 w_row = (LIS_INT *)lis_malloc( (n+1)*sizeof(LIS_INT),"lis_matrix_malloc_vrco::w_row" );
652 if( w_row==NULL )
653 {
654 LIS_SETERR_MEM((n+1)*sizeof(LIS_INT));
655 return LIS_OUT_OF_MEMORY;
656 }
657 w_index = (LIS_INT **)lis_malloc( n*sizeof(LIS_INT *),"lis_matrix_malloc_vrco::w_index" );
658 if( w_index==NULL )
659 {
660 LIS_SETERR_MEM(n*sizeof(LIS_INT *));
661 lis_free2(3,w_row,w_index,w_value);
662 return LIS_OUT_OF_MEMORY;
663 }
664 w_value = (LIS_SCALAR ***)lis_malloc( n*sizeof(LIS_SCALAR **),"lis_matrix_malloc_vrco::w_value" );
665 if( w_value==NULL )
666 {
667 LIS_SETERR_MEM(n*sizeof(LIS_SCALAR **));
668 lis_free2(3,w_row,w_index,w_value);
669 return LIS_OUT_OF_MEMORY;
670 }
671 #pragma omp parallel for private(i)
672 for(i=0;i<n;i++)
673 {
674 w_nnz[i] = 0;
675 w_index[i] = NULL;
676 }
677 *nnz = w_nnz;
678 *row = w_row;
679 *index = w_index;
680 *value = w_value;
681
682 LIS_DEBUG_FUNC_OUT;
683 return LIS_SUCCESS;
684 }
685
686 */
687