attachment:orb01.c of OrbitCirc

   1 // orb01.c  - compute insertion delta V from launch loop, space elevator
   2 
   3 
   4 #define MUKM   398600.4418   // km3/s2
   5 #define RE     6378.0        // km  earth equatorial radius
   6 #define LOOP   80.0          // km  loop altitude
   7 #define PI2    6.28318530718 
   8 #define SDAT   86164.0905    // sidereal day seconds
   9 
  10 #define RB      6600.0       // km  orbit sweep bottom
  11 #define RT     50000.0       // km  orbit sweep top
  12 #define RSTEP    100.0       // km  orbit sweep step
  13 
  14 #define RAB    29700.0       // km  lowest ra search 
  15 #define RAT    53000.0       // km  highest ra search
  16 
  17 #define NLINES    50         // number of lines to print
  18 
  19 #include <stdio.h>
  20 #include <math.h>
  21 
  22 int main() {
  23    double  omegae = PI2 / SDAT  ;  // earth rotation angular velocity
  24    double  rloop  = RE + LOOP   ;  // km launch loop radius
  25    double  rgeo   = pow( (MUKM/(omegae*omegae)), (1.0/3.0)) ;
  26    double  muom   = 2.0 * MUKM / ( omegae * omegae ) ; 
  27    int     iter ;
  28 
  29    double  rdest ;
  30    int lines = 0 ;
  31 
  32    for( rdest = RB ; rdest <= RT ; rdest += RSTEP ) {
  33       if( lines-- == 0 ) {
  34          printf( "#\n# KMradius  Eradius   DVloop   DVelev   RAelev    Vdest\n" ) ;
  35          lines = NLINES ;
  36       }
  37       double vdest = sqrt( MUKM / rdest ) ;    // destination orbit velocity
  38       double rre   = rdest/RE             ;    // destination earth units
  39 
  40       // launch loop transfer orbit
  41       double lva   = sqrt( 2.0*MUKM / ( rdest * ( 1.0 + rdest/rloop ) ) ) ;
  42       double ldv   = fabs( vdest - lva );
  43 
  44       // space elevator transfer orbit
  45       // set up binary search
  46       double rab   = RAB ;
  47       double rat   = RAT ;
  48       double sra         ;
  49       double srp         ;
  50 
  51       for( iter = 20 ; iter-- > 0 ; ) {
  52          sra = 0.5 * (rab + rat) ;
  53          srp = sra / ( ( muom / ( sra*sra*sra ) ) - 1.0 ) ;
  54          if ( srp < rdest ) { rab = sra ; }
  55          else               { rat = sra ; }
  56       }
  57       double svd   = sra*sra*omegae / rdest ;
  58       double sdv   = fabs( vdest - svd )    ;
  59 
  60       printf( "%10.0f%9.5f%9.5f%9.5f%9.2f%9.5f\n",
  61                rdest, rre, ldv, sdv, sra, vdest );
  62    } 
  63    return (0);
  64 }