1 #ifdef PETSC_RCS_HEADER 2 static char vcid[] = "$Id: shellcnv.c,v 1.7 1999/05/04 20:31:58 balay Exp balay $"; 3 #endif 4 5 6 #include "petsc.h" 7 #include "src/mat/matimpl.h" /*I "mat.h" I*/ 8 #include "src/vec/vecimpl.h" 9 10 #undef __FUNC__ 11 #define __FUNC__ "MatConvert_Shell" 12 int MatConvert_Shell(Mat oldmat,MatType newtype, Mat *mat) 13 { 14 Vec in,out; 15 int ierr,i,M,m,size,*rows,start,end; 16 MPI_Comm comm; 17 Scalar *array,zero = 0.0,one = 1.0; 18 19 PetscFunctionBegin; 20 PetscValidHeaderSpecific(oldmat,MAT_COOKIE); 21 PetscValidPointer(mat); 22 23 if (newtype != MATSEQDENSE || newtype != MATMPIDENSE) { 24 SETERRQ(PETSC_ERR_SUP,1,"Can only convert shell matrices to dense currently"); 25 } 26 comm = oldmat->comm; 27 28 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 29 30 ierr = MatGetOwnershipRange(oldmat,&start,&end);CHKERRQ(ierr); 31 ierr = VecCreateMPI(comm,end-start,PETSC_DECIDE,&in);CHKERRQ(ierr); 32 ierr = VecDuplicate(in,&out);CHKERRQ(ierr); 33 ierr = VecGetSize(in,&M);CHKERRQ(ierr); 34 ierr = VecGetLocalSize(in,&m);CHKERRQ(ierr); 35 rows = (int *) PetscMalloc( (m+1)*sizeof(int) );CHKPTRQ(rows); 36 for ( i=0; i<m; i++ ) {rows[i] = start + i;} 37 38 if (size == 1) { 39 ierr = MatCreateSeqDense(comm,M,M,PETSC_NULL,mat);CHKERRQ(ierr); 40 } else { 41 ierr = MatCreateMPIDense(comm,m,M,M,M,PETSC_NULL,mat);CHKERRQ(ierr); 42 /* ierr = MatCreateMPIAIJ(comm,m,m,M,M,0,0,0,0,mat);CHKERRQ(ierr); */ 43 } 44 45 for ( i=0; i<M; i++ ) { 46 ierr = VecSet(&zero,in);CHKERRQ(ierr); 47 ierr = VecSetValues(in,1,&i,&one,INSERT_VALUES);CHKERRQ(ierr); 48 ierr = VecAssemblyBegin(in);CHKERRQ(ierr); 49 ierr = VecAssemblyEnd(in);CHKERRQ(ierr); 50 51 ierr = MatMult(oldmat,in,out);CHKERRQ(ierr); 52 53 ierr = VecGetArray(out,&array);CHKERRQ(ierr); 54 ierr = MatSetValues(*mat,m,rows,1,&i,array,INSERT_VALUES);CHKERRQ(ierr); 55 ierr = VecRestoreArray(out,&array);CHKERRQ(ierr); 56 57 } 58 ierr = PetscFree(rows);CHKERRQ(ierr); 59 ierr = VecDestroy(in);CHKERRQ(ierr); 60 ierr = VecDestroy(out);CHKERRQ(ierr); 61 ierr = MatAssemblyBegin(*mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 62 ierr = MatAssemblyEnd(*mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 63 PetscFunctionReturn(0); 64 } 65 66 67 68