May 23, 2003 Compiled on September 10, 2023 at 8:33am

1 Problem

2 Assumptions

3 Method

4 Analysis

4.1 Common calculations

4.2 Decide how to correct the LEO plane inclination

4.2.1 First scenario. All plane correction at perigee.

4.2.2 Second scenario. All plane correction at apogee.

4.2.3 Third Scenario. Partial plane correction at perigee. Rest at apogee.

4.2.4 Summary of scenarios to correct plane inclination

4.3 Calculate time to move from LEO to GEO

4.4 Decide when to inject to GEO. Calculate lead angle \(\beta \) and rendezvous with ﬁrst satellite

4.4.1 First option. Inject to GEO at \(t=0\)

4.4.2 second option: phase-wait in the LEO orbit until the correct lead angle \(\beta \) reached at one end of lines of nodes.

4.4.3 Third option: phase-wait in the LEO orbit until the correct lead angle \(\beta \) reached at the other end of lines of nodes.

4.4.4 Fourth option. phase-wait in the LEO orbit until smallest diﬀerence to \(\beta \) reached ﬁrst time.

4.4.5 Summary of lead angle \(\beta \) scenarios

4.5 Rendezvous with the second satellite

4.6 Stay locked in with the second satellite.

4.7 reposition to ﬁnal destination.

5 Summary and Results

6 Appendix

2 Assumptions

3 Method

4 Analysis

4.1 Common calculations

4.2 Decide how to correct the LEO plane inclination

4.2.1 First scenario. All plane correction at perigee.

4.2.2 Second scenario. All plane correction at apogee.

4.2.3 Third Scenario. Partial plane correction at perigee. Rest at apogee.

4.2.4 Summary of scenarios to correct plane inclination

4.3 Calculate time to move from LEO to GEO

4.4 Decide when to inject to GEO. Calculate lead angle \(\beta \) and rendezvous with ﬁrst satellite

4.4.1 First option. Inject to GEO at \(t=0\)

4.4.2 second option: phase-wait in the LEO orbit until the correct lead angle \(\beta \) reached at one end of lines of nodes.

4.4.3 Third option: phase-wait in the LEO orbit until the correct lead angle \(\beta \) reached at the other end of lines of nodes.

4.4.4 Fourth option. phase-wait in the LEO orbit until smallest diﬀerence to \(\beta \) reached ﬁrst time.

4.4.5 Summary of lead angle \(\beta \) scenarios

4.5 Rendezvous with the second satellite

4.6 Stay locked in with the second satellite.

4.7 reposition to ﬁnal destination.

5 Summary and Results

6 Appendix

Inject spacecraft directly into circular parking orbit of altitude 100km. At burnout ﬂight path angle is zero and velocity is \(V_{bo}\).

Even though launch at equator (line of nodes), Orbit has an undesired \(15^{0} \) inclination in the geocentric frame and longitude of ascending node \(\Omega =20^{0}\). Need to correct this inclination so that you may rendezvous with a satellite in the GEO (35,860 km) on the equatorial plane with zero inclination.

This target satellite was \(40^{0}\) behind you at the time you entered your parking orbit of 100 km.

After you rendezvous with the ﬁrst satellite, you need to transfer to a second satellite in GEO which was \(10^{0}\) ahead of you at time you entered your parking orbit.

Finally after making one complete orbit with this second satellite you need to take your ﬁnal position in GEO which is another \(5^{0}\) ahead of the second satellite.

Design the above sequence. You may opt to minimize the time to complete or the fuel needed \(\left ( \Delta V\right ) \). Specify your design criteria.

No drag while in LEO orbit. This allows the spacecraft to orbit as many times as needed to improve rendezvous conditions with the ﬁrst target in GEO.

All impulses applied are assumed to have inﬁnitely small time durations.

In addition, all assumptions used to derived Kepler equation apply as well.

Launch site latitude eﬀect on \(\Delta V\) are ignored. In practice, \(\Delta V \) requirement need to be modiﬁed by small magnitude depending on the launch site location on the surface of the earth.

The geometry of the problem is illustrated in the ﬁgure below, which depicts the state of the system at \(t_{0}\). The time the space vehicle enters its parking orbit at burn out \(V_{bo}\).

The analysis part will show an outline of the maneuvers to achieve the goal of the project.

The criteria for selecting the maneuvering sequence is:

The minimization of fuel

Which directly relates to the minimization of the \(\Delta V.\)

Now, we will show each phase of the maneuver, with diﬀerent scenarios to achieve each phase. Before starting, there are common calculations that will be done now that will be shared by many scenarios, so that we do not have to re-calculate these each time.

For the Hohmann transfer from LEO to GEO:

Let \(V_{p}\) be the speed on the Hohmann orbit at the perigee point.

Let \(V_{a}\) be the speed on the Hohmann orbit at the apogee point.

Let \(2a\) be the semi major axis for the Hohmann orbit.

Let \(T_{h}\) be the period of the Hohmann orbit.

For the LEO and GEO orbits:

Let \(V_{leo}\) be the speed on the LEO orbit.

Let \(V_{geo}\) be the speed on the GEO orbit.

Let \(r_{leo}\) be the radius of the LEO orbit.

Let \(r_{geo}\) be the radius of the GEO orbit.

Let \(T_{leo}\) be the period of the LEO orbit

let \(T_{geo}\) be the period of the GEO orbit

Then

\(r_{leo}\) \(=6378.145+100=6478.145\) km

\(r_{geo}\) \(=6378.145+35,860=42,238.145\) km

\(a=\frac {r_{leo}+r_{geo}}{2}=\frac {6478.145+42238.145}{2}=24358.145\) km

\(V_{p}=\) \(\sqrt {\mu (\frac {2}{r_{p}}-\frac {1}{a})}\) but \(r_{p}=r_{leo}\)

so, \(V_{p}=\) \(\sqrt {3.986012\times 10^{5}\left ( \frac {2}{6478.145}-\frac {1}{24358.145}\right ) }\) \(=10.3294\) km/sec \(=37,185.84\) km/hr

\(V_{a}=\) \(\sqrt {\mu (\frac {2}{r_{a}}-\frac {1}{a})}\) but \(r_{a}=r_{geo}\)

so, \(V_{a}=\) \(\sqrt {3.986012\times 10^{5}\left ( \frac {2}{42238.145}-\frac {1}{24358.145}\right ) }=\) \(1.584\ \)km/sec \(=5,702.4\) km/hr

\(V_{leo}=\sqrt {\frac {\mu }{r_{leo}}}=\sqrt {\frac {3.986012\times 10^{5}}{6478.145}}=7.844\ \) km/sec \(=28,238.82\) km/hr

\(V_{geo}=\sqrt {\frac {\mu }{r_{geo}}}=\sqrt {\frac {3.986012\times 10^{5}}{42238.145}}=3.072\ \)km/\(\sec =11,059.089\ \)km/hr

\(2a=r_{leo}+r_{geo}=\left ( 6378.145+100\right ) +\left ( 6378.145+35860\right ) =48,716.29\) km. Hence \(a=24,358.145\ \)km.

\(T_{leo}=2\pi \sqrt {\frac {r_{leo}^{3}}{\mu }}=2\pi \sqrt {\frac {6478.145^{3}}{3.986012\times 10^{5}}}=\) \(5,189\ \) sec = 1 hr 26 minutes and 28 seconds.

\(T_{geo}=2\pi \sqrt {\frac {r_{geo}^{3}}{\mu }}=2\pi \sqrt {\frac {42238.145^{3}}{3.986012\times 10^{5}}}=\) \(86,391\) sec = 23 hrs and 59 minutes and 50 seconds.

\(T_{h}=2\pi \sqrt {\frac {a^{3}}{\mu }}=2\pi \sqrt {\frac {24358.145^{3}}{3.986012\times 10^{5}}}=37,833\ \)sec \(=10\ \)hr 30 min 22 seconds

The problem is that we are given two non-coplanar circular orbits of diﬀerent radices. A LEO orbit that is inclined at an angle to the plane of another, and larger circular GEO orbit.

We wish to transfer from the inclined (\(i=15^{0})\ \)LEO orbit to the equatorial (\(i=0^{0})\ \)GEO orbit.

We must correct the plane inclination to be able to transfer to the desired GEO orbit.

There are 3 possible ways to achieve this^{1}:

- Correct all of the plane inclination before performing a Hohmann transfer from LEO to GEO. In other words, all of the inclination correction is made at the perigee of the Hohmann ellipse where the Hohmann transfer speed is largest. This will turn out to be the most fuel costly maneuver.
- First perform a Hohmann transfer to transfer from the inclined LEO orbit to an inclined GEO orbit, and then apply all of the plane inclination correction at the apogee of the Hohmann elliptical orbit where the ellipse speed will be smallest. This is less costly in \(\Delta V\ \)than the above sequence, and is a common maneuver.
- Apply a small and partial plane inclination correction (say an angle \(\alpha )\) at the perigee of the Hohmann orbit, then apply the remaining inclination correction (angle \(\beta -\alpha \)) at the apogee. Notice that for \(\alpha =0\), this case becomes the same as case (2) above.

Now, we will analyze each case above in details and ﬁnd the \(\Delta V\) for each case and select the maneuver with the smallest \(\Delta V.\)

Move from the initial circular parking orbit (which has \(15^{0}\) degrees inclination) to a new circular orbit of the same radius but on the equatorial (\(0^{0}\) degrees inclination). This requires one impulse to adjust the inclination. This impulse applied at the point where the parking orbit intersects the equator (line of nodes).

Next, and immediately, apply a coplanar Hohmann transfer (2 impulses) to transfer from the LEO orbit to the outer GEO orbit. (We Can combine the inclination correction velocity impulse vector with the ﬁrst Hohmann velocity impulse vector using vector additions.)

This is illustrated in ﬁgure below.

To ﬁnd impulse 1:

\[ \Delta _{1}V=2V_{leo} \sin \frac {15^{0}}{2}\]

hence, \(\Delta _{1}V=2 \left ( 7.84412 \right ) \sin \frac {15^{0}}{2}\ \)

\[ \Delta _{1}V=2.048 km/sec =7,371.81 km/hr \]

To ﬁnd impulse 2:

\[ \Delta _{2}V=\left \vert V_{l}eo-V_{p}\right \vert =\left \vert 7.844-10.3298\right \vert =2.4858 km/sec =8,948.88 km/hr \]

(speed up).

To ﬁnd impulse 3:

\[ \Delta _{3}V=\left \vert V_{g}eo-V_{a}\right \vert =\left \vert 3.072-1.584\right \vert = 1.488 km/sec = 5,356.8 km/hr \]

(speed up).

Hence, total impulses is found by summing the above

\[ \Delta V=2.048+2.4858+1.488=6.0218 km/sec =21,678.48 km/hr \]

Transfer from the initial LEO orbit (which has \(15^{0}\) degrees inclination) to a GEO orbit at \(35860\) km altitude (which still has a \(15^{0}\) degrees inclination). This is achieved using a normal Hohmann transfer (2 impulses). Next, perform an orbit plane inclination correction (one impulse) to move into the equatorial (\(0^{o}\) inclination) GEO circular orbit on which the ﬁrst target is currently orbiting.

This is illustrated in ﬁgure below.

To ﬁnd impulse 1: \[ \Delta _{1}V=V_{leo}-V_{p}=\left \vert 7.844-10.3298\right \vert =\ 2.4858\ \text {km/sec}=8,948.88\ \text {km/hr}\]

To ﬁnd impulse 2: \[ \Delta _{2}V=\left \vert V_{geo}-V_{a}\right \vert =\left \vert 3.072-1.584\right \vert =1.488\ \ \text {km/sec}=5,356.8\text {\ km/hr}\]

(speed up). To ﬁnd impulse 3 \[ \Delta _{3}V = 2V_{geo}\sin \frac {15^{0}}{2} = 0.80195\ \text {km/sec}=2,887.031\ \text {km/hr}\] Hence, total impulses is founding by summing the above total. \[ \Delta V=2.4858+1.488+0.80195=4.77575\ \text {km/sec}=17,192.7\ \text {km/hr}\]

In this scenario, we will apply a partial orbit plane correction at the perigee and the remaining correction at the apogee. See ﬁgure below.

At the perigee of the Hohmann transfer, apply the law of the cosines to obtain \begin {align*} \Delta V_{1} & =\sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }\\ \Delta V_{2} & =\sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }\\ \Delta V_{total} & =\Delta V_{1}+\Delta V_{2} \end {align*}

To ﬁnd the minimum \(\Delta V_{total}\) for a given \(\alpha \), take \(\alpha \) as the independent variable, and minimize \(\Delta V_{total}\) as a function of \(\alpha \). Hence solve \[ \frac {\partial \Delta V_{total}}{\partial \alpha }=0 \] Let\[ f\left ( \alpha \right ) =\Delta V_{total}=\sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }+ \sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }\] Hence for minimum \[ \frac {\partial f}{\partial \alpha }=\frac {2V_{leo}V_{p}\sin \left ( \alpha \right ) }{\sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }}-\frac {2V_{geo}V_{a}\sin \left ( \beta -\alpha \right ) }{\sqrt { V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }}=0 \] Hence \[ \frac {2V_{leo}V_{p}\sin \left ( \alpha \right ) \sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) } -2V_{geo}V_{a}\sin \left ( \beta -\alpha \right ) \sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }}{\sqrt { V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }\sqrt { V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }}=0 \] Then \[ V_{leo}V_{p}\sin \left ( \alpha \right ) \sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }-V_{geo}V_{a}\sin \left ( \beta -\alpha \right ) \sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }=0 \] Let \[ F\left ( \alpha \right ) =V_{leo}V_{p}\sin \left ( \alpha \right ) \sqrt { V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) } -V_{geo}V_{a}\sin \left ( \beta -\alpha \right ) \sqrt { V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }\] This is a non-linear equation in \(\alpha \). we solved for \(\alpha \) using Newton root ﬁnding method. Hence we need to ﬁnd \(F^{^{\prime }}(\alpha )\) as follows \begin {multline*} F^{^{\prime }}\left ( \alpha \right ) =V_{leo}V_{p}\cos \left ( \alpha \right ) \sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }\\ +V_{leo}V_{p}\sin \left ( \alpha \right ) \left ( \frac {-V_{geo}V_{a}\sin \left ( \beta -\alpha \right ) }{\sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }}\right ) \\ -\left [ V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) \sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }+V_{geo}V_{a}\sin \left ( \beta -\alpha \right ) \frac {V_{leo}V_{p}\sin \left ( \alpha \right ) }{\sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }}\right ] \end {multline*}

To solve for \(\alpha ,\) I wrote a MATLAB function that uses Newton root ﬁnding method to ﬁnd the root of \(F\left ( \alpha \right ) \) for a given \(r_{leo},r_{geo},\beta ,\mu \) where \(\beta \) is the total angle (in degrees) of the inclination of the ﬁrst circular orbit relative to the second circular orbit, and \(\mu \) is the gravitational constant.

The function returns back the angle \(\alpha \) for which \(\Delta V_{total}\) is minimum.

For an initial guess for \(\alpha \), and since \(\alpha \) is expected to be small compared to \(\beta \), I selected \(\alpha _{0}=0.1\ \beta \).

Applying Newton iterative root ﬁnding:

\[ \alpha _{i+1}=\alpha _{i}-\frac {F\left ( \alpha _{i}\right ) }{F^{^{\prime }}\left ( \alpha _{i}\right ) }\]

The tricky part in this problem was ﬁnding a good initial guess for the root (common problem with using Newton roots ﬁnding method).

I had to try diﬀerent values for an initial guess before the root was converged to. For example, when I selected \(\alpha _{0}\) to be \(50\%\) of \(\beta \), Newton method did not converge to the root. Selecting \(\alpha _{0}\) to be close to where one expects it to be (which is a small value compared to \(\beta \)) did work and a root was found.

For this design project, we are given that \(\beta =15^{0}.\) Using this matlab
function^{2}
I found that

\[ \alpha =1.28891^{0}\]

is the solution. Hence, this is the angle I will use for the correction to apply at the perigee.

So, at the perigee, Apply a correction of \(\alpha \) angle, and at the apogee, apply a correction of \(15-1.28891=\) \(13.711\).

To ﬁnd the impulse needed to correct inclination at the perigee and transfer from LEO to a Hohmann orbit with a correction of \(\alpha \)

\[ \Delta V_{1}=\sqrt {V_{leo}^{2}+V_{p}^{2}-2V_{leo}V_{p}\cos \left ( \alpha \right ) }=\sqrt {7.844\ ^{2}+10.3294^{2}-2\left ( 7.844\right ) \left ( 10.3294\right ) \cos \left ( 1.28891^{0}\right ) }\]

so, \(\Delta V_{1}=2.4936\ \)km/sec \(=8,977.085\ \)km/hr

To ﬁnd impulse needed to correct inclination at the apogee and transfer from Hohmann to the GEO orbit for an angle \(\left ( \beta -\alpha \right ) :\)

\[ \Delta V_{2}=\sqrt {V_{geo}^{2}+V_{a}^{2}-2V_{geo}V_{a}\cos \left ( \beta -\alpha \right ) }=\sqrt {3.072\ ^{2}+1.584^{2}-2\left ( 3.072\right ) \left ( 1.584\right ) \cos \left ( 13.71109^{0}\right ) }\]

so, \(\Delta V_{2}=1.578\) km/sec \(=5682.388\) km/hr

So, for this scenario, total \(\Delta V\) is given by \[ \Delta V_{1}+\Delta V_{2}=2.4936+1.578\ =4.0716\ \text {km/sec}=14657.76\, \text {km/hr}\]

Comparing the \(\Delta V\) from the above 3 scenarios, we see this:

scenario 1: \(6.0218\) km/sec

scenario 2: \(4.77575\) km/sec

scenario 3: \(4.0716\) km/sec

So, we can see that splitting the plane correction between the perigee and the apogee points leads to a more economical maneuver.

Hence Choose scenario 3 for the next sequence.

The vehicle was at the lines of node at \(t=0\), hence the time to reach GEO orbit is half the period of the Hohmann transfer orbit.

\[ t=\frac {T_{h}}{2}=18,916.77\, \text {sec}=5.25\, \text {hr}\]

The above time is the same regardless if we inject at one end of the lines of nodes or at the other end. Also, this time is independent of what inclination the Hohmann transfer orbit was at the time of injection.

In the previous step, the time it takes to move from LEO to GEO over a Hohmann orbit was found.

Now, ﬁnd when to make this transfer. That is, we need to ﬁnd the time to inject into the transfer orbit such that the overall \(\Delta V\) is minimized.

The injection must occur when the space vehicle is on the line of nodes. Since this line is where the LEO and the desired GEO plane intersects at.

Hence, there are only 2 points on the LEO orbit that we can use to launch to GEO. (Both ends of the lines of nodes, at both sides of earth). Let me call one end of the lines of nodes, the top end, and the other end, the bottom end. Where the top end is that end which the spacecraft was at when it ﬁrst reached LEO, i.e. at time=0.

In addition to the above restriction, if we want to eliminate the need to make any phasing loops when we arrive at the desired GEO orbit in order to rendezvous with the ﬁrst target, then injection must occur only when the correct lead angle \(\beta \) is encountered. This additional synchronization requirement will turn out to be costly in time to achieve. If LEO and GEO orbits had been coplanar, then we can inject from any point on the LEO orbit as long as the lead angle \(\beta \) requirement is met. There will not be an additional requirement of the injection having to be from only two points in the LEO orbit.

We know that the synodic period between LEO and GEO is \(\frac {T_{leo}T_{geo}}{T_{geo-T_{leo}}}=\frac {5189\times 86391}{86391-5189}=5520\) seconds, or 92 minutes. This means the LEO and the GEO objects will be aligned along a radial vector originating from the center of the earth once every 92 minutes.

But due to the restriction that this radial vector be only the lines of nodes of the space vehicle, using this synodic period is not too much help for me here.

So, what options do we have?

These are the options to investigate:

- Inject from the top end of the lines of nodes. Reach GEO orbit and then phase-wait in that orbit to rendezvous with the target.
- phase-wait in the LEO orbit until the correct lead angle \(\beta \ \)with the target is reached. Calculate this for when the spacecraft is on the top end of the lines of nodes.
- The same as above, but for the case when the spacecraft is on the bottom end of the lines of nodes.
- phase-wait in the LEO orbit until the lead angle is closest to \(\beta \) the ﬁrst time this happens. (i.e. within the ﬁrst \(2\pi \))

At the end, select the option which gives the smallest \(\Delta V\) as long as the time cost is reasonable.

we know that at \(t=0\) we have this state as shown in ﬁgure below

From previous calculations, we found the time needed for space vehicle to reach GEO orbit is \(18,916.77\)
sec. Hence, angle that the GEO satellite will travel in this time is found from

\begin {align*} \frac {18916}{\left ( 0.9972696\right ) \left ( 24\right ) \left ( 60\right ) \left ( 60\right ) } & =\frac {x}{2\pi }\\ x & =1.3793765886\ \ \text {rad}\ \\ & =79.03\ ^{0} \end {align*}

Hence When Vehicle reaches GEO, the ﬁrst satellite will be (see diagram)

\[ 180-\left ( 79.03-40\right ) \ =140.9675^{0}\]

behind the space vehicle.

Hence\[ \Delta L=-140.9675^{0}\]

The reason a minus sign is used, is because \(\Delta L\) is measured positive if the change in longitude desired is eastwards, and since our target is behind us (westwards relative to the spacecraft), this change is negative.

\(-140.9675^{0}\) \(=-140.9675\left ( \frac {\pi }{180}\right ) =-2.46\) radians.

For \(n=1\), \(P_{ph}=\frac {-2.46}{1\times \ 6.3}+0.9972696=0.60676\) days.

Hence \(\overset {\circ }{L}=\frac {\Delta L}{n\ P_{ph}}=\frac {-2.46}{1\left ( 0.60676\right ) }=\) \(-4.054\ \)radians/day \(=-4.054\left ( \frac {180}{\pi }\right ) =\) \(-232.2771\) deg/day.

Hence \(\Delta V=5.8\left ( -232.2771\right ) =\) \(-1347.20718\) m/sec \(=-1.347\) km/sec.

For \(n=2\), \(P_{ph}=\frac {-2.46}{2\times \ 6.3}+0.9972696=0.802\ \) days.

\(\overset {\circ }{L} = \frac {-2.46}{2 (0.802)} = -1.5337\) rad/day or \(-87.875\) deg/day

Hence \(\Delta V=5.8\left ( -87.875\right ) =\) \(-509.675\) m/sec \(=-0.5\) km/sec.

Continuing the above process, I obtain this table.

\(n\) | \(\Delta V\)(Km/s) | Period of phase orbit (in days) | Total time in phasing period in days |

\(1\) | \(-1.347\) | \(0.6067934\) | \(0.6067934\) |

\(2\) | \(-0.509\) | \(0.80203\) | \(1.60406\) |

\(3\) | \(-0.314\) | \(0.86711\) | \(2.60133\) |

\(4\) | \(-0.2272\) | \(0.8996\) | \(3.5984\) |

\(5\) | \(-0.17788\) | \(0.91917436\) | \(4.5958718\) |

\(6\) | \(-0.1462\) | \(0.93219\) | \(5.59314\) |

\(12\) | \(-0.07062\) | \(0.96472\) | \(11.57664\) |

\(24\) | \(-0.0347\) | \(0.980999\) | \(23.543976\) |

\(96\) | \(-0.0085\) | \(0.993202\) | \(95.347392\) |

Any one of the above choices for \(n\) will achieve rendezvous with the ﬁrst satellite.

Plotting \(\Delta V\) against the total time in the phasing orbit results in the plot shown in ﬁgure below

Figure above shows that most saving in fuel is made by staying in the phasing orbit for less than 20 days. For more than 20 days, the additional saving in fuel is not justiﬁed by the time wasted in the phasing orbit.

Zooming in the region of interest in the plot shown in above ﬁgure results in ﬁgure below

From the above, it is clear that the slope after 6 days in the phasing orbit is less steep than earlier. The two options I see is to choose \(n=6\) and save some fuel, or choose \(n=1\) (the smallest possible value) and save time.

For \(n=6\) \[ \text {time}=5.59314\ \text {days}\ \ \ \Delta V=-0.1462\ \text {km/sec.}\]

For \(n=1\)\[ \text {time}=0.6067934\ \text {days}\ \ \ \Delta V=-1.347\text {\ km/sec}. \]

Since the fuel saved is so small compared to the initial fuel needed to send the vehicle into GEO orbit, I decided to use the smaller time here.

Looking at above ﬁgure, the problem can be seen as the following: we need to ﬁnd the time it takes for the spacecraft to be at the top end of the lines of nodes when the target is at the correct lead angle \(\beta \). Because in this case, the spacecraft can injects into the Hohmann orbit and will meet the target at the apogee. This would result in the spacecraft not having to do any phase-waiting loops in the GEO orbit. Hence saving \(\Delta V\). A trade is made between time and \(\Delta V.\)

How do we ﬁnd this time value?

First, ﬁnd \(\beta \). To do this, equate travel time for target and spacecraft.

Travel time for spacecraft is half the Hohmann orbit period = \(37,833\ \)sec (see common calculations section for derivation)

Travel time \(x\) for target is found from \begin {align*} \frac {\text {distance travelled in radians}}{2\pi } & =\frac {x}{\text {Period}}\\ x & =\frac {P\left ( \theta -\beta \right ) }{2\pi } \end {align*}

Where period \begin {align*} P & =2\pi \sqrt {\frac {r_{geo}^{3}}{\mu }}\\ & =2\pi \sqrt {\frac {42,238.145^{3}}{3.986012\times 10^{5}}}\\ & =86390\text { sec} \end {align*}

Hence, and since \(\theta =\pi \), when equating times of travel, we get this relation \(\frac {86390 (\pi -\beta )}{2\pi } = 37833 \)

Solve for \(\beta \)

\begin {align*} 37833\left ( 2\pi \right ) & =86390\left ( \pi -\beta \right ) \\ \frac {37833\left ( 2\pi \right ) }{86390} & =\pi -\beta \\ \beta & =\pi -\frac {37833\left ( 2\pi \right ) }{86390}\\ & =0.389980\ \text {rad}\\ & =22.344^{0} \end {align*}

So, now that we know \(\beta \), we need to ﬁnd when will the spacecraft be at the top end of the lines of nodes when the target has this \(\beta \) with it.

The angle that the target will move by for each one full period that the spacecraft makes in the LEO orbit is \(T_{leo}\varpi \), where \(T_{leo}\) is the LEO period found in the common calculations section to be \(5,189\) seconds, and \(\varpi \) is the average angular speed of the GEO target in radians per second which is \(2\pi \) each \(24\) hrs (since GEO).

Exactly, \(\varpi =\frac {2\pi }{0.9972696\ast 24\ast 60\ast 60}=0.00007292115609102 \) rad per second.

Hence, target will travel \(5,189\times 0.00007292115609102=\) \(0.378388\) rad \(=21.68^{0}\) for each one full LEO period.

Now that we know the angle the target will travel for each one full LEO period, we need to ﬁnd how many times we have to do this so that target will end at the correct \(\beta \ \)location starting from t=0.

This becomes a simple counting problem. If we imagine a straight line, starting at t=0, and then we move a stick from its left end to its right end at an equal increments of 21.68 units, we just need to ﬁnd when this stick will land at the correct point (or close enough) on the line where the point \(\beta \) is located.

When the stick reaches the right end of the line, we carry the reminder back to the start of the line and continue the process.

Figure below shows how to do this counting. Notice that the point I am interested in ﬁnding, which is the angle \(\beta ,\) needs to be compensated for by adding the initial \(40^{0}\) to it (since the counting is starting from the epoch). In other words, this is the degrees the spacecraft was ahead of the target when counting starts.

So, the point on the line we are looking for is \(40+\beta =40+22.344=60.344^{0}\), measured from t=0.

Figure below shows how the counting is actually done.

I wrote a small MATLAB^{3} function
called nma_ﬁndPointOnLine.m to do the counting. This function accepts as input the step size,
the length of the line, and the distance we are looking for (60.344 in the above example), and how
close to the target we want to be (the tolerance).

It returns the number of steps needed to achieve the synchronization.

This is an example call

>> stepSize=21.68; lineLength=360; target=60.344; tol=0.1; >> nSteps=nma_findPointOnLine(stepSize,lineLength,target,tol) nSteps = 3274

This is the result of the MATLAB function

tolerance (in degrees) | number of LEO periods needed | time to achieve (seconds) | time in hrs |

1 | 36 | \(5189\times 36=\ \)186804 | 51.89 |

0.1 | 3274 | \(5189\times 3274=16988786\) | 4719 |

It is clear that to achieve synchronization to 0.1 degree is too costly in time.

For the case of 1 degree tolerance, it will take 36 LEO loops to achieve the optimal situation with the target at the correct \(\beta .\)

This means, if we spend this time in LEO, we can inject and will meet the target at the same time when we reach GEO. Hence no additional \(\Delta V\) would be needed in GEO to phase-wait. We have traded time for fuel.

This case is the same as above, except now we want the spacecraft to be at the other end of the lines of nodes at injection. The only diﬀerence is that now t=0 have been shifted to become time after making one half LEO period. We can ﬁnd this shift since we know the angle the target will travel in one half LEO period. We found from above that target will travel \(21.68^{0}\) for one full LEO period, hence it will travel \(10.84^{0}\) per half that period.

So, the only thing we need to do is determine where the point is located that we need to synchronize with, as illustrated by the ﬁgure below.

So, use the same MATLAB function to ﬁnd the number of full LEO rotations needed.

>> stepSize=21.68; lineLength=360; target= 231.544; tol=0.1; >> nSteps=nma_findPointOnLine(stepSize,lineLength,target,tol) nSteps = 3614 >> nSteps=nma_findPointOnLine(stepSize,lineLength,target,1) nSteps = 293

Hence, we see that for a tolerance of 1 degree, we have to wait 293 full LEO loops.

Compare this with the earlier case where we looked for the other end of the lines of nodes, which achieved the same synchronization for only 36 LEO loops. Hence this maneuver will not be accepted.

To solve this problem, I will calculate the lead angle \(\beta \ \)with the ﬁrst target for a number of time increments of \(0.5P\) each, where \(P\) is the LEO period, and for each such time increment, will calculate where the ﬁrst target will be at the end of a Hohmann transfer. Then will calculate the \(\Delta V\) needed to phase-wait in GEO to close this ﬁnal angle gap.

At \(t=0.5T_{leo}\)

Figure below illustrates this case.

Total time from epoch to reach GEO \(=0.5T_{leo}+0.5T_{h}=0.5\left ( 5189+37833\right ) =\) \(21,511\sec \)

Hence, angle that the GEO satellite will travel in this time is found from

\begin {align*} \frac {21511}{\left ( 0.9972696\right ) \left ( 24\right ) \left ( 60\right ) \left ( 60\right ) } & =\frac {x}{2\pi }\\ x & =1.568606\ \text {rad}\ \\ & =89.874\ ^{0}. \end {align*}

Hence When Vehicle reaches GEO, the ﬁrst satellite will be

\[ 89.874-40\ =49.874^{0}\] ahead.

At \(t=T_{leo}\)

Total time from epoch to reach GEO \(=T_{leo}+0.5T_{h}=5189+0.5\left ( 37833\right ) =\) \(24,105.5\sec \)

Hence, angle that the GEO satellite will travel in this time is found from

\begin {align*} \frac {24105.5}{\left ( 0.9972696\right ) \left ( 24\right ) \left ( 60\right ) \left ( 60\right ) } & =\frac {x}{2\pi }\\ x & =1.7578\ \text {rad}\\ \ & =100.71^{0}\ \end {align*}

Hence When Vehicle reaches GEO, the ﬁrst satellite will be\[ 180-\left ( 100.71-40\right ) \ =119.29^{0}\]

behind

At \(t=1.5T_{leo}\)

Total time from epoch to reach GEO \(=1.5T_{leo}+0.5T_{h}=1.5\left ( 5189\right ) +0.5\left ( 37833\right ) =\) \(26,700\sec \)

Hence, angle that the GEO satellite will travel in this time is found from

\begin {align*} \frac {26700}{0.9972696\ast 24\ast 60\ast 60} & =\frac {x}{2\pi }\\ x & =1.946995\ \text {rad}\ \\ & =111.5546\ ^{0} \end {align*}

Hence When Vehicle reaches GEO, the ﬁrst satellite will be

\[ 111.555-40\ =71.555^{0}\] ahead.

Continue this process. The result is illustrated in this table

time | total time to reach GEO (sec) | \(\phi \) target travelled angle | \(\beta \ \)lead angle | \(\beta \ \)sense relative to spacecraft |

0 | \(18,916\ \) | \(79.03^{0}\) | \(140.9675^{0}\) | behind |

\(0.5T_{leo}\) | \(21,511\) | \(89.874\ ^{0}\) | \(49.874^{0}\) | ahead |

\(T_{leo}\) | \(24,105\ \) | \(100.71^{0}\) | \(119.29^{0}\) | behind |

\(1.5T_{leo}\) | \(26,700\) | \(111.5546\ ^{0}\) | \(71.555^{0}\) | ahead |

\(2T_{leo}\) | \(29,294\ \) | \(122.3946^{0}\) | \(97.605^{0}\) | behind |

\(2.5T_{leo}\) | \(31,889\) | \(133.2346^{0}\) | \(93.235^{0}\) | ahead |

\(3T_{leo}\) | \(34,483\ \) | \(144.0746^{0}\) | \(75.9254^{0}\) | behind |

\(3.5T_{leo}\) | \(37,078\) | \(154.9146^{0}\) | \(114.915^{0}\) | ahead |

\(4T_{leo}\) | \(39,672\ \) | \(165.7547^{0}\) | \(54.2453^{0}\) | behind |

\(4.5T_{leo}\) | \(42,267\) | \(176.5947^{0}\) | \(136.595^{0}\) | ahead |

\(5T_{leo}\) | \(44,861\ \) | \(187.4347^{0}\) | \(32.5653^{0}\) | behind |

\(5.5T_{leo}\) | \(47,456\) | \(198.2747^{0}\) | \(158.275^{0}\) | ahead |

\(6T_{leo}\) | \(50,050\) | \(209.1147^{0}\) | \(10.8853^{0}\) | behind |

\(6.5T_{leo}\) | \(52,645\) | \(219.9547^{0}\) | \(179.957^{0}\) | ahead |

\(7T_{leo}\) | \(55,239\) | \(230.7947^{0}\) | \(-10.7947^{0}\) | behind |

From the above table, we see that the closest the target gets to the lines of nodes (within the ﬁrst \(2\pi \)) at the same time as the vehicle is \(10.8853^{0}\) and is achieved after \(6\) periods in the LEO orbit.

What is left to do is to determine is the \(\Delta V\) needed in the GEO orbit to wait-phase so as to close this ﬁnal remaining \(\Delta L=-10.8853^{0}\)

Now, for the \(-10.8853^{0}\ \)case calculate the cost \(\Delta V\) for phase-waiting in the GEO orbit:

Let \(P_{ph}\) be the period of the phasing orbit while in the GEO orbit.

Let \(\omega _{E}\) be the angular rate of axial rotation of the earth, which is \(360.985647\) deg/day \(=360.985647\ \frac {\pi }{180}=6.3\ \)radians/day.

Let \(P_{o}\) be the period of the geosynchronous orbit, which is \(1436.068\) min or \(0.9972696\) days.

Let \(\Delta L\) be the change in longitude desired. Which we found it to be \(-10^{0}\) \(=-10\left ( \frac {\pi }{180}\right ) =-0.174533\) radians.

Let \(n\) be the number of revolutions spent in the phasing orbit.

Let \(\overset {\circ }{L}\) be the drift rate, positive eastwards.

First, ﬁnd the period of the phasing orbit \(P_{ph}\) using the equation

\[ P_{ph}=\frac {\Delta L}{n\ \omega _{E}}+P_{o}\]

for a speciﬁc \(n.\) Then solve for \(\overset {\circ }{L}\) from the equation

\[ \Delta L=\overset {\circ }{L}n\ P_{ph}\]

Then ﬁnd \(\Delta V\) corresponding to this \(\overset {\circ }{L}\) from ﬁgure 7.14, Orbital Mechanics book or by using the relation

\[ \Delta V=5.8\overset {\circ }{L}. \]

Where 5.8 is the slope of the line relating \(\overset {\circ }{L}\ \)to \(\Delta V\). (The above slope is not exact, but it is close enough).

Try the above for a number of diﬀerent values for \(n\).

For \(n=1\), \(P_{ph}=\frac {-0.174533}{1\times \ 6.3}+0.9972696=\) \(0.96957\) days.

Hence \(\overset {\circ }{L}=\frac {\Delta L}{n\ P_{ph}}=\frac {-0.174533}{1\left ( 0.96957\right ) }=\) \(-0.18\ \ \)radians/day \(=-0.18\left ( \frac {180}{\pi }\right ) =\) \(-10.31385\ \) deg/day.

Hence \[ \Delta V=5.8\left ( -10.31385\right ) =-59.82033\ \text {m/sec}=-0.0598\ \text {km/sec}. \]

I do not need to look for n=2 and higher for this case, since the saving in \(\Delta V\) is clearly not worth spending an extra day for each additional increment in \(n\). We see that the \(\Delta V\) is very small already (0.0598 km/sec), and for larger n, it will only get smaller and smaller.

For ﬁrst option, move to GEO without phase-waiting in LEO, and instead phase-wait in GEO, results in \(\Delta V=-1.347\) with time spend \((0.6067934)(24)(60)(60) = 52426\) seconds, or \(14\) hrs and \(33\) minutes.

For the second and third options, the time spend is all in LEO, and zero \(\Delta V\) was needed to phase-wait in GEO.

For 4th option, part of time spend is in LEO, and some part of time spend is in GEO. \(\Delta V\) is not zero, but smaller than ﬁrst option.

Notice that in this table, the \(\Delta V\) cost refer only to the cost of phase wait in GEO to rendezvous with the ﬁrst target.

scenario | \(\Delta V\) cost | LEO loop | total time cost to rendezvous |

ﬁrst. Inject at t=0 | \(1.347\) km/sec | 0 | 52,426 (sec)\(\ \)= 14 hrs 33 minutes |

second. optimal \(\beta \) at top end | 0 km/sec | 36 | 186,804 (sec) = 51.89 hrs = 2 days 3 hrs 53 minutes |

third. optimal \(\beta \) at bottom end | 0 km/sec | 293 | 1,520,377 (sec) |

fourth. smallest \(\beta \) at apogee | \(0.0598\) km/sec | 6 | 133,820 (sec) = 1 day, 13 hrs, 10 minutes, 20 seconds. |

Clearly option 3 is bad. With option 2, we get the same saving for much less time.

Between options 1 and 4, I prefer option 4, since for the cost of about 1.5 days, we reduced \(\Delta V\) from \(1.347\) to \(0.0598\ \)km/sec.

So, the ﬁnal choice is between option 4 and option 2.

With option 2, we have zero \(\Delta V\) but we have to spend about 15 more hours in LEO to save \(0.0598\ \)km/sec. Is this good or not?

Compared to the \(\Delta V\) needed to inject from LEO to GEO which is \(4.0716 \) km/sec, this amount is 1.5% of that. It takes about 5 hrs to go from LEO to GEO. So, should I spend about 3 times as many hours to save 1.5% as many in fuel? (This is interesting that spending more time ﬂying ends up saving fuel! only in space this is possible).

Will not consider option 2 as the time needed is not worth the saving in fuel.

Option 4 is selected

Simply perform an in-orbit repositioning using a phasing orbit to rendezvous with the second satellite. The second satellite is \(50^{0}\) ahead of the ﬁrst satellite (this is given), hence \(\Delta L=50^{0}=0.87267\) rad. Apply the same process of in-orbit repositioning to decide on the procedure to select. This table was generated:

\(n\) | \(\Delta V\)(Km/s) | Period of phase orbit (in days) | Total time in phasing period in days |

\(1\) | \(0.255\) | \(1.13578\) | \(1.13578\) |

\(2\) | \(0.13596\) | \(1.066529\) | \(2.133058\) |

\(3\) | \(0.093\) | \(1.04344\) | \(3.13032\) |

\(4\) | \(0.0703\) | \(1.0318993\) | \(4.1275972\) |

\(5\) | \(0.0566\) | \(1.0249734\) | \(5.124867\) |

\(6\) | \(0.04746\) | \(1.020356\) | \(6.122136\) |

\(12\) | \(0.02396\) | \(1.008812\) | \(12.105744\) |

\(24\) | \(0.012\) | \(1.00304122\) | \(24.0729\) |

\(96\) | \(0.003\) | \(0.998712\) | \(95.87635\) |

Plotting \(\Delta V\) against the total time in the orbit results in

Similarly, at \(n=6\), the fuel saving is best for the time spend in the phasing orbit. This gives \(6.122136\) days in the phasing orbit, and \[ \Delta V=0.04746\ \text {km/sec}. \]

However, since the fuel saved as a percentage of initial fuel is small, it does not seem that spending extra 6 days in orbit is worth it. I will use \(n=1\)\[ \text {Time}=1.3578\text { days \ }\Delta V=0.255\ \text {km/sec}. \]

Also, since from the speciﬁcations, it seems that it says that the second satellite needs to be reached quickly, so I choose \(n=1\), the smallest possible value.

Stay in orbit for one complete orbit revolution. This adds one day to the total time in ﬂight. No additional \(\Delta V\) needed.

Perform an in-orbit repositioning using a phasing orbit to reposition \(\Delta L=5^{0}=0.087266\) rad ahead.

Apply the same process of in-orbit repositioning to decide on the procedure to select. This table was generated:

\(n\) | \(\Delta V\)(Km/s) | Period of phase orbit (in days) | Total time in phasing period in days |

\(1\) | \(0.02868\) | \(1.0111213\) | \(1.0111213\) |

\(2\) | \(0.0144\) | \(1.0041954\) | \(2.0083908\) |

\(3\) | \(0.00965\) | \(1.00188684\) | \(3.00566052\) |

\(4\) | \(0.007\) | \(1.00073253\) | \(4.00293012\) |

\(5\) | \(0.0058\) | \(1.00003994\) | \(5.0001997\) |

\(6\) | \(0.00484\) | \(0.9995782\) | \(5.9974692\) |

\(12\) | \(0.00242\) | \(0.9984239\) | \(11.9810868\) |

\(24\) | \(0.001211\) | \(0.99784675\) | \(23.948322\) |

\(96\) | \(0.0003\) | \(0.9974138\) | \(95.7517248\) |

As before, choose \(n=1\) with \[ \text {Time}=1.0111213\ \text {day\ }\ \Delta V=0.02868\ \text {km/sec.}\] This adds little over one day to the total ﬂight time.

This completes the required sequence.

See ﬁgure below for the ﬁnal decision tree.

This is a summary, in table format, of the entire orbiting sequences using the most \(\Delta V\) optimal maneuvers selected out of the diﬀerent scenarios shown above. This table below is the ﬁnal result of selecting the best option out of each phase. This table only shows the \(\Delta V\) needed starting from LEO orbit. It does not include the \(\Delta V\) needed to reach LEO which I will add next.

sequence: Inject to LEO. wait in LEO for 6 full orbits. Perform Hohmann transfer with partial plane correction at each end. phase-wait in GEO to rendezvous with ﬁrst target. Phase-wait in GEO to rendezvous with second target. Lock into second target. Reposition to ﬁnal location in GEO.

sequence | \(\Delta V\)(km/sec) | duration |

wait in LEO | 0 | 50,050 (sec) |

Hohmann transfer | 4.0716 | 18,916 (sec) |

Rendezvous with 1st target | 0.0598 | 83,770 (sec) |

Rendezvous with 2nd target | 0.255 | 98,131 (sec) |

lock with 2nd target | 0 | 86,400 (sec) |

position to ﬁnal destination | 0.02868 | 87,360 (sec) |

TOTAL | 4.41508 | 424,627 (sec)=117.95 hrs = 4 days 21hrs 57 min |

To reach LEO, we have found that \(V_{leo}=7.844\) km/sec. To be more realistic, we need to account for the gravitational loss and drag. Typical time to reach LEO is about 2 minutes or 120 seconds. Hence additional \(\Delta V\) to account for gravity loss is \(g\Delta t=9.8\left ( 120\right ) =1176\) m/sec \(=1.176\) km/sec.

For the drag eﬀect, it of course depends on the rocket cross sectional area, the drag coeﬃcient and air density. A typical value I have seen in the literature for rockets is to use \(5\%\) of the LEO velocity to account for drag. Hence an additional \(\Delta V\) for drag will be \(5\%\) of \(7.844\) km/sec or \(0.392\) km/sec.

Hence total \(\Delta V\) to reach LEO = \(7.844+1.176+0.392=9.412\) km/sec

Hence the total \(\Delta V\) to achieve ﬁnal position of spacecraft is

\[ 4.41508+9.412=13.827\ \text {km/sec}\]

see ﬁgure below for a graphical display of the deltaV used at each stage.

This is the MATLAB function that solves for \(\alpha \) to ﬁnd what partial correction in inclination angle can be done at the perigee for the Hohmann transfer.

Caller script:

function nma_testfindAlphaForMinDeltaV r0 = 6378.145; r1 = 100+r0; r2 = 35860+r0; beta = 15; mu = 3.986012*10^5; alpha = nma_findAlphaForMinDeltaV(r1,r2,beta,mu)

This is the function that solves for alpha.

function alpha=nma_findAlphaForMinDeltaV(r1,r2,beta,mu) %function alpha=nma_findAlphaForMinDeltaV(r1,r2,beta,mu) % % Finds the minimum alpha (initial inclination correction) % for an orbit relative to a larger circular orbit. % see design note for more details. % %INPUT: % r1: The radius of the smaller circular orbit % r2: the radius of the larger circular orbit % beta: the angle, in degrees, in which the two % circular orbits are non co-planers to % each others. % mu: gravitational constant % % OUTPUT: % alpha: The angle in degrees to use for initial % correction such that minimum delta V is % obtained to move from the smaller circular % orbit to the larger circular orbit % % Author: Nasser Abbasi % May 19,2003 a = (r1+r2)/2; rp = r1; ra = r2; Va=sqrt( mu*(2/ra - 1/a )); Vp=sqrt( mu*(2/rp - 1/a )); Vc1=sqrt( mu/r1 ); Vc2=sqrt( mu/r2 ); beta=beta*pi/180; \end{Verbatim} \newpage\begin{Verbatim} root(1)=0.1*beta; keepLooking = true; i=0; while(keepLooking) i=i+1; root(i+1)=root(i)- ( F(Vc1,Vc2,Vp,Va,root(i),beta)/dF(Vc1,Vc2,Vp,Va,root(i),beta) ); root(i+1) if( abs ( (root(i+1) - root(i)) / root(i+1) ) * 100 < 0.001 ) keepLooking=false; else if( ( root(i+1) * root(i) )<0.0) error('jumped out of root'); end end if(i>100) error('Failed to converge'); end end alpha=root(end); alpha=alpha*180/pi; %%%%%%%%%%%%%%% % % %%%%%%%%%%%%%%%% function v=F(Vc1,Vc2,Vp,Va,alpha,beta) v=Vc1*Vp*sin(alpha)*sqrt(Vc2^2+Va^2-2*Vc2*Va*cos(beta-alpha)) ... - Vc2*Va*sin(beta-alpha)*sqrt(Vc1^2+Vp^2-2*Vc1*Vp*cos(alpha)); %%%%%%%%%%%%%%% % % %%%%%%%%%%%%%%%% function v=dF(Vc1,Vc2,Vp,Va,alpha,beta) v=Vc1*Vp*cos(alpha)*sqrt(Vc2^2+Va^2-2*Vc2*Va*cos(beta-alpha)) ... +Vc1*Vp*sin(alpha)* ... ( -Vc2*Va*sin(beta-alpha)/sqrt(Vc2^2+Va^2-2*Vc2*Va*cos(beta-alpha))) ... - (Vc2*Va*cos(beta-alpha)*sqrt(Vc1^2+Vp^2-2*Vc1*Vp*cos(alpha)) ... + Vc2*Va*sin(beta-alpha)*... ( Vc1*Vp*sin(alpha) / sqrt(Vc1^2+Vp^2-2*Vc1*Vp*cos(alpha)))); \end{Verbatim} \newpage\begin{Verbatim} function nSteps=nma_findPointOnLine(stepSize,lineLength,target,tol) %function nSteps=nma_findPointOnLine(stepSize,lineLength,target,tol) % % Function to find how many steps needed to reach withing tolearance % close to a point on a line by taking fixed number of steps. Line wrapes % around. % %INPUT: % stepSize: the step size % lineLength: The line length % target: The distance from leftend of line that we want to reach % tol: tolerance in abs. value % %OUTPUT: % nSteps: number of steps needed % %Author Nasser Abbasi. May 22, 2003 % currentDist = 0; nSteps = 0; while true currentDist = currentDist+stepSize; nSteps = 1; while(currentDist <= lineLength) if( abs (currentDist-target) < tol ) return; end rem = lineLength-currentDist; if(rem < stepSize) currentDist = stepSize-rem; else currentDist = currentDist+stepSize; end nSteps = nSteps+1; end end