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 #ifdef _OPENMP
43 #include <omp.h>
44 #endif
45 #ifdef USE_MPI
46 #include <mpi.h>
47 #endif
48 #include "lislib.h"
49
50 /************************************************
51 * lis_matrix_init
52 * lis_matrix_check
53 * lis_matrix_create
54 * lis_matrix_destroy
55 * lis_matrix_duplicate
56 * lis_matrix_malloc
57 * lis_matrix_assemble
58 * lis_matrix_is_assembled
59 * lis_matrix_get_range
60 * lis_matrix_get_size
61 * lis_matrix_get_nnz
62 * lis_matrix_get_type
63 * lis_matrix_get_option
64 * lis_matrix_set_type
65 * lis_matrix_set_value
66 * lis_matrix_set_option
67 ************************************************/
68
69 #undef __FUNC__
70 #define __FUNC__ "lis_matrix_init"
lis_matrix_init(LIS_MATRIX * Amat)71 LIS_INT lis_matrix_init(LIS_MATRIX *Amat)
72 {
73
74 LIS_DEBUG_FUNC_IN;
75
76 memset(*Amat,0,sizeof(struct LIS_MATRIX_STRUCT));
77
78 (*Amat)->label = LIS_LABEL_MATRIX;
79 (*Amat)->matrix_type = LIS_MATRIX_CSR;
80 (*Amat)->status = LIS_MATRIX_DECIDING_SIZE;
81 /* (*Amat)->status = LIS_MATRIX_NULL;*/
82 (*Amat)->w_annz = LIS_MATRIX_W_ANNZ;
83 (*Amat)->conv_bnr = 2;
84 (*Amat)->conv_bnc = 2;
85 (*Amat)->is_destroy = LIS_TRUE;
86
87 LIS_DEBUG_FUNC_OUT;
88 return LIS_SUCCESS;
89 }
90
91 #undef __FUNC__
92 #define __FUNC__ "lis_matrix_check"
lis_matrix_check(LIS_MATRIX A,LIS_INT level)93 LIS_INT lis_matrix_check(LIS_MATRIX A, LIS_INT level)
94 {
95 LIS_DEBUG_FUNC_IN;
96
97 switch( level )
98 {
99 case LIS_MATRIX_CHECK_NULL:
100 if( !lis_is_malloc(A) )
101 {
102 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is undefined\n");
103 return LIS_ERR_ILL_ARG;
104 }
105 break;
106 case LIS_MATRIX_CHECK_SIZE:
107 if( !lis_is_malloc(A) )
108 {
109 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is undefined\n");
110 return LIS_ERR_ILL_ARG;
111 }
112 if( A->status==LIS_MATRIX_DECIDING_SIZE )
113 {
114 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix size is undefined\n");
115 return LIS_ERR_ILL_ARG;
116 }
117 break;
118 case LIS_MATRIX_CHECK_TYPE:
119 if( !lis_is_malloc(A) )
120 {
121 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is undefined\n");
122 return LIS_ERR_ILL_ARG;
123 }
124 if( A->status==LIS_MATRIX_DECIDING_SIZE )
125 {
126 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix size is undefined\n");
127 return LIS_ERR_ILL_ARG;
128 }
129 /* only want this condition satisfied for the case of non-zero number of DOFs */
130 /* on the current process */
131 if( A->status==LIS_MATRIX_NULL && A->n>0 )
132 {
133 /* as the status is actually changed from ...NULL to ...ASSEMBLING in lis_matrix_set_value, */
134 /* ...NULL here really means that the components are undefined, not the type */
135 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix components undefined\n");
136 return LIS_ERR_ILL_ARG;
137 }
138 break;
139 case LIS_MATRIX_CHECK_NOT_ASSEMBLED:
140 if( !lis_is_malloc(A) )
141 {
142 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is undefined\n");
143 return LIS_ERR_ILL_ARG;
144 }
145 if( A->status!=LIS_MATRIX_DECIDING_SIZE && A->status!=LIS_MATRIX_NULL && A->status!=LIS_MATRIX_ASSEMBLING )
146 {
147 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A has already been assembled\n");
148 return LIS_ERR_ILL_ARG;
149 }
150 break;
151 case LIS_MATRIX_CHECK_SET:
152 if( !lis_is_malloc(A) )
153 {
154 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is undefined\n");
155 return LIS_ERR_ILL_ARG;
156 }
157 if( A->status==LIS_MATRIX_DECIDING_SIZE )
158 {
159 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix size is undefined\n");
160 return LIS_ERR_ILL_ARG;
161 }
162 if( A->status!=LIS_MATRIX_NULL )
163 {
164 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A has already been assembled\n");
165 return LIS_ERR_ILL_ARG;
166 }
167 break;
168 default:
169 if( !lis_is_malloc(A) )
170 {
171 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is undefined\n");
172 return LIS_ERR_ILL_ARG;
173 }
174 if( A->status==LIS_MATRIX_DECIDING_SIZE )
175 {
176 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix size is undefined\n");
177 return LIS_ERR_ILL_ARG;
178 }
179 /* if( A->status==LIS_MATRIX_NULL )*/
180 if( A->status==LIS_MATRIX_NULL && A->n>0 )
181 {
182 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix type is undefined\n");
183 return LIS_ERR_ILL_ARG;
184 }
185 /* if( A->status<=LIS_MATRIX_ASSEMBLING )*/
186 if( A->status<=LIS_MATRIX_ASSEMBLING && A->n>0 )
187 {
188 LIS_SETERR(LIS_ERR_ILL_ARG,"matrix A is assembling\n");
189 return LIS_ERR_ILL_ARG;
190 }
191 break;
192 }
193 LIS_DEBUG_FUNC_OUT;
194 return LIS_SUCCESS;
195 }
196
197
198 #undef __FUNC__
199 #define __FUNC__ "lis_matrix_create"
lis_matrix_create(LIS_Comm comm,LIS_MATRIX * Amat)200 LIS_INT lis_matrix_create(LIS_Comm comm, LIS_MATRIX *Amat)
201 {
202 LIS_INT nprocs,my_rank;
203 int int_nprocs,int_my_rank;
204
205 LIS_DEBUG_FUNC_IN;
206
207 *Amat = NULL;
208
209
210 *Amat = (LIS_MATRIX)lis_malloc( sizeof(struct LIS_MATRIX_STRUCT),"lis_matrix_create::Amat" );
211 if( NULL==*Amat )
212 {
213 LIS_SETERR_MEM(sizeof(struct LIS_MATRIX_STRUCT));
214 return LIS_OUT_OF_MEMORY;
215 }
216 lis_matrix_init(Amat);
217
218 #ifdef USE_MPI
219 #if _DEBUG
220 printf("c_comm = %d f_comm = %d\n",MPI_COMM_WORLD,comm);
221 #endif
222 MPI_Comm_size(comm,&int_nprocs);
223 MPI_Comm_rank(comm,&int_my_rank);
224 nprocs = int_nprocs;
225 my_rank = int_my_rank;
226 #else
227 nprocs = 1;
228 my_rank = 0;
229 #endif
230 (*Amat)->comm = comm;
231 (*Amat)->nprocs = nprocs;
232 (*Amat)->my_rank = my_rank;
233
234 LIS_DEBUG_FUNC_OUT;
235 return LIS_SUCCESS;
236 }
237
238 #undef __FUNC__
239 #define __FUNC__ "lis_matrix_set_size"
lis_matrix_set_size(LIS_MATRIX Amat,LIS_INT local_n,LIS_INT global_n)240 LIS_INT lis_matrix_set_size(LIS_MATRIX Amat, LIS_INT local_n, LIS_INT global_n)
241 {
242 LIS_INT nprocs,my_rank;
243 int int_nprocs,int_my_rank;
244 LIS_INT is,ie;
245 LIS_INT err;
246 LIS_INT *ranges;
247
248 LIS_DEBUG_FUNC_IN;
249
250 err = lis_matrix_check(Amat,LIS_MATRIX_CHECK_NULL);
251 if( err ) return err;
252
253 if( global_n>0 && local_n>global_n )
254 {
255 LIS_SETERR2(LIS_ERR_ILL_ARG,"local n(=%D) is larger than global n(=%D)\n",local_n,global_n);
256 return LIS_ERR_ILL_ARG;
257 }
258 if( local_n<0 || global_n<0 )
259 {
260 LIS_SETERR2(LIS_ERR_ILL_ARG,"local n(=%D) or global n(=%D) are less than 0\n",local_n,global_n);
261 return LIS_ERR_ILL_ARG;
262 }
263 /* actually want to allow this case for MPI, but not otherwise */
264 #ifndef USE_MPI
265 if( local_n==0 && global_n==0 )
266 {
267 LIS_SETERR2(LIS_ERR_ILL_ARG,"local n(=%D) and global n(=%D) are 0\n",local_n,global_n);
268 return LIS_ERR_ILL_ARG;
269 }
270 #endif
271 #ifdef USE_MPI
272 MPI_Comm_size(Amat->comm,&int_nprocs);
273 nprocs = int_nprocs;
274 /* change the logic to something a little more direct, i.e., that will not */
275 /* throw an error if local_n=0 . . . */
276 if( global_n>0 && global_n<nprocs )
277 {
278 LIS_SETERR2(LIS_ERR_ILL_ARG,"global n(=%D) is smaller than nprocs(=%D)\n",global_n,nprocs);
279 return LIS_ERR_ILL_ARG;
280 }
281 #endif
282
283 err = lis_ranges_create(Amat->comm,&local_n,&global_n,&ranges,&is,&ie,&nprocs,&my_rank);
284 if( err )
285 {
286 return err;
287 }
288
289 Amat->status = LIS_MATRIX_NULL;
290 Amat->ranges = ranges;
291 Amat->n = local_n;
292 Amat->gn = global_n;
293 Amat->np = local_n;
294 Amat->my_rank = my_rank;
295 Amat->nprocs = nprocs;
296 Amat->is = is;
297 Amat->ie = ie;
298
299 LIS_DEBUG_FUNC_OUT;
300 return LIS_SUCCESS;
301 }
302
303 #undef __FUNC__
304 #define __FUNC__ "lis_matrix_LU_create"
lis_matrix_LU_create(LIS_MATRIX A)305 LIS_INT lis_matrix_LU_create(LIS_MATRIX A)
306 {
307 LIS_MATRIX_CORE L,U;
308
309 LIS_DEBUG_FUNC_IN;
310
311 L = NULL;
312 U = NULL;
313
314 L = (LIS_MATRIX_CORE)lis_calloc(sizeof(struct LIS_MATRIX_CORE_STRUCT),"lis_matrix_LU_create::L");
315 if( L==NULL )
316 {
317 LIS_SETERR_MEM(sizeof(struct LIS_MATRIX_CORE_STRUCT));
318 return LIS_OUT_OF_MEMORY;
319 }
320 U = (LIS_MATRIX_CORE)lis_calloc(sizeof(struct LIS_MATRIX_CORE_STRUCT),"lis_matrix_LU_create::U");
321 if( U==NULL )
322 {
323 LIS_SETERR_MEM(sizeof(struct LIS_MATRIX_CORE_STRUCT));
324 lis_free(L);
325 return LIS_OUT_OF_MEMORY;
326 }
327
328 A->L = L;
329 A->U = U;
330
331 LIS_DEBUG_FUNC_OUT;
332 return LIS_SUCCESS;
333 }
334
335 #undef __FUNC__
336 #define __FUNC__ "lis_matrix_LU_destroy"
lis_matrix_LU_destroy(LIS_MATRIX_CORE Amat)337 LIS_INT lis_matrix_LU_destroy(LIS_MATRIX_CORE Amat)
338 {
339 LIS_DEBUG_FUNC_IN;
340
341 if( Amat )
342 {
343 if( Amat->ptr ) lis_free( Amat->ptr );
344 if( Amat->row ) lis_free( Amat->row );
345 if( Amat->col ) lis_free( Amat->col );
346 if( Amat->index) lis_free( Amat->index );
347 if( Amat->bptr ) lis_free( Amat->bptr );
348 if( Amat->bindex ) lis_free( Amat->bindex );
349 if( Amat->value ) lis_free( Amat->value );
350 if( Amat->work ) lis_free( Amat->work );
351 lis_free(Amat);
352 }
353
354 LIS_DEBUG_FUNC_OUT;
355 return LIS_SUCCESS;
356 }
357
358 #undef __FUNC__
359 #define __FUNC__ "lis_matrix_DLU_destroy"
lis_matrix_DLU_destroy(LIS_MATRIX Amat)360 LIS_INT lis_matrix_DLU_destroy(LIS_MATRIX Amat)
361 {
362 LIS_DEBUG_FUNC_IN;
363
364 if( Amat )
365 {
366 if( Amat->D ) lis_matrix_diag_destroy(Amat->D);
367 if( Amat->L ) lis_matrix_LU_destroy(Amat->L);
368 if( Amat->U ) lis_matrix_LU_destroy(Amat->U);
369 Amat->D = NULL;
370 Amat->L = NULL;
371 Amat->U = NULL;
372 }
373
374 LIS_DEBUG_FUNC_OUT;
375 return LIS_SUCCESS;
376 }
377
378 #undef __FUNC__
379 #define __FUNC__ "lis_matrix_storage_destroy"
lis_matrix_storage_destroy(LIS_MATRIX Amat)380 LIS_INT lis_matrix_storage_destroy(LIS_MATRIX Amat)
381 {
382
383 LIS_DEBUG_FUNC_IN;
384
385 if( Amat )
386 {
387 if( Amat->is_destroy )
388 {
389 if( Amat->ptr ) lis_free( Amat->ptr );
390 if( Amat->row ) lis_free( Amat->row );
391 if( Amat->col ) lis_free( Amat->col );
392 if( Amat->index) lis_free( Amat->index );
393 if( Amat->bptr ) lis_free( Amat->bptr );
394 if( Amat->bindex ) lis_free( Amat->bindex );
395 if( Amat->value ) lis_free( Amat->value );
396 if( Amat->work ) lis_free( Amat->work );
397 if( Amat->conv_row ) lis_free( Amat->conv_row );
398 if( Amat->conv_col ) lis_free( Amat->conv_col );
399 if( Amat->w_nnz ) lis_free( Amat->w_nnz );
400 if( Amat->w_row )
401 {
402 #if 0
403 for(i=0;i<Amat->n;i++)
404 {
405 lis_free(Amat->w_index[i]);
406 lis_free(Amat->w_value[i]);
407 }
408 #else
409 lis_free_mat(Amat);
410 #endif
411 lis_free2(3,Amat->w_row,Amat->w_index,Amat->w_value);
412 }
413 }
414 Amat->row = NULL;
415 Amat->col = NULL;
416 Amat->ptr = NULL;
417 Amat->index = NULL;
418 Amat->bptr = NULL;
419 Amat->bindex = NULL;
420 Amat->value = NULL;
421 Amat->work = NULL;
422 Amat->conv_row = NULL;
423 Amat->conv_col = NULL;
424 Amat->w_nnz = NULL;
425 Amat->w_row = NULL;
426 Amat->w_index = NULL;
427 Amat->w_value = NULL;
428 Amat->v_value = NULL;
429 }
430
431 LIS_DEBUG_FUNC_OUT;
432 return LIS_SUCCESS;
433 }
434
435 #undef __FUNC__
436 #define __FUNC__ "lis_matrix_destroy"
lis_matrix_destroy(LIS_MATRIX Amat)437 LIS_INT lis_matrix_destroy(LIS_MATRIX Amat)
438 {
439 LIS_DEBUG_FUNC_IN;
440
441 if( lis_is_malloc(Amat) )
442 {
443 if( !Amat->is_fallocated ) lis_matrix_storage_destroy(Amat);
444 lis_matrix_DLU_destroy(Amat);
445 lis_matrix_diag_destroy(Amat->WD);
446 if( Amat->l2g_map ) lis_free( Amat->l2g_map );
447 if( Amat->commtable ) lis_commtable_destroy( Amat->commtable );
448 if( Amat->ranges ) lis_free( Amat->ranges );
449 lis_free(Amat);
450 }
451
452 LIS_DEBUG_FUNC_OUT;
453 return LIS_SUCCESS;
454 }
455
456 #undef __FUNC__
457 #define __FUNC__ "lis_matrix_duplicate"
lis_matrix_duplicate(LIS_MATRIX Ain,LIS_MATRIX * Aout)458 LIS_INT lis_matrix_duplicate(LIS_MATRIX Ain, LIS_MATRIX *Aout)
459 {
460 LIS_INT err;
461 LIS_INT nprocs;
462 #ifdef USE_MPI
463 LIS_INT i;
464 LIS_INT *ranges;
465 LIS_INT *l2g_map;
466 #endif
467
468 LIS_DEBUG_FUNC_IN;
469
470 err = lis_matrix_check(Ain,LIS_MATRIX_CHECK_ALL);
471 if( err ) return err;
472
473 nprocs = Ain->nprocs;
474 *Aout = NULL;
475
476 *Aout = (LIS_MATRIX)lis_malloc( sizeof(struct LIS_MATRIX_STRUCT),"lis_matrix_duplicate::Aout" );
477 if( NULL==*Aout )
478 {
479 LIS_SETERR_MEM(sizeof(struct LIS_MATRIX_STRUCT));
480 return LIS_OUT_OF_MEMORY;
481 }
482 lis_matrix_init(Aout);
483
484 #ifdef USE_MPI
485 l2g_map = (LIS_INT *)lis_malloc( (Ain->np-Ain->n)*sizeof(LIS_INT),"lis_matrix_duplicate::l2g_map" );
486 if( l2g_map==NULL )
487 {
488 LIS_SETERR_MEM((Ain->np-Ain->n)*sizeof(LIS_INT));
489 lis_matrix_destroy(*Aout);
490 *Aout = NULL;
491 return LIS_OUT_OF_MEMORY;
492 }
493 memcpy(l2g_map,Ain->l2g_map,(Ain->np-Ain->n)*sizeof(LIS_INT));
494 ranges = (LIS_INT *)lis_malloc( (nprocs+1)*sizeof(LIS_INT),"lis_matrix_duplicate::range" );
495 if( ranges==NULL )
496 {
497 LIS_SETERR_MEM((nprocs+1)*sizeof(LIS_INT));
498 lis_free(l2g_map);
499 lis_matrix_destroy(*Aout);
500 *Aout = NULL;
501 return LIS_OUT_OF_MEMORY;
502 }
503 for(i=0;i<nprocs+1;i++) ranges[i] = Ain->ranges[i];
504 (*Aout)->ranges = ranges;
505 (*Aout)->l2g_map = l2g_map;
506 #else
507 (*Aout)->ranges = NULL;
508 #endif
509
510 (*Aout)->status = LIS_MATRIX_NULL;
511 (*Aout)->is_block = Ain->is_block;
512 (*Aout)->n = Ain->n;
513 (*Aout)->gn = Ain->gn;
514 (*Aout)->np = Ain->np;
515 (*Aout)->comm = Ain->comm;
516 (*Aout)->my_rank = Ain->my_rank;
517 (*Aout)->nprocs = Ain->nprocs;
518 (*Aout)->is = Ain->is;
519 (*Aout)->ie = Ain->ie;
520 (*Aout)->origin = Ain->origin;
521 (*Aout)->is_destroy = Ain->is_destroy;
522
523 #ifdef USE_MPI
524 err = lis_commtable_duplicateM(Ain,Aout);
525 if( err )
526 {
527 lis_matrix_destroy(*Aout);
528 *Aout = NULL;
529 return err;
530 }
531 (*Aout)->is_comm = LIS_TRUE;
532 LIS_DEBUG_FUNC_OUT;
533 return LIS_SUCCESS;
534 #else
535 LIS_DEBUG_FUNC_OUT;
536 return LIS_SUCCESS;
537 #endif
538 }
539
540 #undef __FUNC__
541 #define __FUNC__ "lis_matrix_copy_struct"
lis_matrix_copy_struct(LIS_MATRIX Ain,LIS_MATRIX Aout)542 LIS_INT lis_matrix_copy_struct(LIS_MATRIX Ain, LIS_MATRIX Aout)
543 {
544 LIS_DEBUG_FUNC_IN;
545
546 memcpy(Aout,Ain,sizeof(struct LIS_MATRIX_STRUCT));
547
548 LIS_DEBUG_FUNC_OUT;
549 return LIS_SUCCESS;
550 }
551
552 #undef __FUNC__
553 #define __FUNC__ "lis_matrix_get_range"
lis_matrix_get_range(LIS_MATRIX A,LIS_INT * is,LIS_INT * ie)554 LIS_INT lis_matrix_get_range(LIS_MATRIX A, LIS_INT *is, LIS_INT *ie)
555 {
556 LIS_INT err;
557
558 LIS_DEBUG_FUNC_IN;
559
560 err = lis_matrix_check(A,LIS_MATRIX_CHECK_SIZE);
561 if( err ) return err;
562
563 *is = A->is;
564 *ie = A->ie;
565
566 LIS_DEBUG_FUNC_OUT;
567 return LIS_SUCCESS;
568 }
569
570 #undef __FUNC__
571 #define __FUNC__ "lis_matrix_get_size"
lis_matrix_get_size(LIS_MATRIX A,LIS_INT * local_n,LIS_INT * global_n)572 LIS_INT lis_matrix_get_size(LIS_MATRIX A, LIS_INT *local_n, LIS_INT *global_n)
573 {
574 LIS_INT err;
575
576 LIS_DEBUG_FUNC_IN;
577
578 err = lis_matrix_check(A,LIS_MATRIX_CHECK_SIZE);
579 if( err ) return err;
580
581 *local_n = A->n;
582 *global_n = A->gn;
583
584 LIS_DEBUG_FUNC_OUT;
585 return LIS_SUCCESS;
586 }
587
588 #undef __FUNC__
589 #define __FUNC__ "lis_matrix_get_nnz"
lis_matrix_get_nnz(LIS_MATRIX A,LIS_INT * nnz)590 LIS_INT lis_matrix_get_nnz(LIS_MATRIX A, LIS_INT *nnz)
591 {
592 LIS_INT err;
593
594 LIS_DEBUG_FUNC_IN;
595
596 err = lis_matrix_check(A,LIS_MATRIX_CHECK_SIZE);
597 if( err ) return err;
598
599 *nnz = A->nnz;
600
601 LIS_DEBUG_FUNC_OUT;
602 return LIS_SUCCESS;
603 }
604
605 #undef __FUNC__
606 #define __FUNC__ "lis_matrix_assemble"
lis_matrix_assemble(LIS_MATRIX A)607 LIS_INT lis_matrix_assemble(LIS_MATRIX A)
608 {
609 LIS_INT err;
610 LIS_INT matrix_type;
611 LIS_MATRIX B;
612
613 LIS_DEBUG_FUNC_IN;
614
615 err = lis_matrix_check(A,LIS_MATRIX_CHECK_TYPE);
616 if( err ) return err;
617
618 matrix_type = A->matrix_type;
619 if( A->status==LIS_MATRIX_ASSEMBLING )
620 {
621 A->matrix_type = LIS_MATRIX_RCO;
622 A->status = LIS_MATRIX_RCO;
623 #ifdef USE_MPI
624 err = lis_matrix_g2l(A);
625 if( err ) return err;
626 err = lis_commtable_create(A);
627 if( err ) return err;
628 A->is_comm = LIS_TRUE;
629 #endif
630 err = lis_matrix_duplicate(A,&B);
631 if( err ) return err;
632 lis_matrix_set_type(B,matrix_type);
633 err = lis_matrix_convert(A,B);
634 if( err ) return err;
635 lis_matrix_storage_destroy(A);
636 lis_matrix_DLU_destroy(A);
637 lis_matrix_diag_destroy(A->WD);
638 if( A->l2g_map ) lis_free( A->l2g_map );
639 if( A->commtable ) lis_commtable_destroy( A->commtable );
640 if( A->ranges ) lis_free( A->ranges );
641 err = lis_matrix_copy_struct(B,A);
642 if( err ) return err;
643 lis_free(B);
644 if( A->matrix_type==LIS_MATRIX_JAD )
645 {
646 A->work = (LIS_SCALAR *)lis_malloc(A->n*sizeof(LIS_SCALAR),"lis_matrix_assemble::A->work");
647 if( A->work==NULL )
648 {
649 LIS_SETERR_MEM(A->n*sizeof(LIS_SCALAR));
650 return LIS_OUT_OF_MEMORY;
651 }
652 }
653
654 LIS_DEBUG_FUNC_OUT;
655 return LIS_SUCCESS;
656 }
657 else if( A->status<0 && A->n>0 )
658 {
659 A->status = -A->status;
660 A->matrix_type = A->status;
661 if( A->matrix_type==LIS_MATRIX_JAD )
662 {
663 A->work = (LIS_SCALAR *)lis_malloc(A->n*sizeof(LIS_SCALAR),"lis_matrix_assemble::A->work");
664 if( A->work==NULL )
665 {
666 LIS_SETERR_MEM(A->n*sizeof(LIS_SCALAR));
667 return LIS_OUT_OF_MEMORY;
668 }
669 }
670 }
671 #ifdef USE_MPI
672 if( !A->is_pmat )
673 {
674 if( A->l2g_map==NULL )
675 {
676 err = lis_matrix_g2l(A);
677 if( err ) return err;
678 }
679 if( A->commtable==NULL )
680 {
681 err = lis_commtable_create(A);
682 if( err ) return err;
683 }
684 }
685 #endif
686 LIS_DEBUG_FUNC_OUT;
687 return LIS_SUCCESS;
688 }
689
690 #undef __FUNC__
691 #define __FUNC__ "lis_matrix_is_assembled"
lis_matrix_is_assembled(LIS_MATRIX A)692 LIS_INT lis_matrix_is_assembled(LIS_MATRIX A)
693 {
694 if( A->status!=LIS_MATRIX_NULL ) return !LIS_SUCCESS;
695 else return LIS_SUCCESS;
696 }
697
698 #undef __FUNC__
699 #define __FUNC__ "lis_matrix_set_value"
lis_matrix_set_value(LIS_INT flag,LIS_INT i,LIS_INT j,LIS_SCALAR value,LIS_MATRIX A)700 LIS_INT lis_matrix_set_value(LIS_INT flag, LIS_INT i, LIS_INT j, LIS_SCALAR value, LIS_MATRIX A)
701 {
702 LIS_INT n,gn,is,k;
703 LIS_INT err;
704
705 LIS_DEBUG_FUNC_IN;
706
707 err = lis_matrix_check(A,LIS_MATRIX_CHECK_NOT_ASSEMBLED);
708 if( err ) return err;
709
710 n = A->n;
711 gn = A->gn;
712 is = A->is;
713 if( A->origin )
714 {
715 i--;
716 j--;
717 }
718 if( i<0 || j<0 )
719 {
720 k = 0;
721 if( A->origin )
722 {
723 i++;
724 j++;
725 k++;
726 }
727 LIS_SETERR3(LIS_ERR_ILL_ARG,"i(=%D) or j(=%D) are less than %D\n",i,j,k);
728 return LIS_ERR_ILL_ARG;
729 }
730 if( i>=gn || j>=gn )
731 {
732 if( A->origin )
733 {
734 i++;
735 j++;
736 }
737 LIS_SETERR3(LIS_ERR_ILL_ARG,"i(=%D) or j(=%D) are larger than global n=(%D)\n",i,j,gn);
738 return LIS_ERR_ILL_ARG;
739 }
740
741 if( A->status==LIS_MATRIX_NULL )
742 {
743 if( A->w_nnz==NULL )
744 {
745 A->w_nnz = (LIS_INT *)lis_malloc(n*sizeof(LIS_INT),"lis_matrix_set_value::A->w_nnz");
746 if( A->w_nnz==NULL )
747 {
748 LIS_SETERR_MEM(n*sizeof(LIS_INT));
749 return LIS_OUT_OF_MEMORY;
750 }
751 for(k=0;k<n;k++) A->w_nnz[k] = A->w_annz;
752 }
753 err = lis_matrix_malloc_rco(n,A->w_nnz,&A->w_row,&A->w_index,&A->w_value);
754 if( err )
755 {
756 lis_free(A->w_nnz);
757 return err;
758 }
759 A->status = LIS_MATRIX_ASSEMBLING;
760 A->is_copy = LIS_TRUE;
761 }
762 if( A->w_nnz[i-is]==A->w_row[i-is] )
763 {
764 A->w_nnz[i-is] += A->w_annz;
765 err = lis_matrix_realloc_rco(i-is,A->w_nnz[i-is],&A->w_index,&A->w_value);
766 if( err )
767 {
768 for(k=0;k<n;k++)
769 {
770 lis_free(A->w_index[k]);
771 lis_free(A->w_value[k]);
772 }
773 lis_free2(4,A->w_nnz,A->w_row,A->w_index,A->w_value);
774 return err;
775 }
776 }
777 for(k=0;k<A->w_row[i-is];k++)
778 {
779 if( A->w_index[i-is][k]==j ) break;
780 }
781 if( k<A->w_row[i-is] )
782 {
783 if( flag==LIS_INS_VALUE )
784 {
785 A->w_value[i-is][k] = value;
786 }
787 else
788 {
789 A->w_value[i-is][k] += value;
790 }
791 }
792 else
793 {
794 k = A->w_row[i-is]++;
795 A->w_index[i-is][k] = j;
796 A->w_value[i-is][k] = value;
797 }
798
799 LIS_DEBUG_FUNC_OUT;
800 return LIS_SUCCESS;
801 }
802
803 /*NEH support for extended "solve_kernel" workflow*/
804 #undef __FUNC__
805 #define __FUNC__ "lis_matrix_psd_set_value"
lis_matrix_psd_set_value(LIS_INT flag,LIS_INT i,LIS_INT j,LIS_SCALAR value,LIS_MATRIX A)806 LIS_INT lis_matrix_psd_set_value(LIS_INT flag, LIS_INT i, LIS_INT j, LIS_SCALAR value, LIS_MATRIX A)
807 {
808 LIS_INT err;
809
810 LIS_DEBUG_FUNC_IN;
811
812 /* presuming that the matrix is assembled */
813 /* still, if not, default defined below */
814 switch(A->status) {
815 case LIS_MATRIX_CSR:
816 err=lis_matrix_psd_set_value_csr(flag,i,j,value,A);
817 if (err) return err;
818 break;
819 case LIS_MATRIX_CSC:
820 err = LIS_ERR_NOT_IMPLEMENTED;
821 return err;
822 case LIS_MATRIX_MSR:
823 err = LIS_ERR_NOT_IMPLEMENTED;
824 return err;
825 case LIS_MATRIX_DIA:
826 err = LIS_ERR_NOT_IMPLEMENTED;
827 return err;
828 case LIS_MATRIX_ELL:
829 err = LIS_ERR_NOT_IMPLEMENTED;
830 return err;
831 case LIS_MATRIX_JAD:
832 err = LIS_ERR_NOT_IMPLEMENTED;
833 return err;
834 case LIS_MATRIX_BSR:
835 err = LIS_ERR_NOT_IMPLEMENTED;
836 return err;
837 case LIS_MATRIX_BSC:
838 err = LIS_ERR_NOT_IMPLEMENTED;
839 return err;
840 case LIS_MATRIX_VBR:
841 err = LIS_ERR_NOT_IMPLEMENTED;
842 return err;
843 case LIS_MATRIX_COO:
844 err = LIS_ERR_NOT_IMPLEMENTED;
845 return err;
846 case LIS_MATRIX_DNS:
847 err = LIS_ERR_NOT_IMPLEMENTED;
848 return err;
849 default:
850 err = LIS_ERR_NOT_IMPLEMENTED;
851 return err;
852 }
853
854 LIS_DEBUG_FUNC_OUT;
855 return LIS_SUCCESS;
856 }
857
858 #undef __FUNC__
859 #define __FUNC__ "lis_matrix_set_values"
lis_matrix_set_values(LIS_INT flag,LIS_INT n,LIS_SCALAR value[],LIS_MATRIX A)860 LIS_INT lis_matrix_set_values(LIS_INT flag, LIS_INT n, LIS_SCALAR value[], LIS_MATRIX A)
861 {
862 LIS_INT i,j;
863 LIS_DEBUG_FUNC_IN;
864
865 for (i=0;i<n;i++)
866 {
867 for (j=0;j<n;j++)
868 {
869 lis_matrix_set_value(flag, i, j, value[i*n+j], A);
870 }
871 }
872
873 LIS_DEBUG_FUNC_OUT;
874 return LIS_SUCCESS;
875 }
876
877 #undef __FUNC__
878 #define __FUNC__ "lis_matrix_set_type"
lis_matrix_set_type(LIS_MATRIX A,LIS_INT matrix_type)879 LIS_INT lis_matrix_set_type(LIS_MATRIX A, LIS_INT matrix_type)
880 {
881 LIS_INT err;
882
883 LIS_DEBUG_FUNC_IN;
884
885 err = lis_matrix_check(A,LIS_MATRIX_CHECK_NOT_ASSEMBLED);
886 if( err ) return err;
887 if( matrix_type < LIS_MATRIX_CSR || matrix_type > LIS_MATRIX_DNS )
888 {
889 LIS_SETERR2(LIS_ERR_ILL_ARG,"matrix_type is %D (Set between 1 to %D)\n", matrix_type, LIS_MATRIX_DNS);
890 return LIS_ERR_ILL_ARG;
891 }
892
893 A->matrix_type = matrix_type;
894
895 LIS_DEBUG_FUNC_OUT;
896 return LIS_SUCCESS;
897 }
898
899 #undef __FUNC__
900 #define __FUNC__ "lis_matrix_get_type"
lis_matrix_get_type(LIS_MATRIX A,LIS_INT * matrix_type)901 LIS_INT lis_matrix_get_type(LIS_MATRIX A, LIS_INT *matrix_type)
902 {
903 LIS_INT err;
904
905 LIS_DEBUG_FUNC_IN;
906
907 err = lis_matrix_check(A,LIS_MATRIX_CHECK_NULL);
908 if( err ) return err;
909
910 *matrix_type = A->matrix_type;
911
912 LIS_DEBUG_FUNC_OUT;
913 return LIS_SUCCESS;
914 }
915
916 #undef __FUNC__
917 #define __FUNC__ "lis_matrix_set_destroyflag"
lis_matrix_set_destroyflag(LIS_MATRIX A,LIS_INT flag)918 LIS_INT lis_matrix_set_destroyflag(LIS_MATRIX A, LIS_INT flag)
919 {
920 LIS_INT err;
921
922 LIS_DEBUG_FUNC_IN;
923
924 err = lis_matrix_check(A,LIS_MATRIX_CHECK_NULL);
925 if( err ) return err;
926
927 if( flag )
928 {
929 A->is_destroy = LIS_TRUE;
930 }
931 else
932 {
933 A->is_destroy = LIS_FALSE;
934 }
935
936 LIS_DEBUG_FUNC_OUT;
937 return LIS_SUCCESS;
938 }
939
940 #undef __FUNC__
941 #define __FUNC__ "lis_matrix_get_destroyflag"
lis_matrix_get_destroyflag(LIS_MATRIX A,LIS_INT * flag)942 LIS_INT lis_matrix_get_destroyflag(LIS_MATRIX A, LIS_INT *flag)
943 {
944 LIS_INT err;
945
946 LIS_DEBUG_FUNC_IN;
947
948 err = lis_matrix_check(A,LIS_MATRIX_CHECK_NULL);
949 if( err ) return err;
950
951 *flag = A->is_destroy;
952
953 LIS_DEBUG_FUNC_OUT;
954 return LIS_SUCCESS;
955 }
956
957 #undef __FUNC__
958 #define __FUNC__ "lis_matrix_malloc"
lis_matrix_malloc(LIS_MATRIX A,LIS_INT nnz_row,LIS_INT nnz[])959 LIS_INT lis_matrix_malloc(LIS_MATRIX A, LIS_INT nnz_row, LIS_INT nnz[])
960 {
961 LIS_INT n,k;
962 LIS_INT err;
963
964 LIS_DEBUG_FUNC_IN;
965
966 err = lis_matrix_check(A,LIS_MATRIX_CHECK_NOT_ASSEMBLED);
967 if( err ) return err;
968
969 n = A->n;
970 if( A->w_nnz==NULL )
971 {
972 A->w_nnz = (LIS_INT *)lis_malloc(n*sizeof(LIS_INT),"lis_matrix_malloc::A->w_nnz");
973 if( A->w_nnz==NULL )
974 {
975 LIS_SETERR_MEM(n*sizeof(LIS_INT));
976 return LIS_OUT_OF_MEMORY;
977 }
978 }
979 if( nnz==NULL )
980 {
981 A->w_annz = nnz_row;
982 for(k=0;k<n;k++) A->w_nnz[k] = nnz_row;
983 }
984 else
985 {
986 for(k=0;k<n;k++) A->w_nnz[k] = nnz[k];
987 }
988
989 LIS_DEBUG_FUNC_OUT;
990 return LIS_SUCCESS;
991 }
992
993 #undef __FUNC__
994 #define __FUNC__ "lis_matrix_unset"
lis_matrix_unset(LIS_MATRIX A)995 LIS_INT lis_matrix_unset(LIS_MATRIX A)
996 {
997 LIS_INT err;
998
999 LIS_DEBUG_FUNC_IN;
1000
1001 err = lis_matrix_check(A,LIS_MATRIX_CHECK_SIZE);
1002 if( err ) return err;
1003
1004 if( A->is_copy )
1005 {
1006 lis_matrix_storage_destroy(A);
1007 }
1008 A->row = NULL;
1009 A->col = NULL;
1010 A->ptr = NULL;
1011 A->index = NULL;
1012 A->bptr = NULL;
1013 A->bindex = NULL;
1014 A->value = NULL;
1015 A->is_copy = LIS_FALSE;
1016 A->status = LIS_MATRIX_NULL;
1017
1018 LIS_DEBUG_FUNC_OUT;
1019 return LIS_SUCCESS;
1020 }
1021
1022