Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP

Loading…

Boundary Conditions Bug Fixes #164

Closed
wants to merge 2 commits into from

2 participants

@mandli
Owner

This pull request includes a bug fix and additional functionality to the current boundary condition routines:

  • The is_valid routine checked to see if a custom boundary condition was set and then checked if a user defined was provided. The routine mistakenly checked to see if the lower boundary condition was present in both cases.

  • Time dependent aux array support has been added. This required the setting and getting of aux arrays to and from the auxbc array in the solver at every time step. A slight refactoring of the boundary condition routines was also included.

There may be debate about whether the aux arrays should be set at every time step as that theoretically this could lead to overhead. One way around this would be to let the user set a parameter indicating a time dependent aux array but this is currently not included in the pull request.

mandli added some commits
@mandli mandli Bug fix in is_valid 1f22944
@mandli mandli Add set and get auxbc functions, add where needed
This commit seperate and adds functions for setting and getting aux
arrays from and to auxbc arrays.  These calls have also been added to
the stepping routines in clawpack as without them time dependent aux
arrays are not updated properly.  Sharpclaw was not modified.
1989bfe
@ketch ketch referenced this pull request
Merged

Bc bugfixes #168

@ketch
Owner

I rebased this and merged it as pull request #168.

@ketch ketch closed this
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Commits on Apr 30, 2012
  1. @mandli

    Bug fix in is_valid

    mandli authored
Commits on May 1, 2012
  1. @mandli

    Add set and get auxbc functions, add where needed

    mandli authored
    This commit seperate and adds functions for setting and getting aux
    arrays from and to auxbc arrays.  These calls have also been added to
    the stepping routines in clawpack as without them time dependent aux
    arrays are not updated properly.  Sharpclaw was not modified.
This page is out of date. Refresh to see the latest.
View
42 src/petclaw/state.py
@@ -235,22 +235,46 @@ def set_q_from_qbc(self,num_ghost,qbc):
raise NotImplementedError("The case of 3D is not handled in "\
+"this function yet")
- def get_qbc_from_q(self,num_ghost,whichvec,qbc):
+ def set_aux_from_auxbc(self,num_ghost,auxbc):
+ """
+ Set the value of aux using the array auxbc. for PetSolver, this
+ involves setting auxbc as the local vector array then perform
+ a local to global communication.
+ """
+
+ patch = self.patch
+ if patch.num_dim == 1:
+ self.aux = auxbc[:,num_ghost:-num_ghost]
+ elif patch.num_dim == 2:
+ self.aux = auxbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost]
+ elif patch.num_dim == 3:
+ self.aux = auxbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost]
+ else:
+ raise NotImplementedError("The case of 3D is not handled in "\
+ +"this function yet")
+
+
+ def get_qbc_from_q(self,num_ghost,qbc):
"""
Returns q with ghost cells attached. For PetSolver,
this means returning the local vector.
"""
shape = [n + 2*num_ghost for n in self.grid.num_cells]
- if whichvec == 'q':
- self.q_da.globalToLocal(self.gqVec, self.lqVec)
- shape.insert(0,self.num_eqn)
- return self.lqVec.getArray().reshape(shape, order = 'F')
+ self.q_da.globalToLocal(self.gqVec, self.lqVec)
+ shape.insert(0,self.num_eqn)
+ return self.lqVec.getArray().reshape(shape, order = 'F')
- elif whichvec == 'aux':
- self.aux_da.globalToLocal(self.gauxVec, self.lauxVec)
- shape.insert(0,self.num_aux)
- return self.lauxVec.getArray().reshape(shape, order = 'F')
+ def set_auxbc_from_aux(self,num_ghost,auxbc):
+ """
+ Returns aux with ghost cells attached. For PetSolver,
+ this means returning the local vector.
+ """
+ shape = [n + 2*num_ghost for n in self.grid.num_cells]
+
+ self.aux_da.globalToLocal(self.gauxVec, self.lauxVec)
+ shape.insert(0,self.num_aux)
+ return self.lauxVec.getArray().reshape(shape, order = 'F')
def set_num_ghost(self,num_ghost):
r"""
View
12 src/pyclaw/clawpack/clawpack.py
@@ -290,6 +290,8 @@ def step_hyperbolic(self,solution):
grid = state.grid
self.apply_q_bcs(state)
+ if state.num_aux > 0:
+ self.apply_aux_bcs(state)
num_eqn,num_ghost = state.num_eqn,self.num_ghost
@@ -383,6 +385,8 @@ def step_hyperbolic(self,solution):
self.cfl.update_global_max(cfl)
state.set_q_from_qbc(num_ghost,self.qbc)
+ if state.num_aux > 0:
+ state.set_aux_from_auxbc(num_ghost,self.auxbc)
# ============================================================================
@@ -507,6 +511,8 @@ def step_hyperbolic(self,solution):
maxm = max(mx,my)
self.apply_q_bcs(state)
+ if state.num_aux > 0:
+ self.apply_aux_bcs(state)
qold = self.qbc.copy('F')
rpn2 = self.rp.rpn2._cpointer
@@ -534,6 +540,8 @@ def step_hyperbolic(self,solution):
self.cfl.update_global_max(cfl)
state.set_q_from_qbc(self.num_ghost,self.qbc)
+ if state.num_aux > 0:
+ state.set_aux_from_auxbc(self.num_ghost,self.auxbc)
else:
raise NotImplementedError("No python implementation for step_hyperbolic in 2D.")
@@ -649,6 +657,8 @@ def step_hyperbolic(self,solution):
maxm = max(mx,my,mz)
self.apply_q_bcs(state)
+ if state.num_aux > 0:
+ self.apply_aux_bcs(state)
qnew = self.qbc
qold = qnew.copy('F')
@@ -682,6 +692,8 @@ def step_hyperbolic(self,solution):
self.cfl.update_global_max(cfl)
state.set_q_from_qbc(self.num_ghost,self.qbc)
+ if state.num_aux > 0:
+ state.set_aux_from_auxbc(self.num_ghost,self.auxbc)
else:
raise NotImplementedError("No python implementation for step_hyperbolic in 3D.")
View
7 src/pyclaw/sharpclaw/sharpclaw.py
@@ -397,6 +397,8 @@ def dq_hyperbolic(self,state):
import numpy as np
self.apply_q_bcs(state)
+ if state.num_aux > 0:
+ self.apply_aux_bcs(state)
q = self.qbc
grid = state.grid
@@ -470,7 +472,8 @@ def dq_hyperbolic(self,state):
dq[m,LL:UL] = -dtdx[LL:UL]*(amdq[m,LL:UL] + apdq[m,LL-1:UL-1] \
+ apdq2[m,LL:UL] + amdq2[m,LL:UL])
- else: raise Exception('Unrecognized value of solver.kernel_language.')
+ else:
+ raise Exception('Unrecognized value of solver.kernel_language.')
self.cfl.update_global_max(cfl)
return dq[:,self.num_ghost:-self.num_ghost]
@@ -536,6 +539,8 @@ def dq_hyperbolic(self,state):
"""
self.apply_q_bcs(state)
+ if state.num_aux > 0:
+ self.apply_aux_bcs(state)
q = self.qbc
grid = state.grid
View
8 src/pyclaw/solver.py
@@ -238,7 +238,7 @@ def is_valid(self):
self.logger.debug('Lower custom BC function has not been set.')
valid = False
if any([bcmeth == BC.custom for bcmeth in self.bc_upper]):
- if self.user_bc_lower is None:
+ if self.user_bc_upper is None:
self.logger.debug('Upper custom BC function has not been set.')
valid = False
return valid
@@ -330,8 +330,8 @@ def apply_q_bcs(self,state):
"""
import numpy as np
-
- self.qbc = state.get_qbc_from_q(self.num_ghost,'q',self.qbc)
+
+ self.qbc = state.get_qbc_from_q(self.num_ghost,self.qbc)
grid = state.grid
for idim,dim in enumerate(grid.dimensions):
@@ -470,7 +470,7 @@ def apply_aux_bcs(self,state):
import numpy as np
- self.auxbc = state.get_qbc_from_q(self.num_ghost,'aux',self.auxbc)
+ self.auxbc = state.get_auxbc_from_aux(self.num_ghost,self.auxbc)
patch = state.patch
View
51 src/pyclaw/state.py 100644 → 100755
@@ -196,26 +196,59 @@ def set_q_from_qbc(self,num_ghost,qbc):
self.q = qbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost]
else:
raise Exception("Assumption (1 <= num_dim <= 3) violated.")
+
+ def set_aux_from_auxbc(self,num_ghost,auxbc):
+ """
+ Set the value of aux using the array auxbc. for PetSolver, this
+ involves setting auxbc as the local vector array then perform
+ a local to global communication.
+ """
+
+ patch = self.patch
+ if patch.num_dim == 1:
+ self.aux = auxbc[:,num_ghost:-num_ghost]
+ elif patch.num_dim == 2:
+ self.aux = auxbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost]
+ elif patch.num_dim == 3:
+ self.aux = auxbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost]
+ else:
+ raise Exception("Assumption (1 <= num_dim <= 3) violated.")
- def get_qbc_from_q(self,num_ghost,whichvec,qbc):
+
+ def get_qbc_from_q(self,num_ghost,qbc):
"""
Fills in the interior of qbc (local vector) by copying q (global vector) to it.
"""
num_dim = self.patch.num_dim
- if whichvec == 'q':
- q = self.q
- elif whichvec == 'aux':
- q = self.aux
-
if num_dim == 1:
- qbc[:,num_ghost:-num_ghost] = q
+ qbc[:,num_ghost:-num_ghost] = self.q
elif num_dim == 2:
- qbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost] = q
+ qbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost] = self.q
elif num_dim == 3:
- qbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost] = q
+ qbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost] = self.q
+ else:
+ raise Exception("Assumption (1 <= num_dim <= 3) violated.")
return qbc
+
+ def get_auxbc_from_aux(self,num_ghost,auxbc):
+ """
+ Fills in the interior of auxbc (local vector) by copying aux (global vector) to it.
+ """
+ num_dim = self.patch.num_dim
+
+ if num_dim == 1:
+ auxbc[:,num_ghost:-num_ghost] = self.aux
+ elif num_dim == 2:
+ auxbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost] = self.aux
+ elif num_dim == 3:
+ auxbc[:,num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost] = self.aux
+ else:
+ raise Exception("Assumption (1 <= num_dim <= 3) violated.")
+
+ return auxbc
+
# ========== Copy functionality ==========================================
def __copy__(self):
Something went wrong with that request. Please try again.