Skip to content

Commit

Permalink
Merge pull request #337 from donnaaboise/app-sphere-swe
Browse files Browse the repository at this point in the history 
WIP : Add SWE on the sphere example
  • Loading branch information
@scottaiton
scottaiton authored Jul 11, 2024
2 parents 4e6f7f5 + cc84495 commit 084afac
Show file tree
Hide file tree
Showing 118 changed files with 28,573 additions and 1,670 deletions.
4 changes: 3 additions & 1 deletion applications/clawpack/advection/2d/sphere/afterframe.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,9 @@
cv([1 end]) = [];
drawcontourlines(cv);

view(vbot)
% view(vbot)
view(3)
showgridlines

MaxFrames = 64;
NoQuery = 0;
Expand Down
16 changes: 8 additions & 8 deletions applications/clawpack/advection/2d/sphere/fclaw_options.ini
View file
Original file line number Diff line number Diff line change
@@ -1,12 +1,12 @@
[user]
example = 0 # 0 = cubed_sphere; 1 = pillowsphere
example = 1 # 0 = cubed_sphere; 1 = pillowsphere
revs-per-second = 0.5 # Rigid body rotation
claw-version = 4

[clawpatch]
# Grid dimensions
mx = 8 # mx_leaf
my = 8 # my_leaf
mx = 64 # mx_leaf
my = 64 # my_leaf
mbc = 2 # Number of ghost cells

maux = 3 # Auxilliary arrays (needed by Clawpack solvers)
Expand All @@ -26,8 +26,8 @@
[Options]

# Regridding
minlevel = 0 # Minimum level
maxlevel = 5 # Maximum levels of refinement
minlevel = 1 # Minimum level
maxlevel = 1 # Maximum levels of refinement
regrid_interval = 1 # Regrid every 'regrid_interval' time steps.
refine_threshold = 0.25
coarsen_threshold = 0.05
Expand All @@ -37,8 +37,8 @@
# Time stepping
tfinal = 2.0 #

outstyle = 1 # OutputStyle = 1,2,3 (see below, e.g. [OutputStyle_1]
nout = 16 # Used for all three out styles; has different meaning, though.
outstyle = 3 # OutputStyle = 1,2,3 (see below, e.g. [OutputStyle_1]
nout = 1 # Used for all three out styles; has different meaning, though.
nstep = 1 # Only used if outstyle is 3

use_fixed_dt = F # Take a fixed time step
Expand All @@ -49,7 +49,7 @@
subcycle = T # subcycle in time ( = T or F).

# Output
verbosity = essential # verbose
verbosity = production # verbose
output = T # Set this to F for parallel runs. T is dead slow.

# Diagnostics and debugging
Expand Down
2 changes: 1 addition & 1 deletion applications/clawpack/advection/2d/sphere/mapc2m.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@
global notpillowsphere;

map = 'cubedsphere';
% map = 'pillowsphere';
map = 'pillowsphere';

b = getblocknumber();

Expand Down
1 change: 1 addition & 0 deletions applications/clawpack/clawpack.apps
View file
Original file line number Diff line number Diff line change
Expand Up @@ -90,6 +90,7 @@ include applications/clawpack/burgers/2d/pwconst/Makefile.am
## Shallow
include applications/clawpack/shallow/2d/radialdam/Makefile.am
include applications/clawpack/shallow/2d/bump/Makefile.am
include applications/clawpack/shallow/2d/sphere/Makefile.am

## Euler
include applications/clawpack/euler/2d/shockbubble/Makefile.am
Expand Down
1 change: 1 addition & 0 deletions applications/clawpack/clawpack.cmake
View file
Original file line number Diff line number Diff line change
Expand Up @@ -91,6 +91,7 @@ add_subdirectory(clawpack/burgers/2d/pwconst)
## Shallow
add_subdirectory(clawpack/shallow/2d/radialdam)
add_subdirectory(clawpack/shallow/2d/bump)
add_subdirectory(clawpack/shallow/2d/sphere)

## Euler
add_subdirectory(clawpack/euler/2d/shockbubble)
Expand Down
51 changes: 51 additions & 0 deletions applications/clawpack/shallow/2d/rp/rpn2cons_update.f
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,51 @@
c # This routine computes the flux of cell-centered
c #
c # All scaling and potential sign changes are handled by the
c # calling routine.
c

subroutine rpn2_cons_update(meqn,maux, idir, iface, q,
& auxvec_center,
& auxvec_edge,flux)
implicit none

integer meqn,maux,idir, iface
double precision q(meqn), flux(meqn)

c # maux == 0 in this routine.
double precision auxvec_center(maux), auxvec_edge(maux)

double precision grav
common /cparam/ grav

double precision hun, hut, h, un, ut

integer mu, mv

c # f1 = (hu; hu^2 + 0.5*g*h^2; huv)
c # f2 = (hv; huv; hv^2 + 0.5*gh^2)

if (idir .eq. 0) then
mu = 2
mv = 3
else
mu = 3
mv = 2
endif

hun = q(mu)
hut = q(mv)

h = q(1)
un = hun/h
ut = hut/h

flux(1) = hun
flux(mu) = hun**2/h + 0.5*grav*h**2
flux(mv) = un*hut
c !! This should only be set if we are solving for the manifold case.
c flux(4) = un*q(4)

end


73 changes: 73 additions & 0 deletions applications/clawpack/shallow/2d/rp/rpn2cons_update_manifold.f
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,73 @@
c # This routine rotates cell-centered data into frame specified
c # by vectors at the edge. The flux is then computed from this
c # rotated data.
c #
c # All scaling and potential sign changes are handled by the
c # calling routine.
c
c # Note that edge tangents must be normalized.

subroutine rpn2cons_update_manifold(meqn,maux, idir, iface, q,
& auxvec_center,
& auxvec_edge,flux)
implicit none

integer meqn,maux,idir, iface
double precision q(meqn), flux(meqn)
double precision auxvec_center(maux), auxvec_edge(maux)

double precision grav
common /swe_model_parms/ grav

double precision enx, eny, enz
double precision etx, ety, etz
double precision hun, hut, h, un, ut, gamma, f(4)
double precision b, psi
integer ioff

c # f1 = (hu; hu^2 + 0.5*g*h^2; huv)
c # f2 = (hv; huv; hv^2 + 0.5*gh^2)

if (idir .eq. 0) then
ioff = 1
else
ioff = 7
endif

enx = auxvec_edge(ioff+1)
eny = auxvec_edge(ioff+2)
enz = auxvec_edge(ioff+3)
etx = auxvec_edge(ioff+4)
ety = auxvec_edge(ioff+5)
etz = auxvec_edge(ioff+6)

c !! Normalize the edge lengths
gamma = dsqrt(etx**2 + ety**2 + etz**2)
etx = etx / gamma
ety = ety / gamma
etz = etz / gamma

hun = enx*q(2) + eny*q(3) + enz*q(4)
hut = etx*q(2) + ety*q(3) + etz*q(4)

h = q(1)
un = hun/h
ut = hut/h

b = auxvec_edge(18)
psi = grav*h*b


f(1) = hun
f(2) = hun**2/h + 0.5*grav*h**2 + psi
f(3) = un*hut
f(4) = un*q(4)

flux(1) = hun
flux(2) = enx*f(2) + etx*f(3)
flux(3) = eny*f(2) + ety*f(3)
flux(4) = enz*f(2) + etz*f(3)

end


19 changes: 19 additions & 0 deletions applications/clawpack/shallow/2d/rp/rpn2cons_update_zero.f
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,19 @@
subroutine rpn2_cons_update_zero(meqn,maux, idir, iface, q,
& auxvec_center,
& auxvec_edge,flux)

implicit none

integer meqn,maux,idir, iface
double precision q(meqn), flux(meqn)
double precision auxvec_center(maux), auxvec_edge(maux)
integer m

c # f(q) = (n dot u)*q
do m = 1,meqn
c # No flux function available for equations in non-conservative form
flux(m) = 0
enddo

end

73 changes: 73 additions & 0 deletions applications/clawpack/shallow/2d/sphere/1d_latitude/Makefile
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,73 @@

# Makefile for Clawpack code in this directory.
# This version only sets the local files and frequently changed
# options, and then includes the standard makefile pointed to by CLAWMAKE.
CLAWMAKE = $(CLAW)/clawutil/src/Makefile.common

# See the above file for details and a list of make options, or type
# make .help
# at the unix prompt.


# Adjust these variables if desired:
# ----------------------------------

CLAW_PKG = geoclaw # Clawpack package to use
EXE = xgeo # Executable to create
SETRUN_FILE = setrun.py # File containing function to make data
OUTDIR = _output # Directory for output
SETPLOT_FILE = setplot.py # File containing function to set plots
PLOTDIR = _plots # Directory for plots

OVERWRITE ?= True # False ==> make a copy of OUTDIR first
RESTART ?= False # Should = clawdata.restart in setrun

# Environment variable FC should be set to fortran compiler, e.g. gfortran

# Compiler flags can be specified here or set as an environment variable
FFLAGS ?=

# ---------------------------------
# package sources for this program:
# ---------------------------------

GEOLIB = $(CLAW)/geoclaw/src/1d_classic/shallow
include $(GEOLIB)/Makefile.geoclaw_1d_classic

# ---------------------------------------
# package sources specifically to exclude
# (i.e. if a custom replacement source
# under a different name is provided)
# ---------------------------------------

EXCLUDE_MODULES = \

EXCLUDE_SOURCES = \
setprob.f90 \

# ----------------------------------------
# List of custom sources for this program:
# ----------------------------------------

MODULES = \

SOURCES = \
../setprob.f90 \
qinit.f90 \


#-------------------------------------------------------------------
# Include Makefile containing standard definitions and make options:
include $(CLAWMAKE)

# Construct the topography data
.PHONY: topo all

topo:
python makegrid.py

all:
$(MAKE) topo
$(MAKE) .plots
$(MAKE) .htmls

50 changes: 50 additions & 0 deletions applications/clawpack/shallow/2d/sphere/1d_latitude/_reference/claw.data
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,50 @@
########################################################
### DO NOT EDIT THIS FILE: GENERATED AUTOMATICALLY ####
### To modify data, edit setrun.py ####
### and then "make .data" ####
########################################################

1 =: num_dim
-90.0 =: lower
90.0 =: upper
3600 =: num_cells

2 =: num_eqn
2 =: num_waves
2 =: num_aux

0.0 =: t0

1 =: output_style
4 =: num_output_times
4.0 =: tfinal
T =: output_t0

1 =: output_format
1 1 =: iout_q
1 1 =: iout_aux
T =: output_aux_onlyonce

0.01 =: dt_initial
1000000000.0 =: dt_max
1.0 =: cfl_max
0.75 =: cfl_desired
50000 =: steps_max

T =: dt_variable
2 =: order
0 =: verbosity
1 =: source_split
2 =: capa_index
T =: use_fwaves

4 4 =: limiter

2 =: num_ghost
3 =: bc_lower
3 =: bc_upper

F =: restart
'fort.q0006' =: restart_file
0 =: checkpt_style

10 changes: 10 additions & 0 deletions applications/clawpack/shallow/2d/sphere/1d_latitude/_reference/dtopo.data
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,10 @@
########################################################
### DO NOT EDIT THIS FILE: GENERATED AUTOMATICALLY ####
### To modify data, edit setrun.py ####
### and then "make .data" ####
########################################################

0 =: mdtopofiles


1e+99 =: dt_max_dtopo
Loading

0 comments on commit 084afac

Please sign in to comment.