using System;
using System.Collections.Generic;
using System.Linq;


namespace Xirr
{
    public class Program
    {
        private const Double DaysPerYear = 365.0;
        private const int MaxIterations = 100;
        private const double DefaultTolerance = 1E-6;
        private const double DefaultGuess = 0.1;

    private static readonly Func, Double> NewthonsMethod =
        cf => NewtonsMethodImplementation(cf, Xnpv, XnpvPrime);

    private static readonly Func, Double> BisectionMethod =
        cf => BisectionMethodImplementation(cf, Xnpv);

    public static void Main(string[] args)
    {
        RunScenario(new[]
            {
                // this scenario fails with Newton's but succeeds with slower Bisection
                new CashItem(new DateTime(2012, 6, 1), 0.01),
                new CashItem(new DateTime(2012, 7, 23), 3042626.18),
                new CashItem(new DateTime(2012, 11, 7), -491356.62),
                new CashItem(new DateTime(2012, 11, 30), 631579.92),
                new CashItem(new DateTime(2012, 12, 1), 19769.5),
                new CashItem(new DateTime(2013, 1, 16), 1551771.47),
                new CashItem(new DateTime(2013, 2, 8), -304595),
                new CashItem(new DateTime(2013, 3, 26), 3880609.64),
                new CashItem(new DateTime(2013, 3, 31), -4331949.61)
            });
        RunScenario(new[]
            {
                new CashItem(new DateTime(2001, 5, 1), 10000),
                new CashItem(new DateTime(2002, 3, 1), 2000),
                new CashItem(new DateTime(2002, 5, 1), -5500),
                new CashItem(new DateTime(2002, 9, 1), 3000),
                new CashItem(new DateTime(2003, 2, 1), 3500),
                new CashItem(new DateTime(2003, 5, 1), -15000)
            });
    }

    private static void RunScenario(IEnumerable cashFlow)
    {
        try
        {
            try
            {
                var result = CalcXirr(cashFlow, NewthonsMethod);
                Console.WriteLine("XIRR [Newton's] value is {0}", result);
            }
            catch (InvalidOperationException)
            {
                // Failed: try another algorithm
                var result = CalcXirr(cashFlow, BisectionMethod);
                Console.WriteLine("XIRR [Bisection] (Newton's failed) value is {0}", result);
            }
        }
        catch (ArgumentException e)
        {
            Console.WriteLine(e.Message);
        }
        catch (InvalidOperationException exception)
        {
            Console.WriteLine(exception.Message);
        }
    }

    private static double CalcXirr(IEnumerable cashFlow, Func, double> method)
    {
        if (cashFlow.Count(cf => cf.Amount > 0) == 0)
            throw new ArgumentException("Add at least one positive item");

        if (cashFlow.Count(c => c.Amount < 0) == 0)
            throw new ArgumentException("Add at least one negative item");

        var result = method(cashFlow);

        if (Double.IsInfinity(result))
            throw new InvalidOperationException("Could not calculate: Infinity");

        if (Double.IsNaN(result))
            throw new InvalidOperationException("Could not calculate: Not a number");

        return result;
    }

    private static Double NewtonsMethodImplementation(IEnumerable cashFlow,
                                                      Func, Double, Double> f,
                                                      Func, Double, Double> df,
                                                      Double guess = DefaultGuess,
                                                      Double tolerance = DefaultTolerance,
                                                      int maxIterations = MaxIterations)
    {
        var x0 = guess;
        var i = 0;
        Double error;
        do
        {
            var dfx0 = df(cashFlow, x0);
            if (Math.Abs(dfx0 - 0) < Double.Epsilon)
                throw new InvalidOperationException("Could not calculate: No solution found. df(x) = 0");

            var fx0 = f(cashFlow, x0);
            var x1 = x0 - fx0/dfx0;
            error = Math.Abs(x1 - x0);

            x0 = x1;
        } while (error > tolerance && ++i < maxIterations);
        if (i == maxIterations)
            throw new InvalidOperationException("Could not calculate: No solution found. Max iterations reached.");

        return x0;
    }

    internal static Double BisectionMethodImplementation(IEnumerable cashFlow,
                                                         Func, Double, Double> f,
                                                         Double tolerance = DefaultTolerance,
                                                         int maxIterations = MaxIterations)
    {
        // From "Applied Numerical Analysis" by Gerald
        var brackets = Brackets.Find(Xnpv, cashFlow);
        if (Math.Abs(brackets.First - brackets.Second) < Double.Epsilon)
            throw new ArgumentException("Could not calculate: bracket failed");

        Double f3;
        Double result;
        var x1 = brackets.First;
        var x2 = brackets.Second;

        var i = 0;
        do
        {
            var f1 = f(cashFlow, x1);
            var f2 = f(cashFlow, x2);

            if (Math.Abs(f1) < Double.Epsilon && Math.Abs(f2) < Double.Epsilon)
                throw new InvalidOperationException("Could not calculate: No solution found");

            if (f1*f2 > 0)
                throw new ArgumentException("Could not calculate: bracket failed for x1, x2");

            result = (x1 + x2)/2;
            f3 = f(cashFlow, result);

            if (f3*f1 < 0)
                x2 = result;
            else
                x1 = result;
        } while (Math.Abs(x1 - x2)/2 > tolerance && Math.Abs(f3) > Double.Epsilon && ++i < maxIterations);

        if (i == maxIterations)
            throw new InvalidOperationException("Could not calculate: No solution found");

        return result;
    }

    private static Double Xnpv(IEnumerable cashFlow, Double rate)
    {
        if (rate <= -1)
            rate = -1 + 1E-10; // Very funky ... Better check what an IRR <= -100% means

        var startDate = cashFlow.OrderBy(i => i.Date).First().Date;
        return
            (from item in cashFlow
             let days = -(item.Date - startDate).Days
             select item.Amount*Math.Pow(1 + rate, days/DaysPerYear)).Sum();
    }

    private static Double XnpvPrime(IEnumerable cashFlow, Double rate)
    {
        var startDate = cashFlow.OrderBy(i => i.Date).First().Date;
        return (from item in cashFlow
                let daysRatio = -(item.Date - startDate).Days/DaysPerYear
                select item.Amount*daysRatio*Math.Pow(1.0 + rate, daysRatio - 1)).Sum();
    }

    public struct Brackets
    {
        public readonly Double First;
        public readonly Double Second;

        public Brackets(Double first, Double second)
        {
            First = first;
            Second = second;
        }

        internal static Brackets Find(Func, Double, Double> f,
                                      IEnumerable cashFlow,
                                      Double guess = DefaultGuess,
                                      int maxIterations = MaxIterations)
        {
            const Double bracketStep = 0.5;
            var leftBracket = guess - bracketStep;
            var rightBracket = guess + bracketStep;
            var i = 0;
            while (f(cashFlow, leftBracket)*f(cashFlow, rightBracket) > 0 && i++ < maxIterations)
            {
                leftBracket -= bracketStep;
                rightBracket += bracketStep;
            }

            return i >= maxIterations
                       ? new Brackets(0, 0)
                       : new Brackets(leftBracket, rightBracket);
        }
    }

    public struct CashItem
    {
        public DateTime Date;
        public Double Amount;

        public CashItem(DateTime date, Double amount)
        {
            Date = date;
            Amount = amount;
        }
    }
}

SQL> set serverout on
SQL> declare
  2
  3    cDaysPerYear      number := 365.0;
  4    cMaxIterations    number := 100;
  5    cDefaultTolerance number := 1E-6;
  6    cDefaultGuess     number := 0.1;
  7    cDouble_Epsilon   number := 0.00000000001;
  8
  9    type cash_obj is record ( dte date, amount number );
 10    type cash_list is table of cash_obj index by pls_integer;
 11    cash_item cash_list;
 12
 13    type bracket_rec is record ( frst number, second number );
 14
 15    function Xnpv(cashFlow cash_list, prate number) return number is
 16      startdate date := cashFlow(1).dte;
 17      dys number;
 18      tot number := 0;
 19      rate number := prate;
 20    begin
 21        if (rate <= -1) then rate := -1 + 1E-10; end if;
 22
 23        for item in 1 .. cashFlow.count loop
 24             dys := -(cashFlow(item).Dte - startDate);
 25             tot := tot + cashFlow(item).Amount*Power(1 + rate, dys/cDaysPerYear);
 26        end loop;
 27        return tot;
 28    end;
 29
 30    function XnpvPrime(cashFlow cash_list, prate number) return number is
 31      startdate date := cashFlow(1).dte;
 32      daysRatio number;
 33      tot number := 0;
 34      rate number := prate;
 35    begin
 36        for item in 1 .. cashFlow.count loop
 37             daysRatio := -(cashFlow(item).Dte - startDate)/cDaysPerYear;
 38             tot := tot + cashFlow(item).Amount*daysRatio*Power(1 + rate, daysRatio - 1);
 39        end loop;
 40        return tot;
 41    end;
 42
 43    function BracketsFind(cashFlow cash_list,guess number := cDefaultGuess) return bracket_rec is
 44      x bracket_rec;
 45      bracketStep number := 0.5;
 46      leftBracket number := guess - bracketStep;
 47      rightBracket number := guess + bracketStep;
 48      i int := 0;
 49    begin
 50      while Xnpv(cashFlow, leftBracket)*Xnpv(cashFlow, rightBracket) > 0
 51      loop
 52          leftBracket := leftBracket - bracketStep;
 53          rightBracket:= rightBracket + bracketStep;
 54          i := i + 1;
 55      end loop;
 56
 57      if i >= cmaxIterations then
 58         x.frst := 0;
 59         x.second := 0;
 60         return x;
 61      else
 62         x.frst := leftBracket;
 63         x.second := rightBracket;
 64         return x;
 65      end if;
 66    end;
 67
 68    function BisectionMethodImplementation(cashFlow cash_list,tolerance number := cDefaultTolerance, maxIterations number := cMaxIterations) return number is
 69
 70       brackets bracket_rec;
 71       f3 number;
 72       result number;
 73       x1 number;
 74       x2 number;
 75       i int := 0;
 76       f1 number;
 77       f2 number;
 78    begin
 79       brackets := BracketsFind(cashFlow);
 80       if Abs(brackets.Frst - brackets.Second) < cDouble_Epsilon then
 81          raise_application_error(-20000,'Could not calculate: bracket failed');
 82       end if;
 83       x1 := brackets.Frst;
 84       x2 := brackets.Second;
 85
 86       loop
 87
 88          f1 := Xnpv(cashFlow, x1);
 89          f2 := Xnpv(cashFlow, x2);
 90
 91          if Abs(f1) < cDouble_Epsilon and Abs(f2) < cDouble_Epsilon then
 92              raise_application_error(-20000,'Could not calculate: No solution found');
 93          end if;
 94
 95          if (f1*f2 > 0) then
 96              raise_application_error(-20000,'Could not calculate: bracket failed for x1, x2');
 97          end if;
 98
 99          result := (x1 + x2)/2;
100          f3 := Xnpv(cashFlow, result);
101
102          if (f3*f1 < 0) then
103              x2 := result;
104          else
105              x1 := result;
106          end if;
107          i := i + 1;
108          exit when not
109            (Abs(x1 - x2)/2 > tolerance and Abs(f3) > cDouble_Epsilon and i < cmaxIterations);
110       end loop;
111
112       if (i = cmaxIterations) then
113          raise_application_error(-20000,'Could not calculate: No solution found');
114       end if;
115
116       return result;
117    end;
118
119    function NewtonsMethodImplementation(cashFlow cash_list,tolerance number := cDefaultTolerance, maxIterations number := cMaxIterations,guess number := cDefaultGuess) return number is
120        x0 number := guess;
121        i int := 0;
122        error number;
123        fx0 number;
124        x1 number;
125        dfx0 number;
126    begin
127        loop
128            dfx0 := XnpvPrime(cashFlow, x0);
129            if (Abs(dfx0 - 0) < cDouble_Epsilon) then
130                raise_application_error(-20000,'Could not calculate: No solution found. df(x) = 0');
131            end if;
132
133            fx0 := Xnpv(cashFlow, x0);
134            x1 := x0 - fx0/dfx0;
135            error := Abs(x1 - x0);
136
137            x0 := x1;
138          i := i + 1;
139          exit when not (error > tolerance and  i < cmaxIterations );
140        end loop;
141
142        if (i = maxIterations) then
143            raise_application_error(-20000,'Could not calculate: No solution found. Max iterations reached.');
144        end if;
145        return x0;
146      end;
147
148      function CalcXirr(cashFlow cash_list) return number is
149        ok boolean := false;
150      begin
151        for i in 1 .. cashFlow.count loop
152           ok := ok or cashFlow(i).amount > 0;
153        end loop;
154        if not ok then raise_application_error(-20000,'Add at least one positive item');   end if;
155
156        ok := false;
157        for i in 1 .. cashFlow.count loop
158           ok := ok or cashFlow(i).amount < 0;
159        end loop;
160        if not ok then raise_application_error(-20000,'Add at least one negative item');   end if;
161
162  --      return BisectionMethodImplementation(cashFlow);
163        return NewtonsMethodImplementation(cashFlow);
164
165      end;
166
167  begin
168
169    cash_item(1).dte := date '2008-01-01';
170    cash_item(2).dte := date '2008-03-01';
171    cash_item(3).dte := date '2008-10-30';
172    cash_item(4).dte := date '2009-02-15';
173    cash_item(5).dte := date '2009-04-01';
174
175    cash_item(1).Amount := -10000;
176    cash_item(2).Amount := 2750;
177    cash_item(3).Amount := 4250;
178    cash_item(4).Amount := 3250;
179    cash_item(5).Amount := 2750;
180
181    dbms_output.put_line(CalcXirr(cash_item));
182  end;
183  /
.3733625335188315103083924482521883278399

PL/SQL procedure successfully completed.