Actual source code: convert.c
1: #define PETSCMAT_DLL
3: #include src/mat/matimpl.h
7: /*
8: MatConvert_Basic - Converts from any input format to another format. For
9: parallel formats, the new matrix distribution is determined by PETSc.
11: Does not do preallocation so in general will be slow
12: */
13: PetscErrorCode MatConvert_Basic(Mat mat, MatType newtype,MatReuse reuse,Mat *newmat)
14: {
15: Mat M;
16: const PetscScalar *vwork;
17: PetscErrorCode ierr;
18: PetscInt i,nz,m,n,rstart,rend,lm,ln;
19: const PetscInt *cwork;
22: MatGetSize(mat,&m,&n);
23: MatGetLocalSize(mat,&lm,&ln);
25: if (ln == n) ln = PETSC_DECIDE; /* try to preserve column ownership */
26:
27: MatCreate(mat->comm,&M);
28: MatSetSizes(M,lm,ln,m,n);
29: MatSetType(M,newtype);
31: MatGetOwnershipRange(mat,&rstart,&rend);
32: for (i=rstart; i<rend; i++) {
33: MatGetRow(mat,i,&nz,&cwork,&vwork);
34: MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES);
35: MatRestoreRow(mat,i,&nz,&cwork,&vwork);
36: }
37: MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);
38: MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);
40: if (reuse == MAT_REUSE_MATRIX) {
41: MatHeaderReplace(mat,M);
42: } else {
43: *newmat = M;
44: }
45: return(0);
46: }