*23456789a123456789b123456789c123456789d123456789e123456789f123456789g12
*
* 
      program math
      include '/homes/jeffery/jef/syn/const.f'
      dimension echo(4,4),asemi(7)
*
      open(unit=1,file='/homes/jeffery/jef/syn/const.dat',
     &     status='old')
      read(1,const)  ! const is cgs recall
      read(1,astro)
      close(unit=1)
*
      pi=acos(-1.)
      raddeg=180./pi
*     

      epsilon=8.854187817e-12
      xkcon=1./(4.*pi*epsilon)
*
      print*
      q=2.e-6
      xmass=2.e-5
      ri=3.
      vinf=sqrt(2.*xkcon*q**2/(ri*xmass))
      print*,'Velocity at infinity is ',vinf !=34.6170502
*
      print*
      amu=1.66e-27
      xmpenny=2.5e-3
      ampenny=63.5
      acount=xmpenny/(amu*ampenny)
      echarge=1.602e-19
      anpenny=29.
      aproton=anpenny*acount
      charge=aproton*echarge
      qnet=charge*1.e-6
      xkcon=8.99e+9
      force=xkcon*qnet**2/1.**2
      print*,'acount,aproton,charge,qnet,force'
      print*,acount,aproton,charge,qnet,force
* 2.37169149E+22  6.87790528E+23  110184.047  0.110184044  109143304.
*
      print*
      print*,.05/(1.+1836.)
*             2.72182915E-05
*
      print*
      bohr=0.529e-10
      xkcon=8.99e+9
      echarge=1.602e-19
      force=xkcon*echarge**2/bohr**2
      emass=9.11e-31
      acc=force/emass
      vel=sqrt(acc*bohr)
      clight=2.998e8
      beta=vel/clight
      print*,'force,acc,vel,beta'
      print*,force,acc,vel,beta
* 8.2446725E-08  9.05013387E+22  2188040.5  0.00729833404
*
      print*
      volt=1.e+5
      charge=2*1.602e-19
      bfield=.1
      r=1.5/2.
      amu=1.66e-27
      xmass=bfield**2*r**2*charge/(2.*volt)
      print*,'xmass,xmass/amu'
      print*,xmass,xmass/amu
* 9.01125056E-27  5.42846394
*
      print*
      pi=acos(-1.)
      xmu=8.e+22
      rad=3500.e+3
      xi=xmu/(pi*rad**2)
      print*,'xi'
      print*,xi
* 2.0787584E+09
* 
      print*
      p=1250.
      v=115.
      xi=p/v
      r=v/xi
      en=p*3600.
      enkwhr=(p*1.e-3)*1.
      eninf=(v**2/r)*3600.
      print*,'xi,r,en,enkwhr,eninf'
      print*,xi,r,en,enkwhr,eninf
* 10.869565  10.5799999  4500000.  1.25  4500000.
*
      print*
      pi=acos(-1.)
      echarge=1.602e-19
      enkev=12.0
      emass=9.11e-31
      bfield=55.e-6
      en=enkev*1.e+3*echarge
      v=sqrt(2.*en/emass)
      rad=(emass*v)/(echarge*bfield)
      dia=2.*rad
      period=(2.*pi*rad/v)
      print*,'en,v,rad,dia,period'
      print*,en,v,rad,dia,period
* 1.92240014E-15  64964740.  6.71693087  13.4338617  6.49640469E-07
*
      print*
      echarge=1.602e-19
      sig=2.70e-14
      efield=120.
      xnp=620.e+6
      xnm=550.e+6
      xj=sig*efield
      vd=xj/( (xnp+xnm)*echarge )
      print*,'xj,vd'
      print*,xj,vd
*     3.24000007E-12  0.0172860846
*
      print*
      time=1.30e-6
      vend=12.
      v=5.
      r=15.e+3
      tau=time/log(vend/(vend-v))
      c=tau/r
      print*,'tau,c'
      print*,tau,c
* 2.41188945E-06  1.60792629E-10
*
      print*
      echarge=1.602e-19
      emass=9.11e-31
      acc=2.e+12
      vy=5.e+3
      vz=5.e+3
      bfield=400.e-6
      efield=emass*acc/echarge
      speed=sqrt(vy**2+vz**2)
      rad=emass*speed/(echarge*bfield)
      print*,'efield,speed,rad'
      print*,efield,speed,rad
* 11.3732834  7071.06787  0.00010052658
*
      print*
      pi=acos(-1.)
      xturns=200.
      radius=.5
      bfield=4.
      freq=1000./60.
      omega=2.*pi*freq
      emfamp=omega*xturns*pi*radius**2*bfield
      rr=100.e+3
      powermax=emfamp**2/rr
      powerave=powermax*.5
      emfampearth=emfamp*1.e-4/4.
      print*,'emfamp,powermax,powerave,emfampearth'
      print*,emfamp,powermax,powerave,emfampearth
*   65797.3672  43292.9336  21646.4668  1.64493418
*
      print*
      pi=acos(-1.)
      perm=4.*pi*1.e-7
      xL=1.e+6
      vol=125.
      xi=100
      en=.5*xl*xi**2/vol
      turnden=sqrt(xL/(perm*vol))
      bfield=perm*turnden*xi
      print*,'en,turnden,bfield'
      print*,en,turnden,bfield
*  40000000.  79788.4531  10.0265131
*
      print*
      clight=2.99792458e+8     ! in mks
      freq=clight/(1.e+4*6.37e+6)
      print*,'freq,1./freq'
      print*,freq,1./freq
*  0.00470631802  212.480331
*
      print*
      clight=2.99792458e+8     ! in mks
      xlam=632.8e-9
      dxlam=.01e-9
      freq=clight/xlam
      dfreq=freq*dxlam/xlam
      print*,'freq,dfreq'
      print*,freq,dfreq
*  4.73755462E+14  7.48665395E+09
*
*
      print*
*      dimension echo(4,4)
      iecho=4
      ix=4
      cair=30.
      csoil=10.
      cmarble=10.5
      echo(1,1)=80.1
      echo(2,1)=164.3
*23456789a123456789b123456789c123456789d123456789e123456789f123456789g12
      echo(1,2)=80.5 
      echo(2,2)=90.8 
      echo(3,2)=111.0
      echo(4,2)=122.7
      echo(1,3)=80.4
      echo(2,3)=90.2
      echo(3,3)=112.0
      echo(4,3)=121.2
      echo(1,4)=80.6
      echo(2,4)=170.0
      do 410 i=1,ix
      do 420 j=1,iecho
      if(j .gt. 1) then
         decho=echo(j,i)-echo(j-1,i)
       else
         decho=echo(j,i)
      end if
      decho=.5*decho
      if(i .eq. 1  .or.  i .eq. 4) then
          if(j .eq. 1) then
              thick=csoil*decho
            else
              thick=cmarble*decho
          end if
        else
          if(j .eq. 1) then
              thick=csoil*decho
            else if(j .eq. 3) then
              thick=cair*decho
            else
              thick=cmarble*decho
          end if
      end if
      print*,i,j,thick*.01
  420 continue
  410 continue
*
      print*
      pi=acos(-1.)
      xlum=3.86e+26
      astunit=1.496e+11
      flux=xlum/(4.*pi*astunit**2)
      fluxp=flux/40.**2
      print*,'flux,fluxp'
      print*,flux,fluxp
* 1372.5061  0.857816339
*
      print*
      clight=2.99792458e+8     ! in mks
      wave=.589e-6
      xindexair=1.000293
      freq=(clight/xindexair)/wave
      xindex=1.52
      waveglass=wave*(xindexair/xindex) ! HZ-38 air index
      vglass=clight/xindex
      print*,'freq,waveglass,vglass'
      print*,freq,waveglass,vglass
* 5.08836352E+14  3.87613568E-07  197231872.
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      v1=4.
      v2=3
      theta1=40./raddeg
      theta2=raddeg*asin((v2/v1)*sin(theta1))
      print*,'theta2=',theta2  ! = 28.8220367 degrees
*
      print*
      ds=0.10
      xlam=632.8e-9
      ell=42.
      dd=ell*xlam/ds
      print*,'slit separation is ',dd  ! = 0.000265775976
*
*
      print*
      vsound=343. ! at 20 C in 1 atm (HRW-400).
      flow=30.    ! Clark-66
      fhigh=3.e+4 ! Clark-66
      wlow=vsound/fhigh
      whigh=vsound/flow
      fconcerta=440.
      wconcerta=vsound/fconcerta
      theta=asin(2.*wconcerta/4.)
      x=100.*tan(theta)
      theta2=asin(1.*wconcerta/1.)
      x2=100.*tan(theta2)
      print*,'wlow,whigh,wconcerta,x,x2'
      print*,wlow,whigh,wconcerta,x,x2
* 0.0114333332  11.4333334  0.779545426  42.3246841  124.45929
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      wave=.550e-6
      dd=(1.*wave)/sin(26./raddeg)
      theta2=raddeg*asin(2.*wave/dd)
      red=.7e-6
      blue=.4e-6
      thetar1=raddeg*asin(red/dd)
      thetab2=raddeg*asin(2.*blue/dd)
      print*,'dd,theta2,thetar1,thetab2'
      print*,dd,theta2,thetar1,thetab2
* 1.25464453E-06  61.2518501  33.9125481  39.6153793
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      wave=1.2
      theta2=28.
      dd=(2.*wave)/(2.*sin(theta2/raddeg))
      theta1=raddeg*asin(1.*wave/(2.*dd))
      print*,'dd,theta1'
      print*,dd,theta1
*   2.55606532  13.5760489
*
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      wave=.633e-6
      delta=3.
      aa=wave/sin((delta/2.)/raddeg)
      print*,'Width=',aa   ! 2.41815796E-05
*
      print*
      en=7.5e+6*1.602e-19
      hemass=4.0026*1.66e-27
      hplanck=6.63e-34
      wave=hplanck/sqrt(2.*hemass*en)
      print*,'wave=',wave    !  5.24700488E-15 
*
      print*
      hplanck=6.626e-34
      bolt=1.381e-23
      amu=1.661e-27
      hemass=4.0026*amu
      tem=300.
      pres=1.01e+5
      wave=hplanck/sqrt(2.*hemass*1.5*bolt*tem)
      den=pres/(bolt*tem)
      vol=1./den
      dist=vol**(1./3.)
      ratio=den*wave**3
      ratio3=ratio**(1./3.)
      print*,'wave,den,vol,dist'
      print*,wave,den,vol,dist
* 7.28915955E-11  2.43784699E+25  4.10198026E-26  3.44877038E-09
      print*,'ratio,ratio3'
      print*,ratio,ratio3
*  9.44145268E-06  0.0211355183
*
      print*
      y0=1000
      v0=100
      gmoon=1.6
      ymax=y0+v0**2/(2.*gmoon)
      vfin=sqrt(v0**2+2.*(-gmoon)*(0.-y0))
      print*,'ymax,vfin'
      print*,ymax,vfin
*  4125.  114.891251
*
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      gg=9.8
      theta=30.
      xmass=2.
      fn=xmass*gg*cos(theta/raddeg)
      fd=xmass*gg*sin(theta/raddeg)
      aa=gg*sin(theta/raddeg)
      dist=.5*aa*10.**2
      print*,'fn,fd,aa,dist'
      print*,fn,fd,aa,dist
* 16.9740982  9.80000019  4.9000001  245.
*
      print*
      eryd=13.6
      en=-eryd*(1./2.**2-1.)
      hc=12398.
      wave=hc/en
      echarge=1.602e-19
      clight=2.998e+8
      aa=1.00797
      amu=1.66e-27
      amucsq=.9315e+9
      ratio=.5*en/(aa*amucsq)
      enh=ratio*en
      boltev=.8617e-4
      enth=1.5*boltev*300.
      vh=sqrt(2.*enh*echarge/(aa*amu))
      print*,'en,wave,ratio,enh,enth,vh'
      print*,en,wave,ratio,enh,enth,vh
* 10.2000008 1215.49011 5.43174972E-09 5.54038522E-08 0.038776502 3.2571547
*
      print*
      aa=.529
      rr=8e-6
      pp=(4./3.)*(rr/aa)
      print*,'pp=',pp     !
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      aa=5.
      wave1=.50
      del1=2.*asin(wave1/aa)*raddeg
      wave2=.010
      del2=2.*asin(wave2/aa)*raddeg
      print*,'del1,del2'
      print*,del1,del2
* 11.4783401  0.229183257
*
      print*
      v=30.
      a=1.
      x0=-200.
      t1=(v-sqrt(v**2+2.*a*x0) )/a
      t2=(v+sqrt(v**2+2.*a*x0) )/a
      print*,'t1,t2'
      print*,t1,t2
*  7.63932037  52.3606796
*
      print*
      v=50.
      x0=20.
      a=6.
      conv=(1./1609.)*(3600./1.)
      vmph=v*conv
      t1=(v-sqrt(v**2+2.*a*x0))/a
      t2=(v+sqrt(v**2+2.*a*x0))/a
      xpass=.5*a*t2**2
      vpass=a*t2
      vpassmph=vpass*conv
      print*,'conv,vmph,t1,t2,xpass,vpass,vpassmph'
      print*,conv,vmph,t1,t2,xpass,vpass,vpassmph
* 2.2374146  111.870728 -0.390834808  17.0575008  872.875  102.345001 228.988205
*
      print*
      va0=60.
      vb0=-50.
      x0=1000.
      v0=vb0-va0
      acel=-v0**2/(2.*(0.-x0))
      t=-v0/acel
      acela=-va0/t
      acelb=-vb0/t
      acelg=acel/9.8
      va0kmhr=va0*3600./1000.
      print*,'acel,t,acela,acelb,acelg,va0kmhr'
      print*,acel,t,acela,acelb,acelg,va0kmhr
*   6.05000019  18.181818 -3.29999995  2.75  0.617346942  216.
*
      print*
      g=9.8
      acel=g*.05
      conv=1000./3600.
      v=216.*conv
      rmin=v**2/acel
      rad=1000.
      vmax=sqrt(acel*rad)
      vmaxkmhr=vmax/conv
      print*,'v,rmin,vmax,vmaxkmhr'
      print*,v,rmin,vmax,vmaxkmhr
*   60.0000038  7346.93945  22.1359444  79.6893997
*
      print*
      pi=acos(-1.)
      raddeg=180./pi
      vpad=.5
      vx=.2
      ymax=20.
      tcross=ymax/vpad
      xcross=vx*tcross
      v=sqrt(vpad**2+vx**2)
      theta=atan(vpad/vx)*raddeg
      vxmax=.8
      xcross2=.5*vxmax*(vpad/ymax)*tcross**2
      print*,'tcross,xcross,v,theta,xcross2'
      print*,tcross,xcross,v,theta,xcross2
*   40.  8.  0.538516462  68.1985855  16.
*
      t=(50.+sqrt(50**2+4.*3.3*20) )/6.6
      print*,t
*
      print*
      pi=acos(-1.)
      frpm=78.
      g=9.8
      r=0.06
      w=2.*pi*frpm/60.
      xmu=r*w**2/g
      print*,'xmu=',xmu   ! = 0.408480793
*
      print*
      pi=acos(-1.)
      daysec=86400.
      earthmass=5.9742e24
      earthradmean=6378.138e+3
      omega=2.*pi/daysec
      alpha=omega/(pi*1.e7) ! constant acceleration over a year
      xirot=0.4*earthmass*earthradmean**2
      torque=xirot*alpha
      xke=.5*xirot*omega**2
      print*,'omega,alpha,xirot,torque,xke'
      print*,omega,alpha,xirot,torque,xke
* 7.2722054E-05  2.31481479E-12 9.72137211E+37 2.25031757E+26 2.5705724E+29
*
      print*
      epsilon=.28
      power=800.
      tzero=273.15
      th=250.+tzero
      tc=30.+tzero
      epcarnot=1.-(tc/th)
      pextr=power/epsilon
      pwaste=pextr-power
      xnum=pwaste*1.e+6/1.e+4
      print*,'epcarnot,pextr,pwaste,xnum'
      print*,epcarnot,pextr,pwaste,xnum
*   0.420529515  2857.14282  2057.14282  205714.281
*
      print*
      xn=6.
      rgas=8.31
      f=5.
      ss=xn*rgas*(f/2.)*log(10.)
      print*,'Change in entropy is ',ss   ! 287.017242
*
      print*
      pi=acos(-1.)
      gg=6.67407e-11   ! the 2002 value, but not yet confirmed
      p=27.*3600.
      r=100.*1.e+3
      xmass=(1./gg)*(4.*pi**2*r**3)/p**2
      vol=55.*(pi*11.**2)  ! in km**3
      den=xmass/(vol*1.e+9)
      print*,'xmass,vol,den'
      print*,xmass,vol,den
* 6.26089537E+16  20907.2988  2994.5979
*
      print*
      tt=1. ! in years
      cc=1. ! 1 lyr/yr
      dd=8.7
      beta=1./sqrt(1.+(cc*tt/dd)**2)
      print*,'beta'
      print*,beta
*
      print*
      tc=26.85
      tf=1.8*tc+32.
      print*,'tc,tf'
      print*,tc,tf
*
      print*
      radmau=rearthmoon/astunit
      radmea=rearthmoon/reartheq
      print*,'rearthmoon,reartheq'
      print*,rearthmoon,reartheq
      print*,radmau,radmea
*
      print*
      dimoon=2.*asin(1728./384400.)*raddeg
      disun=2.*asin(6.9599e+8/astunit)*raddeg
      print*,'dimoon,disun,astunit'
      print*,dimoon,disun,astunit
      print*,dimoon*60.,disun*60.
*
      print*
      xmoon=.0123
      cm=xmoon*radmea
      print*,'Center of mass location ', cm
*
      print*
      theta=360.
      xlun=29.53059 
      xsid=27.321661
      rlun=theta/xlun
      rsid=theta/xsid
      print*,'Lunar angular speeds'
      print*,'rlun,rsid'
      print*,rlun,rsid
*
      print*
      xjyr=365.25
      daysec=86400.
      yearsec=xjyr*daysec
      pi=acos(-1.)
      print*,'yearsec,yearsec*1.e-7/pi'
      print*,yearsec,yearsec*1.e-7/pi
*          3.15576E+07    1.00451
*
      print*
      pi=acos(-1.)
      radeartheq=6.378136e+6  ! Cox-340
      daysec=86400.
      veq=2.*pi*radeartheq/daysec ! This is, of course relative to Sun.
      print*,'The approximate equatorial speed is ',veq ! 463.831
      acf=veq**2/radeartheq
      print*,'The equatorial centrifugal force per mass is ',acf ! 3.37308E-02
*
      print*
      xlunar=29.53059  ! Cox-16
      xsolar=365.2421897 ! Cox-15
      xlyear=12.*xlunar
      dif=3.*(xsolar-xlyear)
      xlun=12.*12.*xlunar+7.*13.*xlunar
      xsol=19.*xsolar
      print*,'xlyear,dif,xlun,xsol'
      print*,xlyear,dif,xlun,xsol
* 354.367    32.6254    6939.69    6939.60
*
      print*
      rmoon=384400e+3   ! Cox-303
      xmoon=7.3483e+22  ! Cox-305
      grav=6.673e-11    ! Cox-8
      reartheq=6.378136e+6  ! Cox-240
      ftidalmax=2.*grav*xmoon*reartheq/rmoon**3
      print*,'The maximum tidal force per unit mass is ',ftidalmax 
*      
      print*
      print*,100.*29./35.,100.*29./36.
      print*,100.*31./35.,100.*31./36.
*
      print*
      rmars=3397.
      rphobos=9378.
      rdeimos=23459.
      ratp=rphobos/rmars
      ratd=rdeimos/rmars
      print*,'ratp,ratd'
      print*,ratp,ratd
*
      print*
      rjup=71492.
      rio=422000.
      reu=671000.
      rga=1070000. 
      rca=1883000.
      print*,'rio/rjup,reu/rjup,rga/rjup,rca/rjup'
      print*,rio/rjup,reu/rjup,rga/rjup,rca/rjup
*
      gg=6.673e-11
      xmass=5.9737e+24
      rr=6.378136e+6
      vv1=sqrt(gg*xmass/rr)
      vv=sqrt(2.*gg*xmass/rr)
      print*,'vv1,vv'
      print*,vv1,vv
*          7905.61    11180.2
*
      print*
      pi=acos(-1.)
      xmpc=3.0856776e+24  ! 1 Mpc in cm.
      h0=70.
      daysec=86400.
      clight=2.99792485e+5 ! speed of light in km/s
      grav=6.673e-8 ! cgs gravity constant
      h0asec=h0*1.e+5/xmpc
      th0sec=1./h0asec
      th0=th0sec/(365.25*daysec)
      dh0=clight/h0
      rhoc=( 3./(8.*pi) )/(grav*th0sec**2 )
      print*,'h0asec,th0,dh0,rhoc'
      print*,h0asec,th0,dh0,rhoc 
*
      print*
      hplanck=6.6260755e-34
      clight=2.99792485e+8
      hc=hplanck*clight
      print*,'hc'
      print*,hc
*
      print*
      earthradmean=6378.138e+3
      acharon=19.6e+6   ! Cox-305  semi-major axis
      radpluto=1195.e+3
      radcharon=593.e+3  ! Cox-306
      ratio=acharon/earthradmean
      ratiop=acharon/radpluto
      radch=radcharon/earthradmean
      radchp=radcharon/radpluto
      print*,'ratio,ratiop,radch,radchp'
      print*,ratio,ratiop,radch,radchp
*
      print*
      pmass=1.3e+22  ! Cox-295
      cmass=1.62e+21  ! Cox-306
      emass=5.9742e+24 ! Cox-295
      rpe=pmass/emass
      rce=cmass/emass
      rcp=cmass/pmass
      print*,'rpe,rce,rcp'
      print*,rpe,rce,rcp
*
      print*
      pi=acos(-1.)
      vol=4.*pi/3.
      denp=pmass/(vol*radpluto**3)
      denc=cmass/(vol*radcharon**3)
      print*,'denp,denc'
      print*,denp,denc
*
      print*
      solcon=1367.  ! Cox-340 average value
*      aa=.367      ! Cox-299 this is only visual light
      aa=.3         ! Rampino & Caldiera - 85 presumably this
*                   ! the good average for all wavelengths
      sigma=5.67e-8
      tem=( solcon*(1.-aa)/(4.*sigma) )**.25
      print*,'tem'
      print*,tem
*            254.862
*
      print*
      p1=-243.01
      p2=224.68
      p3=583.9214
      pday=(p2*p1)/(p2-p1)
      pres=p3/243.16
      print*,'pday,pres'
      print*,pday,pres
*
      print*
      ymars=1.88071105   ! Cox-294 in Julian years
      rotmars=1.02595675   ! Cox-296
      ymarsd=ymars*365.25
      dmars=ymarsd*rotmars/(ymarsd-rotmars)
      print*,'Martian day ',dmars 
*
      print*
      asemi(1)=39.48168677  ! Pluto
      asemi(2)=39.473       ! 2004 DW
      asemi(3)=532.         ! Sedna VB12 
      asemi(4)=43.377       ! Quaoar 
      asemi(5)=39.485       ! Ixion
      asemi(6)=43.129       ! Varuna
      asemi(7)=47.501       ! 2002 AW197
      iasemi=7
      conv=365.256/365.25  ! conversion from sidereal to Julian years
      do 490 i=1,iasemi
      pp=asemi(i)**1.5*conv
      print*,i,asemi(i),pp
  490 continue
* 
      print*
      print*,'The inner moons of Jupiter.'
      asemi(1)=128.e+3
      asemi(2)=129.e+3
      asemi(3)=181.e+3
      asemi(4)=222.e+3
      do 492 i=1,4
      asemi(i)=asemi(i)/71492.
      print*,i,asemi(i)
  492 continue
*
      print*
      aaa=1.458
      eee=.223
      peri=(1.-eee)*aaa
      print*,'peri'
      print*,peri
*
      print*
      xm=40.e+6
      rho=3000.
      pi=acos(-1.)
      rr=( xm/(4.*(pi/3.)*rho) )**(1./3.)
      print*,'rr'
      print*,rr
*
      print*
      yrven=.61518257
      yrven=yrven*365.25
      print*,'yrven'
      print*,yrven
*
      print*
      pi=acos(-1.)
      radeartheq=6.378136e+6  ! Cox-340
      daysec=86400.
      daysid=.99726968
      veq=2.*pi*radeartheq/daysec ! This is, of course relative to Sun.
      print*,'The synodic equatorial speed is ',veq ! 463.831 m/s
      veq=2.*pi*radeartheq/(daysec*daysid) ! This is relative to fixed stars 
      print*,'The sidereal equatorial speed is ',veq ! 465.101 m/s
      acf=veq**2/radeartheq
      print*,'The equatorial centrifugal force per mass is ',acf !3.39157E-02 
*
      print*
      payus=33000.
      exch=1.35
      paycdn=payus*exch
      print*,'payus,paycdn'
      print*,payus,paycdn
*
      print*
      f1=10.**14.875
      f2=f1/10.
      clight=3.e+10
      xlam1=clight/f1 
      xlam2=clight/f2
      print*,'xlam1,xlam2'
      print*,xlam1,xlam2
*
      print*
      pi=acos(-1.)
      xlum=3.86e26
      solconst=1373.
      dd=sqrt(xlum/(4.*pi*solconst))
      print*,'The calculated Earth-Sun distance is ',dd,
     & ' m'  !  1.49573E+11
*
      print*
      grav=6.6742e-11
      clight=2.99792458e+8
      solmass=1.9891e+30
      rsch=2.*grav*solmass/clight**2
      rsch6AU=(2.*grav*solmass*1.e+6/clight**2)*(1./1.49597870e+11)
      print*,'rsch,rsch6AU'
      print*,rsch,rsch6AU
*      2954.23    1.97478E-02
*
      print*
      ckms=2.99792458e+5
      hubble=71.
      dh=ckms/hubble
      th=3.085678e+19/71. 
      thjyr=th/3.15576e+7
      ph=46.e+3  ! particle horizon in Mlyr 
      ph=ph*(1./3.262) ! particle horizon in Mpc
      print*,'hubble,dh,th,thjyr,ph'
      print*,hubble,dh,th,thjyr,ph
*     71.00000     4222.429    4.3460255E+17  1.3771724E+10   14101.78
*
      print*
      radch=593.
      radea=6378.14
      radpl=1195
      radchea=radch/radea
      radchpl=radch/radpl
      print*,'radchea,radchpl'
      print*,radchea,radchpl
*     9.2973813E-02  0.4962343 
*
      print*
      au=1.4959787066e+13
      rsun=6.95508e+10
      rsunau=rsun/au
      print*,'au= ',rsunau ! = 4.6491837E-03
* 
      print*
      grav=6.6742e-11
      rr=100.*3.0856776e+16
      vv=100*10.e+3
      xmsun=1.9891e+30
      xmm=vv**2*(rr/xmsun)/grav
      print*,'xmm=',xmm   !   2.3243135E+10
*
      print*
      xm=.012300034   ! Moon in Earth masses
      xmi=1./xm
      print*,'xmi=',xmi ! inverse Moon mass
*
      print*
      xmass=100.
      gg=9.8 
      hh=4.e+3
      tt=3600.*10.
      energy=xmass*gg*hh
      power=xmass*gg*hh/tt
      print*,'energy,power'
      print*,energy,power
*            3920000.       108.8889 
*
      print*
      dd=1.83
      ddin=dd*100.*(1./2.54)
      print*,'ddin'
      print*,ddin
*
*
      print*
      pp=0.44401
      pp=pp*24.
      pp2=pp-10.
      pp2=pp2*60.
      print*,'pp,pp-10,pp2'
      print*,pp,10.,pp2
*   10.65624    10.   39.37437
*
      print*
      clight=2.99792458e+8 
      ex=5.4e6
      xmcsq=1.*clight**2
      frac=ex/xmcsq
      xmult=ex*clight**2
      print*,'xmcsq,xmult,frac'
      print*,xmcsq,xmult,frac
*    8.9875515E+16  4.8532778E+23  6.0083105E-11
*
      print*
      rho=19.3   ! density of gold, ordinary terrestrial conditions.
      xm=30.
      vol=xm/rho
      print*,'rho/xm,rho*xm,vol'
      print*,rho/xm,rho*xm,vol 
*    0.6433333       579.0000       1.554404
*
      print*
      theta=.232
      dd=1./theta
      print*,'distance d=',dd   ! 4.310345 
*
      print*
      tt=.387e+9
      yearj=3.15576e+7
      tty=tt/yearj
      print*,'Triton half-life in years ',tty  !  12.26329
*
      print*
      rsun=6.95508e+8      ! Cox-12
      au=1.4959787066e+11  ! Cox-12
      rsunau=rsun/au
      rcorona=30*rsunau
      print*,'rsunau,rcorona'
      print*,rsunau,rcorona
*     4.6491837E-03  0.1394755 
* 
      print*
      print*,111./144.
*
      end
