積分計算使用梯形或辛普森規則

-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

積分計算使用梯形或辛普森規則

回答

相關問題