积分计算使用梯形或辛普森规则

-1

基本上，分配“说找到解决方案f（x）= 0 ....”。其中从x = 0到x-0.45，f（x）=积分（1/sqrt（2 * pi））exp（-t^2/2）dt。对不起，我无法找到如何把积分放在这里。积分计算使用梯形或辛普森规则

我的教授。给了我们一些提示如何开始。

start x 
do i = 1,2,3... 
    fp = 1/sqrt(2*pi)exp(-x^2/2) 
    f = use trap,or simpson's rule to find the integration than subtract 0.45 
    x = x - (f/fp) 
end do

这里是我做过什么

program main 
implicit none 

integer :: n, k, i 
double precision :: h, a, fp, f, x1, x2, pi, blub, integ, e, dx, j, m 

a = 0 
n = 25 
x1 = 0.5 
pi = 4.0d0 * atan(1.0d0) 

do i = 1, n 
    e = exp((-x1)**2.0d0/2.0d0) 
    print*, e 

    fp = (1.0d0/sqrt(2.0d0 * pi))* e 

    !start of the trapezoid 
    dx = (x1-a)/n !find dx by subtracting starting point from endpoint and divide by number of steps taken 
    m = (blub(a) + blub(x1))/2 !find the mean value of the integral 


    j = 0 
     do k=1, n-1 
      h = i 
      j = j + blub(h) !calculate the each step of the integral from one to n and add all together. 
     end do 

    integ = dx*(m+j) !result of the trapezoid method 

    print*, "integ: ", integ 
    f = integ - 0.45 
    a = x1 
    x1 = a - (f/fp) 
    print*, "x1: ",x1 
    print*, "a: ",a 
    print*, "f: ",f 
    print*, "fp: ", fp 
    print*, (x1-a)/n 
end do 

stop 
end 

double precision function blub (x) 
double precision :: x 

blub(x) = (1.0d0/sqrt(2.0d0 * (4.0d0 * atan(1.0d0))))*exp((-x)**2.0d0/2.0d0) 

return 
end

我做了所有的印刷找出我的错误可能。我意识到我的DX正在变得小，基于x1小于n。因为我的整合变成e^-310的数字，对0.45几乎没有影响，因此我的f保持不变，并且x1向上混乱...

如果有人可以向我解释我如何解决这个错误。

编辑：我确定我对梯形部分犯了一个错误，我只是不知道在哪里。 f应该随着每个循环而改变，但dx太小会阻止这种情况发生。

来源

2014-12-01 un known

你会得到哪些错误？ – 2014-12-01 17:43:13

他们都是语法错误。像无效字符一样，将一个变量从int变为real。我一直工作到凌晨5点，并开始失去焦点 – 2014-12-01 17:47:15

因此，您只是在这里发布它，以便某人完成您的硬件，甚至在您发现问题时甚至不给他提示？ – 2014-12-01 18:18:59

其中一个错误是这样的。我不要求它是唯一一个：

double precision function blub (x) 
double precision :: x 

blub(x) = (1.0d0/sqrt(2.0d0 * (4.0d0 * atan(1.0d0))))*exp((-x)**2.0d0/2.0d0) 

return 
end

如果用gfortran -Wall -std=f95 -c fun.f90编译你

fun.f90:5.76: 

lub(x) = (1.0d0/sqrt(2.0d0 * (4.0d0 * atan(1.0d0))))*exp((-x)**2.0d0/2.0d0) 
                      1 
Warning: Obsolescent feature: Statement function at (1) 
fun.f90:1.30: 

double precision function blub (x) 
           1 
Warning: Return value of function 'blub' at (1) not set

这可能提醒你一些奇怪的事情。实际上，您已经在原始函数内创建了一个新的语句函数（一个过时的Fortran函数），而不是给它正确的值。您应该使用：

double precision function blub (x) 
    implicit none 
    double precision, intent(in) :: x 

    blub = (1/sqrt(8 * atan(1.0d0))) * exp((-x)**2/2) 
end function

你可以忘记return末stop你end前放。没有理由使用它们。

我建议您使用更多的调试编译器标志，例如-g -fbacktrace -fcheck=all，并将所有功能和子程序放入模块中或使其成为内部模块（使用contains）。

来源

2014-12-01 19:45:36

积分计算使用梯形或辛普森规则

回答

相关问题