From d25a746067cd174a8a81e1698bfaa95f1576f7c3 Mon Sep 17 00:00:00 2001 From: gabrielequatrana Date: Sun, 7 Jul 2024 23:27:00 +0200 Subject: [PATCH] Some bugs fixed --- base/modules/serial/psb_d_base_vect_mod.F90 | 6 ----- cuda/psb_d_cuda_vect_mod.F90 | 2 +- krylov/psb_dbgmres.f90 | 30 ++++++++++----------- 3 files changed, 15 insertions(+), 23 deletions(-) diff --git a/base/modules/serial/psb_d_base_vect_mod.F90 b/base/modules/serial/psb_d_base_vect_mod.F90 index 8c87d362..78e21cda 100644 --- a/base/modules/serial/psb_d_base_vect_mod.F90 +++ b/base/modules/serial/psb_d_base_vect_mod.F90 @@ -2945,9 +2945,6 @@ contains external :: dgemm integer(psb_ipk_) :: x_n, y_n, lda, ldb - write(*,*) 'CPU PROD' - - if (x%is_dev()) call x%sync() if (y%is_dev()) call y%sync() select type(yy => y) @@ -2982,9 +2979,6 @@ contains external :: dgemm integer(psb_ipk_) :: x_n, y_n, lda, ldb - write(*,*) 'CPU PROD' - - if (x%is_dev()) call x%sync() x_n = x%get_ncols() y_n = size(y,dim=2) diff --git a/cuda/psb_d_cuda_vect_mod.F90 b/cuda/psb_d_cuda_vect_mod.F90 index 97b91d2e..b1e57f7e 100644 --- a/cuda/psb_d_cuda_vect_mod.F90 +++ b/cuda/psb_d_cuda_vect_mod.F90 @@ -1668,7 +1668,7 @@ contains if (yy%is_host()) call yy%sync() if (zz%is_host()) call zz%sync() if (info == 0) & - & info = prodMultiVecDevice('N',nr,yy%get_ncols(),x%get_ncols(),done,x%deviceVect,yy%deviceVect,done,zz%deviceVect) + & info = prodMultiVecDevice('N',nr,yy%get_ncols(),x%get_ncols(),done,x%deviceVect,yy%deviceVect,dzero,zz%deviceVect) call zz%set_dev() class default ! TODO diff --git a/krylov/psb_dbgmres.f90 b/krylov/psb_dbgmres.f90 index 159e8bda..f65e595f 100644 --- a/krylov/psb_dbgmres.f90 +++ b/krylov/psb_dbgmres.f90 @@ -63,12 +63,12 @@ subroutine psb_dbgmres_multivect(a, prec, b, x, eps, desc_a, info, itmax, iter, real(psb_dpk_), allocatable :: aux(:), c(:,:), s(:,:), h(:,:), beta(:,:) type(psb_d_multivect_type), allocatable :: v(:) - type(psb_d_multivect_type) :: w, r0, rm, pd + type(psb_d_multivect_type) :: w, pd integer(psb_ipk_) :: naux, itrace_, n_row, n_col, nrhs, nrep integer(psb_lpk_) :: mglob, n_add, ncv - integer(psb_ipk_) :: i, j, istop_, err_act, idx_i, idx_j, idx + integer(psb_ipk_) :: i, j, k, istop_, err_act, idx_i, idx_j, idx integer(psb_ipk_) :: debug_level, debug_unit type(psb_ctxt_type) :: ctxt @@ -152,13 +152,9 @@ subroutine psb_dbgmres_multivect(a, prec, b, x, eps, desc_a, info, itmax, iter, & beta((nrep+1)*nrhs,nrhs),r0n2(nrhs),rmn2(nrhs),stat=info) if (info == psb_success_) call psb_geall(v,desc_a,info,m=nrep+1,n=nrhs) if (info == psb_success_) call psb_geall(w,desc_a,info,n=nrhs) - if (info == psb_success_) call psb_geall(r0,desc_a,info,n=nrhs) - if (info == psb_success_) call psb_geall(rm,desc_a,info,n=nrhs) if (info == psb_success_) call psb_geall(pd,desc_a,info,n=nrhs) if (info == psb_success_) call psb_geasb(v(1),desc_a,info,mold=x%v,n=nrhs) if (info == psb_success_) call psb_geasb(w,desc_a,info,mold=x%v,n=nrhs) - if (info == psb_success_) call psb_geasb(r0,desc_a,info,mold=x%v,n=nrhs) - if (info == psb_success_) call psb_geasb(rm,desc_a,info,mold=x%v,n=nrhs) if (info == psb_success_) call psb_geasb(pd,desc_a,info,mold=x%v,n=nrhs) if (info /= psb_success_) then info=psb_err_from_subroutine_non_ @@ -173,19 +169,19 @@ subroutine psb_dbgmres_multivect(a, prec, b, x, eps, desc_a, info, itmax, iter, ! Compute norm2 of R(0) if (istop_ == 1) then - call psb_geaxpby(done,b,dzero,r0,desc_a,info) + call psb_geaxpby(done,b,dzero,v(1),desc_a,info) if (info /= psb_success_) then info=psb_err_from_subroutine_non_ call psb_errpush(info,name) goto 9999 end if - call psb_spmm(-done,a,x,done,r0,desc_a,info,work=aux) + call psb_spmm(-done,a,x,done,v(1),desc_a,info,work=aux) if (info /= psb_success_) then info=psb_err_from_subroutine_non_ call psb_errpush(info,name) goto 9999 end if - r0n2 = psb_genrm2(r0,desc_a,info) + r0n2 = psb_genrm2(v(1),desc_a,info) endif if (info /= psb_success_) then info=psb_err_from_subroutine_non_ @@ -197,7 +193,7 @@ subroutine psb_dbgmres_multivect(a, prec, b, x, eps, desc_a, info, itmax, iter, c = dzero s = dzero beta = dzero - itx = dzero + itx = 0 deps = eps n_add = nrhs-1 @@ -370,8 +366,6 @@ subroutine psb_dbgmres_multivect(a, prec, b, x, eps, desc_a, info, itmax, iter, if (info == psb_success_) call psb_gefree(v,desc_a,info) if (info == psb_success_) call psb_gefree(w,desc_a,info) - if (info == psb_success_) call psb_gefree(r0,desc_a,info) - if (info == psb_success_) call psb_gefree(rm,desc_a,info) if (info == psb_success_) call psb_gefree(pd,desc_a,info) if (info == psb_success_) deallocate(aux,c,s,h,beta,r0n2,rmn2,stat=info) if (info /= psb_success_) then @@ -402,6 +396,9 @@ contains ! Allocate output allocate(res(nrhs,nrhs)) + + ! Initialize params + res = dzero ! Start factorization do i=1,nrhs @@ -473,16 +470,17 @@ contains call drotg(rti,rti1,c(idx_col,j),s(idx_col,j)) ! Apply Givens rotation to H - call drot(nrhs,h(back_idx-1,idx_rep:idx_rep+n_add),1,& - & h(back_idx,idx_rep:idx_rep+n_add),& + call drot(nrhs-(i-1),h(back_idx-1,idx_col:idx_rep+n_add),1,& + & h(back_idx,idx_col:idx_rep+n_add),& & 1,c(idx_col,j),s(idx_col,j)) ! Eliminate rotated values h(back_idx,idx_rep:idx_col) = dzero ! Apply Givens rotation to G=E1*Beta - call drot(nrhs,beta(back_idx-1,:),1,beta(back_idx,:),1,& - & c(idx_col,j),s(idx_col,j)) + call drot(j,beta(back_idx-1,nrhs-(j-1):nrhs),& + & 1,beta(back_idx,nrhs-(j-1):nrhs),1,& + & c(idx_col,j),s(idx_col,j)) end do end do