Skip to content

Commit

Permalink
Merge remote-tracking branch 'origin/develop' into dempsey/dust_update
Browse files Browse the repository at this point in the history
  • Loading branch information
@adamdempsey90
adamdempsey90 committed Jan 13, 2025
2 parents 5346481 + c98a205 commit 72984f6
Show file tree
Hide file tree
Showing 9 changed files with 137 additions and 115 deletions.
2 changes: 1 addition & 1 deletion inputs/diffusion/conduction.in
View file
Original file line number Diff line number Diff line change
Expand Up @@ -71,7 +71,7 @@ dfloor = 1.0e-10
siefloor = 1.0e-15

<gas/conductivity>
type = constant
type = conductivity
cond = 1.0e-1

<gravity/uniform>
Expand Down
2 changes: 1 addition & 1 deletion inputs/diffusion/gaussian_bump.in
View file
Original file line number Diff line number Diff line change
Expand Up @@ -80,7 +80,7 @@ type = constant
nu = 5.0e-3

<gas/conductivity>
type = constant
type = conductivity
cond = 5.0e-3

<problem>
Expand Down
10 changes: 5 additions & 5 deletions src/drag/drag.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -110,8 +110,8 @@ TaskStatus DragSource(MeshData<Real> *md, const Real time, const Real dt) {
auto &gas_pkg = pm->packages.Get("gas");
const auto &dp = gas_pkg->template Param<Diffusion::DiffCoeffParams>("visc_params");
const auto &eos_d = gas_pkg->template Param<EOS>("eos_d");
if (dp.type == Diffusion::DiffType::viscosity_const) {
return SelfDragSourceImpl<Diffusion::DiffType::viscosity_const, GEOM>(
if (dp.type == Diffusion::DiffType::viscosity_plaw) {
return SelfDragSourceImpl<Diffusion::DiffType::viscosity_plaw, GEOM>(
md, time, dt, dp, eos_d, gas_self_par, dust_self_par);
} else if (dp.type == Diffusion::DiffType::viscosity_alpha) {
return SelfDragSourceImpl<Diffusion::DiffType::viscosity_alpha, GEOM>(
Expand All @@ -134,13 +134,13 @@ TaskStatus DragSource(MeshData<Real> *md, const Real time, const Real dt) {
drag_pkg->template Param<StoppingTimeParams>("stopping_time_params");
if (gas_self_par.damp_to_visc) {
auto &dp = gas_pkg->template Param<Diffusion::DiffCoeffParams>("visc_params");
if (dp.type == Diffusion::DiffType::viscosity_const) {
if (dp.type == Diffusion::DiffType::viscosity_plaw) {
if (stop_par.model == DragModel::constant) {
return SimpleDragSourceImpl<Diffusion::DiffType::viscosity_const,
return SimpleDragSourceImpl<Diffusion::DiffType::viscosity_plaw,
DragModel::constant, GEOM>(
md, time, dt, dp, eos_d, gas_self_par, dust_self_par, stop_par);
} else if (stop_par.model == DragModel::stokes) {
return SimpleDragSourceImpl<Diffusion::DiffType::viscosity_const,
return SimpleDragSourceImpl<Diffusion::DiffType::viscosity_plaw,
DragModel::stokes, GEOM>(
md, time, dt, dp, eos_d, gas_self_par, dust_self_par, stop_par);
}
Expand Down
30 changes: 15 additions & 15 deletions src/gas/gas.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -436,9 +436,9 @@ Real EstimateTimestepMesh(MeshData<Real> *md) {
const auto do_viscosity = params.template Get<bool>("do_viscosity");
if (do_viscosity) {
auto dp = params.template Get<Diffusion::DiffCoeffParams>("visc_params");
if (dp.type == Diffusion::DiffType::viscosity_const) {
if (dp.type == Diffusion::DiffType::viscosity_plaw) {
visc_dt = Diffusion::EstimateTimestep<GEOM, Fluid::gas,
Diffusion::DiffType::viscosity_const>(
Diffusion::DiffType::viscosity_plaw>(
md, dp, gas_pkg, eos_d, vmesh);
} else if (dp.type == Diffusion::DiffType::viscosity_alpha) {
visc_dt = Diffusion::EstimateTimestep<GEOM, Fluid::gas,
Expand All @@ -451,13 +451,13 @@ Real EstimateTimestepMesh(MeshData<Real> *md) {
const auto do_conduction = params.template Get<bool>("do_conduction");
if (do_conduction) {
auto dp = params.template Get<Diffusion::DiffCoeffParams>("cond_params");
if (dp.type == Diffusion::DiffType::conductivity_const) {
if (dp.type == Diffusion::DiffType::conductivity_plaw) {
cond_dt = Diffusion::EstimateTimestep<GEOM, Fluid::gas,
Diffusion::DiffType::conductivity_const>(
Diffusion::DiffType::conductivity_plaw>(
md, dp, gas_pkg, eos_d, vmesh);
} else if (dp.type == Diffusion::DiffType::thermaldiff_const) {
} else if (dp.type == Diffusion::DiffType::thermaldiff_plaw) {
cond_dt = Diffusion::EstimateTimestep<GEOM, Fluid::gas,
Diffusion::DiffType::thermaldiff_const>(
Diffusion::DiffType::thermaldiff_plaw>(
md, dp, gas_pkg, eos_d, vmesh);
}
}
Expand Down Expand Up @@ -543,10 +543,10 @@ TaskStatus ViscousFlux(MeshData<Real> *md) {

if (dp.type == Diffusion::DiffType::null) {
return TaskStatus::complete;
} else if (dp.type == Diffusion::DiffType::viscosity_const) {
} else if (dp.type == Diffusion::DiffType::viscosity_plaw) {
return Diffusion::MomentumFluxImpl<GEOM, Fluid::gas,
Diffusion::DiffType::viscosity_const>(md, dp, pkg,
vprim, vf);
Diffusion::DiffType::viscosity_plaw>(md, dp, pkg,
vprim, vf);
} else if (dp.type == Diffusion::DiffType::viscosity_alpha) {
return Diffusion::MomentumFluxImpl<GEOM, Fluid::gas,
Diffusion::DiffType::viscosity_alpha>(md, dp, pkg,
Expand Down Expand Up @@ -581,14 +581,14 @@ TaskStatus ThermalFlux(MeshData<Real> *md) {

if (dp.type == Diffusion::DiffType::null) {
return TaskStatus::complete;
} else if (dp.type == Diffusion::DiffType::conductivity_const) {
return Diffusion::ThermalFluxImpl<GEOM, Fluid::gas,
Diffusion::DiffType::conductivity_const>(
md, dp, pkg, vprim, vf);
} else if (dp.type == Diffusion::DiffType::thermaldiff_const) {
} else if (dp.type == Diffusion::DiffType::conductivity_plaw) {
return Diffusion::ThermalFluxImpl<GEOM, Fluid::gas,
Diffusion::DiffType::thermaldiff_const>(md, dp, pkg,
Diffusion::DiffType::conductivity_plaw>(md, dp, pkg,
vprim, vf);
} else if (dp.type == Diffusion::DiffType::thermaldiff_plaw) {
return Diffusion::ThermalFluxImpl<GEOM, Fluid::gas,
Diffusion::DiffType::thermaldiff_plaw>(md, dp, pkg,
vprim, vf);
} else {
PARTHENON_FAIL("Invalid conductivity type");
}
Expand Down
48 changes: 22 additions & 26 deletions src/gas/params.yaml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -129,14 +129,13 @@ gas:
_description: "Viscosity models"
type:
_type: string
_default: constant
_description: "Viscosity model"
alpha:
_type: opt
_description: "Alpha viscosity"
constant:
powerlaw:
_type: opt
_description: "Constant kinematic viscosity"
_description: "Powerlaw kinematic viscosity"
averaging:
_type: string
_default: "arithmetic"
Expand All @@ -154,6 +153,9 @@ gas:
_type: Real
_description: "Value of the kinematic shear viscosity"
_units: "cm^2 s^-1"
r_exp:
_type: Real
_description: "Value of the cylindrical powerlaw exponent"
eta_bulk:
_type: Real
_default: "0.0"
Expand All @@ -164,14 +166,13 @@ gas:
_description: "Conductivity models"
type:
_type: string
_default: constant
_description: "Conductivity model"
constant:
conductivity:
_type: opt
_description: "Constant thermal conductivity"
diffusivity_constant:
_description: "Powerlaw thermal conductivity"
diffusivity:
_type: opt
_description: "Constant thermal diffusivity"
_description: "Powerlaw thermal diffusivity"
averaging:
_type: string
_default: "arithmetic"
Expand All @@ -191,26 +192,21 @@ gas:
_type: Real
_description: "Value of heat conductivity"
_units: "erg cm^-1 s^-1 K^-1"

averaging:
_type: string
_default: "arithmetic"
_description: "Averaging method for the face diffusion coefficient"
arithmetic:
_type: opt
_description: "Simple average (L+R)/2"
harmonic:
_type: opt
_description: "Harmonic average 2*(L*R)/(L+R)"

kappa:
temp_exp:
_type: Real
_description: "Value of thermal diffusivity"
_units: "cm^2 s^-1"
cond:
_description: "Temperature powerlaw exponent"
T_ref:
_type: Real
_description: "Value of heat conductivity"
_units: "erg cm^-1 s^-1 K^-1"
_description: "Reference temperature"
_units: "K"
rho_exp:
_type: Real
_description: "Density powerlaw exponent"
rho_ref:
_type: Real
_description: "Reference density"
_units: "g cm^-3"


opacity/absorption:
_type: node
Expand Down
48 changes: 24 additions & 24 deletions src/pgen/conduction.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -260,67 +260,67 @@ inline void CondBoundary(std::shared_ptr<MeshBlockData<Real>> &mbd, bool coarse)
const auto dcp = pkg->template Param<Diffusion::DiffCoeffParams>("cond_params");

if constexpr (BDY == IndexDomain::inner_x1) {
if (dcp.type == Diffusion::DiffType::conductivity_const) {
if (dcp.type == Diffusion::DiffType::conductivity_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::inner_x1,
Diffusion::DiffType::conductivity_const>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_const) {
Diffusion::DiffType::conductivity_plaw>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::inner_x1,
Diffusion::DiffType::thermaldiff_const>(mbd, coarse);
Diffusion::DiffType::thermaldiff_plaw>(mbd, coarse);
} else {
PARTHENON_FAIL(
"Chosen conductivity type is not compatible with conductivity boundaries");
}
} else if constexpr (BDY == IndexDomain::outer_x1) {
if (dcp.type == Diffusion::DiffType::conductivity_const) {
if (dcp.type == Diffusion::DiffType::conductivity_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::outer_x1,
Diffusion::DiffType::conductivity_const>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_const) {
Diffusion::DiffType::conductivity_plaw>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::outer_x1,
Diffusion::DiffType::thermaldiff_const>(mbd, coarse);
Diffusion::DiffType::thermaldiff_plaw>(mbd, coarse);
} else {
PARTHENON_FAIL(
"Chosen conductivity type is not compatible with conductivity boundaries");
}
} else if constexpr (BDY == IndexDomain::inner_x2) {
if (dcp.type == Diffusion::DiffType::conductivity_const) {
if (dcp.type == Diffusion::DiffType::conductivity_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::inner_x2,
Diffusion::DiffType::conductivity_const>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_const) {
Diffusion::DiffType::conductivity_plaw>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::inner_x2,
Diffusion::DiffType::thermaldiff_const>(mbd, coarse);
Diffusion::DiffType::thermaldiff_plaw>(mbd, coarse);
} else {
PARTHENON_FAIL(
"Chosen conductivity type is not compatible with conductivity boundaries");
}
} else if constexpr (BDY == IndexDomain::outer_x2) {
if (dcp.type == Diffusion::DiffType::conductivity_const) {
if (dcp.type == Diffusion::DiffType::conductivity_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::outer_x2,
Diffusion::DiffType::conductivity_const>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_const) {
Diffusion::DiffType::conductivity_plaw>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::outer_x2,
Diffusion::DiffType::thermaldiff_const>(mbd, coarse);
Diffusion::DiffType::thermaldiff_plaw>(mbd, coarse);
} else {
PARTHENON_FAIL(
"Chosen conductivity type is not compatible with conductivity boundaries");
}
} else if constexpr (BDY == IndexDomain::inner_x3) {
if (dcp.type == Diffusion::DiffType::conductivity_const) {
if (dcp.type == Diffusion::DiffType::conductivity_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::inner_x3,
Diffusion::DiffType::conductivity_const>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_const) {
Diffusion::DiffType::conductivity_plaw>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::inner_x3,
Diffusion::DiffType::thermaldiff_const>(mbd, coarse);
Diffusion::DiffType::thermaldiff_plaw>(mbd, coarse);
} else {
PARTHENON_FAIL(
"Chosen conductivity type is not compatible with conductivity boundaries");
}
} else if constexpr (BDY == IndexDomain::outer_x3) {
if (dcp.type == Diffusion::DiffType::conductivity_const) {
if (dcp.type == Diffusion::DiffType::conductivity_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::outer_x3,
Diffusion::DiffType::conductivity_const>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_const) {
Diffusion::DiffType::conductivity_plaw>(mbd, coarse);
} else if (dcp.type == Diffusion::DiffType::thermaldiff_plaw) {
return CondBoundaryImpl<GEOM, IndexDomain::outer_x3,
Diffusion::DiffType::thermaldiff_const>(mbd, coarse);
Diffusion::DiffType::thermaldiff_plaw>(mbd, coarse);
} else {
PARTHENON_FAIL(
"Chosen conductivity type is not compatible with conductivity boundaries");
Expand Down
6 changes: 3 additions & 3 deletions src/pgen/disk.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -312,11 +312,11 @@ inline void InitDiskParams(MeshBlock *pmb, ParameterInput *pin) {
disk_params.nu0 = disk_params.alpha * disk_params.gamma_gas *
SQR(disk_params.h0 * disk_params.r0 * disk_params.Omega0);
disk_params.nu_indx = 1.5 + disk_params.q;
} else if (vtype == "constant") {
} else if (vtype == "powerlaw") {
disk_params.nu0 = pin->GetReal("gas/viscosity", "nu");
disk_params.nu_indx = 0.0;
disk_params.nu_indx = pin->GetOrAddReal("gas/viscosity", "r_exp", 0.0);
} else {
PARTHENON_FAIL("Disk pgen is only compatible with alpha or constant viscosity");
PARTHENON_FAIL("Disk pgen is only compatible with alpha or powerlaw viscosity");
}
if (pin->DoesParameterExist("problem", "mdot")) {
disk_params.mdot = pin->GetReal("problem", "mdot");
Expand Down
4 changes: 2 additions & 2 deletions src/utils/diffusion/diffusion.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -94,9 +94,9 @@ Real EstimateTimestep(MeshData<Real> *md, DiffCoeffParams &dp, PKG &pkg, const E
// Get the maximum diffusion coefficient (if there's more than one)
DiffusionCoeff<DIFF, GEOM, FLUID_TYPE> diffcoeff;
Real mu = diffcoeff.Get(dp, coords, dens, sie, eos);
if constexpr (DIFF == DiffType::conductivity_const) {
if constexpr (DIFF == DiffType::conductivity_plaw) {
mu /= (dens * eos.SpecificHeatFromDensityInternalEnergy(dens, sie));
} else if constexpr ((DIFF == DiffType::viscosity_const) ||
} else if constexpr ((DIFF == DiffType::viscosity_plaw) ||
(DIFF == DiffType::viscosity_alpha)) {
mu *= (1.0 + (dp.eta > 1.0) * (dp.eta - 1.0)) / dens;
}
Expand Down
Loading

0 comments on commit 72984f6

Please sign in to comment.