diff --git a/base/serial/psb_crwextd.f90 b/base/serial/psb_crwextd.f90 index e19af57a..d801428a 100644 --- a/base/serial/psb_crwextd.f90 +++ b/base/serial/psb_crwextd.f90 @@ -35,17 +35,18 @@ ! ! We have a problem here: 1. How to handle well all the formats? ! 2. What should we do with rowscale? Does it only -! apply when a%fida='COO' ?????? +! apply when a%fmt()='COO' ?????? ! ! subroutine psb_crwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_crwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_crwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr type(psb_cspmat_type), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info type(psb_cspmat_type), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -89,23 +90,20 @@ subroutine psb_crwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue - call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return - end if +9999 call psb_error_handler(err_act) + return end subroutine psb_crwextd subroutine psb_cbase_rwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_cbase_rwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_cbase_rwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr class(psb_c_base_sparse_mat), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info class(psb_c_base_sparse_mat), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -236,12 +234,213 @@ subroutine psb_cbase_rwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_cbase_rwextd + + +subroutine psb_lcrwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_lcrwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + type(psb_lcspmat_type), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + type(psb_lcspmat_type), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + type(psb_lc_coo_sparse_mat) :: actmp + logical rowscale_ + + name='psb_lcrwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (nr > a%get_nrows()) then + select type(aa=> a%a) + type is (psb_lc_csr_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + type is (psb_lc_coo_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + class default + call aa%mv_to_coo(actmp,info) + if (info == psb_success_) then + if (present(b)) then + call psb_rwextd(nr,actmp,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,actmp,info,rowscale=rowscale) + end if + end if + if (info == psb_success_) call aa%mv_from_coo(actmp,info) + end select + end if + if (info /= psb_success_) goto 9999 + call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return + return + +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_lcrwextd +subroutine psb_lcbase_rwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_lcbase_rwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + class(psb_lc_base_sparse_mat), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + class(psb_lc_base_sparse_mat), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb, ma, mb, na, nb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + logical rowscale_ + + name='psb_lcbase_rwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (present(rowscale)) then + rowscale_ = rowscale + else + rowscale_ = .true. end if + + ma = a%get_nrows() + na = a%get_ncols() + + + select type(a) + type is (psb_lc_csr_sparse_mat) + + call psb_ensure_size(nr+1,a%irp,info) + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + + select type (b) + type is (psb_lc_csr_sparse_mat) + call psb_ensure_size(size(a%ja)+nzb,a%ja,info) + call psb_ensure_size(size(a%val)+nzb,a%val,info) + do i=1, min(nr-ma,mb) + a%irp(ma+i+1) = a%irp(ma+i) + b%irp(i+1) - b%irp(i) + ja = a%irp(ma+i) + do jb = b%irp(i), b%irp(i+1)-1 + a%val(ja) = b%val(jb) + a%ja(ja) = b%ja(jb) + ja = ja + 1 + end do + end do + do j=i,nr-ma + a%irp(ma+i+1) = a%irp(ma+i) + end do + class default + + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + end select + call a%set_ncols(max(na,nb)) + + else + + do i=ma+2,nr+1 + a%irp(i) = a%irp(i-1) + end do + + end if + + call a%set_nrows(nr) + + + type is (psb_lc_coo_sparse_mat) + nza = a%get_nzeros() + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + call a%reallocate(nza+nzb) + + select type(b) + type is (psb_lc_coo_sparse_mat) + + if (rowscale_) then + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = ma + b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + else + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + endif + call a%set_nzeros(nza) + + type is (psb_lc_csr_sparse_mat) + + do i=1, min(nr-ma,mb) + do jb = b%irp(i), b%irp(i+1)-1 + nza = nza + 1 + a%val(nza) = b%val(jb) + a%ia(nza) = ma + i + a%ja(nza) = b%ja(jb) + end do + end do + call a%set_nzeros(nza) + + class default + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + + end select + + call a%set_ncols(max(na,nb)) + endif + + call a%set_nrows(nr) + + class default + info = psb_err_unsupported_format_ + ch_err=a%get_fmt() + call psb_errpush(info,name,a_err=ch_err) + goto 9999 + + end select + + call psb_erractionrestore(err_act) return -end subroutine psb_cbase_rwextd +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_lcbase_rwextd diff --git a/base/serial/psb_drwextd.f90 b/base/serial/psb_drwextd.f90 index c2f39473..29898578 100644 --- a/base/serial/psb_drwextd.f90 +++ b/base/serial/psb_drwextd.f90 @@ -35,17 +35,18 @@ ! ! We have a problem here: 1. How to handle well all the formats? ! 2. What should we do with rowscale? Does it only -! apply when a%fida='COO' ?????? +! apply when a%fmt()='COO' ?????? ! ! subroutine psb_drwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_drwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_drwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr type(psb_dspmat_type), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info type(psb_dspmat_type), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -89,23 +90,20 @@ subroutine psb_drwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue - call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return - end if +9999 call psb_error_handler(err_act) + return end subroutine psb_drwextd subroutine psb_dbase_rwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_dbase_rwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_dbase_rwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr class(psb_d_base_sparse_mat), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info class(psb_d_base_sparse_mat), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -174,8 +172,8 @@ subroutine psb_dbase_rwextd(nr,a,info,b,rowscale) nza = a%get_nzeros() if (present(b)) then - mb = b%get_nrows() - nb = b%get_ncols() + mb = b%get_nrows() + nb = b%get_ncols() nzb = b%get_nzeros() call a%reallocate(nza+nzb) @@ -236,12 +234,213 @@ subroutine psb_dbase_rwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_dbase_rwextd + + +subroutine psb_ldrwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_ldrwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + type(psb_ldspmat_type), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + type(psb_ldspmat_type), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + type(psb_ld_coo_sparse_mat) :: actmp + logical rowscale_ + + name='psb_ldrwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (nr > a%get_nrows()) then + select type(aa=> a%a) + type is (psb_ld_csr_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + type is (psb_ld_coo_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + class default + call aa%mv_to_coo(actmp,info) + if (info == psb_success_) then + if (present(b)) then + call psb_rwextd(nr,actmp,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,actmp,info,rowscale=rowscale) + end if + end if + if (info == psb_success_) call aa%mv_from_coo(actmp,info) + end select + end if + if (info /= psb_success_) goto 9999 + call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return + return + +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_ldrwextd +subroutine psb_ldbase_rwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_ldbase_rwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + class(psb_ld_base_sparse_mat), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + class(psb_ld_base_sparse_mat), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb, ma, mb, na, nb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + logical rowscale_ + + name='psb_ldbase_rwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (present(rowscale)) then + rowscale_ = rowscale + else + rowscale_ = .true. end if + + ma = a%get_nrows() + na = a%get_ncols() + + + select type(a) + type is (psb_ld_csr_sparse_mat) + + call psb_ensure_size(nr+1,a%irp,info) + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + + select type (b) + type is (psb_ld_csr_sparse_mat) + call psb_ensure_size(size(a%ja)+nzb,a%ja,info) + call psb_ensure_size(size(a%val)+nzb,a%val,info) + do i=1, min(nr-ma,mb) + a%irp(ma+i+1) = a%irp(ma+i) + b%irp(i+1) - b%irp(i) + ja = a%irp(ma+i) + do jb = b%irp(i), b%irp(i+1)-1 + a%val(ja) = b%val(jb) + a%ja(ja) = b%ja(jb) + ja = ja + 1 + end do + end do + do j=i,nr-ma + a%irp(ma+i+1) = a%irp(ma+i) + end do + class default + + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + end select + call a%set_ncols(max(na,nb)) + + else + + do i=ma+2,nr+1 + a%irp(i) = a%irp(i-1) + end do + + end if + + call a%set_nrows(nr) + + + type is (psb_ld_coo_sparse_mat) + nza = a%get_nzeros() + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + call a%reallocate(nza+nzb) + + select type(b) + type is (psb_ld_coo_sparse_mat) + + if (rowscale_) then + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = ma + b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + else + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + endif + call a%set_nzeros(nza) + + type is (psb_ld_csr_sparse_mat) + + do i=1, min(nr-ma,mb) + do jb = b%irp(i), b%irp(i+1)-1 + nza = nza + 1 + a%val(nza) = b%val(jb) + a%ia(nza) = ma + i + a%ja(nza) = b%ja(jb) + end do + end do + call a%set_nzeros(nza) + + class default + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + + end select + + call a%set_ncols(max(na,nb)) + endif + + call a%set_nrows(nr) + + class default + info = psb_err_unsupported_format_ + ch_err=a%get_fmt() + call psb_errpush(info,name,a_err=ch_err) + goto 9999 + + end select + + call psb_erractionrestore(err_act) return -end subroutine psb_dbase_rwextd +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_ldbase_rwextd diff --git a/base/serial/psb_srwextd.f90 b/base/serial/psb_srwextd.f90 index 4e286d38..6fc53334 100644 --- a/base/serial/psb_srwextd.f90 +++ b/base/serial/psb_srwextd.f90 @@ -35,17 +35,18 @@ ! ! We have a problem here: 1. How to handle well all the formats? ! 2. What should we do with rowscale? Does it only -! apply when a%fida='COO' ?????? +! apply when a%fmt()='COO' ?????? ! ! subroutine psb_srwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_srwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_srwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr type(psb_sspmat_type), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info type(psb_sspmat_type), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -89,23 +90,20 @@ subroutine psb_srwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue - call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return - end if +9999 call psb_error_handler(err_act) + return end subroutine psb_srwextd subroutine psb_sbase_rwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_sbase_rwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_sbase_rwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr class(psb_s_base_sparse_mat), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info class(psb_s_base_sparse_mat), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -236,12 +234,213 @@ subroutine psb_sbase_rwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_sbase_rwextd + + +subroutine psb_lsrwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_lsrwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + type(psb_lsspmat_type), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + type(psb_lsspmat_type), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + type(psb_ls_coo_sparse_mat) :: actmp + logical rowscale_ + + name='psb_lsrwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (nr > a%get_nrows()) then + select type(aa=> a%a) + type is (psb_ls_csr_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + type is (psb_ls_coo_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + class default + call aa%mv_to_coo(actmp,info) + if (info == psb_success_) then + if (present(b)) then + call psb_rwextd(nr,actmp,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,actmp,info,rowscale=rowscale) + end if + end if + if (info == psb_success_) call aa%mv_from_coo(actmp,info) + end select + end if + if (info /= psb_success_) goto 9999 + call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return + return + +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_lsrwextd +subroutine psb_lsbase_rwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_lsbase_rwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + class(psb_ls_base_sparse_mat), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + class(psb_ls_base_sparse_mat), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb, ma, mb, na, nb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + logical rowscale_ + + name='psb_lsbase_rwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (present(rowscale)) then + rowscale_ = rowscale + else + rowscale_ = .true. end if + + ma = a%get_nrows() + na = a%get_ncols() + + + select type(a) + type is (psb_ls_csr_sparse_mat) + + call psb_ensure_size(nr+1,a%irp,info) + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + + select type (b) + type is (psb_ls_csr_sparse_mat) + call psb_ensure_size(size(a%ja)+nzb,a%ja,info) + call psb_ensure_size(size(a%val)+nzb,a%val,info) + do i=1, min(nr-ma,mb) + a%irp(ma+i+1) = a%irp(ma+i) + b%irp(i+1) - b%irp(i) + ja = a%irp(ma+i) + do jb = b%irp(i), b%irp(i+1)-1 + a%val(ja) = b%val(jb) + a%ja(ja) = b%ja(jb) + ja = ja + 1 + end do + end do + do j=i,nr-ma + a%irp(ma+i+1) = a%irp(ma+i) + end do + class default + + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + end select + call a%set_ncols(max(na,nb)) + + else + + do i=ma+2,nr+1 + a%irp(i) = a%irp(i-1) + end do + + end if + + call a%set_nrows(nr) + + + type is (psb_ls_coo_sparse_mat) + nza = a%get_nzeros() + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + call a%reallocate(nza+nzb) + + select type(b) + type is (psb_ls_coo_sparse_mat) + + if (rowscale_) then + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = ma + b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + else + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + endif + call a%set_nzeros(nza) + + type is (psb_ls_csr_sparse_mat) + + do i=1, min(nr-ma,mb) + do jb = b%irp(i), b%irp(i+1)-1 + nza = nza + 1 + a%val(nza) = b%val(jb) + a%ia(nza) = ma + i + a%ja(nza) = b%ja(jb) + end do + end do + call a%set_nzeros(nza) + + class default + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + + end select + + call a%set_ncols(max(na,nb)) + endif + + call a%set_nrows(nr) + + class default + info = psb_err_unsupported_format_ + ch_err=a%get_fmt() + call psb_errpush(info,name,a_err=ch_err) + goto 9999 + + end select + + call psb_erractionrestore(err_act) return -end subroutine psb_sbase_rwextd +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_lsbase_rwextd diff --git a/base/serial/psb_zrwextd.f90 b/base/serial/psb_zrwextd.f90 index a45beedc..949017bc 100644 --- a/base/serial/psb_zrwextd.f90 +++ b/base/serial/psb_zrwextd.f90 @@ -35,17 +35,18 @@ ! ! We have a problem here: 1. How to handle well all the formats? ! 2. What should we do with rowscale? Does it only -! apply when a%fida='COO' ?????? +! apply when a%fmt()='COO' ?????? ! ! subroutine psb_zrwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_zrwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_zrwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr type(psb_zspmat_type), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info type(psb_zspmat_type), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -89,23 +90,20 @@ subroutine psb_zrwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue - call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return - end if +9999 call psb_error_handler(err_act) + return end subroutine psb_zrwextd subroutine psb_zbase_rwextd(nr,a,info,b,rowscale) - use psb_base_mod, psb_protect_name => psb_zbase_rwextd + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_zbase_rwextd implicit none ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) - integer(psb_ipk_), intent(in) :: nr + integer(psb_ipk_), intent(in) :: nr class(psb_z_base_sparse_mat), intent(inout) :: a - integer(psb_ipk_),intent(out) :: info + integer(psb_ipk_),intent(out) :: info class(psb_z_base_sparse_mat), intent(in), optional :: b logical,intent(in), optional :: rowscale @@ -236,12 +234,213 @@ subroutine psb_zbase_rwextd(nr,a,info,b,rowscale) call psb_erractionrestore(err_act) return -9999 continue +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_zbase_rwextd + + +subroutine psb_lzrwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_lzrwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + type(psb_lzspmat_type), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + type(psb_lzspmat_type), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + type(psb_lz_coo_sparse_mat) :: actmp + logical rowscale_ + + name='psb_lzrwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (nr > a%get_nrows()) then + select type(aa=> a%a) + type is (psb_lz_csr_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + type is (psb_lz_coo_sparse_mat) + if (present(b)) then + call psb_rwextd(nr,aa,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,aa,info,rowscale=rowscale) + end if + class default + call aa%mv_to_coo(actmp,info) + if (info == psb_success_) then + if (present(b)) then + call psb_rwextd(nr,actmp,info,b%a,rowscale=rowscale) + else + call psb_rwextd(nr,actmp,info,rowscale=rowscale) + end if + end if + if (info == psb_success_) call aa%mv_from_coo(actmp,info) + end select + end if + if (info /= psb_success_) goto 9999 + call psb_erractionrestore(err_act) - if (err_act == psb_act_abort_) then - call psb_error() - return + return + +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_lzrwextd +subroutine psb_lzbase_rwextd(nr,a,info,b,rowscale) + use psb_mat_mod + use psb_serial_mod, psb_protect_name => psb_lzbase_rwextd + implicit none + + ! Extend matrix A up to NR rows with empty ones (i.e.: all zeroes) + integer(psb_lpk_), intent(in) :: nr + class(psb_lz_base_sparse_mat), intent(inout) :: a + integer(psb_ipk_),intent(out) :: info + class(psb_lz_base_sparse_mat), intent(in), optional :: b + logical,intent(in), optional :: rowscale + + integer(psb_lpk_) :: i,j,ja,jb,nza,nzb, ma, mb, na, nb + integer(psb_ipk_) :: err_act + character(len=20) :: name, ch_err + logical rowscale_ + + name='psb_lzbase_rwextd' + info = psb_success_ + call psb_erractionsave(err_act) + + if (present(rowscale)) then + rowscale_ = rowscale + else + rowscale_ = .true. end if + + ma = a%get_nrows() + na = a%get_ncols() + + + select type(a) + type is (psb_lz_csr_sparse_mat) + + call psb_ensure_size(nr+1,a%irp,info) + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + + select type (b) + type is (psb_lz_csr_sparse_mat) + call psb_ensure_size(size(a%ja)+nzb,a%ja,info) + call psb_ensure_size(size(a%val)+nzb,a%val,info) + do i=1, min(nr-ma,mb) + a%irp(ma+i+1) = a%irp(ma+i) + b%irp(i+1) - b%irp(i) + ja = a%irp(ma+i) + do jb = b%irp(i), b%irp(i+1)-1 + a%val(ja) = b%val(jb) + a%ja(ja) = b%ja(jb) + ja = ja + 1 + end do + end do + do j=i,nr-ma + a%irp(ma+i+1) = a%irp(ma+i) + end do + class default + + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + end select + call a%set_ncols(max(na,nb)) + + else + + do i=ma+2,nr+1 + a%irp(i) = a%irp(i-1) + end do + + end if + + call a%set_nrows(nr) + + + type is (psb_lz_coo_sparse_mat) + nza = a%get_nzeros() + + if (present(b)) then + mb = b%get_nrows() + nb = b%get_ncols() + nzb = b%get_nzeros() + call a%reallocate(nza+nzb) + + select type(b) + type is (psb_lz_coo_sparse_mat) + + if (rowscale_) then + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = ma + b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + else + do j=1,nzb + if ((ma + b%ia(j)) <= nr) then + nza = nza + 1 + a%ia(nza) = b%ia(j) + a%ja(nza) = b%ja(j) + a%val(nza) = b%val(j) + end if + enddo + endif + call a%set_nzeros(nza) + + type is (psb_lz_csr_sparse_mat) + + do i=1, min(nr-ma,mb) + do jb = b%irp(i), b%irp(i+1)-1 + nza = nza + 1 + a%val(nza) = b%val(jb) + a%ia(nza) = ma + i + a%ja(nza) = b%ja(jb) + end do + end do + call a%set_nzeros(nza) + + class default + write(psb_err_unit,*) 'Implement SPGETBLK in RWEXTD!!!!!!!' + + end select + + call a%set_ncols(max(na,nb)) + endif + + call a%set_nrows(nr) + + class default + info = psb_err_unsupported_format_ + ch_err=a%get_fmt() + call psb_errpush(info,name,a_err=ch_err) + goto 9999 + + end select + + call psb_erractionrestore(err_act) return -end subroutine psb_zbase_rwextd +9999 call psb_error_handler(err_act) + + return + +end subroutine psb_lzbase_rwextd