• For Contributors +
• Journal Search +
Journal Search Engine
ISSN : 2093-5587(Print)
ISSN : 2093-1409(Online)
Journal of Astronomy and Space Sciences Vol.34 No.1 pp.55-65
DOI : https://doi.org/10.5140/JASS.2017.34.1.55

# Analysis on the Hour-lines of Yangcheon-cheok

Sang Hyuk Kim1,2, Byeong-Hee Mihn1,2, Yong Sam Lee3
1Korea Astronomy and Space Science Institute, Daejeon 34055, Korea
2Korea University of Science & Technology, Daejeon 34113, Korea
3Department of Astronomy and Space Science, Chungbuk National University, Cheongju 28644, Korea
Corresponding Author:+82-42-865-2046, astro91@kasi.re.kr
September 27, 2016 December 26, 2016 December 28, 2016

## Abstract

We studied the hour-lines of six extant Yangcheon-cheoks in Korea. To find whether Yangcheon-cheok was used in the whole area of Korea, we calculated the length of shadows of Yangcheon-cheok on the Korean Peninsula, Nanjing and Beijing as well as Hanyang (Seoul), according to 24 solar terms. Comparing the length of shadows with hour-lines of those relics, we could find that Yangcheon-cheok was suitable for use at limited times (from 9:00 to 15:00) during the year. Also, this sundial is more appropriate for use at low latitudes than high ones. Among existing relics, that of Seoul Museum of History made with porcelain was much more suitable to use at Hanyang and its higher latitude. Lee’s collection was also suitable to use at Nanjing. It is certain that Yangcheon-cheok was a portable sundial which could be used from nine to fifteen of clock all the year around except for the winter season.

## 1.INTRODUCTION

It has been reported that there were two kinds of a heliotropic sundials in the late period of the Joseon dynasty; one was the Yanggyeong-gyuil-ui (兩景揆日儀) and the other was the Yangcheon-cheok (量天尺) (Kim et al. 2015; Lee et al. 2015). While Yanggyeong-gyuil-ui was introduced in the Gyuil-go (揆 日考, The Study on Solar Motion) and the Uigi-jipseol (儀器輯 說, The Comprehensive Lecture of Astronomical Instruments), Yangcheon-cheok was discussed in Sukjong Sillok (The Veritable Record of the King Sukjong). It has been known that a lot of sundials in the early period of Joseon dynasty were directed with the main axis toward the meridian (Jeon 1963; Lee 1982; Jeon et al. 1984; Kim 2007; Lee & Kim 2011; Mihn et al. 2011; Kim et al. 2012; Mihn et al. 2012, 2013). In the late period of Joseon, a new-fashioned sundial was introduced from the West via China (Lee et al. 2015).

It has been recognized that all of the Yangcheon-cheok have been investigated through the project of cataloguing Korean scientific cultural assets (Jeon et al. 1984). It has been pointed out that various sizes of Yangcheon-cheok all had the same instructions and statistical similarity with respect to hourlines (Kim et al. 2015).

In this study, we analyzed whether Yangcheon-cheok could be used in the whole area of the Korean Peninsula or not. Though the representative place of Korea is Hanyang (Seoul), we performed calculations for a more extensive range - from Nanjing to Beijing. The reason why the latitudinal range covered Nanjing and Beijing was that Yangcheon-cheok was originally introduced by ambassadors of the Qing dynasty. On six Yangcheon-cheok relics (Kim et al. 2015), we calculated the length of shadows by season-changed stylus and compared their lengths with the hour-lines of the six relics. We also analyzed the relics’ time accuracy in terms of solar time measurement divided according to 24 solar terms (氣).

## 2.RELICS

### 2.1.Designation

According to Kim et al. (2015), there are specifications of six Yangcheon-cheoks that remain in Korea. For convenience we numbered these relics. No. 1 and 4 relics were made with wood (W) and porcelain (P), respectively, and are preserved by Seoul Museum of History (SMH). No. 2 and 3 sundials are possessed by Yeoju Museum (YM) and Sungshin Women’s University Museum (SWUM), respectively. No. 5 is in Lee’s private collection (Lee). No. 6 is a reproduction that is kept at National Science Museum (NSM). Except for No. 4, all of the Yangcheon-cheoks were made with wood. These relics were not the same in size; and have different lengths, widths, thickness, interval of hour-lines, and so on. Two relics were distinguished by their materials; SMH(W) and SMH(P).

A previous result that was conducted by Kim et al. (2015) revealed that interval of hour-lines of No. 1 SMH(W) and No. 2 YM had the same proportion to their own stylus even though they had different sizes. No. 3 SWUM relic and No. 4 SMH(P) also had a similar pattern between line-interval and gnomon length only in the summer season.

### 2.2.Size and Position of Stylus

A stylus of Yangcheon-cheok was changed according to the season specified by this sundial itself. Among four kinds of styluses, a 5 pun stylus and a 10 pun stylus were used in the winter season and summer season, respectively. A pun (分) is different to each relics and is corresponding to about 2~3 mm. A 6 pun stylus corresponded to early spring and late autumn (hereafter ‘outer winter season’) and an 8 pun stylus corresponded to late spring and early autumn (hereafter ‘outer summer season’). In the old Korean luni-solar calendar as well as the old Chinese one, the four seasons were subdivided into 24 solar terms (24氣). Table 1 shows the 24 solar terms corresponding to the four season of Yangcheon-cheok and also represent the stylus length within every solar terms. There is the date of the 24 solar terms for the year 2015 in brackets (Ahn & Mihn 2014).

As shown in Fig.1, there are two locations in which to locate the stylus; an upper hole and a lower hole. The 5 pun and the 6 pun stylus are put in the upper hole of this sundial in winter and outer winter season. Meanwhile, the 8 pun and the 10 pun stylus are attached in the lower hole in summer and outer summer season.

## 3.CALCULATION

The Yangcheon-cheok’s instructions say that a stylus can be put in two positions depending on the season and a position is shared in two seasons. Several lines on the lateral side of these sundials are carved to mark the stylus length. Let the longest line interval be the unit of a stylus length. The unit length of No. 1, 2, 3, 4, 5, and 6 artifacts is 32 mm, 16 mm, 14 mm, 16 mm, 20.7 mm, and 30 mm, respectively. If the unit length of each relic was normalized to 1 chi (=10 pun), they could be represented as 5 pun, 6 pun, 8 pun and 10 pun of stylus as in Table 1 (Kim et al. 2015).

The length of shadow by four styli corresponding to 24 solar terms can be calculated. To do so, solar declination is first obtained from the ecliptic longitude at each solar term that increases by 15° from 0° at vernal equinox. With the spherical trigonometry, an altitude of the Sun, h, could be expressed as in the following equation.

(1)
(2)

Here, δ, ε, and λ are the solar declination, the obliquity of the ecliptic, and the ecliptic longitude, respectively. ϕ is observation latitude or altitude of the North Pole. The latitudes of Hanyang, Nanjing, and Beijing are given as 37°.56, 32°.05, and 39°.9, respectively. HA is an hour angle measured from the meridian. The hour-line of Yangcheon-cheok used together two hour angles (that is a plus and minus hour angles) at the same distance from the meridian. Here the unit of time, t, is a true solar time. Because the altitude of the Sun is the function of t, the length of shadow by a stylus, C(t), is represented in the following expression:

$C t = l cot h t$
(3)

Here, l is a stylus length which is shown in Table 1 and changes according to the 24 solar terms.

### 3.2.Hour-lines of Relics

Table 2 shows the normalized length from the upper hole to each hour-line for all relics (Kim et al. 2015). The first hour-line (HL1) corresponded to 11:00 and 13:00 and hourlines every double-hour were drawn after HL1. For 12 solar terms of the winter and outer winter seasons, this sundial user put the long stylus on the upper hole and measures the length of shadow cast by the Sun when this sundial is aimed to (refer to the Table 1). The length from the upper hole to hour-lines at 13:00, 15:00, 17:00, and 19:00 in Table 2 should be identical to the lengths of shadow calculated by Eq. (3). In the same manner, the length from the lower hole to hourlines at 13:00, 15:00, 17:00, and 19:00 which means the lower hole’s value subtracts from every hour-line value in Table 2 should be each calculation of Eq. (3) for 12 solar terms of the summer and outer summer seasons.

In Table 2, each Yangcheon-cheok can be utilized at a minimum of four shadow-length points with respect to a double-hour. The cubic spline interpolation line passing through four shadow-length points of each relic (Press et al. 1992) can be compared with the calculated length of shadow by Eq. (3). Because it can be supposed that the origin of hour in Yangcheon-cheok at the upper or lower hole does not stand for its 12:00 (noon) hour-line, its length of shadow at noon can be an appropriate point close to the induced point which meets the cubic spline curve of four shadow-length points at 12:00. For example, we selected 4.96 pun for SMH(W) as the shadow length at noon in winter season and 1.33 pun as the length in summer season. We also used 8.25 pun and 2.96 pun for SMH(P) in winter and summer season, respectively.

## 4.ANALYSIS

The stylus shadow of Yangcheon-cheok is the same length twice a day before and after noon, because Yangcheon-cheok is a sundial pointing to the Sun. After the Sun rises, the length of shadow shortens until the culmination at noon, and after noon, the length of shadow becomes longer until the setting of the Sun. The length of shadow with respect to time in winter season is longer than that in summer season. The shorter the length of the day is, the earlier the stylus shadow disappeared.

The results of cubic spline curves (or R(t)) of relics and calculated curves (or C(t)) by Eq. (3) at Hanyang are shown in Fig. 2. In this diagram, the abscissa is the time before and after noon and the ordinate is length of shadow by a stylus. This figure is divided into four faces by four seasons which include 6 solar terms, respectively. Figs. 2(a) and 2(b) are the outer summer and summer seasons and Figs. 2(c) and 2(d) are the outer winter and winter seasons. However every season’s graph yielded four calculated curves, not six, because two solar terms are paired as the same solar declination; for example of the pair of ‘Clear and bright’ and ‘White dew’ and the pair of 'Crop rains' and 'End of heat' in the outer summer season. Cubic spline curves of six relics’ hour-line are drawn by a line and point; the hour-line's drift curves of SMH(W), YM, SWUM, SMH(P), Lee's collection, NSM are drawn by a solid line with a solid circle, a dot line with an open inverse triangle, a short dash line with a solid rectangle, a two-dot-and-dash line with an open diamond, a long dash line with a solid triangle, and a one-dot-and-dash line with an open circle, respectively. The bunch of six cubic spline curves is identical to all faces in Fig. 2.

As shown in Fig. 2, the bundle of four computed curves has the smallest binding at noon, which means that the difference of the shadow length in each 6 solar terms is small. The further from noon the time is, the wider the curves’ binding and the bigger the difference of the shadow length is. That is why Eq. (3) is a cotangent function, where it grows gently and extremely before and after 45°, respectively, that is related to the solar altitude. Near to the time at sunrise and sunset, the bundles of computed curves of Figs. 2(a) and 2(c) are spread widely and the bundles of Figs. 2(b) and 2(d) are not. This is because daily changes in the solar declination are much bigger in spring and autumn than in summer and winter. In Eq. (2), if δ is maximum (= ε), then δεsinλ which means δ is a function of sinλ.

If the cubic spline curve of relics intersects with calculated curves at any point, it means that the hour-line of the relic at that instant coincides to the length of shadow by the very same stylus of the solar term. In other words, it is demonstrated which relic is suitable for any solar term by investigating whether the differences between R(t) and C(t) are large or small. For example, if any relic’s cubic spline curve is positioned within the bunch of calculated curves in any time range, then this relic can be used at least in a corresponding solar term. The cubic spline curve cannot satisfy all six solar terms in any season because computed curves for each season greatly change according to time in a day.

For example, when comparing C(t) and R(t) at 09:00 (15:00) in Fig. 2(d), the length of shadow of SMH(P) and SWUM is positioned in the bundle of computed curves. The length of shadow of the other four relics is positioned out of the computed curves as winter begins. While in the case of Fig. 2(b), lengths of shadow of six relics meet with the bundle of four computed curves at 09:00 (15:00). The length of shadow of the relic in Lee’s collection is out of the bundle of computed curves at 07:00 (17:00).

Considering a Yangcheon-cheok’s inflow from China, the region of computed curves needs to be expanded to Beijing and Nanjing in China. It gives a clue on a verification if these relics can be used to the whole of Korea. The results in Beijing and Nanjing are shown in Figs. 3 and 4.

The features of Yangcheon-cheok can be ascertained by the examination of the computed curves and cubic spline curves at three regions in Figs. 2-4. First, this sundial is more applicable at low latitudes than at high latitudes. When the bundle of computed curves on 6 solar terms has a narrow width this means that it is more possible to lessen the value of C(t)-R(t). C(t)-R(t) is a gap of a computed curve by Eq. (3) and drift curve of a relic’s hour-line by the cubic spline interpolation. Second, Yangcheon-cheok is more useful in the summer season than in the winter season. These sundials can be utilized from 9:00 to 15:00 in the winter season, and from 8:00 to 16:00 in the summer season. Because the length of the day in the winter season is short, as the Sun rises late and sets early, the inflection point and asymptotic line of the computed curve is generated earlier from noon in the winter season than in the summer season. Because cubic spline curves change more gently than the bunch of computed curves, computed curves in summer season are closer to cubic spline curves than those of the winter season. Nevertheless in comparison with Figs. 2(d), 3(d), and 4(d), R(t) of all relics at 08:30(15:30) is positioned within the bundle of computed curves in Nanjing, but most of R(t) in Hanyang and Beijing is outside this bundle.

### 4.2.Time Difference Between Two Curves

The credibility of measuring the time of relics can be estimated by comparing the differences of computed curves at 24 solar terms and drift curves for each relic. Figs. 5-7 shows lines of C(t)-R(t) related to Hanyang, Beijing, and Nanjing, applied to 6 relics in relation to 24 solar terms.

In Figs. 5-7, the horizontal axis represents before and after the time of southing, and the vertical axis does C(t)-R(t). If the value of C-R is bigger than zero, calculated shadow length is longer than the time-line (or the drift curve) of the relic. Otherwise, if the value of C-R is smaller than zero, computed length of shadow is shorter than the line of the relic. The ideal result for each Yangcheon-cheok relic is that the value of C-R remains near zero even though time has elapsed. The computed curves of 24 solar terms are distinguished each season by colors. Black, red, blue, and green lines are corresponding to curves for outer summer, summer, outer winter, and winter, respectively. Six solar terms for each season are also represented with four kinds of line; a solid, a short dash, a normal dash, and a long dash line. A normal and a long dash line represent two solar terms to be the same solar declination.

The graph for six relics in each Figs. 5-7 represents the graph box (box-1), where the C-R line of solar terms was analyzed in the range of ±1.0 pun from 13:00 to 15:00 at the upper left corner and expanded in the C-R range of ±1.5 pun to be read conveniently. Also, it represents the box (box- 2) in the range of ±1.0 pun for C-R solar term from 16:30 to 17:30 at the lower right corner. Considering the available time in which Yangcheon-cheok can be used (refer to Section 4.1), we studied the relics for all seasons with the former one (box-1) and studied the relics in all seasons except for winter season with the latter one (box-2). Yangcheon-cheok informs of the time with the line in double-hour, so we made box- 1 in the size of 2 hr. And then, we set the time around 17:00 (7:00), after 2 hr (before 2 hr) as box-2. If C-R solar terms of range of [13:00, 15:00] is in the range of ±1.5 pun at box- 1, Yangcheon-cheok can be evaluated as being usable. Also, if the C-R line of three seasons (outer summer, summer, outer winter) passes 0 at 17:00 (7:00) at box-2, Yangcheon-cheok can be used on the solar term of that line.

Table 3 represents the number of C-R solar term of each relic in box-1 at Hanyang, Beijing, and Nanjing. At box 1 in Figs. 5-7, we score 1 point, if C-R line is in the range of [-1, 1] in 13:00-15:00, while we record 0.5 point, if it is in the range of [-1, 1] at 13:00-14:00 or 14:00-15:00. The former is marked as a square in Table 3 and the latter is marked as a square with parenthesis. For example, ■ is 1 point, (■■) is also 1 point, and ■(■■) is 2 point. In the graph, if C-R line corresponds to a curve that overlaps two solar terms, we also count 2 point for that curve.

Table 4 represents the number of the C-R line of solar terms for each relic passing box-2 at Hanyang, Beijing, and Nanjing. Box-2 is usually useful in the summer season (late spring, summer, and early autumn), so we add 1 point when the C-R line passes 17:00 (7:00). In the same as the case of box-1, if C-R line is a curve for two solar terms, we also count two points.

As shown in Tables 3 and 4, black, red, blue, and green solid squares represent the outer summer, the summer, the outer winter, and the winter season, respectively. According to Figs. 5-7, a total of 24 C-R lines can pass in box-1, and 18 C-R lines can pass in box-2.

In Table 3, SMH(P) has the most points with 16.5 point, and next are SWUM, Lee’s one that obtained 14.0, and 13.0 points, respectively in Hanyang. The average of these three is 14.5 point (60 %), and the average of the six relics is 11.0 point (46 %). In Beijing, SMH(P), Lee’s relic, SWUM obtained 15.5, 11.0, and 9.0, respectively. The average of these three is 11.8 point (49 %), and the average for all the relics is 8.4 point (35 %). In Nanjing, Lee’s relic, SWUM and YM obtained 19.0, 18.5, and 14.5 point, respectively. The average of the three relics is 17.3 point (72 %), and the average of all the relics is 14.7 point (61 %). Comparing Yangcheon-cheok at Beijing, Hanyang, and Nanjing, the average value of all relics is the biggest at Nanjing and smallest at Beijing as we predicted in Section 4.1. As we limited the time from 9:00 to 15:00, SMH(P) is the most useful in Hanyang and Beijing, and Lee’s relic is the most useful at Nanjing in all seasons.

In Table 4, we investigated usable 18 solar terms except the winter season at 17:00 (7:00), and most of the relics are not suitable to use as sundials. However, YM and SWUM are usable at some solar terms in Hanyang and Beijing, and SWUM is usable in only Nanjing.

## 5.CONCLUSIONS

Yangcheon-cheok was a sundial in the Joseon period which kept a standard type with the form of gnomon (Kim et al. 2015). To confirm whether the Yangcheon-cheok was used in Korea or not, we analyzed the possible use of hour-lines carved on six Yangcheon-cheoks. It is presumed that Yangcheon-cheok was used in the whole area of Korea. Considering that this sundial was firstly introduced from China, we calculated the length of shadow over the latitude of Nanjing and Beijing as well as Hanyang. Through comparing the hour-line of relics with the calculated curves of the three regions, we derived the following data.

First, the extant Yangcheon-cheoks can be successfully used during a limited time period from 9:00 to 15:00. It was also learned that the real shadow line by the Sun is more similar to the hour-lines of Yangcheon-cheok in the summer than in the winter.

Second, Yangcheon-cheok is more suitable for use at low latitudes than at high latitudes. The lower the latitude this sundial uses, the narrower the bundle of the shadow lines by the Sun spreads within any season.

Third, there is a useful Yangcheon-cheok for the wide range of latitudes. It was shown that SMH(P) performed best at Hanyang or at higher latitudes, but Lee’s relic makes the best performance when used in Nanjing.

According to the results, though more or less accurate, the six relics can be evaluated as a useful sundial in Korea within the restricted time from 9:00 to 15:00. The merit of Yangcheon-cheok is that it is structurally simple and easy to reproduce (Kim et al. 2015). It is considered that it was a personal and portable clock and that anyone could easily learn the instructions described on the backside. It is also assumed that, in late Joseon period, Yangcheon-cheok was widely used.

## ACKNOWLEDGMENTS

The research was supported by the Korea Astronomy and Space Science Institute under the R&D program supervised by the Ministry of Science, ICT and Future Planning of Korea.

## Figure

The general structure of Yangcheon-cheok (Kim et al. 2015).

The cubic spline curve, R(t) and the computed curve, C(t) at Hanyang (Seoul, 37°.56) by Eq. (3). The hour-line's drift curves of SMH(W), YM, SWUM, SMH(P), Lee's collection, NSM are drawn by a solid line with a solid circle, a dot line with an open inverse triangle, a short dash line with a solid rectangle, a two-dot-and-dash line with an open diamond, a long dash line with a solid triangle, and a one-dot-and-dash line with an open circle, respectively. While four calculated curves of six solar terms are a solid line, a dot line, short dash line, and a two-dot-and-dash line, respectively.

The cubic spline curve, R(t) and the computed curve, C(t) at Beijing (39°.90) by Eq. (3). The description of lines refers to Fig. 2.

The cubic spline curve, R(t) and the computed curve, C(t) at Nanjing (32°.05) by Eq. (3). The description of lines refers to Fig. 2.

Difference between a computed curve, C(t) and a drift curve, R(t) at Hanyang (Seoul, 37°.56). Black, red, blue, and green lines are corresponding to curves for outer summer, summer, outer winter, and winter, respectively. The box at the lower right corner is ranged from 16:30 to 17:30 and between ±1.0 pun for C-R.

Difference between a computed curve, C(t) and a drift curve, R(t) at Beijing (39°.90). The description of lines and a box refers to Fig. 5.

Difference between a computed curve, C(t) and a drift curve, R(t) at Nanjing (32°.05). The description of lines and a box refers to Fig. 5.

## Table

Stylus length and position with respect to the four seasons (Kim et al. 2015)

1The period of the outer summer season covers Mar. 21 to May 5 and Aug. 08 to Sep. 22.
2The period of the summer season spans from May 08 to Sep. 07.
3The period of the outer winter season covers Feb. 04 to Mar. 20 and Sep. 23 to Nov. 07.
4The period of the winter season spans from Nov. 08 to Feb. 03 of the next year.

Normalized hour-line's distance from the origin of the upper hole (unit: pun (分)) (Kim et al. 2015)

*HL1, HL2, and so forth on are short for the first hour-line, the second hour-line, and so on, respectively.

The number of C-R solar term lines included in box-1 (total 24 solar terms)

The number of C-R solar term lines included in box-2 (total 18 solar terms)

## Reference

1. Ahn YS , Mihn BH (2014) The Korean Astronomical Almanac for the Year 2015, Korea Astronomy and Space Science Institute, ; pp.15
2. Jeon S.W (1963) Study on the Clock Manufacture During the Yi dynasty, Hyangto, ; pp.49-114
3. Jeon S.W , Kim S.S , Kim J.H , Nha I.S , Nam C.W (1984) Science and cultural heritage of Korea research report , Korean J. Hist. Sci, Vol.6 ; pp.58-118
4. Kim S.H (2007) A Study on the mechanism and using method of sundials in King Sejong era , Chungbuk Sahak, Vol.19 ; pp.97-123
5. Kim S.H , Mihn B.H , Lee M.S , Lee Y.S (2012) A Study on the structure innovation of astronomical instruments in Joseon dynasty – focusing on the Soganui, the Ilseongjeongsiui, and the Jeokdogyeongwiui , Publ. Korean Astron. Soc, Vol.27 ; pp.61-69https://dx.doi.org/10.5303/PKAS.2012.27.3.061
6. Kim S.H , Mihn B.H , Lee Y.S (2015) A preliminary study on the Yang-cheon-cheok (量天尺) in the late Joseon dynasty , J. Astron. Space Sci, Vol.32 ; pp.395-401https://dx.doi.org/10.5140/JASS.2015.32.4.395
7. Lee E.S (1982) History and principles of the sundial , Dong Bang Hak Chi(東方學志), Vol.33 ; pp.81-156
8. Lee Y.S , Kim S.H (2011) A study for the restoration of the sundials in King Sejong era , J. Astron. Space Sci, Vol.28 ; pp.143-153https://dx.doi.org/10.5140/JASS.2011.28.2.143
9. Lee Y.S , Kim S.H , Mihn B.H (2015) A study of the Yang-gyeong-gyu-il-ui in the Joseon dynasty , J. Astron. Space Sci, Vol.32 ; pp.73-80https://dx.doi.org/10.5140/JASS.2015.32.1.73
10. Mihn B.H , Kim S.H , Lee K.W , Ahn Y.S , Lee Y.S (2011) A development history and structural feature of Sogyupyo in the Joseon dynasty , Publ. Korean Astron. Soc, Vol.26 ; pp.129-140https://dx.doi.org/10.5303/PKAS.2011.26.3.129
11. Mihn B.H , Lee K.W , Kim S.H , Ahn Y.S , Lee Y.S (2012) The study on the structure of Daegyupyo in the early Joseon dynasty , Publ. Korean Astron. Soc, Vol.27 ; pp.29-37https://dx.doi.org/10.5303/PKAS.2012.27.2.029
12. Mihn B.H , Lee K.W , Kim S.H , Lee Y.S (2013) The role of a cross-bar and the enlargement of a Gnomon in Joseon dynasty , Publ. Korean Astron. Soc, Vol.28 ; pp.55-63https://dx.doi.org/10.5303/PKAS.2013.28.3.055
13. Press W.H , Flannery B.P , Teukolsky S.A , Vetterling W.T (1992) Numerical Recipes in Fortran 77 : The Art of Scientific Computing, Cambridge University Press, ; pp.99-110