diff --git a/Makefile b/Makefile
index 8e36d81f8..a0fea5c7c 100644
--- a/Makefile
+++ b/Makefile
@@ -30,7 +30,7 @@ legacy:
 	@echo "=========================================================================================================="
 
 debug:
-	-gcc -g -D DEBUG -Wunused-but-set-variable -Wunused-variable -Wfloat-equal -o bin/vplanet src/*.c -lm -DGITVERSION=\"$(GITVERSION)\"
+	-gcc -g -D DEBUG -o bin/vplanet src/*.c -lm -DGITVERSION=\"$(GITVERSION)\"
 
 debug_no_AE:
 	-gcc -g -o bin/vplanet src/*.c -lm -DGITVERSION=\"$(GITVERSION)\"
@@ -52,7 +52,7 @@ cpp:
 	g++ -o bin/vplanet src/*.c -lm -O3 -fopenmp -fpermissive -w -DGITVERSION=\"$(GITVERSION)\"
 
 warnings:
-	clang -Weverything src/*.c -lm -O3 -DGITVERSION=\"$(GITVERSION)\"
+		-gcc -g -D DEBUG -Wunused-but-set-variable -Wunused-variable -Wfloat-equal -o bin/vplanet src/*.c -lm -DGITVERSION=\"$(GITVERSION)\"
 
 parallel:
 	gcc -o bin/vplanet src/*.c -lm -O3 -fopenmp -DGITVERSION=\"$(GITVERSION)\"
diff --git a/examples/CosmicShoreline/CosmicShoreline.png b/examples/CosmicShoreline/CosmicShoreline.png
new file mode 100644
index 000000000..0c471714e
Binary files /dev/null and b/examples/CosmicShoreline/CosmicShoreline.png differ
diff --git a/examples/CosmicShoreline/README.rst b/examples/CosmicShoreline/README.rst
new file mode 100644
index 000000000..628f8ebb0
--- /dev/null
+++ b/examples/CosmicShoreline/README.rst
@@ -0,0 +1,38 @@
+The "Cosmic Shoreline"
+==========
+
+Overview
+--------
+
+A reproduction of Zahnle & Catling (2018)'s cosmic shoreline,
+a proposed boundary between planet with and without atmosphere.
+
+===================   ============
+**Date**              09/12/2024
+**Author**            Rory Barnes
+**Modules**           Atmesc, STELLAR
+**Approx. runtime**   10 seonds
+===================   ============
+
+
+To run this example
+-------------------
+
+.. code-block:: bash
+
+   python makeplot.py <pdf | png>
+
+
+Expected output
+---------------
+
+.. figure:: CosmicShoreline.png
+   :width: 100%
+   :align: center
+
+The blue line is the cosmic shoreline, which is taken from Fig. 2
+of Zahnle & Catling (2018). The black points are the positions
+of the Solar System planets in this parameter space. The Sun's XUV luminosity
+is assumed to follow the Ribas et al. (2005) model. This figure, while
+not as complete as Zahnle & Catling (2018), is a near exact reproduction
+of their results.
diff --git a/examples/CosmicShoreline/earth.in b/examples/CosmicShoreline/earth.in
new file mode 100644
index 000000000..31db62cd4
--- /dev/null
+++ b/examples/CosmicShoreline/earth.in
@@ -0,0 +1,6 @@
+sName        Earth
+saModules    atmesc
+dMass        3.0018090452e-06
+dRadius      -1.
+dSemi        1.00000321
+saOutputOrder    Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/george.in b/examples/CosmicShoreline/george.in
new file mode 100644
index 000000000..28e6c178e
--- /dev/null
+++ b/examples/CosmicShoreline/george.in
@@ -0,0 +1,6 @@
+sName        George
+saModules    atmesc
+dMass        -14.54
+dRadius      -4.007
+dSemi        19.19203990
+saOutputOrder    Time -fXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/jupiter.in b/examples/CosmicShoreline/jupiter.in
new file mode 100644
index 000000000..24e07eec4
--- /dev/null
+++ b/examples/CosmicShoreline/jupiter.in
@@ -0,0 +1,6 @@
+sName        Jupiter
+saModules    atmesc
+dMass        -317.83
+dRadius      -11.2
+dSemi        5.20880408
+saOutputOrder    Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/makeplot.py b/examples/CosmicShoreline/makeplot.py
new file mode 100644
index 000000000..56419875e
--- /dev/null
+++ b/examples/CosmicShoreline/makeplot.py
@@ -0,0 +1,73 @@
+import os
+import pathlib
+import subprocess
+import sys
+
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+import numpy as np
+import vplot as vpl
+from matplotlib.ticker import MaxNLocator
+from mpl_toolkits.axes_grid1 import make_axes_locatable
+
+import vplanet
+
+path = pathlib.Path(__file__).parents[0].absolute()
+sys.path.insert(1, str(path.parents[0]))
+
+output = vplanet.run(units = False)
+
+# Plot!
+fig = plt.figure(figsize=(8.5, 6))
+
+fxuv_earth = output.log.final.Earth.CumulativeXUVFlux
+
+fxuv = []
+fxuv.append(output.log.final.Mercury.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.Venus.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.Earth.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.Mars.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.Jupiter.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.Saturn.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.George.CumulativeXUVFlux/fxuv_earth)
+fxuv.append(output.log.final.Neptune.CumulativeXUVFlux/fxuv_earth)
+
+escvel = []
+escvel.append(output.log.final.Mercury.EscapeVelocity/1e3)
+escvel.append(output.log.final.Venus.EscapeVelocity/1e3)
+escvel.append(output.log.final.Earth.EscapeVelocity/1e3)
+escvel.append(output.log.final.Mars.EscapeVelocity/1e3)
+escvel.append(output.log.final.Jupiter.EscapeVelocity/1e3)
+escvel.append(output.log.final.Saturn.EscapeVelocity/1e3)
+escvel.append(output.log.final.George.EscapeVelocity/1e3)
+escvel.append(output.log.final.Neptune.EscapeVelocity/1e3)
+
+shorelinex = []
+shorelinex.append(0.2)
+shorelinex.append(60)
+
+shoreliney = []
+shoreliney.append(1e-6)
+shoreliney.append(1e4)
+
+plt.xlabel('Escape Velocity [km/s]')
+plt.ylabel('Normalized Cumulative XUV Flux')
+plt.plot(shorelinex,shoreliney,color=vpl.colors.pale_blue)
+plt.plot(escvel,fxuv,'o',color='k')
+plt.xscale('log')
+plt.yscale('log')
+plt.ylim(1e-6,1e4)
+plt.xlim(0.1,200)
+
+plt.annotate('Mercury',(2.7,9))
+plt.annotate('Venus',(10.2,2.6))
+plt.annotate('Earth',(11.5,0.5))
+plt.annotate('Mars',(5,0.2))
+plt.annotate('Jupiter',(60,0.05))
+plt.annotate('Saturn',(37,0.014))
+plt.annotate('Uranus',(22,0.0034))
+plt.annotate('Neptune',(24,0.0006))
+
+
+# Save figure
+fig.savefig(path / f"CosmicShoreline.png", bbox_inches="tight", dpi=200)
diff --git a/examples/CosmicShoreline/mars.in b/examples/CosmicShoreline/mars.in
new file mode 100644
index 000000000..7a9c9e0d8
--- /dev/null
+++ b/examples/CosmicShoreline/mars.in
@@ -0,0 +1,6 @@
+sName        Mars
+saModules    atmesc
+dMass        -0.107
+dRadius      -0.53202
+dSemi        1.52366290
+saOutputOrder    Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/mercury.in b/examples/CosmicShoreline/mercury.in
new file mode 100644
index 000000000..4debc018d
--- /dev/null
+++ b/examples/CosmicShoreline/mercury.in
@@ -0,0 +1,6 @@
+sName        Mercury
+saModules    atmesc
+dMass        -0.0553
+dRadius      -0.383
+dSemi        0.38709897
+saOutputOrder    Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/neptune.in b/examples/CosmicShoreline/neptune.in
new file mode 100644
index 000000000..8da0ec951
--- /dev/null
+++ b/examples/CosmicShoreline/neptune.in
@@ -0,0 +1,6 @@
+sName        Neptune
+saModules    atmesc
+dMass        -17.15
+dRadius      -3.883
+dSemi        30.07050641
+saOutputOrder    Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/saturn.in b/examples/CosmicShoreline/saturn.in
new file mode 100644
index 000000000..de7e74aeb
--- /dev/null
+++ b/examples/CosmicShoreline/saturn.in
@@ -0,0 +1,6 @@
+sName        Saturn
+saModules    atmesc
+dMass        -95.16
+dRadius      -9.449
+dSemi       9.53999265
+saOutputOrder    Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/sun.in b/examples/CosmicShoreline/sun.in
new file mode 100644
index 000000000..80e98f955
--- /dev/null
+++ b/examples/CosmicShoreline/sun.in
@@ -0,0 +1,18 @@
+# The host star
+sName	        Sun		                    # Body's name
+saModules	    stellar                     # Modules to apply, exact spelling required
+
+# Output
+saOutputOrder Age -Luminosity -LXUVStellar -Radius Temperature
+
+# Physical Parameters
+dMass         1.00                        # Mass, solar masses
+dAge          5e7
+dObliquity    0
+dRotPeriod    -30                         # Rotation period, negative -> days
+dRadGyra      0.5                         # Radius of gyration (moment of inertia constant)
+
+# STELLAR Parameters
+sStellarModel baraffe                     # Stellar evolution model: `baraffe` or `none`
+dSatXUVFrac   1.e-3                       # Saturation level of the XUV luminosity
+dSatXUVTime   -0.1
diff --git a/examples/CosmicShoreline/venus.in b/examples/CosmicShoreline/venus.in
new file mode 100644
index 000000000..5254830d3
--- /dev/null
+++ b/examples/CosmicShoreline/venus.in
@@ -0,0 +1,7 @@
+# Planet a parameters
+sName         Venus                      # Body's name
+saModules     atmesc
+dMass         -0.815
+dRadius       -0.9499
+dSemi         -0.723
+saOutputOrder Time -FXUV -CumulativeXUVFlux
diff --git a/examples/CosmicShoreline/vpl.in b/examples/CosmicShoreline/vpl.in
new file mode 100644
index 000000000..c53e2fbce
--- /dev/null
+++ b/examples/CosmicShoreline/vpl.in
@@ -0,0 +1,25 @@
+sSystemName   solarsystem                   # System Name
+iVerbose      5                             # Verbosity level
+bOverwrite    1                             # Allow file overwrites?
+
+# List of "body files" that contain body-specific parameters
+saBodyFiles   sun.in mercury.in venus.in earth.in mars.in $
+              jupiter.in saturn.in george.in neptune.in
+
+# Input/Output Units
+sUnitMass      solar                        # Options: gram, kg, Earth, Neptune, Jupiter, solar
+sUnitLength    aU                           # Options: cm, m, km, Earth, Jupiter, solar, AU
+sUnitTime      YEARS                        # Options: sec, day, year, Myr, Gyr
+sUnitAngle     d                            # Options: deg, rad
+
+# Input/Output
+bDoLog         1                            # Write a log file?
+iDigits        6                            # Maximum number of digits to right of decimal
+dMinValue      1e-10                        # Minimum value of eccentricity/obliquity
+
+# Evolution Parameters
+bDoForward    1                             # Perform a forward evolution?
+bVarDt        1                             # Use variable timestepping?
+dEta          0.01                           # Coefficient for variable timestepping
+dStopTime     4.5e9                         # Stop time for evolution
+dOutputTime   1e8                           # Output timesteps (assuming in body files)
diff --git a/src/atmesc.c b/src/atmesc.c
index 979afe37e..deee35007 100644
--- a/src/atmesc.c
+++ b/src/atmesc.c
@@ -2001,7 +2001,6 @@ void fnPropsAuxAtmEsc(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
   }
 }
 
-
 /**
 Assigns functions returning the time-derivatives of each variable
 to the magical matrix of function pointers.
diff --git a/src/atmesc.h b/src/atmesc.h
index bb4392110..2eaabd1d1 100644
--- a/src/atmesc.h
+++ b/src/atmesc.h
@@ -171,6 +171,7 @@ void FinalizeUpdateMassAtmEsc(BODY *, UPDATE *, int *, int, int, int);
 #define OUT_HDIFFFLUX 1238 /**< The diffusion limited flux (the true flux in the diffusion regime) */
 #define OUT_HREFODRAGMOD 1239 /**< Multiply by H ref flux to get H flux with drag of oxgyen */
 #define OUT_KTIDE 1240 /**< Gravitational enhancement of mass loss */
+#define OUT_CUMULATIVEFXUV 1250 
 
 void InitializeOutputAtmEsc(OUTPUT *, fnWriteOutput[]);
 void InitializeOutputFunctionAtmEsc(OUTPUT *, int, int);
@@ -212,6 +213,7 @@ void LogBodyAtmEsc(BODY *, CONTROL *, OUTPUT *, SYSTEM *, UPDATE *,
 void fnForceBehaviorAtmEsc(BODY *, MODULE *, EVOLVE *, IO *, SYSTEM *, UPDATE *,
                            fnUpdateVariable ***, int, int);
 void fnPropsAuxAtmEsc(BODY *, EVOLVE *, IO *, UPDATE *, int);
+
 double fdDSurfaceWaterMassDt(BODY *, SYSTEM *, int *);
 double fdDEnvelopeMassDt(BODY *, SYSTEM *, int *);
 double fdDEnvelopeMassDtBondiLimited(BODY *, SYSTEM *, int *);
diff --git a/src/body.c b/src/body.c
index 3c31b2757..59e9b93f3 100644
--- a/src/body.c
+++ b/src/body.c
@@ -1528,4 +1528,11 @@ void fdHabitableZoneKopparapu2013(BODY *body, int iNumBodies,
 double fdEffectiveTemperature(BODY *body,int iBody) {
   double dTeff = pow((body[iBody].dLuminosity/(4*PI*SIGMA*body[iBody].dRadius*body[iBody].dRadius)),0.25);
   return dTeff;
+}
+
+double fdEscapeVelocity(BODY *body,int iBody) {
+  double dEscVel;
+
+  dEscVel = sqrt(2*BIGG*body[iBody].dMass/body[iBody].dRadius);
+  return dEscVel;
 }
\ No newline at end of file
diff --git a/src/body.h b/src/body.h
index f50c57b24..fa3ecf016 100644
--- a/src/body.h
+++ b/src/body.h
@@ -110,6 +110,7 @@ double CalcDynEllipEq(BODY *, int);
 void fdHabitableZoneKopparapu2013(BODY *, int, double *);
 
 double fdEffectiveTemperature(BODY*,int);
+double fdEscapeVelocity(BODY *,int);
 
 // RB: Move
 
diff --git a/src/evolve.c b/src/evolve.c
index 832eb495c..6a2e05488 100644
--- a/src/evolve.c
+++ b/src/evolve.c
@@ -13,6 +13,11 @@
 
 #include "vplanet.h"
 
+double fdTrapezoidalArea(double dY1,double dY2,double dDeltaX) {
+  double dArea = 0.5*dDeltaX*(dY1 + dY2);
+  return dArea;
+}
+
 void PropsAuxGeneral(BODY *body, CONTROL *control) {
   /* Recompute the mean motion, necessary for most modules */
   int iBody; // Dummy counting variable
@@ -50,6 +55,17 @@ void PropertiesAuxiliary(BODY *body, CONTROL *control, SYSTEM *system,
   }
 }
 
+void fvCumulative(BODY *body,CONTROL *control,SYSTEM *system,double dDt) {
+  int iBody;
+
+  if (body[0].bStellar) {
+    for (iBody=1;iBody<control->Evolve.iNumBodies;iBody++) {
+      fvCumulativeXUVFlux(body,&control->Evolve,system,dDt,iBody);
+    }
+  }
+
+}
+
 void CalculateDerivatives(BODY *body, SYSTEM *system, UPDATE *update,
                           fnUpdateVariable ***fnUpdate, int iNumBodies) {
   int iBody, iVar, iEqn;
@@ -634,6 +650,8 @@ void Evolve(BODY *body, CONTROL *control, FILES *files, MODULE *module,
     dDt = control->Evolve.dTimeStep;
   }
 
+  fvCumulative(body,control,system,dDt);
+
   /* Write out initial conditions */
   WriteOutput(body, control, files, output, system, update, fnWrite);
 
@@ -702,6 +720,8 @@ void Evolve(BODY *body, CONTROL *control, FILES *files, MODULE *module,
        was prior to loop. */
     PropertiesAuxiliary(body, control, system, update);
 
+  fvCumulative(body,control,system,dDt);
+
     // If control->Evolve.bFirstStep hasn't been switched off by now, do so.
     if (control->Evolve.bFirstStep) {
       control->Evolve.bFirstStep = 0;
diff --git a/src/evolve.h b/src/evolve.h
index 641e4dedc..6eec7a40b 100644
--- a/src/evolve.h
+++ b/src/evolve.h
@@ -16,6 +16,7 @@
 
 /* @cond DOXYGEN_OVERRIDE */
 
+double fdTrapezoidalArea(double,double,double);
 void PropertiesAuxiliary(BODY *, CONTROL *, SYSTEM *, UPDATE *);
 void fdGetUpdateInfo(BODY *, CONTROL *, SYSTEM *, UPDATE *,
                      fnUpdateVariable ***);
diff --git a/src/output.c b/src/output.c
index 846d5a821..34528cd3d 100644
--- a/src/output.c
+++ b/src/output.c
@@ -82,6 +82,20 @@ void WriteCriticalSemi(BODY *body, CONTROL *control, OUTPUT *output,
   }
 }
 
+void WriteCumulativeXUVFlux(BODY *body, CONTROL *control, OUTPUT *output,
+                       SYSTEM *system, UNITS *units, UPDATE *update, int iBody,
+                       double *dTmp, char **cUnit) {
+  *dTmp = body[iBody].dFXUVCumulative;
+
+    if (output->bDoNeg[iBody]) {
+    *dTmp *= output->dNeg;
+    fvFormattedString(cUnit, output->cNeg);
+  } else {
+    *dTmp *= fdUnitsEnergyFlux(units->iTime,units->iMass,units->iLength);
+    fsUnitsEnergyFlux(units,cUnit);
+  }
+}
+
 /*
  * D
  */
@@ -128,6 +142,23 @@ void WriteDensity(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system,
   }
 }
 
+/*
+  E
+*/
+
+void WriteEscapeVelocity(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system,
+                 UNITS *units, UPDATE *update, int iBody, double *dTmp,
+                 char **cUnit) {
+
+  *dTmp = fdEscapeVelocity(body,iBody);
+  if (output->bDoNeg[iBody]) {
+    *dTmp *= output->dNeg;
+    fvFormattedString(cUnit, output->cNeg);
+  } else {
+    *dTmp *= fdUnitsTime(units->iTime) / fdUnitsLength(units->iLength);
+    fsUnitsVel(units, cUnit);
+  }
+}
 
 /*
  * H
@@ -1281,6 +1312,15 @@ void InitializeOutputGeneral(OUTPUT *output, fnWriteOutput fnWrite[]) {
         "in unstable orbits. This output parameter prints the instantaneous\n"
         "value of that critical distance.");
 
+  fvFormattedString(&output[OUT_CUMULATIVEFXUV].cName, "CumulativeXUVFlux");
+  fvFormattedString(&output[OUT_CUMULATIVEFXUV].cDescr, "Cumulative XUV flux");
+  fvFormattedString(&output[OUT_CUMULATIVEFXUV].cNeg, "W/m^2");
+  output[OUT_CUMULATIVEFXUV].bNeg       = 1;
+  output[OUT_CUMULATIVEFXUV].dNeg       = 1;
+  output[OUT_CUMULATIVEFXUV].iNum       = 1;
+  output[OUT_CUMULATIVEFXUV].iModuleBit = 1;
+  fnWrite[OUT_CUMULATIVEFXUV]           = &WriteCumulativeXUVFlux;
+
   /*
    * D
    */
@@ -1305,6 +1345,19 @@ void InitializeOutputGeneral(OUTPUT *output, fnWriteOutput fnWrite[]) {
   output[OUT_DENSITY].iModuleBit = 1;
   fnWrite[OUT_DENSITY]           = &WriteDensity;
 
+/* 
+  E
+*/
+
+  fvFormattedString(&output[OUT_ESCAPEVELOCITY].cName, "EscapeVelocity");
+  fvFormattedString(&output[OUT_ESCAPEVELOCITY].cDescr, "Escape Velocity");
+  fvFormattedString(&output[OUT_ESCAPEVELOCITY].cNeg, "m/s");
+  output[OUT_ESCAPEVELOCITY].bNeg       = 1;
+  output[OUT_ESCAPEVELOCITY].dNeg       = 1;
+  output[OUT_ESCAPEVELOCITY].iNum       = 1;
+  output[OUT_ESCAPEVELOCITY].iModuleBit = 1;
+  fnWrite[OUT_ESCAPEVELOCITY]           = &WriteEscapeVelocity;
+
   /*
    * H
    */
diff --git a/src/output.h b/src/output.h
index 369fea3a5..f5c5d3a97 100644
--- a/src/output.h
+++ b/src/output.h
@@ -108,6 +108,7 @@
 #define OUT_LOSTENG 691
 #define OUT_LOSTANGMOM 692
 
+#define OUT_ESCAPEVELOCITY 693
 /* @cond DOXYGEN_OVERRIDE */
 
 
diff --git a/src/system.c b/src/system.c
index 21a598107..321a0333d 100644
--- a/src/system.c
+++ b/src/system.c
@@ -1114,4 +1114,19 @@ double fdBondiRadius(BODY *body, int iBody) {
     dBondiRadius = -1;
   }
   return dBondiRadius;
-}
\ No newline at end of file
+}
+
+void fvCumulativeXUVFlux(BODY *body,EVOLVE *evolve,SYSTEM *system,double dDt,int iBody) {
+  if (evolve->bFirstStep) {
+    body[iBody].dFXUVCumulative = 0;
+    if (body[iBody].bCalcFXUV) {
+      body[iBody].dFXUV = fdXUVFlux(body, iBody);
+    }
+  } else {
+    if (body[iBody].bCalcFXUV) {
+      body[iBody].dFXUV = fdXUVFlux(body, iBody);
+    }
+    body[iBody].dFXUVCumulative += fdTrapezoidalArea(body[iBody].dFXUV,body[iBody].dFXUVLast,dDt);
+    body[iBody].dFXUVLast = body[iBody].dFXUV;
+  }
+}
diff --git a/src/system.h b/src/system.h
index c38dabb4a..686bca3d3 100644
--- a/src/system.h
+++ b/src/system.h
@@ -56,4 +56,7 @@ void fdMergePlanet(BODY *, UPDATE *, fnUpdateVariable ***, int);
 double fdBondiRadius(BODY *, int);
 double fdRocheRadius(BODY *, int, int);
 
+void fvCumulativeXUVFlux(BODY *,EVOLVE *,SYSTEM *,double,int);
+
+
 /* @endcond */
diff --git a/src/vplanet.h b/src/vplanet.h
index de3ed3395..ff52fb16a 100644
--- a/src/vplanet.h
+++ b/src/vplanet.h
@@ -262,6 +262,8 @@ struct BODY {
   double dThermTemp; /**< Thermosphere's temperature in Lehmer-Catling model */
   double dAtmGasConst; /**< Atmosphere's gas constant in Lehmer-Catling model */
   double dFXUV;        /**< XUV Flux at planet's atmosphere */
+  double dFXUVLast;    /**< XUV at planet during last time step */
+  double dFXUVCumulative;
   double dJeansTime;   /**< Jeans timescale for atmospheric escape */
   double dFlowTemp;    /**< Temperature of the hydrodynamic flow */
   double dRocheRadius; /**< Radius of the Roche lobe */