fortran question

[Karkouri]

Hello people ,can anyone help me with .

i have the following and want it to put it in fortran so as to find the result.

http://img22.imageshack.us/img22/8291/typos.jpg

and y*=0.55y² + 0.0165y + 0.00012375

 

[Karkouri]

and yb=0.1
yt=0.005

 

[mikrom]

..have the following and want it to put it in fortran..

When you want to to put it in fortran or any other programming language you need first to know how you want to compute it, also the algorithm.
There are several formulas for numerical integration, for example the Newton–Cotes formulas

 

[Karkouri]

ok but i cant learn fortran in one day .
can you just tell me the algorithm for my example?

 

[mikrom]

Hi Karkouri,

What you ask is not a fortran speciifc question. It's more a basic question from numerical mathematics. You should be clear how the numerical integration works and what method you have to choose. Then try to code it and post your code here.

This is a free help forum, but not a free coding service.
You should first show some effort. Try to code self something and if you have a problem post your code here. Then we try to help you.

You can google for infos about numerical integration incl. fortran code. E.g. here is something:

http://oregonstate.edu/instruct/ch490/lessons/lesson11.htm

 

[Karkouri]

mirkom , my problem is not to find the answer using mathematics (i already did this). I just want to solve this problem using fortran with any method. As i 'm not expert in this programming language , i cannot write any algorithm.
i dont even know any commands in fortran.

I m asking if anyone knows how to type my example in fortran so as to get a result.
I didnt ask something so tragic imo , some help maybe by anyone?

 

[mikrom]

Hi Karkouri,

Here is an example for you taken from the pseudo (python) code on wikipedia

http://en.wikipedia.org/wiki/Trapezoid_rule

http://en.wikipedia.org/wiki/Simpson's_rule

integration.f95

[COLOR=#a020f0]program[/color] integration
 [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
 [COLOR=#2e8b57][b]external[/b][/color] f1, f2 [COLOR=#0000ff]! define function to be integrated as external[/color]
[COLOR=#2e8b57][b] real[/b][/color] :: integral
 [COLOR=#0000ff]! integrating he function f1[/color]
 integral[COLOR=#804040][b]=[/b][/color]trapezoid_rule(f1, [COLOR=#ff00ff]0.0[/color], [COLOR=#ff00ff]1.0[/color], [COLOR=#ff00ff]100[/color])
 [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]'Trapezoid rule = '[/color], integral
 integral[COLOR=#804040][b]=[/b][/color]simpson_rule(f1, [COLOR=#ff00ff]0.0[/color], [COLOR=#ff00ff]1.0[/color], [COLOR=#ff00ff]100[/color])
 [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]"Simpson's rule = "[/color], integral

 [COLOR=#0000ff]! integrating the function f2[/color]
 integral[COLOR=#804040][b]=[/b][/color]trapezoid_rule(f2, [COLOR=#ff00ff]0.005[/color], [COLOR=#ff00ff]0.1[/color], [COLOR=#ff00ff]100[/color])
 [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]'Trapezoid rule = '[/color], integral
 integral[COLOR=#804040][b]=[/b][/color]simpson_rule(f2, [COLOR=#ff00ff]0.005[/color], [COLOR=#ff00ff]0.1[/color], [COLOR=#ff00ff]100[/color])
 [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]"Simpson's rule = "[/color], integral

[COLOR=#a020f0]contains[/color]
[COLOR=#0000ff]! ********************** Functions/Procedures *************************[/color]
[COLOR=#2e8b57][b]real[/b][/color] [COLOR=#a020f0]function[/color] trapezoid_rule(f, a, b, n)
  [COLOR=#0000ff]! Approximate the definite integral of f from a to b by the[/color]
  [COLOR=#0000ff]! composite trapezoidal rule, using N subintervals[/color]
  [COLOR=#0000ff]! see: <a href="[URL unfurl="true"]http://en.wikipedia.org/wiki/Trapezoid_rule">http://en.wikipedia.org/wiki/Trapezoid_rule</a>[/URL][/color]
  [COLOR=#0000ff]! function arguments ---------------------[/color]
  [COLOR=#0000ff]! f: real function[/color]
[COLOR=#2e8b57][b]  real[/b][/color] :: f
  [COLOR=#0000ff]! [a, b] : the interval of integration[/color]
[COLOR=#2e8b57][b]  real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) ::  a, b
  [COLOR=#0000ff]! n : number of subintervals used[/color]
  [COLOR=#2e8b57][b]integer[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: n
  [COLOR=#0000ff]! ----------------------------------------[/color]
  [COLOR=#0000ff]! temporary variables[/color]
  [COLOR=#2e8b57][b]integer[/b][/color] :: k
[COLOR=#2e8b57][b]  real[/b][/color] :: s
  s [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0[/color]
  [COLOR=#804040][b]do[/b][/color] k[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n[COLOR=#804040][b]-[/b][/color][COLOR=#ff00ff]1[/color]
    s [COLOR=#804040][b]=[/b][/color] s [COLOR=#804040][b]+[/b][/color] f(a [COLOR=#804040][b]+[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a)[COLOR=#804040][b]*[/b][/color]k[COLOR=#804040][b]/[/b][/color]n)
  [COLOR=#804040][b]end do[/b][/color]  
  trapezoid_rule [COLOR=#804040][b]=[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a) [COLOR=#804040][b]*[/b][/color] (f(a)[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color] [COLOR=#804040][b]+[/b][/color] f(b)[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color] [COLOR=#804040][b]+[/b][/color] s) [COLOR=#804040][b]/[/b][/color] n
[COLOR=#a020f0]end function[/color] trapezoid_rule

[COLOR=#2e8b57][b]real[/b][/color] [COLOR=#a020f0]function[/color] simpson_rule(f, a, b, n)
  [COLOR=#0000ff]! Approximate the definite integral of f from a to b by the[/color]
  [COLOR=#0000ff]! composite Simpson's rule, using N subintervals[/color]
  [COLOR=#0000ff]! see: <a href="[URL unfurl="true"]http://en.wikipedia.org/wiki/Simpson%27s_rule">http://en.wikipedia.org/wiki/Simpson%27s_rule</a>[/URL][/color]
  [COLOR=#0000ff]! function arguments ---------------------[/color]
  [COLOR=#0000ff]! f: real function[/color]
[COLOR=#2e8b57][b]  real[/b][/color] :: f
  [COLOR=#0000ff]! [a, b] : the interval of integration[/color]
[COLOR=#2e8b57][b]  real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) ::  a, b
  [COLOR=#0000ff]! n : number of subintervals used[/color]
  [COLOR=#2e8b57][b]integer[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: n
  [COLOR=#0000ff]! ----------------------------------------[/color]
  [COLOR=#0000ff]! temporary variables[/color]
  [COLOR=#2e8b57][b]integer[/b][/color] :: k
[COLOR=#2e8b57][b]  real[/b][/color] :: s
  s [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0[/color]
  [COLOR=#804040][b]do[/b][/color] k[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n[COLOR=#804040][b]-[/b][/color][COLOR=#ff00ff]1[/color]
    s [COLOR=#804040][b]=[/b][/color] s [COLOR=#804040][b]+[/b][/color] [COLOR=#ff00ff]2[/color][COLOR=#804040][b]**[/b][/color]([COLOR=#008080]mod[/color](k,[COLOR=#ff00ff]2[/color])[COLOR=#804040][b]+[/b][/color][COLOR=#ff00ff]1[/color]) [COLOR=#804040][b]*[/b][/color] f(a [COLOR=#804040][b]+[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a)[COLOR=#804040][b]*[/b][/color]k[COLOR=#804040][b]/[/b][/color]n)
  [COLOR=#804040][b]end do[/b][/color]  
  simpson_rule [COLOR=#804040][b]=[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a) [COLOR=#804040][b]*[/b][/color] (f(a) [COLOR=#804040][b]+[/b][/color] f(b) [COLOR=#804040][b]+[/b][/color] s) [COLOR=#804040][b]/[/b][/color] ([COLOR=#ff00ff]3[/color][COLOR=#804040][b]*[/b][/color]n)
[COLOR=#a020f0]end function[/color] simpson_rule

[COLOR=#a020f0]end program[/color] integration 

[COLOR=#0000ff]! define here your own functions to be integrated [/color]
[COLOR=#2e8b57][b]real[/b][/color] [COLOR=#a020f0]function[/color] f1(x)
  [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
[COLOR=#2e8b57][b]  real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: x
  f1 [COLOR=#804040][b]=[/b][/color] x[COLOR=#804040][b]*[/b][/color]x
[COLOR=#a020f0]end function[/color] f1

[COLOR=#2e8b57][b]real[/b][/color] [COLOR=#a020f0]function[/color] f2(y)
  [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
[COLOR=#2e8b57][b]  real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: y
  f2 [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color][COLOR=#804040][b]/[/b][/color](y[COLOR=#804040][b]-[/b][/color]([COLOR=#ff00ff]0.55[/color][COLOR=#804040][b]*[/b][/color]y[COLOR=#804040][b]**[/b][/color][COLOR=#ff00ff]2[/color] [COLOR=#804040][b]+[/b][/color] [COLOR=#ff00ff]0.0165[/color][COLOR=#804040][b]*[/b][/color]y [COLOR=#804040][b]+[/b][/color] [COLOR=#ff00ff]0.00012375[/color]))
[COLOR=#a020f0]end function[/color] f2

Output:

$ g95 integration.f95 -o integration

$ integration
 Trapezoid rule =  0.33335
 Simpson's rule =  0.3333334
 Trapezoid rule =  3.1299372
 Simpson's rule =  3.1267815

The first funcion f1(x)=x^2 which integrated on [0,1] gives exactly 1/3. You can compare here the accuracy of the 2 methods (trapezoid and Simpson's rule)
The secon function is yours. You can increase the value of n to get better accuracy. The result seems to be 3.1267..

 

[mikrom]

What I posted above compiles fine with g95 but doesn't compile with gfortran because of these errors

$ gfortran integration.f95 -o integration
integration.f95:6.25:

 integral=trapezoid_rule(f1, 0.0, 1.0, 100)
                        1
Error: Type/rank mismatch in argument 'f' at (1)
integration.f95:8.23:

 integral=simpson_rule(f1, 0.0, 1.0, 100)
                      1
Error: Type/rank mismatch in argument 'f' at (1)
integration.f95:12.25:

So, here is another version which uses modules. I created 2 modules: one for the numerical methods and other for the user defined functions.
This version compiles fine with both compilers

integration.f95

[COLOR=#a020f0]module[/color] integration_methods
  [COLOR=#0000ff]! This module contains the integration methods  [/color]
[COLOR=#a020f0]contains[/color]
[COLOR=#2e8b57][b]  real[/b][/color] [COLOR=#a020f0]function[/color] trapezoid_rule(f, a, b, n)
    [COLOR=#0000ff]! Approximate the definite integral of f from a to b by the[/color]
    [COLOR=#0000ff]! composite trapezoidal rule, using N subintervals[/color]
    [COLOR=#0000ff]! see: <a href="[URL unfurl="true"]http://en.wikipedia.org/wiki/Trapezoid_rule">http://en.wikipedia.org/wiki/Trapezoid_rule</a>[/URL][/color]
    [COLOR=#0000ff]! function arguments ---------------------[/color]
    [COLOR=#0000ff]! f: real function[/color]
[COLOR=#2e8b57][b]    real[/b][/color] :: f
    [COLOR=#0000ff]! [a, b] : the interval of integration[/color]
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) ::  a, b
    [COLOR=#0000ff]! n : number of subintervals used[/color]
    [COLOR=#2e8b57][b]integer[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: n
    [COLOR=#0000ff]! ----------------------------------------[/color]
    [COLOR=#0000ff]! temporary variables[/color]
    [COLOR=#2e8b57][b]integer[/b][/color] :: k
[COLOR=#2e8b57][b]    real[/b][/color] :: s
    s [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0[/color]
    [COLOR=#804040][b]do[/b][/color] k[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n[COLOR=#804040][b]-[/b][/color][COLOR=#ff00ff]1[/color]
      s [COLOR=#804040][b]=[/b][/color] s [COLOR=#804040][b]+[/b][/color] f(a [COLOR=#804040][b]+[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a)[COLOR=#804040][b]*[/b][/color]k[COLOR=#804040][b]/[/b][/color]n)
    [COLOR=#804040][b]end do[/b][/color]  
    trapezoid_rule [COLOR=#804040][b]=[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a) [COLOR=#804040][b]*[/b][/color] (f(a)[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color] [COLOR=#804040][b]+[/b][/color] f(b)[COLOR=#804040][b]/[/b][/color][COLOR=#ff00ff]2[/color] [COLOR=#804040][b]+[/b][/color] s) [COLOR=#804040][b]/[/b][/color] n
  [COLOR=#a020f0]end function[/color] trapezoid_rule

[COLOR=#2e8b57][b]  real[/b][/color] [COLOR=#a020f0]function[/color] simpson_rule(f, a, b, n)
    [COLOR=#0000ff]! Approximate the definite integral of f from a to b by the[/color]
    [COLOR=#0000ff]! composite Simpson's rule, using N subintervals[/color]
    [COLOR=#0000ff]! see: <a href="[URL unfurl="true"]http://en.wikipedia.org/wiki/Simpson%27s_rule">http://en.wikipedia.org/wiki/Simpson%27s_rule</a>[/URL][/color]
    [COLOR=#0000ff]! function arguments ---------------------[/color]
    [COLOR=#0000ff]! f: real function[/color]
[COLOR=#2e8b57][b]    real[/b][/color] :: f
    [COLOR=#0000ff]! [a, b] : the interval of integration[/color]
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) ::  a, b
    [COLOR=#0000ff]! n : number of subintervals used[/color]
    [COLOR=#2e8b57][b]integer[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: n
    [COLOR=#0000ff]! ----------------------------------------[/color]
    [COLOR=#0000ff]! temporary variables[/color]
    [COLOR=#2e8b57][b]integer[/b][/color] :: k
[COLOR=#2e8b57][b]    real[/b][/color] :: s
    s [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]0[/color]
    [COLOR=#804040][b]do[/b][/color] k[COLOR=#804040][b]=[/b][/color][COLOR=#ff00ff]1[/color], n[COLOR=#804040][b]-[/b][/color][COLOR=#ff00ff]1[/color]
      s [COLOR=#804040][b]=[/b][/color] s [COLOR=#804040][b]+[/b][/color] [COLOR=#ff00ff]2[/color][COLOR=#804040][b]**[/b][/color]([COLOR=#008080]mod[/color](k,[COLOR=#ff00ff]2[/color])[COLOR=#804040][b]+[/b][/color][COLOR=#ff00ff]1[/color]) [COLOR=#804040][b]*[/b][/color] f(a [COLOR=#804040][b]+[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a)[COLOR=#804040][b]*[/b][/color]k[COLOR=#804040][b]/[/b][/color]n)
    [COLOR=#804040][b]end do[/b][/color]  
    simpson_rule [COLOR=#804040][b]=[/b][/color] (b[COLOR=#804040][b]-[/b][/color]a) [COLOR=#804040][b]*[/b][/color] (f(a) [COLOR=#804040][b]+[/b][/color] f(b) [COLOR=#804040][b]+[/b][/color] s) [COLOR=#804040][b]/[/b][/color] ([COLOR=#ff00ff]3[/color][COLOR=#804040][b]*[/b][/color]n)
  [COLOR=#a020f0]end function[/color] simpson_rule
[COLOR=#a020f0]end module[/color] integration_methods

[COLOR=#a020f0]module[/color] user_functions
  [COLOR=#0000ff]! Define here your own functions to be integrated [/color]
[COLOR=#a020f0]contains[/color]
[COLOR=#2e8b57][b]  real[/b][/color] [COLOR=#a020f0]function[/color] f1(x)
    [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: x
    f1 [COLOR=#804040][b]=[/b][/color] x[COLOR=#804040][b]*[/b][/color]x
  [COLOR=#a020f0]end function[/color] f1

[COLOR=#2e8b57][b]  real[/b][/color] [COLOR=#a020f0]function[/color] f2(y)
    [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
[COLOR=#2e8b57][b]    real[/b][/color], [COLOR=#2e8b57][b]intent[/b][/color]([COLOR=#2e8b57][b]in[/b][/color]) :: y
    f2 [COLOR=#804040][b]=[/b][/color] [COLOR=#ff00ff]1[/color][COLOR=#804040][b]/[/b][/color](y[COLOR=#804040][b]-[/b][/color]([COLOR=#ff00ff]0.55[/color][COLOR=#804040][b]*[/b][/color]y[COLOR=#804040][b]**[/b][/color][COLOR=#ff00ff]2[/color] [COLOR=#804040][b]+[/b][/color] [COLOR=#ff00ff]0.0165[/color][COLOR=#804040][b]*[/b][/color]y [COLOR=#804040][b]+[/b][/color] [COLOR=#ff00ff]0.00012375[/color]))
  [COLOR=#a020f0]end function[/color] f2
[COLOR=#a020f0]end module[/color] user_functions

[COLOR=#a020f0]program[/color] integration
  [COLOR=#a020f0]use[/color] integration_methods
  [COLOR=#a020f0]use[/color] user_functions
  
  [COLOR=#2e8b57][b]implicit[/b][/color] [COLOR=#2e8b57][b]none[/b][/color]
[COLOR=#2e8b57][b]  real[/b][/color] :: integral

  [COLOR=#0000ff]! integrating he function f1[/color]
  integral[COLOR=#804040][b]=[/b][/color]trapezoid_rule(f1, [COLOR=#ff00ff]0.0[/color], [COLOR=#ff00ff]1.0[/color], [COLOR=#ff00ff]1000[/color])
  [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]'Trapezoid rule = '[/color], integral
  integral[COLOR=#804040][b]=[/b][/color]simpson_rule(f1, [COLOR=#ff00ff]0.0[/color], [COLOR=#ff00ff]1.0[/color], [COLOR=#ff00ff]1000[/color])
  [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]"Simpson's rule = "[/color], integral

  [COLOR=#0000ff]! integrating the function f2[/color]
  integral[COLOR=#804040][b]=[/b][/color]trapezoid_rule(f2, [COLOR=#ff00ff]0.005[/color], [COLOR=#ff00ff]0.1[/color], [COLOR=#ff00ff]1000[/color])
  [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]'Trapezoid rule = '[/color], integral
  integral[COLOR=#804040][b]=[/b][/color]simpson_rule(f2, [COLOR=#ff00ff]0.005[/color], [COLOR=#ff00ff]0.1[/color], [COLOR=#ff00ff]1000[/color])
  [COLOR=#804040][b]write[/b][/color] ([COLOR=#804040][b]*[/b][/color],[COLOR=#804040][b]*[/b][/color]) [COLOR=#ff00ff]"Simpson's rule = "[/color], integral
[COLOR=#a020f0]end program[/color] integration

Output:

$ g95 integration.f95 -o integration

$ integration
 Trapezoid rule =  0.33333343
 Simpson's rule =  0.33333328
 Trapezoid rule =  3.1267667
 Simpson's rule =  3.126737

$ gfortran integration.f95 -o integration

$ integration
 Trapezoid rule =   0.33333343    
 Simpson's rule =   0.33333328    
 Trapezoid rule =    3.1267667    
 Simpson's rule =    3.1267371