'********************************************************
'* This program calculates the discrete primitive of a  *
'* user defined real function F(x) in n points by using *
'* a classical Runge-Kutta method of order 4 (better    *
'* than a Simpson method of order 2).                   *
'* ---------------------------------------------------- *
'* SAMPLE RUN:                                          *
'* (Calculate 20 points of the primitive of F(x)=x      *
'*  from x=-5 to x=5 with a finesse of 2. the value of  *
'*  the primitive at x=-5 will be given as 25/2).       *
'*                                                      *
'* Calculate the discrete primitive of a function F(x)  *
'*                                                      *
'*   Value of x begin: -5                               *
'*   Value of x end  : 5                                *
'*   Value of primitive in x begin: 12.5                *
'*   Number of points to calculate: 20                  *
'*   Number of intermediate points: 2                   *
'*                                                      *
'*  x1 =  -4.500000  y1 =  10.125000                    *
'*  x2 =  -4.000000  y2 =   8.000000                    *
'*  x3 =  -3.500000  y3 =   6.125000                    *
'*  x4 =  -3.000000  y4 =   4.500000                    *
'*  x5 =  -2.500000  y5 =   3.125000                    *
'*  x6 =  -2.000000  y6 =   2.000000                    *
'*  x7 =  -1.500000  y7 =   1.125000                    *
'*  x8 =  -1.000000  y8 =   0.500000                    *
'* -----------------------------------                  *
'*  x16 =   3.000000  y16 =   4.500000                  *
'*  x17 =   3.500000  y17 =   6.125000                  *
'*  x18 =   4.000000  y18 =   8.000000                  *
'*  x19 =   4.500000  y19 =  10.125000                  *
'*  x20 =   5.000000  y20 =  12.500000                  *
'*                                                      *
'* (Exact values are yi = xi*xi/2).                     *
'*                                                      *
'*                  Basic Version By J-P Moreau, Paris. *
'*                           (www.jpmoreau.fr)          *
'********************************************************
'PROGRAM Test_Primitive
DEFDBL A-H, O-Z
DEFINT I-N

'integer ifi, n
'double  fa,xa,xb

  CLS
  PRINT
  PRINT " Calculate the discrete primitive of a function F(x)"
  PRINT
  INPUT "   Value of x begin: ", xa
  INPUT "   Value of x end  : ", xb
  INPUT "   Value of primitive in x begin: ", fa
  INPUT "   Number of points to calculate: ", n
  INPUT "   Number of intermediate points: ", ifi
    
  PRINT : PRINT

  GOSUB 1000  'call Primitive(xa,xb,fa,n,ifi)

  PRINT : PRINT

END


'given user defined function F(x)
500 'Function F(xx)
  F = xx
RETURN


'*******************************************************
'* Calculate the discrete primitive of a real function *
'* F(x) by using a Runge-Kutta method of order 4.      *
'* --------------------------------------------------- *
'* INPUTS:                                             *
'*           xa       begin x                          *
'*           xb       end x                            *
'*           fa       value of primitive at x=a        *
'*           n        number of calculated points      *
'*           ifi      number of intermediate points    *
'*                    (finesse)                        *
'*                                                     *
'* OUTPUTS:                                            *
'*           The subroutine displays the n calculated  *
'*           points (xi,yi) of the primitive of the    *
'*           user defined function F(x).               *
'* --------------------------------------------------- *
'* Ref.: "Mathematiques en Turbo Pascal By Alain Rever-*
'*        chon and Marc Ducamp, Armand Colin Editor,   *
'*        Paris, 1990" [BIBLI 03].                     *
'*******************************************************
1000 'Subroutine Primitive(xa, xb, fa, n, ifi)
'h,h2,h6,xly,xlr,x,y,u: doubles
  h = (xb - xa) / ifi / n: h2 = h / 2: h6 = h / 6
  xlr = fa
  xx = xa: GOSUB 500: xly = F             'xly=F(xa)
  FOR i = 1 TO n
    ni = (i - 1) * ifi - 1
    FOR j = 1 TO ifi
      x = xa + h * (ni + j)
      xx = x + h2: GOSUB 500: u = 4# * F  'u=4*F(x+h2)
      x = x + h
      xx = x: GOSUB 500: y = F            'y=F(x)
      xlr = xlr + h6 * (xly + u + y): xly = y
    NEXT j
    F$ = "###.######"                     'format real numbers
    PRINT " x"; i; " = "; : PRINT USING F$; x;
    PRINT " y"; i; " = "; : PRINT USING F$; xlr
  NEXT i
  RETURN
RETURN

' end of file primitiv.bas