Skip to content

Commit

Permalink
add a spherical advection problem. temporarily disable aritficial vis…
Browse files Browse the repository at this point in the history
…cosity for sphericalpolar
  • Loading branch information
zhichen3 committed Aug 1, 2024
1 parent 4a13b12 commit 893c670
Show file tree
Hide file tree
Showing 3 changed files with 120 additions and 4 deletions.
39 changes: 39 additions & 0 deletions pyro/compressible/problems/inputs.spherical_advect
Original file line number Diff line number Diff line change
@@ -0,0 +1,39 @@
[driver]
max_steps = 500
tmax = 1.0

init_tstep_factor = 1.0
fix_dt = 0.005


[compressible]
limiter = 0
cvisc = 0.1
riemann = CGF

[io]
basename = spherical_advect_


[eos]
gamma = 1.4


[mesh]
grid_type = SphericalPolar
nx = 64
ny = 64
xmin = 1.0
xmax = 2.0
ymin = 0.523
ymax = 2.617

xlboundary = periodic
xrboundary = periodic

ylboundary = periodic
yrboundary = periodic


[vis]
dovis = 1
75 changes: 75 additions & 0 deletions pyro/compressible/problems/spherical_advect.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,75 @@
import sys

import numpy as np

from pyro.mesh import patch
from pyro.util import msg


def init_data(my_data, rp):
""" initialize a smooth advection problem for testing convergence """

msg.bold("initializing the spherical advect problem...")

# make sure that we are passed a valid patch object
if not isinstance(my_data, patch.CellCenterData2d):
print("ERROR: patch invalid in advect.py")
print(my_data.__class__)
sys.exit()

# get the density, momenta, and energy as separate variables
dens = my_data.get_var("density")
xmom = my_data.get_var("x-momentum")
ymom = my_data.get_var("y-momentum")
ener = my_data.get_var("energy")

# initialize the components, remember, that ener here is rho*eint
# + 0.5*rho*v**2, where eint is the specific internal energy
# (erg/g)
dens[:, :] = 1.0
xmom[:, :] = 0.0
ymom[:, :] = 0.0

gamma = rp.get_param("eos.gamma")

xmin = rp.get_param("mesh.xmin")
xmax = rp.get_param("mesh.xmax")

ymin = rp.get_param("mesh.ymin")
ymax = rp.get_param("mesh.ymax")

# Actual x- and y- coordinate
myg = my_data.grid

x = myg.scratch_array()
y = myg.scratch_array()

xctr = 0.5*(xmin + xmax) * np.sin((ymin + ymax) * 0.25)
yctr = 0.5*(xmin + xmax) * np.cos((ymin + ymax) * 0.25)

x[:, :] = myg.x2d.v(buf=myg.ng) * np.sin(myg.y2d.v(buf=myg.ng))
y[:, :] = myg.x2d.v(buf=myg.ng) * np.cos(myg.y2d.v(buf=myg.ng))

# Initial density on the very top.

# this is identical to the advection/smooth problem
dens[:, :] = 1.0 + np.exp(-60.0*((x-xctr)**2 +
(y-yctr)**2))

# velocity in theta direction.
u = 0.0
v = 3.0

xmom[:, :] = dens[:, :]*u
ymom[:, :] = dens[:, :]*v

# pressure is constant
p = 1.0
ener[:, :] = p/(gamma - 1.0) + 0.5*(xmom[:, :]**2 + ymom[:, :]**2)/dens[:, :]


def finalize():
""" print out any information to the user at the end of the run """

print("""
""")
10 changes: 6 additions & 4 deletions pyro/compressible/simulation.py
Original file line number Diff line number Diff line change
Expand Up @@ -280,11 +280,13 @@ def evolve(self):
self.solid, self.tc)

# Apply artificial viscosity to fluxes
q = cons_to_prim(self.cc_data.data, gamma, self.ivars, myg)

F_x, F_y = flx.apply_artificial_viscosity(F_x, F_y, q,
self.cc_data, self.rp,
self.ivars)
if myg.coord_type == 0:
q = cons_to_prim(self.cc_data.data, gamma, self.ivars, myg)

F_x, F_y = flx.apply_artificial_viscosity(F_x, F_y, q,
self.cc_data, self.rp,
self.ivars)

old_dens = dens.copy()
old_xmom = xmom.copy()
Expand Down

0 comments on commit 893c670

Please sign in to comment.