change variable names in ppde/spde to make it more readable.
psblas3-type-indexed
Salvatore Filippone 15 years ago
parent f736089c0d
commit d724c14b5a

@ -352,9 +352,9 @@ contains
type(psb_d_csc_sparse_mat) :: acsc type(psb_d_csc_sparse_mat) :: acsc
type(psb_d_coo_sparse_mat) :: acoo type(psb_d_coo_sparse_mat) :: acoo
type(psb_d_csr_sparse_mat) :: acsr type(psb_d_csr_sparse_mat) :: acsr
real(psb_dpk_) :: zt(nb),glob_x,glob_y,glob_z real(psb_dpk_) :: zt(nb),x,y,z
integer :: m,n,nnz,glob_row,nlr,i,ii,ib,k integer :: m,n,nnz,glob_row,nlr,i,ii,ib,k
integer :: x,y,z,ia,indx_owner integer :: ix,iy,iz,ia,indx_owner
integer :: np, iam, nr, nt integer :: np, iam, nr, nt
integer :: element integer :: element
integer, allocatable :: irow(:),icol(:),myidx(:) integer, allocatable :: irow(:),icol(:),myidx(:)
@ -447,20 +447,20 @@ contains
glob_row=myidx(i) glob_row=myidx(i)
! compute gridpoint coordinates ! compute gridpoint coordinates
if (mod(glob_row,(idim*idim)) == 0) then if (mod(glob_row,(idim*idim)) == 0) then
x = glob_row/(idim*idim) ix = glob_row/(idim*idim)
else else
x = glob_row/(idim*idim)+1 ix = glob_row/(idim*idim)+1
endif endif
if (mod((glob_row-(x-1)*idim*idim),idim) == 0) then if (mod((glob_row-(ix-1)*idim*idim),idim) == 0) then
y = (glob_row-(x-1)*idim*idim)/idim iy = (glob_row-(ix-1)*idim*idim)/idim
else else
y = (glob_row-(x-1)*idim*idim)/idim+1 iy = (glob_row-(ix-1)*idim*idim)/idim+1
endif endif
z = glob_row-(x-1)*idim*idim-(y-1)*idim iz = glob_row-(ix-1)*idim*idim-(iy-1)*idim
! glob_x, glob_y, glob_x coordinates ! x, y, x coordinates
glob_x=x*deltah x = ix*deltah
glob_y=y*deltah y = iy*deltah
glob_z=z*deltah z = iz*deltah
! check on boundary points ! check on boundary points
zt(k) = 0.d0 zt(k) = 0.d0
@ -468,99 +468,81 @@ contains
! !
! term depending on (x-1,y,z) ! term depending on (x-1,y,z)
! !
if (x==1) then if (ix==1) then
val(element)=-b1(glob_x,glob_y,glob_z)& val(element)=-b1(x,y,z)-a1(x,y,z)
& -a1(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*&
& deltah) & deltah)
zt(k) = exp(-glob_y**2-glob_z**2)*(-val(element)) zt(k) = exp(-y**2-z**2)*(-val(element))
else else
val(element)=-b1(glob_x,glob_y,glob_z)& val(element)=-b1(x,y,z)-a1(x,y,z)
& -a1(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*&
& deltah) & deltah)
icol(element) = (x-2)*idim*idim+(y-1)*idim+(z) icol(element) = (ix-2)*idim*idim+(iy-1)*idim+(iz)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y-1,z) ! term depending on (x,y-1,z)
if (y==1) then if (iy==1) then
val(element)=-b2(glob_x,glob_y,glob_z)& val(element)=-b2(x,y,z)-a2(x,y,z)
& -a2(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*&
& deltah) & deltah)
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element)) zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
else else
val(element)=-b2(glob_x,glob_y,glob_z)& val(element)=-b2(x,y,z)-a2(x,y,z)
& -a2(glob_x,glob_y,glob_z) val(element) = val(element)/(deltah*deltah)
val(element) = val(element)/(deltah*& icol(element) = (ix-1)*idim*idim+(iy-2)*idim+(iz)
& deltah)
icol(element) = (x-1)*idim*idim+(y-2)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y,z-1) ! term depending on (x,y,z-1)
if (z==1) then if (iz==1) then
val(element)=-b3(glob_x,glob_y,glob_z)& val(element)=-b3(x,y,z)-a3(x,y,z)
& -a3(glob_x,glob_y,glob_z) val(element) = val(element)/(deltah*deltah)
val(element) = val(element)/(deltah*& zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
& deltah)
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))
else else
val(element)=-b3(glob_x,glob_y,glob_z)& val(element)=-b3(x,y,z)-a3(x,y,z)
& -a3(glob_x,glob_y,glob_z) val(element) = val(element)/(deltah*deltah)
val(element) = val(element)/(deltah*& icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz-1)
& deltah)
icol(element) = (x-1)*idim*idim+(y-1)*idim+(z-1)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y,z) ! term depending on (x,y,z)
val(element)=2*b1(glob_x,glob_y,glob_z)& val(element)=2*b1(x,y,z) + 2*b2(x,y,z)&
& +2*b2(glob_x,glob_y,glob_z)& & + 2*b3(x,y,z) + a1(x,y,z)&
& +2*b3(glob_x,glob_y,glob_z)& & + a2(x,y,z) + a3(x,y,z)
& +a1(glob_x,glob_y,glob_z)& val(element) = val(element)/(deltah*deltah)
& +a2(glob_x,glob_y,glob_z)& icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz)
& +a3(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*&
& deltah)
icol(element) = (x-1)*idim*idim+(y-1)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
! term depending on (x,y,z+1) ! term depending on (x,y,z+1)
if (z==idim) then if (iz==idim) then
val(element)=-b1(glob_x,glob_y,glob_z) val(element)=-b1(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))
else else
val(element)=-b1(glob_x,glob_y,glob_z) val(element)=-b1(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz+1)
icol(element) = (x-1)*idim*idim+(y-1)*idim+(z+1)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y+1,z) ! term depending on (x,y+1,z)
if (y==idim) then if (iy==idim) then
val(element)=-b2(glob_x,glob_y,glob_z) val(element)=-b2(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))
else else
val(element)=-b2(glob_x,glob_y,glob_z) val(element)=-b2(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) icol(element) = (ix-1)*idim*idim+(iy)*idim+(iz)
icol(element) = (x-1)*idim*idim+(y)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x+1,y,z) ! term depending on (x+1,y,z)
if (x<idim) then if (ix<idim) then
val(element)=-b3(glob_x,glob_y,glob_z) val(element)=-b3(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) icol(element) = (ix)*idim*idim+(iy-1)*idim+(iz)
icol(element) = (x)*idim*idim+(y-1)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif

@ -351,9 +351,9 @@ contains
type(psb_s_sparse_mat) :: a type(psb_s_sparse_mat) :: a
type(psb_s_coo_sparse_mat) :: acoo type(psb_s_coo_sparse_mat) :: acoo
type(psb_s_csr_sparse_mat) :: acsr type(psb_s_csr_sparse_mat) :: acsr
real(psb_spk_) :: zt(nb),glob_x,glob_y,glob_z real(psb_spk_) :: zt(nb),x,y,z
integer :: m,n,nnz,glob_row,nlr,i,ii,ib,k integer :: m,n,nnz,glob_row,nlr,i,ii,ib,k
integer :: x,y,z,ia,indx_owner integer :: ix,iy,iz,ia,indx_owner
integer :: np, iam, nr, nt integer :: np, iam, nr, nt
integer :: element integer :: element
integer, allocatable :: irow(:),icol(:),myidx(:) integer, allocatable :: irow(:),icol(:),myidx(:)
@ -445,20 +445,20 @@ contains
glob_row=myidx(i) glob_row=myidx(i)
! compute gridpoint coordinates ! compute gridpoint coordinates
if (mod(glob_row,(idim*idim)) == 0) then if (mod(glob_row,(idim*idim)) == 0) then
x = glob_row/(idim*idim) ix = glob_row/(idim*idim)
else else
x = glob_row/(idim*idim)+1 ix = glob_row/(idim*idim)+1
endif endif
if (mod((glob_row-(x-1)*idim*idim),idim) == 0) then if (mod((glob_row-(ix-1)*idim*idim),idim) == 0) then
y = (glob_row-(x-1)*idim*idim)/idim iy = (glob_row-(ix-1)*idim*idim)/idim
else else
y = (glob_row-(x-1)*idim*idim)/idim+1 iy = (glob_row-(ix-1)*idim*idim)/idim+1
endif endif
z = glob_row-(x-1)*idim*idim-(y-1)*idim iz = glob_row-(ix-1)*idim*idim-(iy-1)*idim
! glob_x, glob_y, glob_x coordinates ! x, y, x coordinates
glob_x=x*deltah x = ix*deltah
glob_y=y*deltah y = iy*deltah
glob_z=z*deltah z = iz*deltah
! check on boundary points ! check on boundary points
zt(k) = 0.d0 zt(k) = 0.d0
@ -466,99 +466,81 @@ contains
! !
! term depending on (x-1,y,z) ! term depending on (x-1,y,z)
! !
if (x==1) then if (ix==1) then
val(element)=-b1(glob_x,glob_y,glob_z)& val(element)=-b1(x,y,z)-a1(x,y,z)
& -a1(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*&
& deltah) & deltah)
zt(k) = exp(-glob_y**2-glob_z**2)*(-val(element)) zt(k) = exp(-y**2-z**2)*(-val(element))
else else
val(element)=-b1(glob_x,glob_y,glob_z)& val(element)=-b1(x,y,z)-a1(x,y,z)
& -a1(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*&
& deltah) & deltah)
icol(element) = (x-2)*idim*idim+(y-1)*idim+(z) icol(element) = (ix-2)*idim*idim+(iy-1)*idim+(iz)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y-1,z) ! term depending on (x,y-1,z)
if (y==1) then if (iy==1) then
val(element)=-b2(glob_x,glob_y,glob_z)& val(element)=-b2(x,y,z)-a2(x,y,z)
& -a2(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*&
& deltah) & deltah)
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element)) zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
else else
val(element)=-b2(glob_x,glob_y,glob_z)& val(element)=-b2(x,y,z)-a2(x,y,z)
& -a2(glob_x,glob_y,glob_z) val(element) = val(element)/(deltah*deltah)
val(element) = val(element)/(deltah*& icol(element) = (ix-1)*idim*idim+(iy-2)*idim+(iz)
& deltah)
icol(element) = (x-1)*idim*idim+(y-2)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y,z-1) ! term depending on (x,y,z-1)
if (z==1) then if (iz==1) then
val(element)=-b3(glob_x,glob_y,glob_z)& val(element)=-b3(x,y,z)-a3(x,y,z)
& -a3(glob_x,glob_y,glob_z) val(element) = val(element)/(deltah*deltah)
val(element) = val(element)/(deltah*& zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
& deltah)
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))
else else
val(element)=-b3(glob_x,glob_y,glob_z)& val(element)=-b3(x,y,z)-a3(x,y,z)
& -a3(glob_x,glob_y,glob_z) val(element) = val(element)/(deltah*deltah)
val(element) = val(element)/(deltah*& icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz-1)
& deltah)
icol(element) = (x-1)*idim*idim+(y-1)*idim+(z-1)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y,z) ! term depending on (x,y,z)
val(element)=2*b1(glob_x,glob_y,glob_z)& val(element)=2*b1(x,y,z) + 2*b2(x,y,z)&
& +2*b2(glob_x,glob_y,glob_z)& & + 2*b3(x,y,z) + a1(x,y,z)&
& +2*b3(glob_x,glob_y,glob_z)& & + a2(x,y,z) + a3(x,y,z)
& +a1(glob_x,glob_y,glob_z)& val(element) = val(element)/(deltah*deltah)
& +a2(glob_x,glob_y,glob_z)& icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz)
& +a3(glob_x,glob_y,glob_z)
val(element) = val(element)/(deltah*&
& deltah)
icol(element) = (x-1)*idim*idim+(y-1)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
! term depending on (x,y,z+1) ! term depending on (x,y,z+1)
if (z==idim) then if (iz==idim) then
val(element)=-b1(glob_x,glob_y,glob_z) val(element)=-b1(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))
else else
val(element)=-b1(glob_x,glob_y,glob_z) val(element)=-b1(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz+1)
icol(element) = (x-1)*idim*idim+(y-1)*idim+(z+1)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x,y+1,z) ! term depending on (x,y+1,z)
if (y==idim) then if (iy==idim) then
val(element)=-b2(glob_x,glob_y,glob_z) val(element)=-b2(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))
zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))
else else
val(element)=-b2(glob_x,glob_y,glob_z) val(element)=-b2(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) icol(element) = (ix-1)*idim*idim+(iy)*idim+(iz)
icol(element) = (x-1)*idim*idim+(y)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif
! term depending on (x+1,y,z) ! term depending on (x+1,y,z)
if (x<idim) then if (ix<idim) then
val(element)=-b3(glob_x,glob_y,glob_z) val(element)=-b3(x,y,z)
val(element) = val(element)/(deltah*& val(element) = val(element)/(deltah*deltah)
& deltah) icol(element) = (ix)*idim*idim+(iy-1)*idim+(iz)
icol(element) = (x)*idim*idim+(y-1)*idim+(z)
irow(element) = glob_row irow(element) = glob_row
element = element+1 element = element+1
endif endif

Loading…
Cancel
Save