4. Numerical example

For a nonlinear estimation of CPT functions, some numerical examples are chosen from GNSS {Global navigation satellite systems = GPS (USA) + GLONASS (RU), GALILEO (EU), COM-PASS (CHN)} artificial satellite orbits whose coordinates are downloaded from the internet address ftp://ftp.glonass-iac.ru/MCC/PRODUCTS/17091/final/Sta19426.sp3 [26].

For this purpose, two estimation software have been developed in 64Bit Python (in accordance with the 2.7 and 3.6 version) and 32Bit C++ (Code::Blocks) environments to see the convergence rate of the mathematical model given in Eq. (12) [27, 28]. The computed elements of CPT functions for the selected four satellites R01 (GLONASS), G03 (GPS), E01 (GALILEO) and C06 (COMPASS) are summarized in Table 1 in which they are ordered from the nearest satellite to the farthest one.

The motions of the CPT functions estimated satellites (in Table 1) with respect to earth- (left column of Figure 3) and space-fixed (right column of Figure 3) coordinate systems are demonstrated in Figure 3. Moreover, coherence between the estimated CPT function (black solid line) of the C06 satellite and its data points given (colorful circles) is represented in detail in Figure 4.

We know that accuracy of precise SP3 file coordinates is about σ<sup>0</sup> ¼ �5 cm. If we compare the value with its estimations given in Table 1, we can say that our predicted model is not meet our demands. We should expand the model by raising the degree of polynomial part or/and order of trigonometric part of CPT functions. In fact, we can readily see that the projected model with Eq. (12) will never cover the data. The model is only chosen for this chapter. The more suitable model established on Keplerian orbital elements can be found in the orbit determination literature and in [18–20].

If we chose the threshold as 51e-13, we can see the distinctions of solution convergences between 64Bit and 32Bit running on the estimations of C06 satellite from 1000 (X), 8 (Y), 7 (Z), 1000 (S) in 32Bit C++ and 10 (X), 8 (Y), 6 (Z), 28 (S) in 64Bit Python in Windows. In here, 1000 is the maximum iteration number. If the mathematical model would be more complicated and its data number would be bigger than the number used in the example part, we would see the

Sat. φ X Y Z Unit S

4 5

2.771 0.052 22,965.361 �30�31<sup>0</sup>

�31�57<sup>0</sup>

�1.824 �0.085 21,757.547 �7�51<sup>0</sup>

30�05<sup>0</sup>

11<sup>h</sup>

4 4

�4.673 0.004 24,878.432 �16�23<sup>0</sup>

25�34<sup>0</sup>

14<sup>h</sup>

6 7

266.178 �7.892 34,235.695 67�36<sup>0</sup> 9.42<sup>00</sup> �15�02<sup>0</sup>

23h

11h

4 4 00.32<sup>00</sup>

56.31<sup>00</sup>

47.16"

29.47"

36.74"

55'49.94"

18.57"

04'40.78"

00.02"

58'00.32"

15'43.53"

— —

On Non-Linearity and Convergence in Non-Linear Least Squares

km km/h km deg deg/h h

— —

km km/h km deg. deg./h h

— —

km km/h km deg. deg./h h

— —

km km/h km deg. deg./h h

10 9

http://dx.doi.org/10.5772/intechopen.76313

25,508.091 0.001 8.127 42�41<sup>0</sup> 35.87<sup>00</sup>

�31�58<sup>0</sup>

26,561.324 �0.003 14.211 �72�24<sup>0</sup>

�30�08<sup>0</sup>

29,600.332 �0.000 3.720 8�29<sup>0</sup> 58.23" �25�42<sup>0</sup>

11<sup>h</sup>

32 93

14<sup>h</sup>

28 29

42,175.353 �0.358 227.599 18�32<sup>0</sup>

�14�57<sup>0</sup>

24<sup>h</sup>

27.85"

04'21.41"

17.11"

11h

30 34 13.55<sup>00</sup>

65

19.52"

11.17"

10.68"

00'22.19"

56'44.42"

15'37.46"

4 4

2.420 0.003 11,178.453 �23�03<sup>0</sup>

�31�57<sup>0</sup>

11h

4 4

7.332 �0.045 17,959.166 �48�06<sup>0</sup>

30�04<sup>0</sup>

11<sup>h</sup> 58'3.31"

5 5

1.358 �0.048 26,031.918 85�57<sup>0</sup>

�25�34<sup>0</sup>

228.037 �6.376 24,832.649 59�24<sup>0</sup> 50.73<sup>00</sup>

�15�02<sup>0</sup>

23h

14h

6 7 56.93<sup>00</sup>

44.92<sup>00</sup>

43.68"

� <sup>b</sup>σ<sup>0</sup> 4.898 3.648 3.961 km 0.183

52.50"

22.67"

04'42.81"

14.87"

16.41"

� <sup>b</sup>σ<sup>0</sup> 56.552 38.018 51.901 km 0.491

Table 1. Computed elements of the CPT functions for G03, R01, E01, C06 satellites by IterMAX = 1000 and thres = 5e-12 in loops {iteration numbers of 64Bit Python and 32Bit C++ software in windows are denoted as iter and as iter++ respectively}.

56'22.30"

� <sup>b</sup>σ<sup>0</sup> 1.796 1.428 1.058 km 0.236

� <sup>b</sup>σ<sup>0</sup> 2.902 1.266 2.653 km 0.211

15'47.54"

state more prominently.

iter iter++

aφ bφ cφ dφ eφ Tφ

iter iter++

aφ bφ cφ dφ eφ Tφ

iter iter++

aφ bφ cφ dφ eφ Tφ

iter iter++

aφ bφ cφ dφ eφ Tφ

R01 (RU)

G03 (USA)

E01 (EU)

C06 (CHN) 5 5

�11.860 0.028 25,474.503 �60�54<sup>0</sup>

31�57<sup>0</sup> 53.92<sup>00</sup>

11h

6 6

�17.545 �0.132 24,800.602 66�34<sup>0</sup>

30�04<sup>0</sup>

11h

4 5

1.866 0.052 21,349.011 �35�22<sup>0</sup>

�25�34<sup>0</sup>

�407.495 61.945 41,716.806 �24�21<sup>0</sup>

�15�07<sup>0</sup>

23<sup>h</sup>

14<sup>h</sup>

7 10 04.37<sup>00</sup>

05'44.37"

39.73"

56.03"

15.46<sup>00</sup>

13.06"

36.29"

02.75"

48'48.86"

04'43.81"

58'01.91"

Comparing the solution velocities (from the iteration numbers with respect to 5e � 12 threshold in Table 1) in different platforms, we can say that the solution velocities in 64Bit Python are generally better then 32Bit ones.


A<sup>j</sup> ¼ 1 tj sin d<sup>0</sup> þ f <sup>0</sup> tj

1

cated models, the readers can utilize [19–25] resources.

its solution velocity.

64 Optimization Algorithms - Examples

<sup>δ</sup><sup>b</sup> <sup>¼</sup> <sup>Q</sup> <sup>X</sup><sup>n</sup>

4. Numerical example

the farthest one.

determination literature and in [18–20].

generally better then 32Bit ones.

0 @

j¼1 AT j lj

l<sup>j</sup> ¼ ϕ<sup>j</sup> � a<sup>0</sup> þ b<sup>0</sup> tj þ c<sup>0</sup> sin d<sup>0</sup> þ e<sup>0</sup> tj � � � � h i

� � <sup>c</sup><sup>0</sup> cos <sup>d</sup><sup>0</sup> <sup>þ</sup> <sup>e</sup><sup>0</sup> tj

<sup>A</sup>, <sup>Q</sup> <sup>¼</sup> <sup>X</sup><sup>n</sup>

address ftp://ftp.glonass-iac.ru/MCC/PRODUCTS/17091/final/Sta19426.sp3 [26].

We can use a recursive solution for Eq. (14) instead of the batch solution as Eq. (11) because of

Continuation of the solution of Eq. (15) can be performed according to Eq. (11). The model given by Eq. (12) is a simple model to determine the satellite orbit motions. For more compli-

For a nonlinear estimation of CPT functions, some numerical examples are chosen from GNSS {Global navigation satellite systems = GPS (USA) + GLONASS (RU), GALILEO (EU), COM-PASS (CHN)} artificial satellite orbits whose coordinates are downloaded from the internet

For this purpose, two estimation software have been developed in 64Bit Python (in accordance with the 2.7 and 3.6 version) and 32Bit C++ (Code::Blocks) environments to see the convergence rate of the mathematical model given in Eq. (12) [27, 28]. The computed elements of CPT functions for the selected four satellites R01 (GLONASS), G03 (GPS), E01 (GALILEO) and C06 (COMPASS) are summarized in Table 1 in which they are ordered from the nearest satellite to

The motions of the CPT functions estimated satellites (in Table 1) with respect to earth- (left column of Figure 3) and space-fixed (right column of Figure 3) coordinate systems are demonstrated in Figure 3. Moreover, coherence between the estimated CPT function (black solid line) of the C06 satellite and its data points given (colorful circles) is represented in detail in Figure 4.

We know that accuracy of precise SP3 file coordinates is about σ<sup>0</sup> ¼ �5 cm. If we compare the value with its estimations given in Table 1, we can say that our predicted model is not meet our demands. We should expand the model by raising the degree of polynomial part or/and order of trigonometric part of CPT functions. In fact, we can readily see that the projected model with Eq. (12) will never cover the data. The model is only chosen for this chapter. The more suitable model established on Keplerian orbital elements can be found in the orbit

Comparing the solution velocities (from the iteration numbers with respect to 5e � 12 threshold in Table 1) in different platforms, we can say that the solution velocities in 64Bit Python are

0 @

j¼1 AT <sup>j</sup> A<sup>j</sup>

� � tj <sup>c</sup><sup>0</sup> cos <sup>d</sup><sup>0</sup> <sup>þ</sup> <sup>e</sup><sup>0</sup> tj

1 A �1

<sup>¼</sup> <sup>X</sup><sup>n</sup> j¼1 AT <sup>j</sup> A<sup>j</sup>

0 @

, j∈ f g 1; 2;…; n : (14b)

1 A þ

: (15)

� � � � , (14a)

Table 1. Computed elements of the CPT functions for G03, R01, E01, C06 satellites by IterMAX = 1000 and thres = 5e-12 in loops {iteration numbers of 64Bit Python and 32Bit C++ software in windows are denoted as iter and as iter++ respectively}.

If we chose the threshold as 51e-13, we can see the distinctions of solution convergences between 64Bit and 32Bit running on the estimations of C06 satellite from 1000 (X), 8 (Y), 7 (Z), 1000 (S) in 32Bit C++ and 10 (X), 8 (Y), 6 (Z), 28 (S) in 64Bit Python in Windows. In here, 1000 is the maximum iteration number. If the mathematical model would be more complicated and its data number would be bigger than the number used in the example part, we would see the state more prominently.

Figure 3. Earth (left) and space (right) fixed orbital traces (see the appendix) with time tags of R01, G03, E02, C06 satellites and the motion of X-coordinate axis shown as GR (XGR(t0) position on the intersection Greenwich meridian and equator) symbol at t0 (=2017 April 01 00:00:00).

Approximated values and the loop element for the unknowns are computed following the order in all solutions of the satellites when running in 32Bit C++ and 64Bit Python platforms for

Figure 4. Temporal changing of space fixed coordinates of C06. The circles and solid lines represent the data points and

x ¼ ½ � a<sup>0</sup> b<sup>0</sup> c<sup>0</sup> d<sup>0</sup> e<sup>0</sup>

Maximum iteration number and threshold loop elements are iterMAX ¼ 1000 and thres ¼ 5e � 12

n o, j <sup>∈</sup>f g <sup>1</sup>; <sup>2</sup>; …; <sup>n</sup> <sup>¼</sup> <sup>96</sup>

On Non-Linearity and Convergence in Non-Linear Least Squares

http://dx.doi.org/10.5772/intechopen.76313

67

T

� � � arcsin <sup>ϕ</sup>1=c<sup>0</sup> � � � � <sup>=</sup>f g <sup>t</sup><sup>2</sup> � <sup>t</sup><sup>1</sup>

a<sup>0</sup> ¼ 0:0, b<sup>0</sup> ¼ 0:0, c<sup>0</sup> ¼ max ϕ<sup>j</sup>

estimated functions under LS respectively for X (red), Y (green), Z (blue), S (cyan).

� �, e<sup>0</sup> <sup>¼</sup> arcsin <sup>ϕ</sup>2=c<sup>0</sup>

the estimations in Table 1.

to break the iteration loop.

d<sup>0</sup> ¼ arcsin ϕ1=c<sup>0</sup>

Figure 4. Temporal changing of space fixed coordinates of C06. The circles and solid lines represent the data points and estimated functions under LS respectively for X (red), Y (green), Z (blue), S (cyan).

Approximated values and the loop element for the unknowns are computed following the order in all solutions of the satellites when running in 32Bit C++ and 64Bit Python platforms for the estimations in Table 1.

$$a\_0 = 0.0, b\_0 = 0.0, c\_0 = \max\{\phi\_j\}, \quad j \in \{ \ 1, 2, \dots, n = 96 \}$$

$$d\_0 = \arcsin(\phi\_1/c\_0), e\_0 = \left\{ \arcsin(\phi\_2/c\_0) - \arcsin(\phi\_1/c\_0) \right\} / \left\{ t\_2 - t\_1 \right\}$$

$$\mathbf{x} = \begin{bmatrix} a\_0 & b\_0 & c\_0 & d\_0 & e\_0 \end{bmatrix}^T$$

Maximum iteration number and threshold loop elements are iterMAX ¼ 1000 and thres ¼ 5e � 12 to break the iteration loop.

Figure 3. Earth (left) and space (right) fixed orbital traces (see the appendix) with time tags of R01, G03, E02, C06 satellites and the motion of X-coordinate axis shown as GR (XGR(t0) position on the intersection Greenwich meridian and equator)

symbol at t0 (=2017 April 01 00:00:00).

66 Optimization Algorithms - Examples

Taking the approximated values as <sup>x</sup> <sup>¼</sup> ½ � 0 0 thres 0 0 <sup>T</sup> <sup>≈</sup> <sup>0</sup><sup>T</sup> and the same loop elements given above, the iteration numbers are handled as 138 (X), 178 (Y), 78 (Z), 16 (S) in 64Bit Python. For <sup>x</sup> <sup>¼</sup> ½ � 0 0 <sup>c</sup><sup>0</sup> 0 0 <sup>T</sup> they are 33 (X), 235 (Y), 96 (Z), 11 (S) in 64Bit Python. Different approximated value selections cause different iteration numbers (namely convergence rate).

#### 4.1. An expanded model example by an auxiliary cosine wave

Since the estimated standard deviations <sup>b</sup>σ<sup>0</sup> = {�56.552, �38.018, �51.901, �0.491} (for X, Y, Z, S in Table 1) of the CPT functions for the coordinates of the C06 satellite are not statically equal to their expected values (σ<sup>0</sup> ¼ �5 cm), the CPT model should be expanded. As an example, three more unknowns are added to the model given in Eq. (12)

$$\boldsymbol{a}\phi\_{j} = \boldsymbol{a}\_{\phi} + \boldsymbol{b}\_{\phi}\ \boldsymbol{t}\_{j} + \boldsymbol{c}\_{\phi}\ \sin\left(\boldsymbol{d}\_{\phi} + \boldsymbol{e}\_{\phi}\ \boldsymbol{t}\_{j}\right) + \boldsymbol{f}\_{\phi}\ \cos\left(\mathbf{g}\_{\phi} + \boldsymbol{h}\_{\phi}\ \boldsymbol{t}\_{j}\right)$$

The added terms represent an amplitude f <sup>ϕ</sup>, an initial phase g<sup>ϕ</sup> and a frequency h<sup>ϕ</sup> of a new wave carried by first (sine) wave. After the first estimation with respect to Eq. (12),

(Table 2) in kilometers for the C06 satellite. We can develop the last model more by means of

17.24"

� <sup>b</sup>σ<sup>0</sup> 0.178 0.137 0.191 km 0.003

Sat. φ X Y Z Unit S

149.296 0.224 66.183 11�15<sup>0</sup> 24.23<sup>00</sup>

30�04<sup>0</sup> 30.44<sup>00</sup> 24,808.034 30�28<sup>0</sup> 55.61<sup>00</sup>

15�02<sup>0</sup>

Table 2. The results by the expanded model for the C06 satellite in Python 3.6 in windows.

iter 27 158 53 — 19

182.012 0.060 91.527 3�15<sup>0</sup> 58.76<sup>00</sup> 30�04<sup>0</sup> 02.09<sup>00</sup> 34,217.994 �22�21<sup>0</sup>

�15�02<sup>0</sup>

14.60<sup>00</sup>

km km/h km deg deg/h km deg. deg./h

http://dx.doi.org/10.5772/intechopen.76313

On Non-Linearity and Convergence in Non-Linear Least Squares

42,169.956 0.055 225.161 19�47<sup>0</sup> 37.54<sup>00</sup>

�15�02<sup>0</sup>

0.838 8�56<sup>0</sup> 45.50<sup>00</sup> �29�40<sup>0</sup>

43.26<sup>00</sup>

69

58.52"

19.27"

As another example, the expanded model has been trained on the coordinate components of R01 GLONASS satellite by means of Python 3.6 software on Windows. We can also see the improvements upon the downs of the standard deviations (its iteration numbers) from <sup>b</sup>σ<sup>0</sup> <sup>=</sup> {�2.902 (5), �1.266 (4), �2.653 (4), �0.211 (10)} (Table 1) to <sup>b</sup>σ<sup>0</sup> <sup>=</sup> {�0.226 (76), �0.448 (32),

Condition numbers computed from a rate of maximum and minimum eigenvalues or a rate of singular values under LS are very effective tools for determining the consistency as well. Therefore, a larger condition number can cause larger iteration number (related to convergence rate), We can see those states from the CPT estimations of the C06 satellite with the iteration (iter) and condition numbers (cond). These are given for X, Y, Z, S as iter = {7, 6, 6, 28} and cond = {3.4e + 13, 2.0e + 12, 2.8e + 12, 1.3e + 09} (Table 1), and as iter = {27, 158, 53, 19} and cond = {9.5e +

In this chapter, the least squares (LS) estimations of the artificial satellite orbital movements by a combination of polynomial and trigonometric (CPT) functions have been given after a general overview has been made on the hard and soft computations. In practice, the orbital motions are modeled on Keplerian orbital elements. In contrary to this, the coordinate components have been selected for this chapter due to the nonlinear relations of the components and the unknowns which are the elements of CPT functions. The relations cause inconsistencies in the LS solutions. The inconsistencies result from the two injectivity defects, c-defects and i-defects. We can readily see the defects from the differences of the convergence rates (in other words the iteration numbers) in different computer platforms and architectures as shown in the chapter. The defects are not fully removed as long as not change the mathematic models. However, we can surpass the effects of those defects in part by means of the pseudoinverse based on the eigendecomposition or the singular value decomposition (SVD) as in here. The

the same manner if we want.

aφ bφ cφ dφ eφ fφ gφ hφ 223.514 1.515 42,064.675 �25�07<sup>0</sup>

�15�02<sup>0</sup>

110.710 5�08<sup>0</sup> 57.82" 30�10<sup>0</sup>

54.78<sup>00</sup>

24.75"

58.43"

C06 (CHN)

�0.201 (43), �0.037 (57)} in kilometers for the R01 satellite.

14, 2.3e + 13, 3.0e + 13, 2.6e + 10} (Table 2).

5. Conclusions

we can choose the approximate values of the new parameters as

$$f\_0 = \max(abs(\hat{\boldsymbol{\varepsilon}})) $$

$$g\_0 = \arcsin\left(\hat{\boldsymbol{\varepsilon}}\_1/f\_0\right) $$

$$h\_0 = \left\{\arcsin\left(\hat{\boldsymbol{\varepsilon}}\_2/f\_0\right) - \arcsin\left(\hat{\boldsymbol{\varepsilon}}\_1/f\_0\right)\right\}/\left\{t\_2 - t\_1\right\}.$$

from <sup>b</sup>ε<sup>j</sup> <sup>¼</sup> <sup>ϕ</sup><sup>j</sup> � <sup>b</sup><sup>a</sup> <sup>þ</sup> <sup>b</sup>b tj <sup>þ</sup>b<sup>c</sup> sin <sup>b</sup><sup>d</sup> <sup>þ</sup>be tj h i n o � � , <sup>j</sup>∈f g <sup>1</sup>; <sup>2</sup>;…; <sup>n</sup> <sup>¼</sup> <sup>96</sup> . The approximate value vector of the expanded model by a new wave is:

$$\mathbf{x} = \begin{bmatrix} a\_0 & b\_0 & c\_0 & d\_0 & e\_0 & f\_0 & g\_0 & h\_0 \end{bmatrix}^T = \begin{bmatrix} \hat{a} & \hat{b} & \hat{c} & \hat{d} & \hat{e} & f\_0 & g\_0 & h\_0 \end{bmatrix}^T$$

The approximate values are substituted in the following linearized model as initial values for the loop in LS estimation.

$$\mathbf{1}\_{j} = \begin{bmatrix} \phi\_{j} - \begin{pmatrix} a\_{0} + b\_{0} \ t\_{j} + c\_{0} & \sin\left(d\_{0} + f\_{0} \ t\_{j}\right) + e\_{0}\cos\left(h\_{0} + g\_{0} \ t\_{j}\right) \end{pmatrix} \end{bmatrix}$$

$$\mathbf{A}\_{j}^{T} = \begin{bmatrix} 1 \\ t\_{j} \\ \sin\left(d\_{0} + e\_{0} \ t\_{j}\right) \\\ c\_{0}\cos\left(d\_{0} + e\_{0} \ t\_{j}\right) \\\ t\_{j}c\_{0}\cos\left(d\_{0} + e\_{0} \ t\_{j}\right) \\\ \cos\left(g\_{0} + h\_{0} \ t\_{j}\right) \\\ -f\_{0}\sin\left(g\_{0} + h\_{0} \ t\_{j}\right) \\\ -t\_{j}\ f\_{0}\sin\left(g\_{0} + h\_{0} \ t\_{j}\right) \end{bmatrix}$$

After the evaluation, the improved solution for the C06 satellite is represented in Table 2. We can readily see the improvements upon the downs of the standard deviations from <sup>b</sup>σ<sup>0</sup> <sup>=</sup> {�56.552, �38.018, �51.901, �0.491} (Table 1) into <sup>b</sup>σ<sup>0</sup> <sup>=</sup> {�0.178, �0.137, �0.191, �0.003}


Table 2. The results by the expanded model for the C06 satellite in Python 3.6 in windows.

(Table 2) in kilometers for the C06 satellite. We can develop the last model more by means of the same manner if we want.

As another example, the expanded model has been trained on the coordinate components of R01 GLONASS satellite by means of Python 3.6 software on Windows. We can also see the improvements upon the downs of the standard deviations (its iteration numbers) from <sup>b</sup>σ<sup>0</sup> <sup>=</sup> {�2.902 (5), �1.266 (4), �2.653 (4), �0.211 (10)} (Table 1) to <sup>b</sup>σ<sup>0</sup> <sup>=</sup> {�0.226 (76), �0.448 (32), �0.201 (43), �0.037 (57)} in kilometers for the R01 satellite.

Condition numbers computed from a rate of maximum and minimum eigenvalues or a rate of singular values under LS are very effective tools for determining the consistency as well. Therefore, a larger condition number can cause larger iteration number (related to convergence rate), We can see those states from the CPT estimations of the C06 satellite with the iteration (iter) and condition numbers (cond). These are given for X, Y, Z, S as iter = {7, 6, 6, 28} and cond = {3.4e + 13, 2.0e + 12, 2.8e + 12, 1.3e + 09} (Table 1), and as iter = {27, 158, 53, 19} and cond = {9.5e + 14, 2.3e + 13, 3.0e + 13, 2.6e + 10} (Table 2).
