Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP
Browse files

Merge branch 'optional-args'

  • Loading branch information...
commit a4c4d87905f35e0170975e4ea314d31524a933a8 2 parents eeea741 + 29a0a6a
@ketch ketch authored
View
26 src/fortran/1d/classic/step1.f
@@ -34,18 +34,24 @@ subroutine step1(maxmx,meqn,mwaves,mbc,maux,mx, q,aux,dx,dt,
c --------------------------------------------------------------------
c
implicit double precision (a-h,o-z)
- dimension q(meqn,1-mbc:maxmx+mbc)
+ double precision q(meqn,1-mbc:maxmx+mbc)
+ double precision aux(maux,1-mbc:maxmx+mbc)
+ double precision f(meqn,1-mbc:maxmx+mbc)
+ double precision s(mwaves,1-mbc:maxmx+mbc)
+ double precision wave(meqn, mwaves,1-mbc:maxmx+mbc)
+ double precision amdq(meqn,1-mbc:maxmx+mbc)
+ double precision apdq(meqn,1-mbc:maxmx+mbc)
+ double precision dtdx(1-mbc:maxmx+mbc)
+ integer method(7),mthlim(mwaves)
+ logical limit
+
cf2py intent(in,out) q
cf2py intent(out) cfl
- dimension aux(maux,1-mbc:maxmx+mbc)
- dimension f(meqn,1-mbc:maxmx+mbc)
- dimension s(mwaves,1-mbc:maxmx+mbc)
- dimension wave(meqn, mwaves,1-mbc:maxmx+mbc)
- dimension amdq(meqn,1-mbc:maxmx+mbc)
- dimension apdq(meqn,1-mbc:maxmx+mbc)
- dimension dtdx(1-mbc:maxmx+mbc)
- dimension method(7),mthlim(mwaves)
- logical limit
+cf2py intent(in) meqn
+cf2py intent(in) mbc
+cf2py intent(in) maxmx
+cf2py optional f, amdq, apdq, dtdx, s, wave
+
c
c # check if any limiters are used:
limit = .false.
View
36 src/fortran/2d/classic/dimsp2.f
@@ -22,25 +22,27 @@ subroutine dimsp2(maxm,maxmx,maxmy,meqn,mwaves,maux,mbc,mx,my,
c # conditions are handled properly.
c
implicit double precision (a-h,o-z)
-cf2py intent(in,out) cfl
external rpn2,rpt2
- dimension qold(meqn, 1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
- dimension qnew(meqn, 1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
-cf2py intent(in,out) qnew
- dimension q1d(meqn, 1-mbc:maxm+mbc)
- dimension cflv(4)
- dimension qadd(meqn, 1-mbc:maxm+mbc)
- dimension fadd(meqn, 1-mbc:maxm+mbc)
- dimension gadd(meqn, 2, 1-mbc:maxm+mbc)
- dimension aux(maux,1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
- dimension aux1(maux,1-mbc:maxm+mbc)
- dimension aux2(maux,1-mbc:maxm+mbc)
- dimension aux3(maux,1-mbc:maxm+mbc)
+ double precision qold(meqn, 1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision qnew(meqn, 1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision q1d(meqn, 1-mbc:maxm+mbc)
+ double precision cflv(4)
+ double precision qadd(meqn, 1-mbc:maxm+mbc)
+ double precision fadd(meqn, 1-mbc:maxm+mbc)
+ double precision gadd(meqn, 2, 1-mbc:maxm+mbc)
+ double precision aux(maux,1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision aux1(maux,1-mbc:maxm+mbc)
+ double precision aux2(maux,1-mbc:maxm+mbc)
+ double precision aux3(maux,1-mbc:maxm+mbc)
+
+ double precision dtdx1d(1-mbc:maxm+mbc)
+ double precision dtdy1d(1-mbc:maxm+mbc)
+ integer method(7),mthlim(mwaves)
+ double precision work(mwork)
- dimension dtdx1d(1-mbc:maxm+mbc)
- dimension dtdy1d(1-mbc:maxm+mbc)
- dimension method(7),mthlim(mwaves)
- dimension work(mwork)
+cf2py intent(in,out) cfl
+cf2py intent(in,out) qnew
+cf2py optional q1d, qadd, fadd, gadd, dtdx1d, dtdy1d
c
c # If method(3) = -1, take a full time step in x.
View
32 src/fortran/2d/classic/step2.f
@@ -18,23 +18,25 @@ subroutine step2(maxm,maxmx,maxmy,meqn,mwaves,maux,mbc,mx,my,
c
implicit double precision (a-h,o-z)
external rpn2,rpt2
- dimension qold(meqn,1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
- dimension qnew(meqn,1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision qold(meqn,1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision qnew(meqn,1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision q1d(meqn,1-mbc:maxm+mbc)
+ double precision qadd(meqn,1-mbc:maxm+mbc)
+ double precision fadd(meqn,1-mbc:maxm+mbc)
+ double precision gadd(meqn, 2, 1-mbc:maxm+mbc)
+ double precision aux(maux, 1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
+ double precision aux1(maux, 1-mbc:maxm+mbc)
+ double precision aux2(maux, 1-mbc:maxm+mbc)
+ double precision aux3(maux, 1-mbc:maxm+mbc)
+
+ double precision dtdx1d(1-mbc:maxm+mbc)
+ double precision dtdy1d(1-mbc:maxm+mbc)
+ integer method(7),mthlim(mwaves)
+ double precision work(mwork)
+
cf2py intent(in,out) cfl
cf2py intent(in,out) qnew
- dimension q1d(meqn,1-mbc:maxm+mbc)
- dimension qadd(meqn,1-mbc:maxm+mbc)
- dimension fadd(meqn,1-mbc:maxm+mbc)
- dimension gadd(meqn, 2, 1-mbc:maxm+mbc)
- dimension aux(maux, 1-mbc:maxmx+mbc, 1-mbc:maxmy+mbc)
- dimension aux1(maux, 1-mbc:maxm+mbc)
- dimension aux2(maux, 1-mbc:maxm+mbc)
- dimension aux3(maux, 1-mbc:maxm+mbc)
-
- dimension dtdx1d(1-mbc:maxm+mbc)
- dimension dtdy1d(1-mbc:maxm+mbc)
- dimension method(7),mthlim(mwaves)
- dimension work(mwork)
+cf2py optional q1d, qadd, fadd, gadd, dtdx1d, dtdy1d
common /comxyt/ dtcom,dxcom,dycom,tcom,icom,jcom
c
View
28 src/pyclaw/clawpack.py
@@ -309,12 +309,6 @@ def set_fortran_parameters(self,solutions):
import numpy as np
state = solutions['n'].state
- grid = state.grid
-
- # Number of equations
- meqn,maux,mwaves,mbc = state.meqn,state.maux,self.mwaves,self.mbc
- mx = grid.ng[0]
-
self.method =np.ones(7, dtype=int)
self.method[0] = self.dt_variable
@@ -326,19 +320,13 @@ def set_fortran_parameters(self,solutions):
self.method[5] = 0
else:
self.method[5] = 1
- self.method[6] = maux # aux
+ self.method[6] = state.maux # aux
if self.fwave:
import classic1fw as classic1
else:
import classic1
-
state.set_cparam(classic1)
- self.f = np.empty( (meqn,mx+2*mbc) )
- self.wave = np.empty( (meqn,mwaves,mx+2*mbc) )
- self.s = np.empty( (mwaves,mx+2*mbc) )
- self.amdq = np.empty( (meqn,mx+2*mbc) )
- self.apdq = np.empty( (meqn,mx+2*mbc) )
def teardown(self):
@@ -424,7 +412,7 @@ def homogeneous_step(self,solutions):
if(maux == 0): aux = np.empty( (maux,mx+2*mbc) )
- q,self.cfl = classic1.step1(mx,mbc,mx,q,aux,dx,dt,self.method,self.mthlim,self.f,self.wave,self.s,self.amdq,self.apdq,dtdx)
+ q,self.cfl = classic1.step1(mx,mbc,mx,q,aux,dx,dt,self.method,self.mthlim)
elif(self.kernel_language == 'Python'):
@@ -616,12 +604,6 @@ def set_fortran_parameters(self,solutions):
# These work arrays really ought to live inside a fortran module
# as is done for sharpclaw
- self.qadd = np.empty((meqn,maxm+2*mbc))
- self.fadd = np.empty((meqn,maxm+2*mbc))
- self.gadd = np.empty((meqn,2,maxm+2*mbc))
- self.q1d = np.empty((meqn,maxm+2*mbc))
- self.dtdx1d = np.empty((maxm+2*mbc))
- self.dtdy1d = np.empty((maxm+2*mbc))
self.aux1 = np.empty((maux,maxm+2*mbc))
self.aux2 = np.empty((maux,maxm+2*mbc))
self.aux3 = np.empty((maux,maxm+2*mbc))
@@ -709,13 +691,11 @@ def homogeneous_step(self,solutions):
if self.dim_split:
q, cfl = classic2.dimsp2(maxm,mx,my,mbc,mx,my, \
qold,qnew,aux,dx,dy,dt,self.method,self.mthlim,self.cfl,self.cflv, \
- self.qadd,self.fadd,self.gadd,self.q1d,self.dtdx1d,\
- self.dtdy1d,self.aux1,self.aux2,self.aux3,self.work)
+ self.aux1,self.aux2,self.aux3,self.work)
else:
q, cfl = classic2.step2(maxm,mx,my,mbc,mx,my, \
qold,qnew,aux,dx,dy,dt,self.method,self.mthlim,self.cfl, \
- self.qadd,self.fadd,self.gadd,self.q1d,self.dtdx1d,\
- self.dtdy1d,self.aux1,self.aux2,self.aux3,self.work)
+ self.aux1,self.aux2,self.aux3,self.work)
self.cfl = cfl
state.q=q[:,mbc:mx+mbc,mbc:my+mbc]
Please sign in to comment.
Something went wrong with that request. Please try again.