jclcomplex.pas

East make Tray Icon in delphi

{**************************************************************************************************}
{                                                                                                  }
{ Project JEDI Code Library (JCL)                                                                  }
{                                                                                                  }
{ The contents of this file are subject to the Mozilla Public License Version 1.1 (the "License"); }
{ you may not use this file except in compliance with the License. You may obtain a copy of the    }
{ License at http://www.mozilla.org/MPL/                                                           }
{                                                                                                  }
{ Software distributed under the License is distributed on an "AS IS" basis, WITHOUT WARRANTY OF   }
{ ANY KIND, either express or implied. See the License for the specific language governing rights  }
{ and limitations under the License.                                                               }
{                                                                                                  }
{ The Original Code is JclComplex.pas.                                                             }
{                                                                                                  }
{ The Initial Developer of the Original Code is Alexei Koudinov. Portions created by               }
{ Alexei Koudinov are Copyright (C) of Alexei Koudinov. All Rights Reserved.                       }                         
{                                                                                                  }
{ Contributor(s):                                                                                  }
{   Marcel van Brakel                                                                              }
{   Alexei Koudinov                                                                                }
{   Robert Marquardt (marquardt)                                                                   }
{   Robert Rossmair (rrossmair)                                                                    }
{   Matthias Thoma  (mthoma)                                                                       }
{   Petr Vones (pvones)                                                                            }        
{                                                                                                  }
{**************************************************************************************************}
{                                                                                                  }
{ Class for working with complex numbers.                                                          }
{                                                                                                  }
{**************************************************************************************************}

// Last modified: $Date: 2005/03/08 16:10:07 $
// For history see end of file

unit JclComplex;

{$I jcl.inc}

interface

uses
  SysUtils,
  JclBase, JclMath, JclResources, JclStrings;

const
  TComplex_VERSION = 5.01;

type
  TComplexKind = (crRectangular, crPolar);

  TCoords = record
    X: Float;     // rectangular real
    Y: Float;     // rectangular imaginary
    R: Float;     // polar 1
    Theta: Float; // polar 2
  end;

  TRectCoord = record
    X: Float;
    Y: Float;
  end;

  TJclComplex = class(TObject)
  private   {z = x + yi}
    FCoord: TCoords;
    FFracLen: Byte;
    function MiscalcSingle(const X: Float): Float;
    procedure MiscalcComplex; // eliminates miscalculation
    procedure FillCoords(const ComplexType: TComplexKind);
    function GetRectangularString: string;
    function GetPolarString: string;
    procedure SetRectangularString(StrToParse: string);
    procedure SetPolarString(StrToParse: string);
    procedure SetFracLen(const X: Byte);
    function GetRadius: Float;
    function GetAngle: Float;
    function NormalizeAngle(Value: Float): Float;
  protected
    function Assign(const Coord: TCoords; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CoreAdd(const First, Second: TRectCoord): TRectCoord;
    function CoreDiv(const First, Second: TRectCoord): TRectCoord;
    function CoreMul(const First, Second: TRectCoord): TRectCoord;
    function CoreSub(const First, Second: TRectCoord): TRectCoord;
    function CoreLn (const LnValue: TRectCoord): TRectCoord;
    function CoreExp(const ExpValue: TRectCoord): TRectCoord;
    function CorePwr(First, Second, Polar: TRectCoord): TRectCoord;
    function CoreIntPwr(First: TRectCoord; const Polar: TRectCoord; const Pwr: Integer): TRectCoord;
    function CoreRealPwr(First: TRectCoord; const Polar: TRectCoord; const Pwr: Float): TRectCoord;
    function CoreRoot(First: TRectCoord; const Polar: TRectCoord; const K, N: Word): TRectCoord;
    function CoreCos(const Value: TRectCoord): TRectCoord;
    function CoreSin(const Value: TRectCoord): TRectCoord;
    function CoreTan(const Value: TRectCoord): TRectCoord;
    function CoreCot(const Value: TRectCoord): TRectCoord;
    function CoreSec(const Value: TRectCoord): TRectCoord;
    function CoreCsc(const Value: TRectCoord): TRectCoord;
    function CoreCosH(const Value: TRectCoord): TRectCoord;
    function CoreSinH(const Value: TRectCoord): TRectCoord;
    function CoreTanH(const Value: TRectCoord): TRectCoord;
    function CoreCotH(const Value: TRectCoord): TRectCoord;
    function CoreSecH(const Value: TRectCoord): TRectCoord;
    function CoreCscH(const Value: TRectCoord): TRectCoord;
    function CoreI0(const Value: TRectCoord): TRectCoord;
    function CoreJ0(const Value: TRectCoord): TRectCoord;
    function CoreApproxLnGamma(const Value: TRectCoord): TRectCoord;
    function CoreLnGamma(Value: TRectCoord): TRectCoord;
    function CoreGamma(const Value: TRectCoord): TRectCoord;
  public
    //----------- constructors
    constructor Create; overload;
    constructor Create(const X, Y: Float; const ComplexType: TComplexKind = crRectangular); overload;

    //----------- complex numbers assignment routines
    function Assign(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function AssignZero: TJclComplex;
    function AssignOne: TJclComplex;
    function Duplicate: TJclComplex;

    //----------- arithmetics -- modify the object itself
    function CAdd(const AddValue: TJclComplex): TJclComplex; overload;
    function CAdd(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CDiv(const DivValue: TJclComplex): TJclComplex; overload;
    function CDiv(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CMul(const MulValue: TJclComplex): TJclComplex; overload;
    function CMul(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CSub(const SubValue: TJclComplex): TJclComplex; overload;
    function CSub(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CNeg: TJclComplex;
    function CConjugate: TJclComplex;

    //----------- arithmetics -- creates new resulting object
    function CNewAdd(const AddValue: TJclComplex): TJclComplex; overload;
    function CNewAdd(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CNewDiv(const DivValue: TJclComplex): TJclComplex; overload;
    function CNewDiv(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CNewMul(const MulValue: TJclComplex): TJclComplex; overload;
    function CNewMul(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CNewSub(const SubValue: TJclComplex): TJclComplex; overload;
    function CNewSub(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CNewNeg: TJclComplex;
    function CNewConjugate: TJclComplex;

    //----------- natural log and exponential functions
    function CLn: TJclComplex;
    function CNewLn: TJclComplex;
    function CExp: TJclComplex;
    function CNewExp: TJclComplex;
    function CPwr(const PwrValue: TJclComplex): TJclComplex; overload;
    function CPwr(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CNewPwr(PwrValue: TJclComplex): TJclComplex; overload;
    function CNewPwr(const X, Y: Float; const ComplexType: TComplexKind = crRectangular): TJclComplex; overload;
    function CIntPwr(const Pwr: Integer): TJclComplex; overload;
    function CNewIntPwr(const Pwr: Integer): TJclComplex; overload;
    function CRealPwr(const Pwr: Float): TJclComplex; overload;
    function CNewRealPwr(const Pwr: Float): TJclComplex; overload;
    function CRoot(const K, N: Word): TJclComplex; overload;
    function CNewRoot(const K, N: Word): TJclComplex; overload;
    function CSqrt: TJclComplex; overload;
    function CNewSqrt: TJclComplex; overload;

    //----------- trigonometric functions
    function CCos: TJclComplex;
    function CNewCos: TJclComplex;
    function CSin: TJclComplex;
    function CNewSin: TJclComplex;
    function CTan: TJclComplex;
    function CNewTan: TJclComplex;
    function CCot: TJclComplex;
    function CNewCot: TJclComplex;
    function CSec: TJclComplex;
    function CNewSec: TJclComplex;
    function CCsc: TJclComplex;
    function CNewCsc: TJclComplex;

    //----------- complex hyperbolic functions
    function CCosH: TJclComplex;
    function CNewCosH: TJclComplex;
    function CSinH: TJclComplex;
    function CNewSinH: TJclComplex;
    function CTanH: TJclComplex;
    function CNewTanH: TJclComplex;
    function CCotH: TJclComplex;
    function CNewCotH: TJclComplex;
    function CSecH: TJclComplex;
    function CNewSecH: TJclComplex;
    function CCscH: TJclComplex;
    function CNewCscH: TJclComplex;

    //----------- complex Bessel functions of order zero
    function CI0: TJclComplex;
    function CNewI0: TJclComplex;
    function CJ0: TJclComplex;
    function CNewJ0: TJclComplex;

    function CApproxLnGamma: TJclComplex;
    function CNewApproxLnGamma: TJclComplex;
    function CLnGamma: TJclComplex;
    function CNewLnGamma: TJclComplex;
    function CGamma: TJclComplex;
    function CNewGamma: TJclComplex;

    //----------- miscellaneous routines
    function AbsoluteValue: Float; overload;
    function AbsoluteValue(const Coord: TRectCoord): Float; overload;
    function AbsoluteValueSqr: Float; overload;
    function AbsoluteValueSqr(const Coord: TRectCoord): Float; overload;
    function FormatExtended(const X: Float): string;

    property FracLength: Byte read FFracLen write SetFracLen default 8;

    //----------- getting different parts of the number
    property RealPart: Float read FCoord.X;
    property ImaginaryPart: Float read FCoord.Y;
    property Radius: Float read GetRadius;
    property Angle: Float read GetAngle;

    //----------- format output
    property AsString: string read GetRectangularString write SetRectangularString;
    property AsPolarString: string read GetPolarString write SetPolarString;
  end;

var
  ComplexPrecision: Float = 1E-14;

const
  MaxTerm: Byte = 35;
  EpsilonSqr: Float = 1E-20;

implementation

const
  MaxFracLen = 18;
  RectOne: TRectCoord = (X: 1.0; Y: 0.0);
  RectZero: TRectCoord = (X: 0.0; Y: 0.0);
  RectInfinity: TRectCoord = (X: Infinity; Y: Infinity);

function Coordinates(const cX, cY: Float; CoordType: TComplexKind): TCoords;
begin
  case CoordType of
    crRectangular:
      begin
        Result.X := cX;
        Result.Y := cY;
        Result.R := 0.0;
        Result.Theta := 0.0;
      end;
    crPolar:
      begin
        Result.X := 0.0;
        Result.Y := 0.0;
        Result.R := cX;
        Result.Theta := cY;
      end;
  end;
end;

function RectCoord(X, Y: Float): TRectCoord; overload;
begin
  Result.X := X;
  Result.Y := Y;
end;

function RectCoord(Value: TJclComplex): TRectCoord; overload;
begin
  Result.X := Value.FCoord.X;
  Result.Y := Value.FCoord.Y;
end;

//=== { TJclComplex } ========================================================

constructor TJclComplex.Create;
begin
  inherited Create;
  AssignZero;
  FFracLen := MaxFracLen;
end;

constructor TJclComplex.Create(const X, Y: Float; const ComplexType: TComplexKind);
begin
  inherited Create;
  Assign(X, Y, ComplexType);
  FFracLen := MaxFracLen;
end;

procedure TJclComplex.FillCoords(const ComplexType: TComplexKind);
begin
  MiscalcComplex;
  case ComplexType of
    crPolar:
      begin
        FCoord.X := FCoord.R * Cos(FCoord.Theta);
        FCoord.Y := FCoord.R * Sin(FCoord.Theta);
      end;
    crRectangular:
      begin
        if FCoord.X = 0.0 then
        begin
          FCoord.R := Abs(FCoord.Y);
          FCoord.Theta := PiOn2 * Sgn(FCoord.Y);
        end
        else
        begin
          FCoord.R := AbsoluteValue;
          FCoord.Theta := System.ArcTan(FCoord.Y / FCoord.X);
          if FCoord.X < 0.0 then
            FCoord.Theta := FCoord.Theta + Pi * Sgn(FCoord.Y);
        end;
      end;
  end;
  MiscalcComplex;
end;

function TJclComplex.MiscalcSingle(const X: Float): Float;
begin
  Result := X;
  if Abs(Result) < ComplexPrecision then
    Result := 0.0;
end;

procedure TJclComplex.MiscalcComplex; // eliminates miscalculation
begin
  FCoord.X := MiscalcSingle(FCoord.X);
  FCoord.Y := MiscalcSingle(FCoord.Y);
  FCoord.R := MiscalcSingle(FCoord.R);
  if FCoord.R = 0.0 then
    FCoord.Theta := 0.0
  else
    FCoord.Theta := MiscalcSingle(FCoord.Theta);
end;

function TJclComplex.Assign(const X, Y: Float; const ComplexType: TComplexKind): TJclComplex;
begin
  Result := Assign(Coordinates(X, Y, ComplexType), ComplexType);
end;

function TJclComplex.Assign(const Coord: TCoords; const ComplexType: TComplexKind): TJclComplex;
begin
  FCoord := Coord;
  FillCoords(ComplexType);
  MiscalcComplex;
  Result := Self;
end;

function TJclComplex.AssignZero: TJclComplex;
begin
  Result := Assign(0.0, 0.0, crRectangular);
end;

function TJclComplex.AssignOne: TJclComplex;
begin
  Result := Assign(1.0, 0.0, crRectangular);
end;

function TJclComplex.GetRectangularString: string;
begin
  MiscalcComplex;
  if (FCoord.X = 0.0) and (FCoord.Y = 0.0) then
    Result := '0'
  else
  if FCoord.X <> 0.0 then
  begin
    Result := FormatExtended(FCoord.X);
    if FCoord.Y > 0.0 then
      Result := Result + '+'
    else
    if FCoord.Y < 0.0 then
      Result := Result + '-';
    if FCoord.Y <> 0.0 then
      Result := Result + FormatExtended(Abs(FCoord.Y)) + 'i';
  end
  else
    Result := FormatExtended(FCoord.Y) + 'i';
end;

function TJclComplex.GetPolarString: string;
begin
  FillCoords(crRectangular);
  Result := FormatExtended(FCoord.R) + '*CIS(' + FormatExtended(FCoord.Theta) + ')';
end;

procedure TJclComplex.SetRectangularString(StrToParse: string);
var
  SignPos: Integer;
  RealPart, ImagPart: Float;
begin
  StrToParse := StrRemoveChars(StrToParse, [' ']);
  SignPos := StrFind('+', StrToParse, 2);
  if SignPos = 0 then
    SignPos := StrFind('-', StrToParse, 2);
  if SignPos > 0 then
  begin
    try
      RealPart := StrToFloat(Copy(StrToParse, 1, SignPos - 1));
    except
      raise EJclMathError.CreateRes(@RsComplexInvalidString);