24小时热门版块排行榜

返回列表

当前主题已经存档。

【有奖交流】积极回复本帖子，参与交流，就有机会分得作者 futurehero 的 100 个金币

futurehero

新虫 (正式写手)

应助: 0 (幼儿园)
金币: 22
散金: 464
红花: 1
帖子: 350
在线: 47.7小时
虫号: 250574
注册: 2006-05-13
性别: GG
专业: 凝聚态物性 II ：电子结构

[交流] 【求助】诚求助并行计算中的一个问题，谢谢大家！

这个程序大致如下，我想用并行计算解决“do循环中包含积分”的计算问题。如果将并行的线程设为1，则和串行计算的结果一直，但是线程大于1时，计算的结果就不正确了。小弟刚刚接触并行计算两天，做了很多种方案都发现线程大于1时计算的结果都不对，请诸位大侠帮助，我想得到争取的并行计算结果。如果将程序中的call omp_set_num_threads(4 )线程改为call omp_set_num_threads(1 )就和串行结果是一致的。(可能您只需要看一下并行区就可以发现问题所在了)
程序如下：

program main
implicit none

integer:: Num, i, p, j, k ,var, r ,Y
real(kind=8):: D, epsl21, epsl22, T, eta, Gama, GamaL, GamaR, epsl0, VL, VR

parameter(D=1.0_8, Num=10001, Gama=0.01_8, GamaL=0.005_8, GamaR=0.005_8, epsl0=-0.1_8, VL=0.0_8, VR= 0.0_8)

real(kind=8):: pi, Dc, Hc, fL(Num), fR(Num), f0(Num), epsl(Num), epsl11, epsl12, dIdV(Num)
complex(kind=8):: sigma0, FG11(Num), FG12(Num), FG21(Num), FG22(Num)

complex(kind=8):: G011(Num), G012(Num), G021(Num), G022(Num), G11(Num), G12(Num), G21(Num), G22(Num)
complex(kind=8):: DDOS(Num), DOS, D11, D12, D21, D22
complex(kind=8):: n11, n12, n21, n22

complex(kind=8):: delta1211(Num), delta2111(Num), delta2211(Num)
complex(kind=8):: delta1112(Num), delta2112(Num), delta2212(Num)

complex(kind=8):: delta1121(Num), delta1221(Num), delta2221(Num)
complex(kind=8):: delta1122(Num), delta1222(Num), delta2122(Num)

complex(kind=8):: FZB11(Num), FZB12(Num), FZB21(Num), FZB22(Num)

complex(kind=8):: B1211(Num), B2111(Num), B2211(Num)
complex(kind=8):: FB1211(Num),  FB2111(Num),  FB2211(Num)

complex(kind=8):: B1112(Num), B2112(Num), B2212(Num)
complex(kind=8)::  FB1112(Num),  FB2112(Num),  FB2212(Num)

complex(kind=8)::  B1121(Num), B1221(Num), B2221(Num)
complex(kind=8)::  FB1121(Num),  FB1221(Num),  FB2221(Num)

complex(kind=8):: B1222(Num), B2122(Num), B1122(Num)
complex(kind=8)::  FB1222(Num),  FB2122(Num),  FB1122(Num)

real(kind=8):: m, Field, SOC, Tk

Dc=2*D / (Num-1)
eta =2.0*Dc
pi=4.0_8*datan(1.0_8)

m=0.00_8
T=0.0_8
SOC=0.00000_8
Field=0.00000_8
epsl11=epsl0 - SOC/2.0_8 - m*Field - Field
epsl12=epsl0 + SOC/2.0_8 + m*Field - Field
epsl21=epsl0 + SOC/2.0_8 - m*Field + Field
epsl22=epsl0 - SOC/2.0_8 + m*Field + Field
do i=1, Num, 1
epsl(i)=-D + (i-1)  *Dc
end do
do i=1, Num, 1
if(  epsl(i)          <  VL    )  fL(i)=1.0_8
if(  epsl(i)          >  VL    )  fL(i)=0.0_8
if(  dabs(epsl(i)-VL)  <  0.01_8*Dc ) fL(i)=0.5_8
if(  epsl(i)          <  VR    )  fR(i)=1.0_8
if(  epsl(i)          >  VR    )  fR(i)=0.0_8
if(  dabs(epsl(i)-VR)  <  0.01_8*Dc ) fR(i)=0.5_8
end do
f0=( fL*GamaL  +  fR*GamaR ) /Gama
delta1211=dcmplx(  epsl12 - epsl11,    eta  ) - epsl
delta2111=dcmplx(  epsl21 - epsl11,    eta  ) - epsl
delta2211=dcmplx(  epsl22 - epsl11,    eta  ) - epsl

delta1112=dcmplx(  epsl11 - epsl12,    eta  ) - epsl
delta2112=dcmplx(  epsl21 - epsl12,    eta  ) - epsl
delta2212=dcmplx(  epsl22 - epsl12,    eta  ) - epsl
delta1121=dcmplx(  epsl11 - epsl21,    eta  ) - epsl
delta1221=dcmplx(  epsl12 - epsl21,    eta  ) - epsl
delta2221=dcmplx(  epsl22 - epsl21,    eta  ) - epsl
delta1122=dcmplx(  epsl11 - epsl22,    eta  ) - epsl
delta1222=dcmplx(  epsl12 - epsl22,    eta  ) - epsl
delta2122=dcmplx(  epsl21 - epsl22,    eta  ) - epsl
!================================
sigma0=dcmplx(0.0_8, -Gama)
G011  =epsl-epsl11-sigma0
G012  =epsl-epsl12-sigma0
G021  =epsl-epsl21-sigma0
G022  =epsl-epsl22-sigma0
!********************************************
call omp_set_num_threads(4 )
!$OMP parallel
!$OMP do private (Y)
do Y=1, Num, 1
FB1211=f0/( epsl(Y) + delta1211  )
FB2111=f0/( epsl(Y) + delta2111  )
FB2211=f0/( epsl(Y) + delta2211  )
FB1112=f0/( epsl(Y) + delta1112  )
FB2112=f0/( epsl(Y) + delta2112  )
FB2212=f0/( epsl(Y) + delta2212  )
FB1121=f0/( epsl(Y) + delta1121  )
FB1221=f0/( epsl(Y) + delta1221  )
FB2221=f0/( epsl(Y) + delta2221  )
FB1122=f0/( epsl(Y) + delta1122  )
FB1222=f0/( epsl(Y) + delta1222  )
FB2122=f0/( epsl(Y) + delta2122  )
call INTEGRAL(Num, Dc,  FB1211, B1211(Y))
call INTEGRAL(Num, Dc,  FB2111, B2111(Y))
call INTEGRAL(Num, Dc,  FB2211, B2211(Y))
call INTEGRAL(Num, Dc,  FB1112, B1112(Y))
call INTEGRAL(Num, Dc,  FB2112, B2112(Y))
call  INTEGRAL(Num, Dc,  FB2212, B2212(Y))
call  INTEGRAL(Num, Dc,  FB1121, B1121(Y))
call  INTEGRAL(Num, Dc,  FB1221, B1221(Y))
call  INTEGRAL(Num, Dc,  FB2221, B2221(Y))
call  INTEGRAL(Num, Dc,  FB1122, B1122(Y))
call  INTEGRAL(Num, Dc,  FB1222, B1222(Y))
call  INTEGRAL(Num, Dc,  FB2122, B2122(Y))
end do
!$OMP end do
!$OMP end parallel
!******************************************
B1211=B1211*Gama/pi
B2111=B2111*Gama/pi
B2211=B2211*Gama/pi

B1112=B1112*Gama/pi
B2112=B2112*Gama/pi
B2212=B2212*Gama/pi

B1121=B1121*Gama/pi
B1221=B1221*Gama/pi
B2221=B2221*Gama/pi

B1122=B1122*Gama/pi
B1222=B1222*Gama/pi
B2122=B2122*Gama/pi
G11= ( 1.0_8 - 3.0_8*0.245_8 )/( G011 - B1211-B2111-B2211 )
G12= ( 1.0_8 - 3.0_8*0.245_8 )/( G012 - B1112-B2112-B2212 )
G21= ( 1.0_8 - 3.0_8*0.245_8 )/( G021 - B1121-B1221-B2221 )
G22= ( 1.0_8 - 3.0_8*0.245_8 )/( G022 - B1122-B1222-B2122 )

FG11=-f0*G11/pi
FG12=-f0*G12/pi
FG21=-f0*G21/pi
FG22=-f0*G22/pi
call INTEGRAL(Num, Dc, FG11, n11)
call INTEGRAL(Num, Dc, FG12, n12)
call INTEGRAL(Num, Dc, FG21, n21)
call INTEGRAL(Num, Dc, FG22, n22)

write(*,*) "n11=", dimag(n11), "n12=", dimag(n12)
write(*,*) "n21=", dimag(n21), "n22=", dimag(n22)

stop
end program main

subroutine INTEGRAL(Num, Dc, Y, answer)
implicit none
integer::  Num
real(kind=8):: Dc
complex(kind=8):: Y(Num), answer
integer:: i
complex:: s1

s1=(Y(1)+ Y(Num))/2
do i=2, Num-1, 1
s1=s1+ Y(i)
end do
answer=s1*Dc

return
end subroutine INTEGRAL

stop
end program main

[ Last edited by futurehero on 2010-3-4 at 23:09 ]

回复此楼

» 猜你喜欢

1楼 2010-03-04 22:14:05

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

futurehero

新虫 (正式写手)

应助: 0 (幼儿园)
金币: 22
散金: 464
红花: 1
帖子: 350
在线: 47.7小时
虫号: 250574
注册: 2006-05-13
性别: GG
专业: 凝聚态物性 II ：电子结构

已经解决了!错误原因在于不同的线程会同时将中间变量写到同一个内存单元FB1211,造成结果错误.只需用将再为每个线程定义不同的变量就解决问题了.

赞一下

回复此楼

2楼2010-03-20 16:08:46

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

相关版块跳转我要订阅楼主 futurehero 的主题更新

返回列表