|  | 
 
| 我用fortran编写了一段程序用于初始以及边界条件的赋予,但在实际运行过程中出现了错误,求解器不进行求解计算即壁钟时间与分析时间一直为0,而不报错。代码如下 MODULE Mandrel datas !芯轴边界条件数据
 IMPLICIT NONE
 SAVE
 real*8,dimension (320600)::InitdarhoLis,InitdaThetaCorLis,InitdazC
 $    orLis,InitdanxLis,InitdanyLis,InitdafxLis,InitdafyLis
 real*8,dimension (7,200,229)::tdarhoLis,tdazCorLis,tdaThetaCorLis,
 $    tdanxLis,tdanyLis,tdafxLis,tdafyLis
 END MODULE Mandrel datas
 
 subroutine usinc(f,n,ndeg,iflag)
 use Mandrel datas
 implicit real*8(a-h,o-z)
 integer n,iflag,ndeg
 real*8 f
 integer u,v
 real*8,dimension (80000,3)::cor1,cor2
 real*8,dimension (100,53)::RhoCorLis,ThetaCorLis,TempLis
 real*8,dimension (6300)::IniTempLis,IniRhoCorLis,IniThetaCorLis
 real*8,dimension (3)::cor
 integer Flag,Numrho,Numtheta,TotaNodeNum,i,j,NodeID,tdNodeID
 $  ,rFlag,tFlag
 real*8:: rho,theta
 
 
 
 dimension f(*)
 Numrho=53
 Numtheta=100
 TotaNodeNum=5300
 
 
 if(flag==0) then
 !!温度数据导入与处理
 open (unit=10,file='2-IniTemp.txt')
 open (unit=20,file='2-RhoCor.txt')
 open (unit=30,file='2-ThetaCor.txt')
 do i=1,5300
 read(10,*) IniTempLis(i)
 read(20,*) IniRhoCorLis(i)
 read(30,*) IniThetaCorLis(i)
 enddo
 close (unit=10)
 close (unit=20)
 close (unit=30)
 
 do i=1,100,1
 do j=1,53,1
 NodeID=(i-1)*53+j
 TempLis(i,j)=IniTempLis(NodeID)
 RhoCorLis(i,j)=IniRhoCorLis(NodeID)
 ThetaCorLis(i,j)=IniThetaCorLis(NodeID)
 enddo
 enddo
 
 !!辊芯非均匀力载荷导入与处理
 open (unit=11,file='tdarho.txt')
 open (unit=12,file='tdatheta.txt')
 open (unit=13,file='tdaz.txt')
 open (unit=14,file='tdanx.txt')
 open (unit=15,file='tdany.txt')
 open (unit=16,file='tdafx.txt')
 open (unit=17,file='tdafy.txt')
 
 do i=1,320600
 read(11,*) InitdarhoLis(i)
 read(12,*) InitdaThetaCorLis(i)
 read(13,*) InitdazCorLis(i)
 read(14,*) InitdanxLis(i)
 read(15,*) InitdanyLis(i)
 read(16,*) InitdafxLis(i)
 read(17,*) InitdafyLis(i)
 enddo
 
 close (unit=11)
 close (unit=12)
 close (unit=13)
 close (unit=14)
 close (unit=15)
 close (unit=16)
 close (unit=17)
 
 do i=1,229,1
 do j=1,200,1
 do k=1,7,1
 tdNodeID=(i-1)*1400+(j-1)*7+k
 tdarhoLis(k,j,i)=InitdarhoLis(tdNodeID)
 tdazCorLis(k,j,i)=InitdazCorLis(tdNodeID)
 tdaThetaCorLis(k,j,i)=InitdaThetaCorLis(tdNodeID)
 tdanxLis(k,j,i)=InitdanxLis(tdNodeID)
 tdanyLis(k,j,i)=InitdanyLis(tdNodeID)
 tdafxLis(k,j,i)=InitdafxLis(tdNodeID)
 tdafyLis(k,j,i)=InitdafyLis(tdNodeID)
 enddo
 enddo
 enddo
 
 Flag=1
 endif
 
 
 
 CALL NODVAR(0,n,cor,NQNCOMP,NQDATATYPE)
 
 rho= sqrt(cor(1)**2 + cor(2)**2)
 
 if(cor(2)>=0)then
 theta = acos(cor(1)/rho)
 elseif(cor(2)<0.and.cor(1)<0)then
 theta = acos(abs(cor(1))/rho)+3.14
 elseif(cor(2)<0.and.cor(1)>0)then
 theta = acos(-(cor(1))/rho)+3.14
 endif
 
 do j=1,53
 if (rho>=RhoCorLis(1,j).and.rho<=RhoCorLis(1,j+1)) then
 
 if (abs(rho-RhoCorLis(1,j))<=abs(rho-RhoCorLis(1,j+1)))then
 rFlag=j
 else
 rFlag=j+1
 endif
 elseif (rho<=RhoCorLis(1,1)) then
 rFlag=1
 elseif (rho>=RhoCorLis(1,Numrho)) then
 rFlag=Numrho
 endif
 enddo
 
 do i=1,100
 if (theta>=ThetaCorLis(i,1).and.theta<=ThetaCorLis(i+1,1)) then
 if (abs(theta-ThetaCorLis(i,1))<=abs(theta-ThetaCorLis(i+1,i)
 $    ))then
 tFlag=i
 else
 tFlag=i+1
 endif
 elseif (theta<=ThetaCorLis(1,1)) then
 tFlag=1
 elseif (theta>=ThetaCorLis(Numtheta,1)) then
 tFlag=Numtheta
 endif
 enddo
 
 if (cor(3)>=-2719.and.cor(3)<=-2294.or.cor(3)>=-2134.and.cor
 $    (3)<=-1674.or.cor(3)>=-1504.and.cor(3)<=-1044.or.cor(3)>
 $    =-884.and.cor(3)<=-649) then
 iflag=3
 f(1)=TempLis(tFlag,rFlag)
 else
 iflag=3
 f(1)=50
 endif
 
 return
 end
 
 
 
 subroutine forcdt(u,v,a,dp,du,time,dtime,ndeg,node,
 *                  ug,xord,ncrd,iacflg,inc,ipass)
 
 use Mandrel datas
 
 #ifdef _IMPLICITNONE
 implicit none
 #else
 implicit logical (a-z)
 #endif
 c     ** Start of generated type statements **
 real*8 a, dp, dtime, du
 integer iacflg, inc, ipass, ncrd, ndeg, node
 real*8 time, u, ug, v, xord
 real*8 rho1,theta1
 
 c     ** End of generated type statements **
 dimension u(ndeg),v(ndeg),a(ndeg),dp(ndeg),du(ndeg),ug(ndeg),
 *          xord(ncrd)
 integer i,j,k,zFlag,Numtdarho,Numtdatheta,Numtdaz,rFlag,tFlag
 
 include 'BCLABEL'
 Numtdarho=7
 Numtdatheta=200
 Numtdaz=229
 
 IF (bcname=='apply6') then
 
 rho1= sqrt(xord(1)**2 + xord(2)**2)
 
 if(xord(2)>=0)then
 theta1 = acos(xord(1)/rho1)
 elseif(xord(2)<0.and.xord(1)<0)then
 theta1 = acos(abs(xord(1))/rho1)+3.14
 elseif(xord(2)<0.and.xord(1)>0)then
 theta1 = acos(-(xord(1))/rho1)+3.14
 endif
 
 do i=1,229
 if (xord(3)>=tdazCorLis(1,1,i).and.xord(3)<=tdazCorLis(1,1,i+1))
 $ then
 if (abs(xord(3)-tdazCorLis(1,1,i))<=abs(xord(3)-tdazCorLis(1,
 $    1,i+1))) then
 zFlag=i
 else
 zFlag=i+1
 endif
 elseif (xord(3)<=tdazCorLis(1,1,1)) then
 zFlag=1
 elseif (xord(3)>=tdazCorLis(1,1,Numtdaz)) then
 zFlag=Numtdaz
 endif
 enddo
 
 do j=1,200
 if (theta1>=tdaThetaCorLis(1,j,1).and.theta1<=tdaThetaCorLis(1,j
 $    +1,1)) then
 if (abs(theta1-tdaThetaCorLis(1,j,1))<=abs(theta1-tdaThetaCor
 $    Lis(1,j+1,1)))then
 tFlag=j
 else
 tFlag=j+1
 endif
 elseif (theta1<=tdaThetaCorLis(1,1,1)) then
 tFlag=1
 elseif (theta1>=tdaThetaCorLis(1,Numtdatheta,1)) then
 tFlag=Numtdatheta
 endif
 enddo
 
 do k=1,7
 if (rho1>=tdarhoLis(k,1,1).and.rho1<=tdarhoLis(k+1,1,1)) then
 if (abs(rho1-tdarhoLis(k,1,1))<=abs(rho1-tdarhoLis(k+1,1,1)))
 $    then
 rFlag=k
 else
 rFlag=k+1
 endif
 elseif (rho1<=tdarhoLis(1,1,1)) then
 rFlag=1
 elseif (rho1>=tdarhoLis(Numtdarho,1,1)) then
 rFlag=Numtdarho
 endif
 enddo
 
 dp(1) = tdanxLis(rFlag,tFlag,zFlag)+tdafxLis(rFlag,tFlag,zFlag)
 dp(2) = tdanyLis(rFlag,tFlag,zFlag)+tdafyLis(rFlag,tFlag,zFlag)
 endif
 return
 end
 
 
 | 
 |