Friday, March 9, 2018

DELPHI - How compute integral - (sin(x))^2 dx from 0 to PI

For integration is used class TIntegral which use for calculation trapezoidal rule.
procedure TForm1.BIntegralClick(Sender: TObject);
var
  i : integer;
  iX, iY : double;
  pIntegral : TIntegral;
  pSerie : TLineSeries;
begin
    pSerie := TLineSeries( Chart1.Series[0] );
    pSerie.Clear;

    { -- prepare integral class - interval 0..PI }

    pIntegral := TIntegral.Create;
    pIntegral.SetInterval( 0, pi );

    for i := 0 to pIntegral.pX.Count - 1 do
      begin
        { calc Y value }
        pIntegral.pY.AddValue( sqr( sin( pIntegral.pX.GetValue( i ) ) ) );

        { -> to graph }
        iX := pIntegral.pX.GetValue( i );
        iY := pIntegral.pY.GetValue( i );
        pSerie.AddXY( iX, iY );
      end;

    ShowMessage( FloatToStr( pIntegral.Calc ) );        // 1.57

    pIntegral.Free;
end;
Here is TIntegral.calc() method:
{ ---------------------------------------------------------------------------
  Method computes integral value (trapezoidal rule).  
  -------------------------------------------------------------------------- }
function TIntegral.Calc( _iCheckSumIntegralValue : double = 999 ) : double;
var
  i, iDeleteXFrom : integer;
  iValue : double;
begin
  result := 0;

  { -- checks }

  if pX.Count < 2 then exit;
  if pX.Count <> pY.Count then exit;

  { -- delete values }

  pIntegral.Clear;
  pCumulation.Clear;

  { -- trapezoidal rule }

  iDeleteXFrom := -1;

  for i := 1 to pX.Count - 1 do
    begin
      { -- compute integral value in interval }

      iValue := (
                  ( pX.GetValueCheckNull( i ) - pX.GetValueCheckNull( i - 1 ) )
                  *
                  ( ( pY.GetValueCheckNull( i ) + pY.GetValueCheckNull( i - 1 ) ) / 2 )
                );

       result := result + iValue;

       { -- check -> is greater than parameter ? }

       if _iCheckSumIntegralValue <> 999 then
       if result > _iCheckSumIntegralValue then
         begin
           { -- subtract }
           result := result - iValue;

           { -- index for removing }
           iDeleteXFrom := i;

           break;
         end;

       { -- record values }

       if bLogIntegral then pIntegral.AddValue( iValue );
       if bLogCumulation then pCumulation.AddValue( result );
     end;

  { -- delete items over test value (if exists) }

  if iDeleteXFrom <> -1 then
    for i := pX.Count - 1 downto iDeleteXFrom do
      begin
        pX.Delete( i );
        pY.Delete( i );
      end;

end;
Output:














And check (same result):

No comments:

Post a Comment