Skip to content

Commit c8ecaa1

Browse files
zingalezingale
authored
Merge branch 'development' into ppm_one_sided_second
2 parents 68702e3 + 2fb7228 commit c8ecaa1

28 files changed

+1067
-246
lines changed

Diff for: ‎CHANGES.md

+21
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,24 @@
1+
# 24.11
2+
3+
* a new well-balanced method was added to the CTU PPM solver. This
4+
does the characteristic projection only on the perturbed pressure
5+
and then adds back in the hydrostatic pressure. It can be enabled
6+
via `castro.ppm_well_balanced` (#2945)
7+
8+
* fixed a bug in the div{U} calculation for artificial viscosity on
9+
symmetry boundaries (#2983)
10+
11+
* more development on 2D spherical geometry (#2973, #2975, #2981)
12+
13+
* updates to the massive star plotting scripts (#2979)
14+
15+
* add some new checks to prevent running unsupported combinations
16+
of solvers (#2978)
17+
18+
* documentation updates (#2977)
19+
20+
* `flame_wave` can now be run in 1D (#2976)
21+
122
# 24.10
223

324
* update initial model for `subchandra` when doing ASE NSE (#2970)

Diff for: ‎Exec/radiation_tests/Rad2Tshock/ci-benchmarks/Rad2TShock-1d.out

+26-26
Original file line numberDiff line numberDiff line change
@@ -1,40 +1,40 @@
11
plotfile = plt_00025
2-
time = 0.00026083835470476599
2+
time = 0.00026084061812361401
33
variables minimum value maximum value
4-
density 5.4596904436e-13 1.2720098871e-12
5-
xmom 1.2327065163e-07 1.4587103381e-07
4+
density 5.4596904436e-13 1.2717694091e-12
5+
xmom 1.2328047524e-07 1.4591226172e-07
66
ymom 0 0
77
zmom 0 0
8-
rho_E 0.021940622286 0.039199658391
9-
rho_e 0.0068091599527 0.032333529333
10-
Temp 100.00003116 217.3564598
11-
rho_X 5.4596904436e-13 1.2720098871e-12
12-
rad 7.5662679486e-07 1.4083425328e-05
13-
pressure 0.0045394399687 0.021555686223
14-
kineng 0.0064467244856 0.015131462333
15-
soundspeed 117717.6283 173551.23637
8+
rho_E 0.021940622285 0.039194728059
9+
rho_e 0.0068091599514 0.032333527949
10+
Temp 100.00003115 217.29660609
11+
rho_X 5.4596904436e-13 1.2717694091e-12
12+
rad 7.5662675907e-07 1.4083418723e-05
13+
pressure 0.0045394399678 0.0215556853
14+
kineng 0.0064473434191 0.015131462334
15+
soundspeed 117717.62829 173527.33922
1616
Gamma_1 1.6666666667 1.6666666667
17-
MachNumber 0.60689740202 1.9999997205
18-
uplusc 269483.97988 362284.15407
19-
uminusc -66810.704217 117717.5954
20-
entropy 2458725989.8 2507243428
17+
MachNumber 0.6068973935 1.9999997207
18+
uplusc 269518.17809 362247.33487
19+
uminusc -66804.446805 117717.59541
20+
entropy 2458725989.8 2507184972.2
2121
magvort 0 0
22-
divu -20286.061564 245.50876324
23-
eint_E 12471696011 27108028479
24-
eint_e 12471696011 27108028479
25-
logden -12.26283198 -11.895509513
26-
StateErr_0 5.4596904436e-13 1.2720098871e-12
27-
StateErr_1 100.00003116 217.3564598
22+
divu -20283.836079 245.94831837
23+
eint_E 12471696009 27100563708
24+
eint_e 12471696009 27100563708
25+
logden -12.26283198 -11.895591626
26+
StateErr_0 5.4596904436e-13 1.2717694091e-12
27+
StateErr_1 100.00003115 217.29660609
2828
StateErr_2 1 1
2929
X(X) 1 1
3030
abar 1 1
31-
x_velocity 102258.95554 235435.2237
31+
x_velocity 102299.57948 235435.22371
3232
y_velocity 0 0
3333
z_velocity 0 0
34-
magvel 102258.95554 235435.2237
35-
radvel -235435.2237 114038.34336
36-
circvel 0 0.005524271728
37-
magmom 1.2327065163e-07 1.4587103381e-07
34+
magvel 102299.57948 235435.22371
35+
radvel -235435.22371 114042.93363
36+
circvel 0 0.0047841596539
37+
magmom 1.2328047524e-07 1.4591226172e-07
3838
angular_momentum_x 0 0
3939
angular_momentum_y 0 0
4040
angular_momentum_z -0 -0

Diff for: ‎Exec/science/xrb_spherical/GNUmakefile

+40
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,40 @@
1+
PRECISION = DOUBLE
2+
PROFILE = FALSE
3+
4+
DEBUG = FALSE
5+
6+
DIM = 2
7+
8+
COMP = gnu
9+
10+
USE_MPI = TRUE
11+
12+
USE_GRAV = TRUE
13+
USE_REACT = FALSE
14+
15+
USE_ROTATION = FALSE
16+
USE_DIFFUSION = FALSE
17+
18+
# define the location of the CASTRO top directory
19+
CASTRO_HOME ?= ../../..
20+
21+
USE_JACOBIAN_CACHING = TRUE
22+
USE_MODEL_PARSER = TRUE
23+
NUM_MODELS := 2
24+
25+
# This sets the EOS directory in $(MICROPHYSICS_HOME)/eos
26+
EOS_DIR := helmholtz
27+
28+
# This sets the network directory in $(MICROPHYSICS_HOME)/networks
29+
NETWORK_DIR := subch_base
30+
31+
INTEGRATOR_DIR := VODE
32+
33+
CONDUCTIVITY_DIR := stellar
34+
35+
PROBLEM_DIR ?= ./
36+
37+
Bpack := $(PROBLEM_DIR)/Make.package
38+
Blocs := $(PROBLEM_DIR)
39+
40+
include $(CASTRO_HOME)/Exec/Make.Castro

Diff for: ‎Exec/science/xrb_spherical/Make.package

+2
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,2 @@
1+
CEXE_headers += initial_model.H
2+

Diff for: ‎Exec/science/xrb_spherical/README.md

+6
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,6 @@
1+
# xrb_spherical
2+
3+
This is the full-star XRB flame setup based on flame_wave.
4+
This setup uses a spherical 2D geometry to model XRB flame
5+
on a spherical shell with initial temperature perturbation
6+
on the north pole.

Diff for: ‎Exec/science/xrb_spherical/_prob_params

+68
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,68 @@
1+
2+
dtemp real 3.81e8_rt y
3+
4+
theta_half_max real 1.745e-2_rt y
5+
6+
theta_half_width real 4.9e-3_rt y
7+
8+
# cutoff mass fraction of the first species for refinement
9+
X_min real 1.e-4_rt y
10+
11+
# do we dynamically refine based on density? or based on height?
12+
tag_by_density integer 1 y
13+
14+
# used for tagging if tag_by_density = 1
15+
cutoff_density real 500.e0_rt y
16+
17+
# used if we are refining based on height rather than density
18+
refine_height real 3600 y
19+
20+
T_hi real 5.e8_rt y
21+
22+
T_star real 1.e8_rt y
23+
24+
T_lo real 5.e7_rt y
25+
26+
dens_base real 2.e6_rt y
27+
28+
H_star real 500.e0_rt y
29+
30+
atm_delta real 25.e0_rt y
31+
32+
fuel1_name string "helium-4" y
33+
34+
fuel2_name string "" y
35+
36+
fuel3_name string "" y
37+
38+
fuel4_name string "" y
39+
40+
ash1_name string "iron-56" y
41+
42+
ash2_name string "" y
43+
44+
ash3_name string "" y
45+
46+
fuel1_frac real 1.0_rt y
47+
48+
fuel2_frac real 0.0_rt y
49+
50+
fuel3_frac real 0.0_rt y
51+
52+
fuel4_frac real 0.0_rt y
53+
54+
ash1_frac real 1.0_rt y
55+
56+
ash2_frac real 0.0_rt y
57+
58+
ash3_frac real 0.0_rt y
59+
60+
low_density_cutoff real 1.e-4_rt y
61+
62+
smallx real 1.e-10_rt y
63+
64+
r_refine_distance real 1.e30_rt y
65+
66+
max_hse_tagging_level integer 2 y
67+
68+
max_base_tagging_level integer 1 y

Diff for: ‎Exec/science/xrb_spherical/analysis/r_profile.py

+56
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,56 @@
1+
#!/usr/bin/env python3
2+
3+
# Spherical R profile at different theta
4+
5+
import os
6+
import sys
7+
import yt
8+
import matplotlib.pyplot as plt
9+
import numpy as np
10+
from functools import reduce
11+
import itertools
12+
13+
import matplotlib.ticker as ptick
14+
from yt.frontends.boxlib.api import CastroDataset
15+
from yt.units import cm
16+
17+
18+
plotfile = sys.argv[1]
19+
ds = CastroDataset(plotfile)
20+
21+
rmin = ds.domain_left_edge[0]
22+
rmax = rmin + 5000.0*cm
23+
#rmax = ds.domain_right_edge[0]
24+
print(ds.domain_left_edge[1])
25+
fig, _ax = plt.subplots(2,2)
26+
27+
axes = list(itertools.chain(*_ax))
28+
29+
fig.set_size_inches(7.0, 8.0)
30+
31+
fields = ["Temp", "density", "x_velocity", "y_velocity"]
32+
nice_names = [r"$T$ (K)", r"$\rho$ (g/${cm}^3$)", r"$u$ (cm/s)", r"$v$ (cm/s)"]
33+
34+
# 4 rays at different theta values
35+
thetal = ds.domain_left_edge[1]
36+
thetar = ds.domain_right_edge[1]
37+
thetas = [thetal, 0.25*thetar, 0.5*thetar, 0.75*thetar]
38+
39+
for i, f in enumerate(fields):
40+
41+
for theta in thetas:
42+
# simply go from (rmin, theta) -> (rmax, theta). Doesn't need to convert to physical R-Z
43+
ray = ds.ray((rmin, theta, 0*cm), (rmax, theta, 0*cm))
44+
isrt = np.argsort(ray["t"])
45+
axes[i].plot(ray['r'][isrt], ray[f][isrt], label=r"$\theta$ = {:.4f}".format(float(theta)))
46+
47+
axes[i].set_xlabel(r"$r$ (cm)")
48+
axes[i].set_ylabel(nice_names[i])
49+
axes[i].set_yscale("symlog")
50+
51+
if i == 0:
52+
axes[0].legend(frameon=False, loc="lower left")
53+
54+
#fig.set_size_inches(10.0, 9.0)
55+
plt.tight_layout()
56+
plt.savefig("{}_profiles.png".format(os.path.basename(plotfile)))

Diff for: ‎Exec/science/xrb_spherical/analysis/slice.py

+94
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,94 @@
1+
#!/usr/bin/env python3
2+
3+
import sys
4+
import os
5+
import yt
6+
import numpy as np
7+
import matplotlib.pyplot as plt
8+
from yt.frontends.boxlib.api import CastroDataset
9+
10+
from yt.units import cm
11+
12+
"""
13+
Given a plot file and field name, it gives slice plots at the top,
14+
middle, and bottom of the domain (shell).
15+
"""
16+
17+
def slice(fname:str, field:str,
18+
loc: str = "top", width_factor: float = 3.0) -> None:
19+
"""
20+
A slice plot of the dataset for Spherical2D geometry.
21+
22+
Parameter
23+
=======================
24+
fname: plot file name
25+
field: field parameter
26+
loc: location on the domain. {top, mid, bot}
27+
"""
28+
29+
ds = CastroDataset(fname)
30+
currentTime = ds.current_time.in_units("s")
31+
print(f"Current time of this plot file is {currentTime} s")
32+
33+
# Some geometry properties
34+
rr = ds.domain_right_edge[0].in_units("km")
35+
rl = ds.domain_left_edge[0].in_units("km")
36+
dr = rr - rl
37+
r_center = 0.5 * (rr + rl)
38+
39+
thetar = ds.domain_right_edge[1]
40+
thetal = ds.domain_left_edge[1]
41+
dtheta = thetar - thetal
42+
theta_center = 0.5 * (thetar + thetal)
43+
44+
# Domain width of the slice plot
45+
width = width_factor * dr
46+
box_widths = (width, width)
47+
48+
loc = loc.lower()
49+
loc_options = ["top", "mid", "bot"]
50+
51+
if loc not in loc_options:
52+
raise Exception("loc parameter must be top, mid or bot")
53+
54+
# Centers for the Top, Mid and Bot panels
55+
centers = {"top":(r_center*np.sin(thetal)+0.5*width, r_center*np.cos(thetal)),
56+
"mid":(r_center*np.sin(theta_center), r_center*np.cos(theta_center)),
57+
"bot":(r_center*np.sin(thetar)+0.5*width, r_center*np.cos(thetar))}
58+
59+
# Note we can also set center during SlicePlot, however then we would enter in [r_center, theta_center, 0]
60+
# rather than the physical R-Z coordinate if we do it via sp.set_center
61+
62+
sp = yt.SlicePlot(ds, 'phi', field, width=box_widths)
63+
sp.set_center(centers[loc])
64+
65+
sp.set_cmap(field, "viridis")
66+
if field in ["x_velocity", "y_velocity", "z_velocity"]:
67+
sp.set_cmap(field, "coolwarm")
68+
elif field == "Temp":
69+
sp.set_zlim(f, 5.e7, 2.5e9)
70+
sp.set_cmap(f, "magma_r")
71+
elif field == "enuc":
72+
sp.set_zlim(f, 1.e18, 1.e20)
73+
elif field == "density":
74+
sp.set_zlim(f, 1.e-3, 5.e8)
75+
76+
# sp.annotate_text((0.05, 0.05), f"{currentTime.in_cgs():8.5f} s")
77+
sp.save(f"{ds}_{loc}")
78+
79+
if __name__ == "__main__":
80+
81+
if len(sys.argv) < 3:
82+
raise Exception("Please enter parameters in order of: fname field_name width_factor[optional] loc[optional]")
83+
84+
fname = sys.argv[1]
85+
field = sys.argv[2]
86+
loc = "top"
87+
width_factor = 3.0
88+
89+
if len(sys.argv) == 4:
90+
width_factor = float(sys.argv[3])
91+
elif len(sys.argv) > 4:
92+
loc = sys.argv[4]
93+
94+
slice(fname, field, loc=loc, width_factor=width_factor)

Diff for: ‎Exec/science/xrb_spherical/initial_model.H

+1
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
../flame_wave/initial_model.H

0 commit comments

Comments
 (0)