如何在 do 循环中区分结果和前一个结果?
How can I make the difference between a result and the preceding one in do loops?
本子程序用来判断复合梯形所以
我想在最终结果(积分)和前一个结果(积分-1)之间进行抽象(差异),并将差异用作重复我的间隔数的限制。
Subroutine Trapezoid(a,b,n,integration)
real,external :: f
real :: h,a,b,summ,p
real,intent(out) :: integration
integer :: n
integer :: i,j
!Here as we have the whole equation is (h/2)*[f(a)+f(b)+2*sum(Xi)
!So we calculate the first part (h/2)*[f(a)+f(b) and then calculate the anoter part
do i=1,n
n=2**i !Double the number of interval
h=(b-a)/n !Calculate the delta X
p=(h/2.)*(f(a)+f(b))
summ=0
do j=1,n-1
summ=summ+h*f(a+j*h) !h/2 *2* sum[f(Xi)
enddo
if(n == 256) then !put a limit for the number of interval
Stop
end if
integration = p + summ !Here the sum the both parts
print*,n,' ',integration
enddo
end Subroutine
所以我想确定差异而不是限制是 250,当这个差异小于 10*-8 时,停止
我尝试了很多,但我没有得到我想要的。
我会像下面那样做(很快就拼凑起来)。请注意,默认种类实数 1e-8 是不切实际的期望精度 - 因此容差较低。如果你想要更高的精度,你需要使用更高精度的 real。
另请注意,我已将其变成一个完整的程序。有问题请自己做。纯粹自私地说,您将更有可能获得有用的答案。
无论如何这是代码
Program integ
Implicit None
Real, Parameter :: pi = 3.1415927
Real :: value, delta
Integer :: n_used
Intrinsic :: sin
Call Trapezoid( sin, 0.0, pi / 2.0, 20, n_used, value, delta )
Write( *, * ) 'final result', value, ' with ', 2 ** n_used, ' intervals'
Contains
Subroutine Trapezoid(f,a,b,n_max,n_used,integration,delta)
Implicit None
Real, Parameter :: tol = 1e-4
Interface
Function f( x ) Result( r )
Real :: r
Real, Intent( In ) :: x
End Function f
End Interface
Real , Intent( In ) :: a
Real , Intent( In ) :: b
Integer, Intent( In ) :: n_max
Integer, Intent( Out ) :: n_used
Real , Intent( Out ) :: integration
Real , Intent( Out ) :: delta
Real :: h,summ,p
Real :: integration_old
Integer :: n
Integer :: i,j
!Here as we have the whole equation is (h/2)*[f(a)+f(b)+2*sum(Xi)
!So we calculate the first part (h/2)*[f(a)+f(b) and then calculate the anoter part
delta = - Huge( delta )
integration_old = Huge( integration_old )
Do i=1,n_max
n=2**i !Double the number of interval
h=(b-a)/n !Calculate the delta X
p=(h/2.)*(f(a)+f(b))
summ=0
Do j=1,n-1
summ=summ+h*f(a+j*h) !h/2 *2* sum[f(Xi)
Enddo
integration = p + summ !Here the sum the both parts
If( i /= 1 ) Then
delta = integration - integration_old
Write( *, * ) n,' ',integration , delta
If( Abs( delta ) < tol ) Exit
End If
integration_old = integration
Enddo
n_used = i
End Subroutine Trapezoid
End Program
ian@eris:~/work/stack$ gfortran --version
GNU Fortran (Ubuntu 7.4.0-1ubuntu1~18.04.1) 7.4.0
Copyright (C) 2017 Free Software Foundation, Inc.
This is free software; see the source for copying conditions. There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
ian@eris:~/work/stack$ gfortran -Wall -Wextra -fcheck=all -O -std=f2008 integ.f90
ian@eris:~/work/stack$ ./a.out
4 0.987115800 3.90563607E-02
8 0.996785223 9.66942310E-03
16 0.999196708 2.41148472E-03
32 0.999799252 6.02543354E-04
64 0.999949872 1.50620937E-04
128 0.999987483 3.76105309E-05
final result 0.999987483 with 128 intervals
本子程序用来判断复合梯形所以
我想在最终结果(积分)和前一个结果(积分-1)之间进行抽象(差异),并将差异用作重复我的间隔数的限制。
Subroutine Trapezoid(a,b,n,integration)
real,external :: f
real :: h,a,b,summ,p
real,intent(out) :: integration
integer :: n
integer :: i,j
!Here as we have the whole equation is (h/2)*[f(a)+f(b)+2*sum(Xi)
!So we calculate the first part (h/2)*[f(a)+f(b) and then calculate the anoter part
do i=1,n
n=2**i !Double the number of interval
h=(b-a)/n !Calculate the delta X
p=(h/2.)*(f(a)+f(b))
summ=0
do j=1,n-1
summ=summ+h*f(a+j*h) !h/2 *2* sum[f(Xi)
enddo
if(n == 256) then !put a limit for the number of interval
Stop
end if
integration = p + summ !Here the sum the both parts
print*,n,' ',integration
enddo
end Subroutine
所以我想确定差异而不是限制是 250,当这个差异小于 10*-8 时,停止 我尝试了很多,但我没有得到我想要的。
我会像下面那样做(很快就拼凑起来)。请注意,默认种类实数 1e-8 是不切实际的期望精度 - 因此容差较低。如果你想要更高的精度,你需要使用更高精度的 real。
另请注意,我已将其变成一个完整的程序。有问题请自己做。纯粹自私地说,您将更有可能获得有用的答案。
无论如何这是代码
Program integ
Implicit None
Real, Parameter :: pi = 3.1415927
Real :: value, delta
Integer :: n_used
Intrinsic :: sin
Call Trapezoid( sin, 0.0, pi / 2.0, 20, n_used, value, delta )
Write( *, * ) 'final result', value, ' with ', 2 ** n_used, ' intervals'
Contains
Subroutine Trapezoid(f,a,b,n_max,n_used,integration,delta)
Implicit None
Real, Parameter :: tol = 1e-4
Interface
Function f( x ) Result( r )
Real :: r
Real, Intent( In ) :: x
End Function f
End Interface
Real , Intent( In ) :: a
Real , Intent( In ) :: b
Integer, Intent( In ) :: n_max
Integer, Intent( Out ) :: n_used
Real , Intent( Out ) :: integration
Real , Intent( Out ) :: delta
Real :: h,summ,p
Real :: integration_old
Integer :: n
Integer :: i,j
!Here as we have the whole equation is (h/2)*[f(a)+f(b)+2*sum(Xi)
!So we calculate the first part (h/2)*[f(a)+f(b) and then calculate the anoter part
delta = - Huge( delta )
integration_old = Huge( integration_old )
Do i=1,n_max
n=2**i !Double the number of interval
h=(b-a)/n !Calculate the delta X
p=(h/2.)*(f(a)+f(b))
summ=0
Do j=1,n-1
summ=summ+h*f(a+j*h) !h/2 *2* sum[f(Xi)
Enddo
integration = p + summ !Here the sum the both parts
If( i /= 1 ) Then
delta = integration - integration_old
Write( *, * ) n,' ',integration , delta
If( Abs( delta ) < tol ) Exit
End If
integration_old = integration
Enddo
n_used = i
End Subroutine Trapezoid
End Program
ian@eris:~/work/stack$ gfortran --version
GNU Fortran (Ubuntu 7.4.0-1ubuntu1~18.04.1) 7.4.0
Copyright (C) 2017 Free Software Foundation, Inc.
This is free software; see the source for copying conditions. There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
ian@eris:~/work/stack$ gfortran -Wall -Wextra -fcheck=all -O -std=f2008 integ.f90
ian@eris:~/work/stack$ ./a.out
4 0.987115800 3.90563607E-02
8 0.996785223 9.66942310E-03
16 0.999196708 2.41148472E-03
32 0.999799252 6.02543354E-04
64 0.999949872 1.50620937E-04
128 0.999987483 3.76105309E-05
final result 0.999987483 with 128 intervals