子例程仅适用于某些值

Question

PROGRAM po
IMPLICIT NONE
INTEGER                                  :: i, j, nc, nd,ok,iter
REAL                                     :: alph, bet, chi
REAL, DIMENSION(:), ALLOCATABLE          :: uexact, x,u,up2, xd , uexactd
REAL                                     :: E, k, Lc, hc, eps, h, Ld
INTEGER, DIMENSION(7)                    :: valnc = (/ 10, 50, 100, 500, 1000, 5000, 10000/)


Do iter=1,7
  nc = valnc(iter)

!--------------------------------------------------------------------[Given DATA]-!
E=50. ; k=125. ; hc=0.01 ; eps=0.01 ; Lc=1 ;

h = Lc/nc ; chi=sqrt((E*hc)/k) ; alph= -(1/h**2) ; bet=(2/h**2)+(k/(E*hc)) ;

Ld=0.2*Lc ; nd=0.2*nc;

!----------------------------------------------------------------------------- 
 OPEN(UNIT=123,FILE="uetexact2.out",ACTION="write",STATUS='old')

 ALLOCATE(x(0:nc), uexact(0:nc), xd(nc:nc+nd), uexactd(nc:nc+nd))

!------------------------------------------------------[ANALYTICAL RESOLUTION]-!
 DO i = 0, nc
   x(i) = h*i-Lc/2
   uexact(i) = eps*chi*((sinh(x(i)/chi))/(cosh(Lc/(2*chi))))
 END DO

 DO i = nc, nc+nd
  xd(i) = (h*i)+(Ld/2)
  uexactd(i) = (eps*xd(i))+eps*(chi*tanh(Lc/(2*chi))-(Lc/2))
 END DO
!-------------------------------------------------------[NUMERICAL SOLUTION]-!
 CALL resolution(0,nc,bet,2*alph,alph,2*alph,-2*eps/h,2*eps/h,u)
 CALL resolution(nc,nc+nd,-2.0,1.0,1.0,2.0,-u(nc),-2*eps*h,up2)

 WRITE(123,fmt='(3E15.6)') xd, uexactd, up2
 CALL SYSTEM('gnuplot graphic.plt')

  END DO


CONTAINS
SUBROUTINE Resolution (n1,n2,a,b1,b,b2,c1,c2,u1)
INTEGER, INTENT(IN) :: n1,n2
REAL, INTENT(IN):: a,b1,b,b2,c1,c2
REAL, DIMENSION (:),ALLOCATABLE, INTENT(OUT) :: u1
REAL, DIMENSION(:), ALLOCATABLE          :: Ap, Ae, Aw, bh, Lw, Lp, Ue, y
INTEGER :: i

ALLOCATE(Ap(n1:n2), Ae(n1:n2), Aw(n1:n2),bh(n1:n2))
ALLOCATE(Lw(n1:n2), Lp(n1:n2), Ue(n1:n2), y(n1:n2),u1(n1:n2))


Aw=0; Ap=0; Ae=0; bh= 0 ; Lw = 0 ; Lp = 0 ; Ue = 0 ; y=0;

DO i = n1,n2
  Ae(i) = b
  Ap(i) = a
  Aw(i) = b
END DO
Ae(n1)=b1
Aw(n2)=b2

bh(n1)=c1
bh(n2)=c2

Lp(n1) = Ap(n1)
Ue(n1) = Ae(n1)/Lp(n1)

DO i = n1+1, n2-1
  Lw(i) = Aw(i)
  Lp(i) = Ap(i) - Lw(i)*Ue(i-1)
  Ue(i) = Ae(i)/Lp(i)
END DO

Lw(n2) = Aw(n2)
Lp(n2) = Ap(n2) - Lw(n2)*Ue(n2-1)

y(n1) = bh(n1)/Lp(n1)

DO i = n1+1, n2
  y(i) = (bh(i) - Lw(i)*y(i-1)) / Lp(i)
END DO

u1(n2) = y(n2)
DO i = n2-1, n1, -1
  u1(i) = y(i) - Ue(i)*u1(i+1)
END DO

DEALLOCATE(Ap, Ae, Aw,bh,Lw, Lp, Ue, y)

END SUBROUTINE

END PROGRAM po

我正在对边缘连续性问题进行数值/分析比较。为此，我为我的数值分辨率创建了一个子程序。该子例程是（LU）h = b分解/分辨率。它并不总是能获得任何价值的问题。

对于第一个CALL（u的第一个分辨率）的任何nc值，但对于第二个CALL（具有连续性边缘问题的up2的分辨率），它的所有nc值都可以正常工作，但并非总是如此。它仅适用于nc = 10,100,1000。...

例如，这是我对于uexactd / up2 的nc = 10的图形>

对于nc = 50（无效）

对于nc = 100，可再次使用

我无法确定是子程序问题还是直接使用CALL System编写的问题。我无法弄清楚为什么它仅适用于某些值。

程序po隐式无整数:: i，j，nc，nd，ok，iter REAL ::阿尔法，赌注，chi REAL，DIMENSION（:)，可分配:: uexact，...

Answer 1

我不确定您到底要做什么，但是您的文件o / p格式错误。