Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Fortran bisection method

  1. Apr 18, 2013 #1

    fluidistic

    User Avatar
    Gold Member

    I'm rusty with fortran and programming in general. I can't see my "error" in a code that I wrote from scratch.
    Basically I wanted to get some fun and solve for a temperature in thermodynamics where I must get "T_f" which appear in a transcendental equation: ##A\ln \left ( \frac{T_f^2}{T_1T_2} \right )+2B[2T_f-(T_1T_2)]=0##.
    Where T_1 is worth 100, T_2 is worth 200, A=2 and B=0.005.
    Here is my code:
    Code (Text):
    program rootfinderbisection
    implicit none

    real(4)::a,b,xn,fa,fb,fn,f,e,x_n,t_1,t_2


    write(*,*)"Enter the temperatures t_1 and t_2"
    read(*,*)t_1, t_2
    write(*,*)"Enter your interval guess, [a,b]"
    read(*,*)a,b
    write(*,*)"Enter the precision you want"
    read(*,*)e
    fa=f(a)
    fb=f(b)

    if (abs(fa)<e) then
    write(*,*) "The solution is",a
    stop
    else if (abs(fb)<e) then
    write(*,*)"The solution is",b
    stop
    end if

    if (f(a)*f(b)>0) then
    write(*,*)"There's an even number of root(s) in that interval"
    stop
    endif

    xn=0.5*(a+b)
    fn=f(xn)

    do while (abs(fn)>e)
    if (fa*fn<0) then
    b=xn
    fb=fn
    else
    a=xn
    fa=fn
    endif

    xn=0.5*(a+b)
    fn=f(xn)
    end do
    write(*,*)"The solution is",x_n


    end program
    that I compile along with
    Code (Text):
    function f(x)
    implicit none

    real(4)::f,x, t_1,t_2
    f=2.*log((x**2.)/(t_1*t_2))+2.*0.005*(2.*x-(t_1+t_2))

    end function
    The problem is that for the values I said above, I get the message "There's an even number of root(s) in that interval". Thus apparently the condition f(a)*f(b)>0 is satisfied. However when I enter the numbers in my pocket calculator I get f(a)=-2.3862... and f(b)=2.3862... where I truncated both values. Thus the condition is not satisfied.
    I really don't see what's wrong. Any idea?

    Edit: I tried with many different values for a and b. I even changed t_1 and t_2, same message.

    Edit2: I've just tested the values of fa and fb. It thinks they are worth infinity. Going to check this out.
    Edit3: No idea why it thinks they are worth infinity. o_0.
    Edit4: ok it does not plug in t_1 and t_2 inside my function.
     
    Last edited: Apr 18, 2013
  2. jcsd
  3. Apr 19, 2013 #2

    SteamKing

    User Avatar
    Staff Emeritus
    Science Advisor
    Homework Helper

    Try including t_1 and t_2 as arguments in the function f and see if this resolves the calculation errors.
     
  4. Apr 19, 2013 #3
    Your function "f(x)" does not agree with your formula at the top of your post. Is the last term in the bracket term of your equation "T_1 * T_2" or "T_1 + T_2".
     
  5. Apr 19, 2013 #4

    jtbell

    User Avatar

    Staff: Mentor

    t_1 and t_2 in the function are completely separate variables from t_1 and t_2 in the main program. If you want t_1 and t_2 in the function to have the same values as in main, you need to make them arguments to the function along with x.
     
  6. Apr 19, 2013 #5

    fluidistic

    User Avatar
    Gold Member

    Thanks all for the tips. To TheoMcCloskey, you're right, I made a latex typo here for T_f. The one in my program is ok.
    So it got rid of the problem! However I now get "The solution is 3.02481484E-39" which is basically 0. I'll investigate and post if stuck.

    Edit: Worked, problem solved. I had to change the line "write(*,*)"The solution is",x_n" to "write(*,*)"The solution is",xn".
    Thanks a lot guys!
     
    Last edited: Apr 19, 2013
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook




Similar Discussions: Fortran bisection method
  1. Bisection Method (Replies: 2)

  2. Bisection method by c (Replies: 10)

Loading...