[lmi] An experimental change for withdrawal solves

[Greg Chicares]

Consider this testcase:
  File | New | Illustration
  solve for withdrawal
    from retirement
    to maturity
For all other parameters, just accept the default values, which are
fairly vanilla: male nonsmoker age 45, $1M specified amount, annual
premium $20K, and a generic default product.

This solve succeeds. Withdrawals reduce the specified amount to the
minimum, $50000, by the maturity duration, when the cash surrender
value equals $50001: it endows for the specified amount. This seems to
be the best answer possible. If you turn off the solve and paste
  0,retirement;56163.67
into the "Withdrawal" field, the ultimate CSV is 49999. But with a
withdrawal of
  0,retirement;56163.66
the ultimate CSV is 50001. A one-cent difference for thirty-five years
leads to a two-dollar difference because of interest and mortality.
The choice between these two withdrawal amounts is not arbitrary: lmi
chooses the one that fully achieves the target, because a CSV of 49999
fails to endow [0].

Now change the solve target from endowment to $0 CSV at maturity. A
messagebox warns:
  Solution not found: using zero instead.
However, a withdrawal of
  0,retirement;56387.68
produces a $1 CSV at maturity. Why wasn't that solution found? [1]

Using cvs HEAD as of 20090130T1135Z, write 'idiosyncrasyT' in the
"Comments" field. That's an expeditious trick we sometimes use to add
a hidden feature; this one writes per-iteration solve details to a
'trace.txt' file in the current directory. Rerun the illustration and
examine that file. It says only this:

  iteration 0 iterand         0    value 6305652.52
  iteration 1 iterand 999999999.99 value   21257.15

The solve algorithm is guaranteed to succeed as long as we know one
input value that's too large and one that's too small:
  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/zero.hpp?annotate=1.16
| Precondition: the input bounds include or bracket a root.
Here, it tries the smallest possible withdrawal (zero), then another
that's implausibly large, but both yield a positive CSV at the target
duration. That's because the generic default policy form, like many
real-world forms, limits withdrawals to CSV less a provision for the
current year's deductions. As long as continuing premiums cover prior
deductions, withdrawals won't lapse this contract.

Then how did the first testcase succeed? Examining its per-iteration
details [2], we see that it satisfies the initial precondition above:

  iteration 0 iterand         0    value 5305652.52
  iteration 1 iterand 999999999.99 value  -28742.85

where "value" is the 'objective function'--what we're trying to reduce
to zero--in this case the difference between actual and target CSV in
the target year. To elaborate:

                      actual  specified     target     objective
    withdrawal           CSV     amount        CSV      function
          0.00    6305652.52    1000000 == 1000000    5305652.52
  999999999.99      21257.15      50000 ==   50000     -28742.85

This is a solve for "endowment", which means that CSV should at least
equal specified amount at the target duration. Of course, withdrawals
affect specified amount, so the target moves:
  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/zero.hpp?annotate=1.16
| coefficients are likely to change with each iteration
What's most important, though, is that the objective function can be
negative, simply because the target happens to be big enough (whereas
in the solve-for-$0 case it isn't). And the target is big enough only
because this generic product has a $50000 minimum specified amount.

The algorithm has begun well enough, by satisfying its precondition,
but it doesn't progress as rapidly as we'd like:

  iteration 0 iterand         0    value 5305652.52
  iteration 1 iterand 999999999.99 value  -28742.85
  iteration 2 iterand 994611788.58 value  -28742.85
  iteration 3 iterand 497305894.29 value  -28742.85
  iteration 4 iterand 248652947.14 value  -28742.85
  iteration 5 iterand 124326473.57 value  -28742.85
  ...
  iteration 14 iterand   242825.14 value  -28741.69
  iteration 15 iterand   121412.57 value  -28741.69
  iteration 16 iterand    60706.28 value  -28742.02

Let's make sure we understand this before we try to fix it [3]. Given
the outcomes of the first two trials, what value should be tried next?
An insurance expert attacking the problem manually might reason thus:

  The billion-dollar annual withdrawal must have been limited, so that
  trial means nothing. With no withdrawal, there's 5305652.52 of extra
  cash value that we want to get out over 35 years, so let's guess:
    151590.07 = 5305652.52 / 35

Without such expert insight, the algorithm examines the first two
trials and "guesses" that the solution is close to a billion, to put
it anthropomorphically. Really, it just interpolates linearly:

  regula falsi: (a * f(b) - b * f(a)) / (f(b) - f(a))
  (0 * -28742.85 - 999999999.99 * 5305652.52) / (-28742.85 - 5305652.52)
  994611788.58

  secant method: b - f(b) * (b - a) / (f(b) - f(a))
  999999999.99 - -28742.85 * (999999999.99 - 0) / (-28742.85 - 5305652.52)
  994611788.58

What's worrisome is that the next fourteen iterations use bisection.
The withdrawal marches down by a factor of 2^14, while the objective
function's value changes by less than a hundredth of a percent.

We could say speed doesn't matter much. The algorithm's author says
  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/zero.hpp?annotate=1.16
| it's "never much slower than bisection"
and bisection is guaranteed to work. In fact, pure bisection would
always take about three dozen steps to find a root to the nearest cent
in the interval [999999999.99, 0.00] (see why? hint: calculate 2^36),
as long as the initial iterands bound a root. But they don't for the
solve-for-$0 testcase.

We could say that's a silly testcase: in the last year, $20000 is
paid, and a somewhat larger amount is withdrawn. Withdrawals usually
bear a fee such as two percent, so it's better not to pay the $20000
and incur a $400 fee to get it back the same day. On the other hand,
it's a simple testcase that might be attempted by a naive user, who
would find the failure discouraging.

Instead, let's look at the cause. We're seeing one-sided convergence
in the endowment solve, and failure in the other testcase, because of
the withdrawal limit. If the billion-dollar withdrawal weren't limited
at all, the objective function would assume a large negative value.

                         gGhH
                        gGhH        -----objective functions-----
                       gGhH
                     gGhH           h: original    solve for $0
    h h h h h h h  gGhH             H: no-WD-limit solve for $0
    -------------------------
    g g g g g g gGH                 g: original    solve to endow
            GH                      G: no-WD-limit solve to endow
        GH
    GH

Probe the behavior of the original objective functions (g and h) for
withdrawals a bit bigger than the solutions given, and you'll see that
the ordinate fluctuates in a way that the ASCII chart above cannot
depict. Those fluctuations can easily trap other algorithms: imagine
applying Newton's method to approach the root of
  f(x) =      2 + sin(x), x < 0
       = -1 + x * 1.05^x, 0 <= x
from the left. Resist the impulse to twiddle the algorithm; see the
"Rationale for choice of algorithm" here:
  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/zero.hpp?annotate=1.16
Instead, recognize that the original objective functions are difficult
to work with, and seek an easier problem whose solution we can use.

We already modify the objective function anyway [4]:

  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/ihs_avsolve.cpp?annotate=1.40
| Therefore, we take certain steps to make the
| objective function more well-behaved as its value approaches zero
| from either direction:
|
| 1. we prevent the policy from lapsing during a solve [...]

so let's try this change:

  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/ihs_avmly.cpp?r1=1.112&r2=1.113

-        if(MaxWD < RequestedWD)
+        if(!Solving && MaxWD < RequestedWD)
             {
             NetWD = MaxWD;
             }

Rerunning the testcases above, we see that the endowment solve finds
the same answer with one-third less work, and the solve for $0 finds
an answer that looks optimal where none was found before. This code
seems ready for serious testing.

---------

[0] "choice...not arbitrary...a CSV of 49999 fails to endow"

It's not as simple as picking the higher CSV in every case. A solve
for tax basis chooses the highest CSV that doesn't exceed the basis,
which is the lower of the final bracketing CSV amounts. Otherwise, a
full surrender wouldn't be tax free.

[1] "Why wasn't that solution found?"

It's quite instructive to try solving for $0 CSV by hand. Apparently
the (original) objective function has no root, even if we let CSV go
negative. What the 20090202T0247Z solve code actually finds is a
minimum that happens to be extremely close to zero.

[2] "per-iteration details"

Due to its size, a complete trace is presented as an appendix below.

[3] "make sure we understand this before we try to fix it"

Moments after I wrote that, google broke. They tried to fix it before
understanding it:

http://news.bbc.co.uk/2/hi/technology/7862840.stm
| "There was a fault. We don't know the nature of it yet. Everything
| has been solved. We are still making initial enquiries," a Google
| spokesperson told BBC News.

and guess what ensued:

http://blogs.zdnet.com/hardware/?p=3420
| it seems that Google has removed the filtering of bad sites. Sites
| that were on StopBadware.org block list are now not being flagged as
| potentially dangerous.

http://digital-lifestyles.info/2009/01/31/googles-gmail-now-freaking-out-mis-marking-mail-as-spam/
| As soon as Google fixed their misreporting "This Site May Harm Your
| Computer" search problem, we had reports of another problem with a
| Google service.

That's why we need to test this work carefully. Solves are subtle.

[4] "We already modify the objective function anyway"

Some other illustration systems use a separate "targeting engine" for
solves. I know one that doesn't produce the same illustration if a
solve result is plugged back in as regular (non-solve) input. We avoid
that problem by making a separate final pass...

  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/ihs_avsolve.cpp?annotate=1.40
| The account and ledger values set as a side effect of solving aren't
| necessarily what we need [...] Therefore, the final solve parameters
| are stored now, and values are regenerated downstream.

...and rounding each iterand carefully:

  http://cvs.savannah.gnu.org/viewvc/lmi/lmi/zero.hpp?annotate=1.16
| The original algorithm returns unrounded values. For life insurance
| illustrations, rounded values are more often wanted, so iterands
| are rounded before each function evaluation. Rounding downstream,
| outside the root-finding algorithm, would not be appropriate: for a
| hypothetical unrounded return value r, f(r) and f(round(r)) might
| easily have different signs.

---------

Appendix: per-iteration details.

20090130T1135Z code, solve for $0:

iteration 0 iterand 0 value 6305652.52
iteration 1 iterand 999999999.99 value 21257.15

20090130T1135Z code, solve for endowment:

iteration 0 iterand 0 value 5305652.52
iteration 1 iterand 999999999.99 value -28742.85
iteration 2 iterand 994611788.58 value -28742.85
iteration 3 iterand 497305894.29 value -28742.85
iteration 4 iterand 248652947.14 value -28742.85
iteration 5 iterand 124326473.57 value -28742.85
iteration 6 iterand 62163236.78 value -28742.85
iteration 7 iterand 31081618.39 value -28742.85
iteration 8 iterand 15540809.2 value -28742.85
iteration 9 iterand 7770404.6 value -28742.85
iteration 10 iterand 3885202.3 value -28742.85
iteration 11 iterand 1942601.15 value -28742.85
iteration 12 iterand 971300.58 value -28742.85
iteration 13 iterand 485650.29 value -28741.69
iteration 14 iterand 242825.14 value -28741.69
iteration 15 iterand 121412.57 value -28741.69
iteration 16 iterand 60706.28 value -28742.02
iteration 17 iterand 30353.14 value 2871686.830000002
iteration 18 iterand 60405.49 value -28742.32
iteration 19 iterand 45379.32 value 1202327.590000001
iteration 20 iterand 60054.67 value -28741.15000000001
iteration 21 iterand 52717 value 385131.2900000006
iteration 22 iterand 59545.11 value -28741.4
iteration 23 iterand 56131.06 value 4451.869999999319
iteration 24 iterand 56588.95 value -45328.72
iteration 25 iterand 56172.01 value -1859.590000000018
iteration 26 iterand 56159.94 value 756.530000000057
iteration 27 iterand 56163.43 value 53.8900000005306
iteration 28 iterand 56163.69 value -5.54999999990105
iteration 29 iterand 56163.67 value -1.209999999780848
iteration 30 iterand 56163.66 value 1.03999999955704

20090202T0247Z code, solve for $0:

iteration 0 iterand 0 value 6305652.52
iteration 1 iterand 999999999.99 value -999369555.86
iteration 2 iterand 6270068.57 value -5639624.43999999
iteration 3 iterand 3309833.15 value -2679389.02
iteration 4 iterand 1429452.74 value -799008.61
iteration 5 iterand 820773.4 value -799861.359999999
iteration 6 iterand 410386.7 value -389027.2399999994
iteration 7 iterand 205193.35 value -185474.4799999999
iteration 8 iterand 102596.68 value -84807.22
iteration 9 iterand 51298.34 value 593218.200000002
iteration 10 iterand 96180.3 value -116609.6199999995
iteration 11 iterand 73739.32 value -57161.5899999996
iteration 12 iterand 62518.83 value -50880.5699999994
iteration 13 iterand 56908.58 value -35223.83000000104
iteration 14 iterand 54103.46 value 280552.2500000005
iteration 15 iterand 56595.68 value -39882.08999999994
iteration 16 iterand 55349.57 value 141620.4400000004
iteration 17 iterand 56321.87 value 14687.20999999867
iteration 18 iterand 56395.57 value -1751.289999999282
iteration 19 iterand 56387.72 value -8.97999999922604
iteration 20 iterand 56387.68 value 0.64000000057149
iteration 21 iterand 56387.69 value -1.879999999618121

20090202T0247Z code, solve for endowment:

iteration 0 iterand 0 value 5305652.52
iteration 1 iterand 999999999.99 value -999419555.86
iteration 2 iterand 5280700.11 value -4700255.97999999
iteration 3 iterand 2800101.54 value -2219657.41
iteration 4 iterand 1235151.92 value -654707.79
iteration 5 iterand 733537.9 value -761848.409999999
iteration 6 iterand 366768.95 value -392034.3899999998
iteration 7 iterand 183384.48 value -226164.0299999995
iteration 8 iterand 91692.24 value -121681.8299999994
iteration 9 iterand 45846.12 value 1150368.730000003
iteration 10 iterand 87306.69 value -113708.5399999995
iteration 11 iterand 66576.4 value -102573.9099999986
iteration 12 iterand 56211.26 value -10621.00999999985
iteration 13 iterand 55104.29 value 118969.3500000002
iteration 14 iterand 56120.53 value 5625.82999999989
iteration 15 iterand 56151.95 value 2030.790000000205
iteration 16 iterand 56166.85 value -711.759999999638
iteration 17 iterand 56162.98 value 150.2400000017806
iteration 18 iterand 56163.65 value 3.429999999971187
iteration 19 iterand 56163.67 value -1.209999999780848
iteration 20 iterand 56163.66 value 1.03999999955704