From 9d081d43a38fc721f6ed01d5d9d11a6dc2a1b10a Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Mon, 25 Nov 2013 19:52:35 +0000 Subject: [PATCH] psblas3: base/serial/impl/psb_c_coo_impl.f90 base/serial/impl/psb_c_csc_impl.f90 base/serial/impl/psb_c_csr_impl.f90 base/serial/impl/psb_d_coo_impl.f90 base/serial/impl/psb_d_csc_impl.f90 base/serial/impl/psb_d_csr_impl.f90 base/serial/impl/psb_s_coo_impl.f90 base/serial/impl/psb_s_csc_impl.f90 base/serial/impl/psb_s_csr_impl.f90 base/serial/impl/psb_z_coo_impl.f90 base/serial/impl/psb_z_csc_impl.f90 base/serial/impl/psb_z_csr_impl.f90 Do not call coo%fix when not necessary --- base/serial/impl/psb_c_coo_impl.f90 | 13 ++-- base/serial/impl/psb_c_csc_impl.f90 | 3 +- base/serial/impl/psb_c_csr_impl.f90 | 97 ++++++++++++++++++++++++++--- base/serial/impl/psb_d_coo_impl.f90 | 24 ++++--- base/serial/impl/psb_d_csc_impl.f90 | 3 +- base/serial/impl/psb_d_csr_impl.f90 | 97 ++++++++++++++++++++++++++--- base/serial/impl/psb_s_coo_impl.f90 | 13 ++-- base/serial/impl/psb_s_csc_impl.f90 | 3 +- base/serial/impl/psb_s_csr_impl.f90 | 97 ++++++++++++++++++++++++++--- base/serial/impl/psb_z_coo_impl.f90 | 13 ++-- base/serial/impl/psb_z_csc_impl.f90 | 3 +- base/serial/impl/psb_z_csr_impl.f90 | 97 ++++++++++++++++++++++++++--- 12 files changed, 396 insertions(+), 67 deletions(-) diff --git a/base/serial/impl/psb_c_coo_impl.f90 b/base/serial/impl/psb_c_coo_impl.f90 index 6d757d01..1c7b83c2 100644 --- a/base/serial/impl/psb_c_coo_impl.f90 +++ b/base/serial/impl/psb_c_coo_impl.f90 @@ -398,6 +398,8 @@ subroutine psb_c_coo_allocate_mnnz(m,n,a,nz) call a%set_triangle(.false.) call a%set_unit(.false.) call a%set_dupl(psb_dupl_def_) + ! An empty matrix is sorted! + call a%set_sorted(.true.) end if if (info /= psb_success_) goto 9999 call psb_erractionrestore(err_act) @@ -2974,7 +2976,8 @@ subroutine psb_c_cp_coo_to_coo(a,b,info) b%ja(1:nz) = a%ja(1:nz) b%val(1:nz) = a%val(1:nz) - call b%fix(info) + + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3020,7 +3023,7 @@ subroutine psb_c_cp_coo_from_coo(a,b,info) a%ja(1:nz) = b%ja(1:nz) a%val(1:nz) = b%val(1:nz) - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 @@ -3132,14 +3135,13 @@ subroutine psb_c_mv_coo_to_coo(a,b,info) info = psb_success_ b%psb_c_base_sparse_mat = a%psb_c_base_sparse_mat call b%set_nzeros(a%get_nzeros()) - call b%reallocate(a%get_nzeros()) call move_alloc(a%ia, b%ia) call move_alloc(a%ja, b%ja) call move_alloc(a%val, b%val) call a%free() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3177,13 +3179,12 @@ subroutine psb_c_mv_coo_from_coo(a,b,info) info = psb_success_ a%psb_c_base_sparse_mat = b%psb_c_base_sparse_mat call a%set_nzeros(b%get_nzeros()) - call a%reallocate(b%get_nzeros()) call move_alloc(b%ia , a%ia ) call move_alloc(b%ja , a%ja ) call move_alloc(b%val, a%val ) call b%free() - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 diff --git a/base/serial/impl/psb_c_csc_impl.f90 b/base/serial/impl/psb_c_csc_impl.f90 index 3279608d..b37cb50b 100644 --- a/base/serial/impl/psb_c_csc_impl.f90 +++ b/base/serial/impl/psb_c_csc_impl.f90 @@ -2273,7 +2273,8 @@ subroutine psb_c_cp_csc_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes + ! We need to make a copy because mv_from will have to + ! sort in column-major order. call tmp%cp_from_coo(b,info) if (info == psb_success_) call a%mv_from_coo(tmp,info) diff --git a/base/serial/impl/psb_c_csr_impl.f90 b/base/serial/impl/psb_c_csr_impl.f90 index 010696d7..aa1d1374 100644 --- a/base/serial/impl/psb_c_csr_impl.f90 +++ b/base/serial/impl/psb_c_csr_impl.f90 @@ -2802,9 +2802,90 @@ subroutine psb_c_cp_csr_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes - call tmp%cp_from_coo(b,info) - if (info == psb_success_) call a%mv_from_coo(tmp,info) + + if (.not.b%is_sorted()) then + ! This is to have fix_coo called behind the scenes + call tmp%cp_from_coo(b,info) + if (info /= psb_success_) return + + nr = tmp%get_nrows() + nc = tmp%get_ncols() + nza = tmp%get_nzeros() + + a%psb_c_base_sparse_mat = tmp%psb_c_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call move_alloc(tmp%ia,itemp) + call move_alloc(tmp%ja,a%ja) + call move_alloc(tmp%val,a%val) + call psb_realloc(max(nr+1,nc+1),a%irp,info) + call tmp%free() + else + + if (info /= psb_success_) return + + nr = b%get_nrows() + nc = b%get_ncols() + nza = b%get_nzeros() + + a%psb_c_base_sparse_mat = b%psb_c_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call psb_safe_ab_cpy(b%ia,itemp,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%ja,a%ja,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%val,a%val,info) + if (info /= psb_success_) call psb_realloc(max(nr+1,nc+1),a%irp,info) + + endif + + if (nza <= 0) then + a%irp(:) = 1 + else + a%irp(1) = 1 + if (nr < itemp(nza)) then + write(debug_unit,*) trim(name),': RWSHR=.false. : ',& + &nr,itemp(nza),' Expect trouble!' + info = 12 + end if + + j = 1 + i = 1 + irw = itemp(j) + + outer: do + inner: do + if (i >= irw) exit inner + if (i>nr) then + write(debug_unit,*) trim(name),& + & 'Strange situation: i>nr ',i,nr,j,nza,irw + exit outer + end if + a%irp(i+1) = a%irp(i) + i = i + 1 + end do inner + j = j + 1 + if (j > nza) exit + if (itemp(j) /= irw) then + a%irp(i+1) = j + irw = itemp(j) + i = i + 1 + endif + if (i>nr) exit + enddo outer + ! + ! Cleanup empty rows at the end + ! + if (j /= (nza+1)) then + write(debug_unit,*) trim(name),': Problem from loop :',j,nza + info = 13 + endif + do + if (i>nr) exit + a%irp(i+1) = j + i = i + 1 + end do + + endif end subroutine psb_c_cp_csr_from_coo @@ -2845,8 +2926,8 @@ subroutine psb_c_cp_csr_to_coo(a,b,info) end do end do call b%set_nzeros(a%get_nzeros()) - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_c_cp_csr_to_coo @@ -2888,8 +2969,8 @@ subroutine psb_c_mv_csr_to_coo(a,b,info) end do end do call a%free() - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_c_mv_csr_to_coo @@ -2920,7 +3001,7 @@ subroutine psb_c_mv_csr_from_coo(a,b,info) debug_level = psb_get_debug_level() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) return nr = b%get_nrows() diff --git a/base/serial/impl/psb_d_coo_impl.f90 b/base/serial/impl/psb_d_coo_impl.f90 index b4ef5f7a..cb5fe775 100644 --- a/base/serial/impl/psb_d_coo_impl.f90 +++ b/base/serial/impl/psb_d_coo_impl.f90 @@ -398,6 +398,8 @@ subroutine psb_d_coo_allocate_mnnz(m,n,a,nz) call a%set_triangle(.false.) call a%set_unit(.false.) call a%set_dupl(psb_dupl_def_) + ! An empty matrix is sorted! + call a%set_sorted(.true.) end if if (info /= psb_success_) goto 9999 call psb_erractionrestore(err_act) @@ -2725,7 +2727,7 @@ contains subroutine d_coo_srch_upd(nz,ia,ja,val,a,& & imin,imax,jmin,jmax,info,gtl) - + use psb_const_mod use psb_realloc_mod use psb_string_mod @@ -2974,7 +2976,8 @@ subroutine psb_d_cp_coo_to_coo(a,b,info) b%ja(1:nz) = a%ja(1:nz) b%val(1:nz) = a%val(1:nz) - call b%fix(info) + + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3020,7 +3023,7 @@ subroutine psb_d_cp_coo_from_coo(a,b,info) a%ja(1:nz) = b%ja(1:nz) a%val(1:nz) = b%val(1:nz) - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 @@ -3132,14 +3135,13 @@ subroutine psb_d_mv_coo_to_coo(a,b,info) info = psb_success_ b%psb_d_base_sparse_mat = a%psb_d_base_sparse_mat call b%set_nzeros(a%get_nzeros()) - call b%reallocate(a%get_nzeros()) call move_alloc(a%ia, b%ia) call move_alloc(a%ja, b%ja) call move_alloc(a%val, b%val) call a%free() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3177,13 +3179,12 @@ subroutine psb_d_mv_coo_from_coo(a,b,info) info = psb_success_ a%psb_d_base_sparse_mat = b%psb_d_base_sparse_mat call a%set_nzeros(b%get_nzeros()) - call a%reallocate(b%get_nzeros()) call move_alloc(b%ia , a%ia ) call move_alloc(b%ja , a%ja ) call move_alloc(b%val, a%val ) call b%free() - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 @@ -3412,7 +3413,6 @@ subroutine psb_d_fix_coo_inner(nzin,dupl,ia,ja,val,nzout,info,idir) use psb_d_base_mat_mod, psb_protect_name => psb_d_fix_coo_inner use psb_string_mod use psb_ip_reord_mod - use psb_sort_mod implicit none integer(psb_ipk_), intent(in) :: nzin, dupl @@ -3423,12 +3423,10 @@ subroutine psb_d_fix_coo_inner(nzin,dupl,ia,ja,val,nzout,info,idir) !locals integer(psb_ipk_), allocatable :: iaux(:) integer(psb_ipk_) :: nza, nzl,iret,idir_, dupl_ - integer(psb_ipk_) :: i,j, irw, icl, err_act, ixp,ki,kx + integer(psb_ipk_) :: i,j, irw, icl, err_act integer(psb_ipk_) :: debug_level, debug_unit integer(psb_ipk_) :: ierr(5) character(len=20) :: name = 'psb_fixcoo' - real(psb_dpk_), allocatable :: vtx(:) - integer(psb_ipk_), allocatable :: itx(:), jtx(:) info = psb_success_ @@ -3458,7 +3456,6 @@ subroutine psb_d_fix_coo_inner(nzin,dupl,ia,ja,val,nzout,info,idir) case(0) ! Row major order - call msort_up(nzin,ia(1:),iaux(1:),iret) if (iret == 0) & & call psb_ip_reord(nzin,val,ia,ja,iaux) @@ -3473,7 +3470,7 @@ subroutine psb_d_fix_coo_inner(nzin,dupl,ia,ja,val,nzout,info,idir) call msort_up(nzl,ja(i:),iaux(1:),iret) if (iret == 0) & & call psb_ip_reord(nzl,val(i:i+nzl-1),& - & ja(i:i+nzl-1),iaux) + & ia(i:i+nzl-1),ja(i:i+nzl-1),iaux) i = j enddo @@ -3538,6 +3535,7 @@ subroutine psb_d_fix_coo_inner(nzin,dupl,ia,ja,val,nzout,info,idir) info =-7 end select + if(debug_level >= psb_debug_serial_)& & write(debug_unit,*) trim(name),': end second loop' diff --git a/base/serial/impl/psb_d_csc_impl.f90 b/base/serial/impl/psb_d_csc_impl.f90 index 5583c7ed..b6d958f7 100644 --- a/base/serial/impl/psb_d_csc_impl.f90 +++ b/base/serial/impl/psb_d_csc_impl.f90 @@ -2273,7 +2273,8 @@ subroutine psb_d_cp_csc_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes + ! We need to make a copy because mv_from will have to + ! sort in column-major order. call tmp%cp_from_coo(b,info) if (info == psb_success_) call a%mv_from_coo(tmp,info) diff --git a/base/serial/impl/psb_d_csr_impl.f90 b/base/serial/impl/psb_d_csr_impl.f90 index 96dee37d..feed133c 100644 --- a/base/serial/impl/psb_d_csr_impl.f90 +++ b/base/serial/impl/psb_d_csr_impl.f90 @@ -2802,9 +2802,90 @@ subroutine psb_d_cp_csr_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes - call tmp%cp_from_coo(b,info) - if (info == psb_success_) call a%mv_from_coo(tmp,info) + + if (.not.b%is_sorted()) then + ! This is to have fix_coo called behind the scenes + call tmp%cp_from_coo(b,info) + if (info /= psb_success_) return + + nr = tmp%get_nrows() + nc = tmp%get_ncols() + nza = tmp%get_nzeros() + + a%psb_d_base_sparse_mat = tmp%psb_d_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call move_alloc(tmp%ia,itemp) + call move_alloc(tmp%ja,a%ja) + call move_alloc(tmp%val,a%val) + call psb_realloc(max(nr+1,nc+1),a%irp,info) + call tmp%free() + else + + if (info /= psb_success_) return + + nr = b%get_nrows() + nc = b%get_ncols() + nza = b%get_nzeros() + + a%psb_d_base_sparse_mat = b%psb_d_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call psb_safe_ab_cpy(b%ia,itemp,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%ja,a%ja,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%val,a%val,info) + if (info /= psb_success_) call psb_realloc(max(nr+1,nc+1),a%irp,info) + + endif + + if (nza <= 0) then + a%irp(:) = 1 + else + a%irp(1) = 1 + if (nr < itemp(nza)) then + write(debug_unit,*) trim(name),': RWSHR=.false. : ',& + &nr,itemp(nza),' Expect trouble!' + info = 12 + end if + + j = 1 + i = 1 + irw = itemp(j) + + outer: do + inner: do + if (i >= irw) exit inner + if (i>nr) then + write(debug_unit,*) trim(name),& + & 'Strange situation: i>nr ',i,nr,j,nza,irw + exit outer + end if + a%irp(i+1) = a%irp(i) + i = i + 1 + end do inner + j = j + 1 + if (j > nza) exit + if (itemp(j) /= irw) then + a%irp(i+1) = j + irw = itemp(j) + i = i + 1 + endif + if (i>nr) exit + enddo outer + ! + ! Cleanup empty rows at the end + ! + if (j /= (nza+1)) then + write(debug_unit,*) trim(name),': Problem from loop :',j,nza + info = 13 + endif + do + if (i>nr) exit + a%irp(i+1) = j + i = i + 1 + end do + + endif end subroutine psb_d_cp_csr_from_coo @@ -2845,8 +2926,8 @@ subroutine psb_d_cp_csr_to_coo(a,b,info) end do end do call b%set_nzeros(a%get_nzeros()) - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_d_cp_csr_to_coo @@ -2888,8 +2969,8 @@ subroutine psb_d_mv_csr_to_coo(a,b,info) end do end do call a%free() - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_d_mv_csr_to_coo @@ -2920,7 +3001,7 @@ subroutine psb_d_mv_csr_from_coo(a,b,info) debug_level = psb_get_debug_level() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) return nr = b%get_nrows() diff --git a/base/serial/impl/psb_s_coo_impl.f90 b/base/serial/impl/psb_s_coo_impl.f90 index 8c84691a..e2e81854 100644 --- a/base/serial/impl/psb_s_coo_impl.f90 +++ b/base/serial/impl/psb_s_coo_impl.f90 @@ -398,6 +398,8 @@ subroutine psb_s_coo_allocate_mnnz(m,n,a,nz) call a%set_triangle(.false.) call a%set_unit(.false.) call a%set_dupl(psb_dupl_def_) + ! An empty matrix is sorted! + call a%set_sorted(.true.) end if if (info /= psb_success_) goto 9999 call psb_erractionrestore(err_act) @@ -2974,7 +2976,8 @@ subroutine psb_s_cp_coo_to_coo(a,b,info) b%ja(1:nz) = a%ja(1:nz) b%val(1:nz) = a%val(1:nz) - call b%fix(info) + + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3020,7 +3023,7 @@ subroutine psb_s_cp_coo_from_coo(a,b,info) a%ja(1:nz) = b%ja(1:nz) a%val(1:nz) = b%val(1:nz) - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 @@ -3132,14 +3135,13 @@ subroutine psb_s_mv_coo_to_coo(a,b,info) info = psb_success_ b%psb_s_base_sparse_mat = a%psb_s_base_sparse_mat call b%set_nzeros(a%get_nzeros()) - call b%reallocate(a%get_nzeros()) call move_alloc(a%ia, b%ia) call move_alloc(a%ja, b%ja) call move_alloc(a%val, b%val) call a%free() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3177,13 +3179,12 @@ subroutine psb_s_mv_coo_from_coo(a,b,info) info = psb_success_ a%psb_s_base_sparse_mat = b%psb_s_base_sparse_mat call a%set_nzeros(b%get_nzeros()) - call a%reallocate(b%get_nzeros()) call move_alloc(b%ia , a%ia ) call move_alloc(b%ja , a%ja ) call move_alloc(b%val, a%val ) call b%free() - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 diff --git a/base/serial/impl/psb_s_csc_impl.f90 b/base/serial/impl/psb_s_csc_impl.f90 index 44aa91e7..59cce665 100644 --- a/base/serial/impl/psb_s_csc_impl.f90 +++ b/base/serial/impl/psb_s_csc_impl.f90 @@ -2273,7 +2273,8 @@ subroutine psb_s_cp_csc_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes + ! We need to make a copy because mv_from will have to + ! sort in column-major order. call tmp%cp_from_coo(b,info) if (info == psb_success_) call a%mv_from_coo(tmp,info) diff --git a/base/serial/impl/psb_s_csr_impl.f90 b/base/serial/impl/psb_s_csr_impl.f90 index 4d281cb3..b6d6ed43 100644 --- a/base/serial/impl/psb_s_csr_impl.f90 +++ b/base/serial/impl/psb_s_csr_impl.f90 @@ -2802,9 +2802,90 @@ subroutine psb_s_cp_csr_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes - call tmp%cp_from_coo(b,info) - if (info == psb_success_) call a%mv_from_coo(tmp,info) + + if (.not.b%is_sorted()) then + ! This is to have fix_coo called behind the scenes + call tmp%cp_from_coo(b,info) + if (info /= psb_success_) return + + nr = tmp%get_nrows() + nc = tmp%get_ncols() + nza = tmp%get_nzeros() + + a%psb_s_base_sparse_mat = tmp%psb_s_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call move_alloc(tmp%ia,itemp) + call move_alloc(tmp%ja,a%ja) + call move_alloc(tmp%val,a%val) + call psb_realloc(max(nr+1,nc+1),a%irp,info) + call tmp%free() + else + + if (info /= psb_success_) return + + nr = b%get_nrows() + nc = b%get_ncols() + nza = b%get_nzeros() + + a%psb_s_base_sparse_mat = b%psb_s_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call psb_safe_ab_cpy(b%ia,itemp,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%ja,a%ja,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%val,a%val,info) + if (info /= psb_success_) call psb_realloc(max(nr+1,nc+1),a%irp,info) + + endif + + if (nza <= 0) then + a%irp(:) = 1 + else + a%irp(1) = 1 + if (nr < itemp(nza)) then + write(debug_unit,*) trim(name),': RWSHR=.false. : ',& + &nr,itemp(nza),' Expect trouble!' + info = 12 + end if + + j = 1 + i = 1 + irw = itemp(j) + + outer: do + inner: do + if (i >= irw) exit inner + if (i>nr) then + write(debug_unit,*) trim(name),& + & 'Strange situation: i>nr ',i,nr,j,nza,irw + exit outer + end if + a%irp(i+1) = a%irp(i) + i = i + 1 + end do inner + j = j + 1 + if (j > nza) exit + if (itemp(j) /= irw) then + a%irp(i+1) = j + irw = itemp(j) + i = i + 1 + endif + if (i>nr) exit + enddo outer + ! + ! Cleanup empty rows at the end + ! + if (j /= (nza+1)) then + write(debug_unit,*) trim(name),': Problem from loop :',j,nza + info = 13 + endif + do + if (i>nr) exit + a%irp(i+1) = j + i = i + 1 + end do + + endif end subroutine psb_s_cp_csr_from_coo @@ -2845,8 +2926,8 @@ subroutine psb_s_cp_csr_to_coo(a,b,info) end do end do call b%set_nzeros(a%get_nzeros()) - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_s_cp_csr_to_coo @@ -2888,8 +2969,8 @@ subroutine psb_s_mv_csr_to_coo(a,b,info) end do end do call a%free() - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_s_mv_csr_to_coo @@ -2920,7 +3001,7 @@ subroutine psb_s_mv_csr_from_coo(a,b,info) debug_level = psb_get_debug_level() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) return nr = b%get_nrows() diff --git a/base/serial/impl/psb_z_coo_impl.f90 b/base/serial/impl/psb_z_coo_impl.f90 index c14ade32..d52fcca5 100644 --- a/base/serial/impl/psb_z_coo_impl.f90 +++ b/base/serial/impl/psb_z_coo_impl.f90 @@ -398,6 +398,8 @@ subroutine psb_z_coo_allocate_mnnz(m,n,a,nz) call a%set_triangle(.false.) call a%set_unit(.false.) call a%set_dupl(psb_dupl_def_) + ! An empty matrix is sorted! + call a%set_sorted(.true.) end if if (info /= psb_success_) goto 9999 call psb_erractionrestore(err_act) @@ -2974,7 +2976,8 @@ subroutine psb_z_cp_coo_to_coo(a,b,info) b%ja(1:nz) = a%ja(1:nz) b%val(1:nz) = a%val(1:nz) - call b%fix(info) + + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3020,7 +3023,7 @@ subroutine psb_z_cp_coo_from_coo(a,b,info) a%ja(1:nz) = b%ja(1:nz) a%val(1:nz) = b%val(1:nz) - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 @@ -3132,14 +3135,13 @@ subroutine psb_z_mv_coo_to_coo(a,b,info) info = psb_success_ b%psb_z_base_sparse_mat = a%psb_z_base_sparse_mat call b%set_nzeros(a%get_nzeros()) - call b%reallocate(a%get_nzeros()) call move_alloc(a%ia, b%ia) call move_alloc(a%ja, b%ja) call move_alloc(a%val, b%val) call a%free() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) goto 9999 @@ -3177,13 +3179,12 @@ subroutine psb_z_mv_coo_from_coo(a,b,info) info = psb_success_ a%psb_z_base_sparse_mat = b%psb_z_base_sparse_mat call a%set_nzeros(b%get_nzeros()) - call a%reallocate(b%get_nzeros()) call move_alloc(b%ia , a%ia ) call move_alloc(b%ja , a%ja ) call move_alloc(b%val, a%val ) call b%free() - call a%fix(info) + if (.not.a%is_sorted()) call a%fix(info) if (info /= psb_success_) goto 9999 diff --git a/base/serial/impl/psb_z_csc_impl.f90 b/base/serial/impl/psb_z_csc_impl.f90 index f6157df3..78a1f439 100644 --- a/base/serial/impl/psb_z_csc_impl.f90 +++ b/base/serial/impl/psb_z_csc_impl.f90 @@ -2273,7 +2273,8 @@ subroutine psb_z_cp_csc_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes + ! We need to make a copy because mv_from will have to + ! sort in column-major order. call tmp%cp_from_coo(b,info) if (info == psb_success_) call a%mv_from_coo(tmp,info) diff --git a/base/serial/impl/psb_z_csr_impl.f90 b/base/serial/impl/psb_z_csr_impl.f90 index 5c6a1a59..482ef039 100644 --- a/base/serial/impl/psb_z_csr_impl.f90 +++ b/base/serial/impl/psb_z_csr_impl.f90 @@ -2802,9 +2802,90 @@ subroutine psb_z_cp_csr_from_coo(a,b,info) character(len=20) :: name info = psb_success_ - ! This is to have fix_coo called behind the scenes - call tmp%cp_from_coo(b,info) - if (info == psb_success_) call a%mv_from_coo(tmp,info) + + if (.not.b%is_sorted()) then + ! This is to have fix_coo called behind the scenes + call tmp%cp_from_coo(b,info) + if (info /= psb_success_) return + + nr = tmp%get_nrows() + nc = tmp%get_ncols() + nza = tmp%get_nzeros() + + a%psb_z_base_sparse_mat = tmp%psb_z_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call move_alloc(tmp%ia,itemp) + call move_alloc(tmp%ja,a%ja) + call move_alloc(tmp%val,a%val) + call psb_realloc(max(nr+1,nc+1),a%irp,info) + call tmp%free() + else + + if (info /= psb_success_) return + + nr = b%get_nrows() + nc = b%get_ncols() + nza = b%get_nzeros() + + a%psb_z_base_sparse_mat = b%psb_z_base_sparse_mat + + ! Dirty trick: call move_alloc to have the new data allocated just once. + call psb_safe_ab_cpy(b%ia,itemp,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%ja,a%ja,info) + if (info /= psb_success_) call psb_safe_ab_cpy(b%val,a%val,info) + if (info /= psb_success_) call psb_realloc(max(nr+1,nc+1),a%irp,info) + + endif + + if (nza <= 0) then + a%irp(:) = 1 + else + a%irp(1) = 1 + if (nr < itemp(nza)) then + write(debug_unit,*) trim(name),': RWSHR=.false. : ',& + &nr,itemp(nza),' Expect trouble!' + info = 12 + end if + + j = 1 + i = 1 + irw = itemp(j) + + outer: do + inner: do + if (i >= irw) exit inner + if (i>nr) then + write(debug_unit,*) trim(name),& + & 'Strange situation: i>nr ',i,nr,j,nza,irw + exit outer + end if + a%irp(i+1) = a%irp(i) + i = i + 1 + end do inner + j = j + 1 + if (j > nza) exit + if (itemp(j) /= irw) then + a%irp(i+1) = j + irw = itemp(j) + i = i + 1 + endif + if (i>nr) exit + enddo outer + ! + ! Cleanup empty rows at the end + ! + if (j /= (nza+1)) then + write(debug_unit,*) trim(name),': Problem from loop :',j,nza + info = 13 + endif + do + if (i>nr) exit + a%irp(i+1) = j + i = i + 1 + end do + + endif end subroutine psb_z_cp_csr_from_coo @@ -2845,8 +2926,8 @@ subroutine psb_z_cp_csr_to_coo(a,b,info) end do end do call b%set_nzeros(a%get_nzeros()) - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_z_cp_csr_to_coo @@ -2888,8 +2969,8 @@ subroutine psb_z_mv_csr_to_coo(a,b,info) end do end do call a%free() - call b%fix(info) - + call b%set_sorted() + call b%set_asb() end subroutine psb_z_mv_csr_to_coo @@ -2920,7 +3001,7 @@ subroutine psb_z_mv_csr_from_coo(a,b,info) debug_level = psb_get_debug_level() - call b%fix(info) + if (.not.b%is_sorted()) call b%fix(info) if (info /= psb_success_) return nr = b%get_nrows()