Fortran bisection method

  • Fortran
  • Thread starter fluidistic
  • Start date
  • #1
fluidistic
Gold Member
3,744
128
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:
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:
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:

Answers and Replies

  • #2
SteamKing
Staff Emeritus
Science Advisor
Homework Helper
12,796
1,668
Try including t_1 and t_2 as arguments in the function f and see if this resolves the calculation errors.
 
  • #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".
 
  • #4
jtbell
Mentor
15,722
3,858
Edit4: ok it does not plug in t_1 and t_2 inside my function.
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.
 
  • #5
fluidistic
Gold Member
3,744
128
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:

Related Threads on Fortran bisection method

Replies
17
Views
19K
  • Last Post
Replies
2
Views
4K
  • Last Post
Replies
10
Views
43K
Replies
6
Views
2K
  • Last Post
Replies
1
Views
11K
  • Last Post
Replies
5
Views
7K
Replies
4
Views
3K
Replies
11
Views
778
Replies
1
Views
6K
Top