Skip to content

Commit

Permalink
Merge pull request #287 from VirtualPlanetaryLaboratory/ConstXUV
Browse files Browse the repository at this point in the history 
Fixed bug in which a body running only atmesc did not lose water if d…
  • Loading branch information
@RoryBarnes
RoryBarnes authored Apr 13, 2024
2 parents f5d1893 + 4507311 commit db1a296
Show file tree
Hide file tree
Showing 24 changed files with 1,354 additions and 13 deletions.
29 changes: 18 additions & 11 deletions src/atmesc.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4138,18 +4138,25 @@ int fbDoesWaterEscape(BODY *body, EVOLVE *evolve, IO *io, int iBody) {
// spectrum! The Kopparapu+14 limit is for a single star only. This
// approximation for a binary is only valid if the two stars have
// similar spectral types, or if body zero dominates the flux.
if (fdInstellation(body, iBody) < fdHZRG14(body, iBody)) {
if (body[iBody].dRGDuration == 0.) {
body[iBody].dRGDuration = body[iBody].dAge;
if (io->iVerbose > VERBPROG && !io->baEnterHZMessage[iBody]) {
printf("%s enters the habitable zone at %.2lf Myr.\n",
body[iBody].cName, evolve->dTime / (YEARSEC * 1e6));
io->baEnterHZMessage[iBody] = 1;
double dInstellation = fdInstellation(body, iBody);
if (dInstellation == -1 && body[iBody].bCalcFXUV == 0) {
// Constant XUV flux, so set water to escape
return 1;
} else {
double dRunawayGreenhouseFlux = fdHZRG14(body, iBody);
if (dInstellation < dRunawayGreenhouseFlux) {
if (body[iBody].dRGDuration == 0.) {
body[iBody].dRGDuration = body[iBody].dAge;
if (io->iVerbose > VERBPROG && !io->baEnterHZMessage[iBody]) {
printf("%s enters the habitable zone at %.2lf Myr.\n",
body[iBody].cName, evolve->dTime / (YEARSEC * 1e6));
io->baEnterHZMessage[iBody] = 1;
}
}
// Only stop water loss if user requested, which is default
if (body[iBody].bStopWaterLossInHZ) {
return 0;
}
}
// Only stop water loss if user requested, which is default
if (body[iBody].bStopWaterLossInHZ) {
return 0;
}
}

Expand Down
119 changes: 119 additions & 0 deletions ...esc/WaterELimConstXUVLB15NoO2SinkBolmont16/test_WaterELimConstXUVLB15NoO2SinkBolmont16.py
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,119 @@
from benchmark import Benchmark, benchmark
import astropy.units as u

@benchmark(
{
"log.initial.system.Age": {"value": 0.000000, "unit": u.sec},
"log.initial.system.Time": {"value": 0.000000, "unit": u.sec},
"log.initial.system.TotAngMom": {"value": 4.416946e+33, "unit": (u.kg * u.m ** 2) / u.sec},
"log.initial.system.TotEnergy": {"value": -2.237790e+32, "unit": u.Joule},
"log.initial.system.PotEnergy": {"value": -2.239397e+32, "unit": u.Joule},
"log.initial.system.KinEnergy": {"value": 1.606047e+29, "unit": u.Joule},
"log.initial.system.DeltaTime": {"value": 0.000000, "unit": u.sec},
"log.initial.earth.Mass": {"value": 5.972186e+24, "unit": u.kg},
"log.initial.earth.Radius": {"value": 6.378100e+06, "unit": u.m},
"log.initial.earth.RadGyra": {"value": 0.500000},
"log.initial.earth.BodyType": {"value": 0.000000},
"log.initial.earth.Density": {"value": 5495.038549, "unit": u.kg / u.m ** 3},
"log.initial.earth.HZLimitDryRunaway": {"value": -1.000000, "unit": u.m},
"log.initial.earth.HZLimRecVenus": {"value": -1.000000},
"log.initial.earth.HZLimRunaway": {"value": -1.000000},
"log.initial.earth.HZLimMoistGreenhouse": {"value": -1.000000},
"log.initial.earth.HZLimMaxGreenhouse": {"value": -1.000000},
"log.initial.earth.HZLimEarlyMars": {"value": -1.000000},
"log.initial.earth.Instellation": {"value": -1.000000, "unit": u.kg / u.sec ** 3},
"log.initial.earth.MeanMotion": {"value": -1.000000, "unit": 1 / u.sec},
"log.initial.earth.OrbPeriod": {"value": -1.000000, "unit": u.sec},
"log.initial.earth.SemiMajorAxis": {"value": -1.000000, "unit": u.m},
"log.initial.earth.LXUVTot": {"value": -1.000000, "unit": u.kg / u.sec ** 3},
"log.initial.earth.SurfWaterMass": {"value": 3.000000, "unit": u.TO},
"log.initial.earth.EnvelopeMass": {"value": 0.000000, "unit": u.kg},
"log.initial.earth.OxygenMass": {"value": 0.000000, "unit": u.bar},
"log.initial.earth.RGLimit": {"value": 0.000000, "unit": u.au},
"log.initial.earth.XO": {"value": 0.333333},
"log.initial.earth.EtaO": {"value": 0.682959},
"log.initial.earth.PlanetRadius": {"value": 6.378100e+06, "unit": u.m},
"log.initial.earth.OxygenMantleMass": {"value": 0.000000, "unit": u.kg},
"log.initial.earth.RadXUV": {"value": -1.000000, "unit": u.m},
"log.initial.earth.RadSolid": {"value": -1.000000, "unit": u.m},
"log.initial.earth.PresXUV": {"value": 5.000000},
"log.initial.earth.ScaleHeight": {"value": -1.000000, "unit": u.m},
"log.initial.earth.ThermTemp": {"value": 400.000000, "unit": u.K},
"log.initial.earth.AtmGasConst": {"value": 4124.000000},
"log.initial.earth.PresSurf": {"value": -1.000000, "unit": u.Pa},
"log.initial.earth.DEnvMassDt": {"value": 0.000000, "unit": u.kg / u.sec},
"log.initial.earth.FXUV": {"value": 100.000000, "unit": u.W / u.m ** 2},
"log.initial.earth.AtmXAbsEffH2O": {"value": 0.010704},
"log.initial.earth.RocheRadius": {"value": 1.037254e+11, "unit": u.m},
"log.initial.earth.BondiRadius": {"value": 1.249016e+08, "unit": u.m},
"log.initial.earth.HEscapeRegime": {"value": 8.000000},
"log.initial.earth.RRCriticalFlux": {"value": 53.697959, "unit": u.W / u.m ** 2},
"log.initial.earth.CrossoverMass": {"value": 8.022481e-26, "unit": u.kg},
"log.initial.earth.WaterEscapeRegime": {"value": 3.000000},
"log.initial.earth.FXUVCRITDRAG": {"value": 4.977066, "unit": u.W / u.m ** 2},
"log.initial.earth.HREFFLUX": {"value": 2.578763e+18, "unit": 1 / u.m ** 2 / u.sec},
"log.initial.earth.XO2": {"value": 0.000000},
"log.initial.earth.XH2O": {"value": 1.000000},
"log.initial.earth.HDiffFlux": {"value": 1.264874e+17, "unit": 1 / u.m ** 2 / u.sec},
"log.initial.earth.HRefODragMod": {"value": 0.154711},
"log.initial.earth.KTide": {"value": 0.999908},
"log.initial.earth.RGDuration": {"value": 0.00000e+00, "unit": u.yr},
"log.final.system.Age": {"value": 3.155760e+13, "unit": u.sec},
"log.final.system.Time": {"value": 3.155760e+13, "unit": u.sec},
"log.final.system.TotAngMom": {"value": 4.416946e+33, "unit": (u.kg * u.m ** 2) / u.sec},
"log.final.system.TotEnergy": {"value": -2.237790e+32, "unit": u.Joule},
"log.final.system.PotEnergy": {"value": -2.239397e+32, "unit": u.Joule},
"log.final.system.KinEnergy": {"value": 1.606047e+29, "unit": u.Joule},
"log.final.system.DeltaTime": {"value": 3.155760e+13, "unit": u.sec},
"log.final.earth.Mass": {"value": 5.972186e+24, "unit": u.kg},
"log.final.earth.Radius": {"value": 6.378100e+06, "unit": u.m},
"log.final.earth.RadGyra": {"value": 0.500000},
"log.final.earth.BodyType": {"value": 0.000000},
"log.final.earth.Density": {"value": 5495.038549, "unit": u.kg / u.m ** 3},
"log.final.earth.HZLimitDryRunaway": {"value": -1.000000, "unit": u.m},
"log.final.earth.HZLimRecVenus": {"value": -1.000000},
"log.final.earth.HZLimRunaway": {"value": -1.000000},
"log.final.earth.HZLimMoistGreenhouse": {"value": -1.000000},
"log.final.earth.HZLimMaxGreenhouse": {"value": -1.000000},
"log.final.earth.HZLimEarlyMars": {"value": -1.000000},
"log.final.earth.Instellation": {"value": -1.000000, "unit": u.kg / u.sec ** 3},
"log.final.earth.MeanMotion": {"value": -1.000000, "unit": 1 / u.sec},
"log.final.earth.OrbPeriod": {"value": -1.000000, "unit": u.sec},
"log.final.earth.SemiMajorAxis": {"value": -1.000000, "unit": u.m},
"log.final.earth.LXUVTot": {"value": -1.000000, "unit": u.kg / u.sec ** 3},
"log.final.earth.SurfWaterMass": {"value": 2.930804, "unit": u.TO},
"log.final.earth.EnvelopeMass": {"value": 0.000000, "unit": u.kg},
"log.final.earth.OxygenMass": {"value": 5.195716, "unit": u.bar},
"log.final.earth.RGLimit": {"value": 0.000000, "unit": u.au},
"log.final.earth.XO": {"value": 0.334993},
"log.final.earth.EtaO": {"value": 0.682959},
"log.final.earth.PlanetRadius": {"value": 6.378100e+06, "unit": u.m},
"log.final.earth.OxygenMantleMass": {"value": 0.000000, "unit": u.kg},
"log.final.earth.RadXUV": {"value": -1.000000, "unit": u.m},
"log.final.earth.RadSolid": {"value": -1.000000, "unit": u.m},
"log.final.earth.PresXUV": {"value": 5.000000},
"log.final.earth.ScaleHeight": {"value": -1.000000, "unit": u.m},
"log.final.earth.ThermTemp": {"value": 400.000000, "unit": u.K},
"log.final.earth.AtmGasConst": {"value": 4124.000000},
"log.final.earth.PresSurf": {"value": -1.000000, "unit": u.Pa},
"log.final.earth.DEnvMassDt": {"value": 0.000000, "unit": u.kg / u.sec},
"log.final.earth.FXUV": {"value": 100.000000, "unit": u.W / u.m ** 2},
"log.final.earth.AtmXAbsEffH2O": {"value": 0.010704},
"log.final.earth.RocheRadius": {"value": 1.037254e+11, "unit": u.m},
"log.final.earth.BondiRadius": {"value": 1.249016e+08, "unit": u.m},
"log.final.earth.HEscapeRegime": {"value": 8.000000},
"log.final.earth.RRCriticalFlux": {"value": 53.697959, "unit": u.W / u.m ** 2},
"log.final.earth.CrossoverMass": {"value": 8.022481e-26, "unit": u.kg},
"log.final.earth.WaterEscapeRegime": {"value": 3.000000},
"log.final.earth.FXUVCRITDRAG": {"value": 4.964677, "unit": u.W / u.m ** 2},
"log.final.earth.HREFFLUX": {"value": 2.578763e+18, "unit": 1 / u.m ** 2 / u.sec},
"log.final.earth.XO2": {"value": 0.003729},
"log.final.earth.XH2O": {"value": 0.996271},
"log.final.earth.HDiffFlux": {"value": 1.261726e+17, "unit": 1 / u.m ** 2 / u.sec},
"log.final.earth.HRefODragMod": {"value": 0.153738},
"log.final.earth.KTide": {"value": 0.999908},
"log.final.earth.RGDuration": {"value": 0.00000e+00, "unit": u.yr},
}
)
class Test_WaterELimConstXUVLB15NoO2SinkBolmont16(Benchmark):
pass
4 changes: 2 additions & 2 deletions tests/Atmesc/WaterELimConstXUVLB15NoO2SinkBolmont16/vpl.in
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,5 +20,5 @@ dMinValue 1e-10 # Minimum value of eccentricity/obli
bDoForward 1 # Perform a forward evolution?
bVarDt 1 # Use variable timestepping?
dEta 0.1 # Coefficient for variable timestepping
dStopTime 1e8 # Stop time for evolution
dOutputTime 1e8 # Output timesteps (assuming in body files)
dStopTime 1e6 # Stop time for evolution
dOutputTime 1e6 # Output timesteps (assuming in body files)
31 changes: 31 additions & 0 deletions tests/Atmesc/WaterELimConstXUVLB15NoO2SinkConstXAbsEffH2O/earth.in
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,31 @@
# Planet a parameters
sName earth # Body's name
saModules atmesc # Modules to apply, exact spelling required

# Physical Properties
dMass -1 # Mass, negative -> Earth masses
dRadius -1 # Radius, negative -> Earth radii
dRotPeriod -1 # Rotation period, negative -> days
dObliquity 23.5 # Retrograde rotation
dRadGyra 0.5 # Radius of gyration (moment of inertia constant)

# ATMESC Properties
dFXUV -100 # Incident XUV flux (constant)
dXFrac 1.0 # X-Ray/XUV absorption radius (fraction of planet radius)
dSurfWaterMass -3.0 # Initial surface water (Earth oceans)
dEnvelopeMass 0 # Initial envelope mass (Earth masses)
bHaltSurfaceDesiccated 0 # Halt when dry?
bHaltEnvelopeGone 0 # Halt when evaporated?
dMinSurfWaterMass -1.e-5 # Planet is desiccated when water content drops below this (Earth oceans)
sWaterLossModel lb15
sPlanetRadiusModel none
bInstantO2Sink 0
sAtmXAbsEffH2OModel none
dAtmXAbsEffH2O 0.1

# Orbital Properties
dSemi -1 # Semi-major axis, negative -> AU
dEcc 0.0167 # Eccentricity

# Output
saOutputOrder Time -SurfWaterMass -RGLimit -OxygenMass
Loading

0 comments on commit db1a296

Please sign in to comment.