[cig-commits] commit: modified InterpField2D for two different BC's

[Mercurial]

changeset:   129:22fca481af77
user:        Marc Spiegelman <mspieg at ldeo.columbia.edu>
date:        Mon Mar 17 10:57:19 2008 -0400
files:       src/FieldFunctions.c src/SolitaryWave2DSL.c src/SolitaryWave2DSL.h
description:
modified InterpField2D for two different BC's
M   src/FieldFunctions.c
	changed InterpFields2D to InterpFields2DPeriodic and added
	InterpFields2DzInflow
	modified BiCubicInterp for efficiency
M   src/SolitaryWave2DSL.c
	check zBoundaryFlag for characteristic callback
M   src/SolitaryWave2DSL.h
	added InterpFields2DPeriodic and InterpFields2DzInflow


diff -r a9b249909823 -r 22fca481af77 src/FieldFunctions.c
--- a/src/FieldFunctions.c	Sun Mar 16 23:54:19 2008 -0400
+++ b/src/FieldFunctions.c	Mon Mar 17 10:57:19 2008 -0400
@@ -404,8 +404,10 @@ PetscErrorCode InterpVelocity2D(void *f,
 
 /*---------------------------------------------------------------------*/
 #undef __FUNCT__
-#define __FUNCT__ "InterpFields2D"
-PetscErrorCode InterpFields2D(void *f, PetscReal ij_real[], PetscInt numComp, PetscInt components[], PetscReal field[], void *ctx)
+#define __FUNCT__ "InterpFields2DPeriodic"
+PetscErrorCode InterpFields2DPeriodic(void *f, PetscReal ij_real[], PetscInt numComp, PetscInt components[], PetscReal field[], void *ctx)
+
+/*  does semi-lagrangian interpolation assuming doubly periodic in x and z */
 /*---------------------------------------------------------------------*/
 {
   AppCtx   *user = (AppCtx *) ctx;
@@ -416,6 +418,7 @@ PetscErrorCode InterpFields2D(void *f, P
   int ierr;
   
   /* check for out of bound interpolation points */
+
   if ( ir < 0 || ir > ni-1 || jr < 0 || jr> nj-1 ) {  
     ierr =  PetscPrintf(PETSC_COMM_WORLD,"Warning: InterpFields2D: ir=%g, jr=%g\n",ir,jr); CHKERRQ(ierr);  
   } 
@@ -423,6 +426,38 @@ PetscErrorCode InterpFields2D(void *f, P
  return 0;
 }
 
+/*---------------------------------------------------------------------*/
+#undef __FUNCT__
+#define __FUNCT__ "InterpFields2DzInflow"
+PetscErrorCode InterpFields2DzInflow(void *f, PetscReal ij_real[], PetscInt numComp, PetscInt components[], PetscReal field[], void *ctx)
+
+/*   does semi-lagrangian interpolation assuming top boundary is inflow with phi=1, C=0 */
+/*   assumes either periodic or reflection in x (should never leave the x boundary...fragile) */
+
+/*---------------------------------------------------------------------*/
+{
+  AppCtx   *user = (AppCtx *) ctx;
+  Field **x = (Field **) f;
+  int       ni=user->grid->niFine, nj=user->grid->njFine;
+  PetscReal ir=ij_real[0], jr=ij_real[1];
+
+  int ierr;
+  
+  /* check for out of bound interpolation points */
+
+  if ( ir < 0 || ir > ni-1 || jr < 0  ) {  
+    ierr =  PetscPrintf(PETSC_COMM_WORLD,"Warning: InterpFields2D: ir=%g, jr=%g\n",ir,jr); CHKERRQ(ierr);  
+  } 
+
+  /* check for inflow boundary on top */
+  if ( jr >= nj-1 ) {
+    field[0] = 0.;
+    field[1] = 1.;
+  } else {
+    BiCubicInterp((TwoVector **) x, ir, jr, (TwoVector *) field);
+  }
+ return 0;
+}
 /*---------------------------------------------------------------------*/
 #undef __FUNCT__
 #define __FUNCT__ "BiCubicInterp"
@@ -453,8 +488,8 @@ PetscReal CubicInterp(PetscReal x, Petsc
 /*---------------------------------------------------------------------*/
 {
   PetscReal  sxth=0.16666666666667, retval;
-  retval = - y_1*(x-1.0)*(x-2.0)*(x-3.0)*sxth + y_2*(x-0.0)*(x-2.0)*(x-3.0)*0.5
-           - y_3*(x-0.0)*(x-1.0)*(x-3.0)*0.5  + y_4*(x-0.0)*(x-1.0)*(x-2.0)*sxth;
+  retval = - y_1*(x-1.0)*(x-2.0)*(x-3.0)*sxth + y_2*x*(x-2.0)*(x-3.0)*0.5
+           - y_3*      x*(x-1.0)*(x-3.0)*0.5  + y_4*x*(x-1.0)*(x-2.0)*sxth;
   return retval;
 }
 
diff -r a9b249909823 -r 22fca481af77 src/SolitaryWave2DSL.c
--- a/src/SolitaryWave2DSL.c	Sun Mar 16 23:54:19 2008 -0400
+++ b/src/SolitaryWave2DSL.c	Mon Mar 17 10:57:19 2008 -0400
@@ -195,8 +195,12 @@ int main(int argc,char **argv)
   ierr = CharacteristicCreate(comm, &(tsCtx.c));CHKERRQ(ierr);
   c = tsCtx.c; Ncomponents = 2; components[0] = 0; components[1] = 1;
   ierr = CharacteristicSetVelocityInterpolationLocal(c, grid.daC,user[fineLevel].Velocity, user[fineLevel].Velocity, Ncomponents, components, InterpVelocity2D, user);CHKERRQ(ierr);
-  ierr = CharacteristicSetFieldInterpolationLocal(c, grid.daC, user[fineLevel].Xold, Ncomponents, components, InterpFields2D, user);CHKERRQ(ierr);CHKERRQ(ierr);
-  
+  if ( param->zBoundaryFlag == 0 ) { /* use periodic BCs */
+    ierr = CharacteristicSetFieldInterpolationLocal(c, grid.daC, user[fineLevel].Xold, Ncomponents, components, InterpFields2DPeriodic, user);CHKERRQ(ierr);CHKERRQ(ierr);
+  } else if ( param->zBoundaryFlag == 1 ) { /* use inflow vertical BC's */
+    ierr = CharacteristicSetFieldInterpolationLocal(c, grid.daC, user[fineLevel].Xold, Ncomponents, components, InterpFields2DzInflow, user);CHKERRQ(ierr);CHKERRQ(ierr);
+  }
+
   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      Initialize the solution
      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
diff -r a9b249909823 -r 22fca481af77 src/SolitaryWave2DSL.h
--- a/src/SolitaryWave2DSL.h	Sun Mar 16 23:54:19 2008 -0400
+++ b/src/SolitaryWave2DSL.h	Mon Mar 17 10:57:19 2008 -0400
@@ -118,7 +118,8 @@ PetscErrorCode KSPMonitorDumpAll(KSP, Pe
 
 /* interpolation */
 PetscErrorCode InterpVelocity2D(void *, PetscReal [], PetscInt, PetscInt [], PetscScalar [], void *);
-PetscErrorCode InterpFields2D(void *, PetscReal [], PetscInt, PetscInt [], PetscReal [], void *);
+PetscErrorCode InterpFields2DPeriodic(void *, PetscReal [], PetscInt, PetscInt [], PetscReal [], void *);
+PetscErrorCode InterpFields2DzInflow(void *, PetscReal [], PetscInt, PetscInt [], PetscReal [], void *);
 PetscErrorCode BiCubicInterp(TwoVector **, PetscScalar, PetscScalar, TwoVector*);
 PetscReal CubicInterp(PetscScalar, PetscScalar, PetscScalar, PetscScalar, PetscScalar);