1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
|
<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
<title>Continued Fraction Evaluation</title>
<link rel="stylesheet" href="../../math.css" type="text/css">
<meta name="generator" content="DocBook XSL Stylesheets V1.79.1">
<link rel="home" href="../../index.html" title="Math Toolkit 4.2.1">
<link rel="up" href="../internals.html" title="Internal tools">
<link rel="prev" href="agm.html" title="Arithmetic-Geometric Mean">
<link rel="next" href="simple_continued_fraction.html" title="Simple Continued Fractions">
<meta name="viewport" content="width=device-width, initial-scale=1">
</head>
<body bgcolor="white" text="black" link="#0000FF" vlink="#840084" alink="#0000FF">
<table cellpadding="2" width="100%"><tr>
<td valign="top"><img alt="Boost C++ Libraries" width="277" height="86" src="../../../../../../boost.png"></td>
<td align="center"><a href="../../../../../../index.html">Home</a></td>
<td align="center"><a href="../../../../../../libs/libraries.htm">Libraries</a></td>
<td align="center"><a href="http://www.boost.org/users/people.html">People</a></td>
<td align="center"><a href="http://www.boost.org/users/faq.html">FAQ</a></td>
<td align="center"><a href="../../../../../../more/index.htm">More</a></td>
</tr></table>
<hr>
<div class="spirit-nav">
<a accesskey="p" href="agm.html"><img src="../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../internals.html"><img src="../../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../../index.html"><img src="../../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="simple_continued_fraction.html"><img src="../../../../../../doc/src/images/next.png" alt="Next"></a>
</div>
<div class="section">
<div class="titlepage"><div><div><h3 class="title">
<a name="math_toolkit.internals.cf"></a><a class="link" href="cf.html" title="Continued Fraction Evaluation">Continued Fraction Evaluation</a>
</h3></div></div></div>
<h5>
<a name="math_toolkit.internals.cf.h0"></a>
<span class="phrase"><a name="math_toolkit.internals.cf.synopsis"></a></span><a class="link" href="cf.html#math_toolkit.internals.cf.synopsis">Synopsis</a>
</h5>
<pre class="programlisting"><span class="preprocessor">#include</span> <span class="special"><</span><span class="identifier">boost</span><span class="special">/</span><span class="identifier">math</span><span class="special">/</span><span class="identifier">tools</span><span class="special">/</span><span class="identifier">fraction</span><span class="special">.</span><span class="identifier">hpp</span><span class="special">></span>
</pre>
<pre class="programlisting"><span class="keyword">namespace</span> <span class="identifier">boost</span><span class="special">{</span> <span class="keyword">namespace</span> <span class="identifier">math</span><span class="special">{</span> <span class="keyword">namespace</span> <span class="identifier">tools</span><span class="special">{</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">,</span> <span class="keyword">class</span> <span class="identifier">U</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_b</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">const</span> <span class="identifier">U</span><span class="special">&</span> <span class="identifier">tolerance</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">uintmax_t</span><span class="special">&</span> <span class="identifier">max_terms</span><span class="special">)</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">,</span> <span class="keyword">class</span> <span class="identifier">U</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_b</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">const</span> <span class="identifier">U</span><span class="special">&</span> <span class="identifier">tolerance</span><span class="special">)</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">,</span> <span class="keyword">class</span> <span class="identifier">U</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_a</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">const</span> <span class="identifier">U</span><span class="special">&</span> <span class="identifier">tolerance</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">uintmax_t</span><span class="special">&</span> <span class="identifier">max_terms</span><span class="special">)</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">,</span> <span class="keyword">class</span> <span class="identifier">U</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_a</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">const</span> <span class="identifier">U</span><span class="special">&</span> <span class="identifier">tolerance</span><span class="special">)</span>
<span class="comment">//</span>
<span class="comment">// These interfaces are present for legacy reasons, and are now deprecated:</span>
<span class="comment">//</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_b</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">int</span> <span class="identifier">bits</span><span class="special">);</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_b</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">int</span> <span class="identifier">bits</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">uintmax_t</span><span class="special">&</span> <span class="identifier">max_terms</span><span class="special">);</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_a</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">int</span> <span class="identifier">bits</span><span class="special">);</span>
<span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">Gen</span><span class="special">></span>
<span class="keyword">typename</span> <span class="identifier">detail</span><span class="special">::</span><span class="identifier">fraction_traits</span><span class="special"><</span><span class="identifier">Gen</span><span class="special">>::</span><span class="identifier">result_type</span>
<span class="identifier">continued_fraction_a</span><span class="special">(</span><span class="identifier">Gen</span><span class="special">&</span> <span class="identifier">g</span><span class="special">,</span> <span class="keyword">int</span> <span class="identifier">bits</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">uintmax_t</span><span class="special">&</span> <span class="identifier">max_terms</span><span class="special">);</span>
<span class="special">}}}</span> <span class="comment">// namespaces</span>
</pre>
<h5>
<a name="math_toolkit.internals.cf.h1"></a>
<span class="phrase"><a name="math_toolkit.internals.cf.description"></a></span><a class="link" href="cf.html#math_toolkit.internals.cf.description">Description</a>
</h5>
<p>
<a href="http://en.wikipedia.org/wiki/Continued_fraction" target="_top">Continued fractions
are a common method of approximation. </a> These functions all evaluate
the continued fraction described by the <span class="emphasis"><em>generator</em></span> type
argument. The functions with an "_a" suffix evaluate the fraction:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/fraction2.svg"></span>
</p></blockquote></div>
<p>
and those with a "_b" suffix evaluate the fraction:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/fraction1.svg"></span>
</p></blockquote></div>
<p>
This latter form is somewhat more natural in that it corresponds with the
usual definition of a continued fraction, but note that the first <span class="emphasis"><em>a</em></span>
value returned by the generator is discarded. Further, often the first <span class="emphasis"><em>a</em></span>
and <span class="emphasis"><em>b</em></span> values in a continued fraction have different
defining equations to the remaining terms, which may make the "_a"
suffixed form more appropriate.
</p>
<p>
The generator type should be a function object which supports the following
operations:
</p>
<div class="informaltable"><table class="table">
<colgroup>
<col>
<col>
</colgroup>
<thead><tr>
<th>
<p>
Expression
</p>
</th>
<th>
<p>
Description
</p>
</th>
</tr></thead>
<tbody>
<tr>
<td>
<p>
Gen::result_type
</p>
</td>
<td>
<p>
The type that is the result of invoking operator(). This can be
either an arithmetic or complex type, or a std::pair<> of
arithmetic or complex types.
</p>
</td>
</tr>
<tr>
<td>
<p>
g()
</p>
</td>
<td>
<p>
Returns an object of type Gen::result_type.
</p>
<p>
Each time this operator is called then the next pair of <span class="emphasis"><em>a</em></span>
and <span class="emphasis"><em>b</em></span> values is returned. Or, if result_type
is an arithmetic type, then the next <span class="emphasis"><em>b</em></span> value
is returned and all the <span class="emphasis"><em>a</em></span> values are assumed
to 1.
</p>
</td>
</tr>
</tbody>
</table></div>
<p>
In all the continued fraction evaluation functions the <span class="emphasis"><em>tolerance</em></span>
parameter is the precision desired in the result, evaluation of the fraction
will continue until the last term evaluated leaves the relative error in
the result less than <span class="emphasis"><em>tolerance</em></span>. The deprecated interfaces
take a number of digits precision here, internally they just convert this
to a tolerance and forward call.
</p>
<p>
If the optional <span class="emphasis"><em>max_terms</em></span> parameter is specified then
no more than <span class="emphasis"><em>max_terms</em></span> calls to the generator will be
made, and on output, <span class="emphasis"><em>max_terms</em></span> will be set to actual
number of calls made. This facility is particularly useful when profiling
a continued fraction for convergence.
</p>
<h5>
<a name="math_toolkit.internals.cf.h2"></a>
<span class="phrase"><a name="math_toolkit.internals.cf.implementation"></a></span><a class="link" href="cf.html#math_toolkit.internals.cf.implementation">Implementation</a>
</h5>
<p>
Internally these algorithms all use the modified Lentz algorithm: refer to
Numeric Recipes in C++, W. H. Press et all, chapter 5, (especially 5.2 Evaluation
of continued fractions, p 175 - 179) for more information, also Lentz, W.J.
1976, Applied Optics, vol. 15, pp. 668-671.
</p>
<h5>
<a name="math_toolkit.internals.cf.h3"></a>
<span class="phrase"><a name="math_toolkit.internals.cf.examples"></a></span><a class="link" href="cf.html#math_toolkit.internals.cf.examples">Examples</a>
</h5>
<p>
All of these examples are in <a href="../../../../example/continued_fractions.cpp" target="_top">continued_fractions.cpp</a>.
</p>
<p>
The <a href="http://en.wikipedia.org/wiki/Golden_ratio" target="_top">golden ratio phi
= 1.618033989...</a> can be computed from the simplest continued fraction
of all:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/fraction3.svg"></span>
</p></blockquote></div>
<p>
We begin by defining a generator function:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="keyword">struct</span> <span class="identifier">golden_ratio_fraction</span>
<span class="special">{</span>
<span class="keyword">typedef</span> <span class="identifier">T</span> <span class="identifier">result_type</span><span class="special">;</span>
<span class="identifier">result_type</span> <span class="keyword">operator</span><span class="special">()()</span>
<span class="special">{</span>
<span class="keyword">return</span> <span class="number">1</span><span class="special">;</span>
<span class="special">}</span>
<span class="special">};</span>
</pre>
<p>
The golden ratio can then be computed to double precision using:
</p>
<pre class="programlisting"><span class="identifier">golden_ratio_fraction</span><span class="special"><</span><span class="keyword">double</span><span class="special">></span> <span class="identifier">func</span><span class="special">;</span>
<span class="keyword">double</span> <span class="identifier">gr</span> <span class="special">=</span> <span class="identifier">continued_fraction_a</span><span class="special">(</span>
<span class="identifier">func</span><span class="special">,</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special"><</span><span class="keyword">double</span><span class="special">>::</span><span class="identifier">epsilon</span><span class="special">());</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">cout</span> <span class="special"><<</span> <span class="string">"The golden ratio is: "</span> <span class="special"><<</span> <span class="identifier">gr</span> <span class="special"><<</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">endl</span><span class="special">;</span>
</pre>
<p>
It's more usual though to have to define both the <span class="emphasis"><em>a</em></span>'s
and the <span class="emphasis"><em>b</em></span>'s when evaluating special functions by continued
fractions, for example the tan function is defined by:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/fraction4.svg"></span>
</p></blockquote></div>
<p>
So its generator object would look like:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="keyword">struct</span> <span class="identifier">tan_fraction</span>
<span class="special">{</span>
<span class="keyword">private</span><span class="special">:</span>
<span class="identifier">T</span> <span class="identifier">a</span><span class="special">,</span> <span class="identifier">b</span><span class="special">;</span>
<span class="keyword">public</span><span class="special">:</span>
<span class="identifier">tan_fraction</span><span class="special">(</span><span class="identifier">T</span> <span class="identifier">v</span><span class="special">)</span>
<span class="special">:</span> <span class="identifier">a</span><span class="special">(-</span><span class="identifier">v</span> <span class="special">*</span> <span class="identifier">v</span><span class="special">),</span> <span class="identifier">b</span><span class="special">(-</span><span class="number">1</span><span class="special">)</span>
<span class="special">{}</span>
<span class="keyword">typedef</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">pair</span><span class="special"><</span><span class="identifier">T</span><span class="special">,</span> <span class="identifier">T</span><span class="special">></span> <span class="identifier">result_type</span><span class="special">;</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">pair</span><span class="special"><</span><span class="identifier">T</span><span class="special">,</span> <span class="identifier">T</span><span class="special">></span> <span class="keyword">operator</span><span class="special">()()</span>
<span class="special">{</span>
<span class="identifier">b</span> <span class="special">+=</span> <span class="number">2</span><span class="special">;</span>
<span class="keyword">return</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">make_pair</span><span class="special">(</span><span class="identifier">a</span><span class="special">,</span> <span class="identifier">b</span><span class="special">);</span>
<span class="special">}</span>
<span class="special">};</span>
</pre>
<p>
Notice that if the continuant is subtracted from the <span class="emphasis"><em>b</em></span>
terms, as is the case here, then all the <span class="emphasis"><em>a</em></span> terms returned
by the generator will be negative. The tangent function can now be evaluated
using:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="identifier">T</span> <span class="identifier">tan</span><span class="special">(</span><span class="identifier">T</span> <span class="identifier">a</span><span class="special">)</span>
<span class="special">{</span>
<span class="identifier">tan_fraction</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="identifier">fract</span><span class="special">(</span><span class="identifier">a</span><span class="special">);</span>
<span class="keyword">return</span> <span class="identifier">a</span> <span class="special">/</span> <span class="identifier">continued_fraction_b</span><span class="special">(</span><span class="identifier">fract</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special"><</span><span class="identifier">T</span><span class="special">>::</span><span class="identifier">epsilon</span><span class="special">());</span>
<span class="special">}</span>
</pre>
<p>
Notice that this time we're using the "_b" suffixed version to
evaluate the fraction: we're removing the leading <span class="emphasis"><em>a</em></span>
term during fraction evaluation as it's different from all the others.
</p>
<p>
Now we'll look at a couple of complex number examples, starting with the
exponential integral which can be calculated via:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/expint_n_3.svg"></span>
</p></blockquote></div>
<p>
So our functor looks like this:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="keyword">struct</span> <span class="identifier">expint_fraction</span>
<span class="special">{</span>
<span class="keyword">typedef</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">pair</span><span class="special"><</span><span class="identifier">T</span><span class="special">,</span> <span class="identifier">T</span><span class="special">></span> <span class="identifier">result_type</span><span class="special">;</span>
<span class="identifier">expint_fraction</span><span class="special">(</span><span class="keyword">unsigned</span> <span class="identifier">n_</span><span class="special">,</span> <span class="identifier">T</span> <span class="identifier">z_</span><span class="special">)</span> <span class="special">:</span> <span class="identifier">b</span><span class="special">(</span><span class="identifier">z_</span> <span class="special">+</span> <span class="identifier">T</span><span class="special">(</span><span class="identifier">n_</span><span class="special">)),</span> <span class="identifier">i</span><span class="special">(-</span><span class="number">1</span><span class="special">),</span> <span class="identifier">n</span><span class="special">(</span><span class="identifier">n_</span><span class="special">)</span> <span class="special">{}</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">pair</span><span class="special"><</span><span class="identifier">T</span><span class="special">,</span> <span class="identifier">T</span><span class="special">></span> <span class="keyword">operator</span><span class="special">()()</span>
<span class="special">{</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">pair</span><span class="special"><</span><span class="identifier">T</span><span class="special">,</span> <span class="identifier">T</span><span class="special">></span> <span class="identifier">result</span> <span class="special">=</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">make_pair</span><span class="special">(-</span><span class="keyword">static_cast</span><span class="special"><</span><span class="identifier">T</span><span class="special">>((</span><span class="identifier">i</span> <span class="special">+</span> <span class="number">1</span><span class="special">)</span> <span class="special">*</span> <span class="special">(</span><span class="identifier">n</span> <span class="special">+</span> <span class="identifier">i</span><span class="special">)),</span> <span class="identifier">b</span><span class="special">);</span>
<span class="identifier">b</span> <span class="special">+=</span> <span class="number">2</span><span class="special">;</span>
<span class="special">++</span><span class="identifier">i</span><span class="special">;</span>
<span class="keyword">return</span> <span class="identifier">result</span><span class="special">;</span>
<span class="special">}</span>
<span class="keyword">private</span><span class="special">:</span>
<span class="identifier">T</span> <span class="identifier">b</span><span class="special">;</span>
<span class="keyword">int</span> <span class="identifier">i</span><span class="special">;</span>
<span class="keyword">unsigned</span> <span class="identifier">n</span><span class="special">;</span>
<span class="special">};</span>
</pre>
<p>
We can finish the example by wrapping everything up in a function:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="keyword">inline</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="identifier">expint_as_fraction</span><span class="special">(</span><span class="keyword">unsigned</span> <span class="identifier">n</span><span class="special">,</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="keyword">const</span><span class="special">&</span> <span class="identifier">z</span><span class="special">)</span>
<span class="special">{</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">uintmax_t</span> <span class="identifier">max_iter</span> <span class="special">=</span> <span class="number">1000</span><span class="special">;</span>
<span class="identifier">expint_fraction</span><span class="special"><</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="special">></span> <span class="identifier">f</span><span class="special">(</span><span class="identifier">n</span><span class="special">,</span> <span class="identifier">z</span><span class="special">);</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="identifier">result</span> <span class="special">=</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">tools</span><span class="special">::</span><span class="identifier">continued_fraction_b</span><span class="special">(</span>
<span class="identifier">f</span><span class="special">,</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">>(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special"><</span><span class="identifier">T</span><span class="special">>::</span><span class="identifier">epsilon</span><span class="special">()),</span>
<span class="identifier">max_iter</span><span class="special">);</span>
<span class="identifier">result</span> <span class="special">=</span> <span class="identifier">exp</span><span class="special">(-</span><span class="identifier">z</span><span class="special">)</span> <span class="special">/</span> <span class="identifier">result</span><span class="special">;</span>
<span class="keyword">return</span> <span class="identifier">result</span><span class="special">;</span>
<span class="special">}</span>
</pre>
<p>
Notice how the termination condition is still expressed as a complex number,
albeit one with zero imaginary part.
</p>
<p>
Our final example will use <code class="literal">continued_fraction_a</code>, in fact
there is only one special function in our code which uses that variant, and
it's the upper incomplete gamma function (Q), which can be calculated via:
</p>
<div class="blockquote"><blockquote class="blockquote"><p>
<span class="inlinemediaobject"><img src="../../../equations/igamma9.svg"></span>
</p></blockquote></div>
<p>
In this case the first couple of terms are different from the rest, so our
fraction will start with the first "regular" a term:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="keyword">struct</span> <span class="identifier">upper_incomplete_gamma_fract</span>
<span class="special">{</span>
<span class="keyword">private</span><span class="special">:</span>
<span class="keyword">typedef</span> <span class="keyword">typename</span> <span class="identifier">T</span><span class="special">::</span><span class="identifier">value_type</span> <span class="identifier">scalar_type</span><span class="special">;</span>
<span class="identifier">T</span> <span class="identifier">z</span><span class="special">,</span> <span class="identifier">a</span><span class="special">;</span>
<span class="keyword">int</span> <span class="identifier">k</span><span class="special">;</span>
<span class="keyword">public</span><span class="special">:</span>
<span class="keyword">typedef</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">pair</span><span class="special"><</span><span class="identifier">T</span><span class="special">,</span> <span class="identifier">T</span><span class="special">></span> <span class="identifier">result_type</span><span class="special">;</span>
<span class="identifier">upper_incomplete_gamma_fract</span><span class="special">(</span><span class="identifier">T</span> <span class="identifier">a1</span><span class="special">,</span> <span class="identifier">T</span> <span class="identifier">z1</span><span class="special">)</span>
<span class="special">:</span> <span class="identifier">z</span><span class="special">(</span><span class="identifier">z1</span> <span class="special">-</span> <span class="identifier">a1</span> <span class="special">+</span> <span class="identifier">scalar_type</span><span class="special">(</span><span class="number">1</span><span class="special">)),</span> <span class="identifier">a</span><span class="special">(</span><span class="identifier">a1</span><span class="special">),</span> <span class="identifier">k</span><span class="special">(</span><span class="number">0</span><span class="special">)</span>
<span class="special">{</span>
<span class="special">}</span>
<span class="identifier">result_type</span> <span class="keyword">operator</span><span class="special">()()</span>
<span class="special">{</span>
<span class="special">++</span><span class="identifier">k</span><span class="special">;</span>
<span class="identifier">z</span> <span class="special">+=</span> <span class="identifier">scalar_type</span><span class="special">(</span><span class="number">2</span><span class="special">);</span>
<span class="keyword">return</span> <span class="identifier">result_type</span><span class="special">(</span><span class="identifier">scalar_type</span><span class="special">(</span><span class="identifier">k</span><span class="special">)</span> <span class="special">*</span> <span class="special">(</span><span class="identifier">a</span> <span class="special">-</span> <span class="identifier">scalar_type</span><span class="special">(</span><span class="identifier">k</span><span class="special">)),</span> <span class="identifier">z</span><span class="special">);</span>
<span class="special">}</span>
<span class="special">};</span>
</pre>
<p>
So now we can implement Q, this time using <code class="literal">continued_fraction_a</code>:
</p>
<pre class="programlisting"><span class="keyword">template</span> <span class="special"><</span><span class="keyword">class</span> <span class="identifier">T</span><span class="special">></span>
<span class="keyword">inline</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="identifier">gamma_Q_as_fraction</span><span class="special">(</span><span class="keyword">const</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">>&</span> <span class="identifier">a</span><span class="special">,</span> <span class="keyword">const</span> <span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">>&</span> <span class="identifier">z</span><span class="special">)</span>
<span class="special">{</span>
<span class="identifier">upper_incomplete_gamma_fract</span><span class="special"><</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="special">></span> <span class="identifier">f</span><span class="special">(</span><span class="identifier">a</span><span class="special">,</span> <span class="identifier">z</span><span class="special">);</span>
<span class="identifier">std</span><span class="special">::</span><span class="identifier">complex</span><span class="special"><</span><span class="identifier">T</span><span class="special">></span> <span class="identifier">eps</span><span class="special">(</span><span class="identifier">std</span><span class="special">::</span><span class="identifier">numeric_limits</span><span class="special"><</span><span class="identifier">T</span><span class="special">>::</span><span class="identifier">epsilon</span><span class="special">());</span>
<span class="keyword">return</span> <span class="identifier">pow</span><span class="special">(</span><span class="identifier">z</span><span class="special">,</span> <span class="identifier">a</span><span class="special">)</span> <span class="special">/</span> <span class="special">(</span><span class="identifier">exp</span><span class="special">(</span><span class="identifier">z</span><span class="special">)</span> <span class="special">*(</span><span class="identifier">z</span> <span class="special">-</span> <span class="identifier">a</span> <span class="special">+</span> <span class="identifier">T</span><span class="special">(</span><span class="number">1</span><span class="special">)</span> <span class="special">+</span> <span class="identifier">boost</span><span class="special">::</span><span class="identifier">math</span><span class="special">::</span><span class="identifier">tools</span><span class="special">::</span><span class="identifier">continued_fraction_a</span><span class="special">(</span><span class="identifier">f</span><span class="special">,</span> <span class="identifier">eps</span><span class="special">)));</span>
<span class="special">}</span>
</pre>
</div>
<div class="copyright-footer">Copyright © 2006-2021 Nikhar Agrawal, Anton Bikineev, Matthew Borland,
Paul A. Bristow, Marco Guazzone, Christopher Kormanyos, Hubert Holin, Bruno
Lalande, John Maddock, Evan Miller, Jeremy Murphy, Matthew Pulver, Johan Råde,
Gautam Sewani, Benjamin Sobotta, Nicholas Thompson, Thijs van den Berg, Daryle
Walker and Xiaogang Zhang<p>
Distributed under the Boost Software License, Version 1.0. (See accompanying
file LICENSE_1_0.txt or copy at <a href="http://www.boost.org/LICENSE_1_0.txt" target="_top">http://www.boost.org/LICENSE_1_0.txt</a>)
</p>
</div>
<hr>
<div class="spirit-nav">
<a accesskey="p" href="agm.html"><img src="../../../../../../doc/src/images/prev.png" alt="Prev"></a><a accesskey="u" href="../internals.html"><img src="../../../../../../doc/src/images/up.png" alt="Up"></a><a accesskey="h" href="../../index.html"><img src="../../../../../../doc/src/images/home.png" alt="Home"></a><a accesskey="n" href="simple_continued_fraction.html"><img src="../../../../../../doc/src/images/next.png" alt="Next"></a>
</div>
</body>
</html>
|
