[cig-commits] r6950 - mc/3D/CitcomCU/branches/inflow-bcs/src
tan2 at geodynamics.org
tan2 at geodynamics.org
Wed May 23 16:10:45 PDT 2007
Author: tan2
Date: 2007-05-23 16:10:45 -0700 (Wed, 23 May 2007)
New Revision: 6950
Modified:
mc/3D/CitcomCU/branches/inflow-bcs/src/Boundary_conditions.c
Log:
Added two functions for inflow/outflow temperature and velocity BC for side (currently x-direction) boundary
Modified: mc/3D/CitcomCU/branches/inflow-bcs/src/Boundary_conditions.c
===================================================================
--- mc/3D/CitcomCU/branches/inflow-bcs/src/Boundary_conditions.c 2007-05-23 22:48:38 UTC (rev 6949)
+++ mc/3D/CitcomCU/branches/inflow-bcs/src/Boundary_conditions.c 2007-05-23 23:10:45 UTC (rev 6950)
@@ -129,6 +129,8 @@
if(E->mesh.periodic_x || E->mesh.periodic_y)
velocity_apply_periodic_bcs(E);
+ else if(E->mesh.sidevbc)
+ velocity_side_bc(E);
else
velocity_refl_vert_bc(E); /* default */
@@ -189,6 +191,8 @@
if(E->mesh.periodic_x || E->mesh.periodic_y)
temperature_apply_periodic_bcs(E);
+ else if(E->mesh.sidevbc)
+ temperature_side_bc(E);
else
temperature_refl_vert_bc(E); /* default */
@@ -358,6 +362,153 @@
return;
}
+
+void velocity_side_bc(struct All_variables *E)
+{
+ int i, j, ii, jj;
+ int node1, node2;
+ int level, nox, noy, noz;
+
+ /* except one side with XOZ and y=0, all others are not reflecting BC */
+ /* for two YOZ planes if 3-D, or two OZ side walls for 2-D */
+
+ if(E->parallel.me_loc[1] == 0 || E->parallel.me_loc[1] == E->parallel.nprocx - 1)
+ for(j = 1; j <= E->lmesh.noy; j++)
+ for(i = 1; i <= E->lmesh.noz; i++)
+ {
+ node1 = i + (j - 1) * E->lmesh.noz * E->lmesh.nox;
+ node2 = node1 + (E->lmesh.nox - 1) * E->lmesh.noz;
+
+ ii = i + E->lmesh.nzs - 1;
+ if(E->parallel.me_loc[1] == 0)
+ {
+ /* Set Vx to a prescribed velocity profile */
+ E->VB[1][node1] = E->segment.Vx[i];
+ if((ii != 1) && (ii != E->mesh.noz))
+ E->VB[3][node1] = 0.0;
+ }
+ if(E->parallel.me_loc[1] == E->parallel.nprocx - 1)
+ {
+ /* Set Vx to a prescribed velocity profile */
+ E->VB[1][node2] = E->segment.Vx[i];
+ if((ii != 1) && (ii != E->mesh.noz))
+ E->VB[3][node2] = 0.0;
+ }
+ } /* end loop for i and j */
+
+ /* for two XOZ planes if 3-D */
+
+ if(E->parallel.me_loc[2] == 0 || E->parallel.me_loc[2] == E->parallel.nprocy - 1)
+ for(j = 1; j <= E->lmesh.nox; j++)
+ for(i = 1; i <= E->lmesh.noz; i++)
+ {
+ node1 = i + (j - 1) * E->lmesh.noz;
+ node2 = node1 + (E->lmesh.noy - 1) * E->lmesh.noz * E->lmesh.nox;
+ ii = i + E->lmesh.nzs - 1;
+
+ if(E->parallel.me_loc[2] == 0)
+ {
+ E->VB[2][node1] = 0.0;
+ if((ii != 1) && (ii != E->mesh.noz))
+ E->VB[3][node1] = 0.0;
+ }
+ if(E->parallel.me_loc[2] == E->parallel.nprocy - 1)
+ {
+ E->VB[2][node2] = 0.0;
+ if((ii != 1) && (ii != E->mesh.noz))
+ E->VB[3][node2] = 0.0;
+ }
+ } /* end of loop i & j */
+
+ /* all vbc's apply at all levels */
+ for(level = E->mesh.levmax; level >= E->mesh.levmin; level--)
+ {
+ nox = E->lmesh.NOX[level];
+ noz = E->lmesh.NOZ[level];
+ noy = E->lmesh.NOY[level];
+
+ if(E->parallel.me_loc[1] == 0 || E->parallel.me_loc[1] == E->parallel.nprocx - 1)
+ for(j = 1; j <= noy; j++)
+ for(i = 1; i <= noz; i++)
+ {
+ node1 = i + (j - 1) * noz * nox;
+ node2 = node1 + (nox - 1) * noz;
+ ii = i + E->lmesh.NZS[level] - 1;
+ if(E->parallel.me_loc[1] == 0)
+ {
+ E->NODE[level][node1] = E->NODE[level][node1] & (~SBX);
+ E->NODE[level][node1] = E->NODE[level][node1] | (VBX);
+ if((ii != 1) && (ii != E->mesh.NOZ[level]))
+ {
+ /* The next two lines differ from those in velocity_refl_vert_bc() */
+ E->NODE[level][node1] = E->NODE[level][node1] & (~SBY);
+ E->NODE[level][node1] = E->NODE[level][node1] | VBY;
+ E->NODE[level][node1] = E->NODE[level][node1] & (~VBZ);
+ E->NODE[level][node1] = E->NODE[level][node1] | SBZ;
+ }
+ }
+ if(E->parallel.me_loc[1] == E->parallel.nprocx - 1)
+ {
+ E->NODE[level][node2] = E->NODE[level][node2] & (~SBX);
+ E->NODE[level][node2] = E->NODE[level][node2] | (VBX);
+ if((ii != 1) && (ii != E->mesh.NOZ[level]))
+ {
+ /* The next two lines differ from those in velocity_refl_vert_bc() */
+ E->NODE[level][node2] = E->NODE[level][node2] & (~SBY);
+ E->NODE[level][node2] = E->NODE[level][node2] | VBY;
+ E->NODE[level][node2] = E->NODE[level][node2] & (~VBZ);
+ E->NODE[level][node2] = E->NODE[level][node2] | SBZ;
+ }
+ }
+ } /* end for loop i & j */
+
+ if(E->parallel.me_loc[2] == 0 || E->parallel.me_loc[2] == E->parallel.nprocy - 1)
+ for(j = 1; j <= nox; j++)
+ for(i = 1; i <= noz; i++)
+ {
+ node1 = i + (j - 1) * noz;
+ node2 = node1 + (noy - 1) * noz * nox;
+ ii = i + E->lmesh.NZS[level] - 1;
+ jj = j + E->lmesh.NXS[level] - 1;
+ if(E->parallel.me_loc[2] == 0)
+ {
+ E->NODE[level][node1] = E->NODE[level][node1] | VBY;
+ E->NODE[level][node1] = E->NODE[level][node1] & (~SBY);
+ if((ii != 1) && (ii != E->mesh.NOZ[level]))
+ {
+ E->NODE[level][node1] = E->NODE[level][node1] & (~VBZ);
+ E->NODE[level][node1] = E->NODE[level][node1] | SBZ;
+ }
+ if((jj != 1) && (jj != E->mesh.NOX[level]) && (ii != 1) && (ii != E->mesh.NOZ[level]))
+ {
+ E->NODE[level][node1] = E->NODE[level][node1] & (~VBX);
+ E->NODE[level][node1] = E->NODE[level][node1] | SBX;
+ }
+ }
+ if(E->parallel.me_loc[2] == E->parallel.nprocy - 1)
+ {
+ E->NODE[level][node2] = E->NODE[level][node2] | VBY;
+ E->NODE[level][node2] = E->NODE[level][node2] & (~SBY);
+ if((ii != 1) && (ii != E->mesh.NOZ[level]))
+ {
+ E->NODE[level][node2] = E->NODE[level][node2] & (~VBZ);
+ E->NODE[level][node2] = E->NODE[level][node2] | SBZ;
+ }
+ if((jj != 1) && (jj != E->mesh.NOX[level]) && (ii != 1) && (ii != E->mesh.NOZ[level]))
+ {
+ E->NODE[level][node2] = E->NODE[level][node2] & (~VBX);
+ E->NODE[level][node2] = E->NODE[level][node2] | SBX;
+ }
+ }
+
+ } /* end for loop i & j */
+ } /* end for loop level */
+
+
+ return;
+}
+
+
void temperature_refl_vert_bc(struct All_variables *E)
{
int i, j;
@@ -411,6 +562,60 @@
}
+void temperature_side_bc(struct All_variables *E)
+{
+ int i, j;
+ int node1, node2;
+ //const int dims = E->mesh.nsd;
+
+ /* Temps and bc-values at top level only */
+ /* fixed temperature at x=0 */
+
+ if(E->parallel.me_loc[1] == 0 || E->parallel.me_loc[1] == E->parallel.nprocx - 1)
+ for(j = 1; j <= E->lmesh.noy; j++)
+ for(i = 1; i <= E->lmesh.noz; i++)
+ {
+ node1 = i + (j - 1) * E->lmesh.noz * E->lmesh.nox;
+ node2 = node1 + (E->lmesh.nox - 1) * E->lmesh.noz;
+ if(E->parallel.me_loc[1] == 0)
+ {
+ /* Set T(x=0) to a prescribed profile */
+ E->node[node1] = E->node[node1] & (~FBX);
+ E->node[node1] = E->node[node1] | TBX;
+ E->TB[1][node1] = E->segment.Tz[i];
+ }
+ if(E->parallel.me_loc[1] == E->parallel.nprocx - 1)
+ {
+ E->node[node2] = E->node[node2] & (~TBX);
+ E->node[node2] = E->node[node2] | FBX;
+ E->TB[1][node2] = 0.0;
+ }
+ } /* end for loop i & j */
+
+ if(E->parallel.me_loc[2] == 0 || E->parallel.me_loc[2] == E->parallel.nprocy - 1)
+ for(j = 1; j <= E->lmesh.nox; j++)
+ for(i = 1; i <= E->lmesh.noz; i++)
+ {
+ node1 = i + (j - 1) * E->lmesh.noz;
+ node2 = node1 + (E->lmesh.noy - 1) * E->lmesh.noz * E->lmesh.nox;
+ if(E->parallel.me_loc[2] == 0)
+ {
+ E->node[node1] = E->node[node1] & (~TBY);
+ E->node[node1] = E->node[node1] | FBY;
+ E->TB[2][node1] = 0.0;
+ }
+ if(E->parallel.me_loc[2] == E->parallel.nprocy - 1)
+ {
+ E->node[node2] = E->node[node2] & (~TBY);
+ E->node[node2] = E->node[node2] | FBY;
+ E->TB[2][node2] = 0.0;
+ }
+ } /* end loop for i and j */
+
+ return;
+}
+
+
void temperature_imposed_botm_bcs(struct All_variables *E, float *BC[], int dirn)
{
//int i, j, node, rowl;
More information about the cig-commits
mailing list