[cig-commits] r15108 - in mc/3D/CitcomS/trunk: CitcomS/Coupler CitcomS/Solver module/Exchanger
tan2 at geodynamics.org
tan2 at geodynamics.org
Tue Jun 2 15:56:48 PDT 2009
Author: tan2
Date: 2009-06-02 15:56:46 -0700 (Tue, 02 Jun 2009)
New Revision: 15108
Added:
mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.cc
mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.h
mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.cc
mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.h
mc/3D/CitcomS/trunk/module/Exchanger/PInlet.cc
mc/3D/CitcomS/trunk/module/Exchanger/PInlet.h
mc/3D/CitcomS/trunk/module/Exchanger/PInterior.cc
mc/3D/CitcomS/trunk/module/Exchanger/PInterior.h
mc/3D/CitcomS/trunk/module/Exchanger/POutlet.cc
mc/3D/CitcomS/trunk/module/Exchanger/POutlet.h
Modified:
mc/3D/CitcomS/trunk/CitcomS/Coupler/ContainingCoupler.py
mc/3D/CitcomS/trunk/CitcomS/Coupler/Coupler.py
mc/3D/CitcomS/trunk/CitcomS/Coupler/EmbeddedCoupler.py
mc/3D/CitcomS/trunk/CitcomS/Coupler/Inlet.py
mc/3D/CitcomS/trunk/CitcomS/Coupler/Outlet.py
mc/3D/CitcomS/trunk/CitcomS/Solver/CoupledSolver.py
mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.cc
mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.h
mc/3D/CitcomS/trunk/module/Exchanger/BoundaryVTInlet.cc
mc/3D/CitcomS/trunk/module/Exchanger/CitcomInterpolator.cc
mc/3D/CitcomS/trunk/module/Exchanger/Makefile.am
mc/3D/CitcomS/trunk/module/Exchanger/bindings.cc
mc/3D/CitcomS/trunk/module/Exchanger/exchangers.cc
mc/3D/CitcomS/trunk/module/Exchanger/exchangers.h
mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.cc
mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.h
Log:
Added parameters 'amending_outflow' and 'exchange_pressure' to help the convergence of esolver
When 'amending_outflow' is set to true, the imposed velocity BC will be
amended slightly to be divergence-free. (The divergence, e.g. outflow, is
caused by the combination of solver inaccuracy and interpolation inaccuracy.)
When 'exchange_pressure' is set to true, the initial pressure (at element
level) of the embedded solver is taken from the pressure solution of the
containing solver.
Modified: mc/3D/CitcomS/trunk/CitcomS/Coupler/ContainingCoupler.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Coupler/ContainingCoupler.py 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/CitcomS/Coupler/ContainingCoupler.py 2009-06-02 22:56:46 UTC (rev 15108)
@@ -38,12 +38,6 @@
def initialize(self, solver):
Coupler.initialize(self, solver)
-
- # allocate space for exchanger objects
- self.remoteBdryList = range(self.remoteSize)
- self.sourceList = range(self.remoteSize)
- self.outletList = range(self.remoteSize)
-
return
@@ -65,8 +59,10 @@
self.interior, self.myBBox = createInterior(self.remoteBBox,
self.all_variables)
- # an empty boundary object, which will be filled by a remote boundary obj.
+ self.remoteBdryList = range(self.remoteSize)
for i in range(self.remoteSize):
+ # an empty boundary object for remote velocity nodes
+ # will be filled by a remote boundary obj.
self.remoteBdryList[i] = createEmptyBoundary()
return
@@ -84,6 +80,7 @@
def createSource(self):
# the source obj's will send boundary conditions to a remote sink
from ExchangerLib import CitcomSource_create
+ self.sourceList = range(self.remoteSize)
for i, comm, b in zip(range(self.remoteSize),
self.srcCommList,
self.remoteBdryList):
@@ -98,6 +95,24 @@
self.myBBox,
self.all_variables)
+ if self.inventory.exchange_pressure:
+ from ExchangerLib import createEmptyPInterior
+ self.pinterior = range(self.remoteSize)
+
+ for i in range(self.remoteSize):
+ self.pinterior[i] = createEmptyPInterior()
+
+ self.psourceList = range(self.remoteSize)
+ for i, comm, b in zip(range(self.remoteSize),
+ self.srcCommList,
+ self.pinterior):
+ # sink is always in the last rank of a communicator
+ sinkRank = comm.size - 1
+ self.psourceList[i] = CitcomSource_create(comm.handle(),
+ sinkRank,
+ b,
+ self.myBBox,
+ self.all_variables)
return
@@ -117,9 +132,18 @@
# boundary conditions will be sent by SVTOutlet, which sends
# stress, velocity, and temperature
import Outlet
+ self.outletList = range(self.remoteSize)
for i, src in zip(range(self.remoteSize),
self.sourceList):
self.outletList[i] = Outlet.SVTOutlet(src, self.all_variables)
+
+ if self.inventory.exchange_pressure:
+ self.poutletList = range(self.remoteSize)
+ for i, src in zip(range(self.remoteSize),
+ self.psourceList):
+
+ self.poutletList[i] = Outlet.POutlet(src, self.all_variables)
+
return
@@ -180,6 +204,10 @@
# send computed velocity to ECPLR for its BCs
for outlet in self.outletList:
outlet.send()
+
+ if self.inventory.exchange_pressure:
+ for outlet in self.poutletList:
+ outlet.send()
return
Modified: mc/3D/CitcomS/trunk/CitcomS/Coupler/Coupler.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Coupler/Coupler.py 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/CitcomS/Coupler/Coupler.py 2009-06-02 22:56:46 UTC (rev 15108)
@@ -151,11 +151,17 @@
# updating the temperature field in the containing solver or not
two_way_communication = prop.bool("two_way_communication", default=True)
- # insuring consistent inititial temperature fields at the overlapping
+
+ # ensuring consistent inititial temperature fields at the overlapping
# domain or not
exchange_initial_temperature = prop.bool("exchange_initial_temperature",
default=True)
+ # pressure from the csolver can be used as a good initial guess
+ # for pressure of the esolver
+ exchange_pressure = prop.bool("exchange_pressure",
+ default=True)
+
# if si_unit is True, quantities exchanged are in SI units
si_unit = prop.bool("si_unit", default=False)
Modified: mc/3D/CitcomS/trunk/CitcomS/Coupler/EmbeddedCoupler.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Coupler/EmbeddedCoupler.py 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/CitcomS/Coupler/EmbeddedCoupler.py 2009-06-02 22:56:46 UTC (rev 15108)
@@ -52,12 +52,6 @@
# otherwise, we have to stop
assert solver.inventory.bc.inventory.side_sbcs == True, \
'Error: esolver.bc.side_sbcs must be on!'
-
- # allocate space for exchanger objects
- self.remoteIntrList = range(self.remoteSize)
- self.sourceList = range(self.remoteSize)
- self.outletList = range(self.remoteSize)
-
return
@@ -84,6 +78,7 @@
# an empty interior object, which will be filled by a remote interior obj.
if inv.two_way_communication:
+ self.remoteIntrList = range(self.remoteSize)
for i in range(self.remoteSize):
self.remoteIntrList[i] = createEmptyInterior()
@@ -107,12 +102,21 @@
self.sink = Sink_create(self.sinkComm.handle(),
self.remoteSize,
self.boundary)
+
+ if self.inventory.exchange_pressure:
+ from ExchangerLib import createPInterior
+ self.pinterior, _ = createPInterior(self.myBBox,
+ self.all_variables)
+ self.psink = Sink_create(self.sinkComm.handle(),
+ self.remoteSize,
+ self.pinterior)
return
def createSource(self):
# the source obj's will send interior temperature to a remote sink
from ExchangerLib import CitcomSource_create
+ self.sourceList = range(self.remoteSize)
for i, comm, b in zip(range(self.remoteSize),
self.srcCommList,
self.remoteIntrList):
@@ -134,15 +138,27 @@
# boundary conditions will be recv. by SVTInlet, which receives
# stress, velocity, and temperature
import Inlet
- self.inlet = Inlet.SVTInlet(self.boundary,
- self.sink,
- self.all_variables)
+ if self.inventory.amending_outflow:
+ self.inlet = Inlet.BoundarySVTInlet(self.boundary,
+ self.sink,
+ self.all_variables,
+ self.communicator.handle())
+ else:
+ self.inlet = Inlet.SVTInlet(self.boundary,
+ self.sink,
+ self.all_variables)
+
+ if self.inventory.exchange_pressure:
+ self.pinlet = Inlet.PInlet(self.pinterior,
+ self.psink,
+ self.all_variables)
return
def createII(self):
# interior temperature will be sent by TOutlet
import Outlet
+ self.outletList = range(self.remoteSize)
for i, src in zip(range(self.remoteSize),
self.sourceList):
self.outletList[i] = Outlet.TOutlet(src, self.all_variables)
@@ -182,11 +198,18 @@
def preVSolverRun(self):
# apply bc before solving the velocity
+
+ # recv'ing vbc after sync'd steps
if self.toApplyBC:
self.inlet.recv()
self.toApplyBC = False
+ if self.inventory.exchange_pressure:
+ self.pinlet.recv()
+ self.pinlet.impose()
+
+ # update vbc every time step
self.inlet.impose()
# applyBC only when previous step is sync'd
@@ -261,6 +284,8 @@
# excluding nodes in bottom boundary?
exclude_bottom = prop.bool("exclude_bottom", default=False)
+ # amending the received vbc to be divergence-free
+ amending_outflow = prop.bool("amending_outflow", default=False)
Modified: mc/3D/CitcomS/trunk/CitcomS/Coupler/Inlet.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Coupler/Inlet.py 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/CitcomS/Coupler/Inlet.py 2009-06-02 22:56:46 UTC (rev 15108)
@@ -66,7 +66,18 @@
+class BoundarySVTInlet(Inlet):
+ def __init__(self, mesh, sink, all_variables, comm):
+ import ExchangerLib
+ self._handle = ExchangerLib.BoundarySVTInlet_create(mesh,
+ sink,
+ all_variables,
+ comm)
+ return
+
+
+
class TInlet(Inlet):
def __init__(self, mesh, sink, all_variables):
@@ -77,6 +88,16 @@
return
+class PInlet(Inlet):
+
+ def __init__(self, mesh, sink, all_variables):
+ import ExchangerLib
+ self._handle = ExchangerLib.PInlet_create(mesh,
+ sink,
+ all_variables)
+ return
+
+
class SInlet(Inlet):
def __init__(self, mesh, sink, all_variables):
Modified: mc/3D/CitcomS/trunk/CitcomS/Coupler/Outlet.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Coupler/Outlet.py 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/CitcomS/Coupler/Outlet.py 2009-06-02 22:56:46 UTC (rev 15108)
@@ -66,6 +66,17 @@
+class POutlet(Outlet):
+
+ def __init__(self, source, all_variables):
+ import ExchangerLib
+ self._handle = ExchangerLib.POutlet_create(source,
+ all_variables)
+ return
+
+
+
+
class VTOutlet(Outlet):
def __init__(self, source, all_variables):
Modified: mc/3D/CitcomS/trunk/CitcomS/Solver/CoupledSolver.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Solver/CoupledSolver.py 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/CitcomS/Solver/CoupledSolver.py 2009-06-02 22:56:46 UTC (rev 15108)
@@ -68,7 +68,6 @@
# insure consistent temperature fields across solvers
self.coupler.exchangeTemperature()
-
self.solveVelocities()
return
@@ -76,9 +75,8 @@
def solveVelocities(self):
# sync boundary conditions before/after vsolver
- vsolver = self.inventory.vsolver
self.coupler.preVSolverRun()
- vsolver.run()
+ self.inventory.vsolver.run()
self.coupler.postVSolverRun()
return
Modified: mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.cc 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -38,6 +38,7 @@
#include "AreaWeightedNormal.h"
extern "C" {
+ // for definition of SIDE_NORTH etc.
#include "element_definitions.h"
}
@@ -51,13 +52,20 @@
const Boundary& boundary,
const Sink& sink,
const All_variables* E) :
- size_(boundary.size()),
- toleranceOutflow_(E->control.accuracy),
- nwght(size_, 0)
+ nwght(boundary.size(), 0)
{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
computeWeightedNormal(boundary, E);
- double total_area = computeTotalArea(comm, sink);
- normalize(total_area);
+ double inv_total_area = 1 / computeTotalArea(comm, sink);
+
+ // normalize
+ for(int n=0; n<sink.size(); n++) {
+ int i = sink.meshNode(n);
+ for(int j=0; j<DIM; j++)
+ nwght[j][i] *= inv_total_area;
+ }
}
@@ -67,22 +75,23 @@
void AreaWeightedNormal::imposeConstraint(Velo& V,
const MPI_Comm& comm,
- const Sink& sink) const
+ const Sink& sink,
+ const All_variables* E) const
{
- journal::info_t info("CitcomS-AreaWeightedNormal-outflow");
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
double outflow = computeOutflow(V, comm, sink);
- info << journal::at(__HERE__)
- << "Net outflow: "
- << outflow << journal::endl;
- if (std::abs(outflow) > toleranceOutflow_) {
- reduceOutflow(V, outflow, sink);
+ reduceOutflow(V, outflow, sink);
- outflow = computeOutflow(V, comm, sink);
- info << journal::at(__HERE__)
- << "Net outflow after correction (SHOULD BE ZERO !): "
- << outflow << journal::endl;
+ double new_outflow = computeOutflow(V, comm, sink);
+
+ if(E->parallel.me == 0) {
+ fprintf(stderr, "Net outflow amended from %e to %e\n",
+ outflow, new_outflow);
+ fprintf(E->fp, "Net outflow amended from %e to %e\n",
+ outflow, new_outflow);
}
}
@@ -94,24 +103,32 @@
{
const int nodes_per_element = 8;
const int vpoints_per_side = 4;
+
+ // converting normal vector [-1, 0, 1] to side index,
const int side_normal[DIM][DIM] = {{SIDE_NORTH, -1, SIDE_SOUTH},
{SIDE_WEST, -1, SIDE_EAST},
{SIDE_BOTTOM, -1, SIDE_TOP}};
+ /* For each node belong to boundary elements, check all 6 faces,
+ * if the node is on the face and the face is part of the boundary,
+ * compute the surface area by summing the 2D determinants. */
for (int m=1; m<=E->sphere.caps_per_proc; m++)
for (int es=1; es<=E->boundary.nel; es++) {
int el = E->boundary.element[m][es];
for(int n=1; n<=nodes_per_element; n++) {
int node = E->ien[m][el].node[n];
int bnode = boundary.bnode(node);
+
+ // skipping nodes not on the boundary
if(bnode < 0) continue;
for(int j=0; j<DIM; j++) {
- int normal = boundary.normal(j,bnode);
+ // normal is either -1, 1, or 0
+ int normal = boundary.normal(j, bnode);
if(normal) {
int side = side_normal[j][normal+1];
for(int k=1; k<=vpoints_per_side; k++)
- nwght[j][bnode] += normal * E->boundary.det[m][side][k][es];
+ nwght[j][bnode] += normal * E->boundary.det[m][side][k][es] / vpoints_per_side;
}
} // end of loop over dim
} // end of loop over element nodes
@@ -131,15 +148,32 @@
Exchanger::util::gatherSum(comm, total_area);
- return total_area;
-}
+ journal::info_t info("CitcomS-AreaWeightedNormal-area");
+ info << journal::at(__HERE__)
+ << "Total surface area: "
+ << total_area << journal::endl;
+ /** debug **
+ for(int j=0; j<DIM; j++) {
+ double in, out;
+ in = out = 0;
+ for(int n=0; n<sink.size(); n++) {
+ int i = sink.meshNode(n);
+ if(nwght[j][i] > 0) out += nwght[j][i];
+ else in += nwght[j][i];
+ }
-void AreaWeightedNormal::normalize(double total_area)
-{
- for(int i=0; i<size_; i++)
- for(int j=0; j<DIM; j++)
- nwght[j][i] /= total_area;
+ Exchanger::util::gatherSum(comm, in);
+ Exchanger::util::gatherSum(comm, out);
+
+ journal::info_t info("CitcomS-AreaWeightedNormal-area");
+ info << journal::at(__HERE__)
+ << "Partial surface area: " << j << " "
+ << in << " " << out << journal::endl;
+ }
+ */
+
+ return total_area;
}
@@ -147,6 +181,7 @@
const MPI_Comm& comm,
const Sink& sink) const
{
+ /* integrate dot(V,n) over the surface */
double outflow = 0;
for(int n=0; n<sink.size(); n++) {
int i = sink.meshNode(n);
Modified: mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.h 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/AreaWeightedNormal.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -40,28 +40,27 @@
class AreaWeightedNormal {
- const int size_;
- const double toleranceOutflow_;
Exchanger::Array2D<double,Exchanger::DIM> nwght;
public:
AreaWeightedNormal(const MPI_Comm& comm,
const Boundary& boundary,
const Exchanger::Sink& sink,
- const All_variables* E);
+ const All_variables* E);
~AreaWeightedNormal();
typedef Exchanger::Array2D<double,Exchanger::DIM> Velo;
void imposeConstraint(Velo& V,
const MPI_Comm& comm,
- const Exchanger::Sink& sink) const;
+ const Exchanger::Sink& sink,
+ const All_variables* E) const;
private:
void computeWeightedNormal(const Boundary& boundary,
const All_variables* E);
double computeTotalArea(const MPI_Comm& comm,
- const Exchanger::Sink& sink) const;
+ const Exchanger::Sink& sink) const;
void normalize(double total_area);
double computeOutflow(const Velo& V,
const MPI_Comm& comm,
Added: mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.cc (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,219 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#include "config.h"
+#include "journal/diagnostics.h"
+#include "global_defs.h"
+#include "Boundary.h"
+#include "Convertor.h"
+#include "Exchanger/Sink.h"
+#include "BaseSVTInlet.h"
+
+extern "C" {
+
+#include "element_definitions.h"
+
+ void check_bc_consistency(const All_variables *E);
+ void construct_id(const All_variables *E);
+}
+
+using Exchanger::Array2D;
+using Exchanger::DIM;
+using Exchanger::STRESS_DIM;
+using Exchanger::Sink;
+
+
+const unsigned vbcFlag[] = {VBX, VBY, VBZ};
+const unsigned sbcFlag[] = {SBX, SBY, SBZ};
+
+
+BaseSVTInlet::BaseSVTInlet(const Boundary& boundary,
+ const Sink& sink,
+ All_variables* e) :
+ Inlet(boundary, sink),
+ E(e),
+ s(sink.size()),
+ s_old(sink.size()),
+ v(sink.size()),
+ v_old(sink.size()),
+ t(sink.size()),
+ t_old(sink.size())
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ // set CitcomS boundary flags
+ setVBCFlag();
+ setTBCFlag();
+
+ check_bc_consistency(E);
+}
+
+
+BaseSVTInlet::~BaseSVTInlet()
+{}
+
+
+void BaseSVTInlet::impose()
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ // impose normal velocity and shear stress as BC for momentum eqn.
+ imposeSV();
+
+ // impose temperature BC for energy eqn.
+ imposeT();
+}
+
+
+// private functions
+
+void BaseSVTInlet::setVBCFlag()
+{
+ // BC: normal velocity and shear traction
+
+ const Boundary& boundary = dynamic_cast<const Boundary&>(mesh);
+ const int m = 1;
+
+ for(int i=0; i<boundary.size(); i++) {
+ int n = boundary.nodeID(i);
+ for(int d=0; d<DIM; ++d)
+ if(boundary.normal(d,i)) {
+ E->node[m][n] = E->node[m][n] | vbcFlag[d];
+ E->node[m][n] = E->node[m][n] & (~sbcFlag[d]);
+ } else {
+ E->node[m][n] = E->node[m][n] | sbcFlag[d];
+ E->node[m][n] = E->node[m][n] & (~vbcFlag[d]);
+ }
+ }
+
+ // reconstruct ID array to reflect changes in VBC
+ construct_id(E);
+}
+
+
+void BaseSVTInlet::setTBCFlag()
+{
+ const int m = 1;
+ for(int i=0; i<mesh.size(); i++) {
+ int n = mesh.nodeID(i);
+ E->node[m][n] = E->node[m][n] | TBX;
+ E->node[m][n] = E->node[m][n] | TBY;
+ E->node[m][n] = E->node[m][n] | TBZ;
+ E->node[m][n] = E->node[m][n] & (~FBX);
+ E->node[m][n] = E->node[m][n] & (~FBY);
+ E->node[m][n] = E->node[m][n] & (~FBZ);
+ }
+}
+
+
+void BaseSVTInlet::imposeSV()
+{
+ const int sidelow[3] = {SIDE_NORTH, SIDE_WEST, SIDE_BOTTOM};
+ const int sidehigh[3] = {SIDE_SOUTH, SIDE_EAST, SIDE_TOP};
+
+ const Boundary& boundary = dynamic_cast<const Boundary&>(mesh);
+
+ double N1, N2;
+ getTimeFactors(N1, N2);
+
+ const int m = 1;
+ for(int i=0; i<sink.size(); i++) {
+ int j = sink.meshNode(i);
+ int n = boundary.nodeID(j);
+ int q = E->sbc.node[m][n];
+
+ for(int d=0; d<DIM; d++)
+ if(E->node[m][n] & vbcFlag[d])
+ E->sphere.cap[m].VB[d+1][n] = N1 * v_old[d][i] + N2 * v[d][i];
+
+ if(E->node[m][n] & (SBX | SBY | SBZ))
+ for(int d=0; d<DIM; d++) {
+ int p;
+ if(boundary.normal(d,j) == -1)
+ p = sidelow[d];
+ else if(boundary.normal(d,j) == 1)
+ p = sidehigh[d];
+ else
+ continue;
+
+ for(int k=0; k<DIM; k++) {
+ E->sbc.SB[m][p][k+1][q] =
+ boundary.normal(d,j) *
+ ( N1 * side_tractions(s_old, i, d, k) +
+ N2 * side_tractions(s, i, d, k) );
+ }
+ }
+ }
+}
+
+
+void BaseSVTInlet::imposeT()
+{
+ journal::debug_t debugBC("CitcomS-imposeT");
+ debugBC << journal::at(__HERE__);
+
+ double N1, N2;
+ getTimeFactors(N1, N2);
+
+ const int m = 1;
+ for(int i=0; i<sink.size(); i++) {
+ int j = sink.meshNode(i);
+ int n = mesh.nodeID(j);
+
+ for(int d=0; d<DIM; d++)
+ E->sphere.cap[m].TB[d+1][n] = N1 * t_old[0][i] + N2 * t[0][i];
+
+ debugBC << E->sphere.cap[m].TB[1][n] << " "
+ << E->sphere.cap[m].TB[2][n] << " "
+ << E->sphere.cap[m].TB[3][n] << journal::newline;
+
+ }
+ debugBC << journal::endl;
+
+ (E->temperatures_conform_bcs)(E);
+}
+
+
+double BaseSVTInlet::side_tractions(const Array2D<double,STRESS_DIM>& stress,
+ int node, int normal_dir, int dim) const
+{
+ const int stress_index[3][3] =
+ { {0, 3, 4},
+ {3, 1, 5},
+ {4, 5, 2} };
+
+ return stress[ stress_index[normal_dir][dim] ][node];
+}
+
+
+// version
+// $Id: BaseSVTInlet.cc 7643 2007-07-11 20:17:32Z tan2 $
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.h (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/BaseSVTInlet.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,82 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+// Role:
+// impose normal velocity and shear traction as velocity b.c.,
+// also impose temperature as temperature b.c.
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#if !defined(pyCitcomSExchanger_BaseSVTInlet_h)
+#define pyCitcomSExchanger_BaseSVTInlet_h
+
+#include "Exchanger/Array2D.h"
+#include "Exchanger/DIM.h"
+#include "Exchanger/Inlet.h"
+
+struct All_variables;
+class Boundary;
+
+
+class BaseSVTInlet : public Exchanger::Inlet{
+protected:
+ All_variables* E;
+ Exchanger::Array2D<double,Exchanger::STRESS_DIM> s;
+ Exchanger::Array2D<double,Exchanger::STRESS_DIM> s_old;
+ Exchanger::Array2D<double,Exchanger::DIM> v;
+ Exchanger::Array2D<double,Exchanger::DIM> v_old;
+ Exchanger::Array2D<double,1> t;
+ Exchanger::Array2D<double,1> t_old;
+
+public:
+ BaseSVTInlet(const Boundary& boundary,
+ const Exchanger::Sink& sink,
+ All_variables* E);
+
+ virtual ~BaseSVTInlet();
+
+ virtual void impose();
+
+ //virtual void recv() should be implemented by child class
+
+protected:
+ void setVBCFlag();
+ void setTBCFlag();
+
+ void imposeSV();
+ void imposeT();
+
+ double side_tractions(const Exchanger::Array2D<double,Exchanger::STRESS_DIM>& stress,
+ int node, int normal_dir, int dim) const;
+};
+
+
+#endif
+
+// version
+// $Id: BaseSVTInlet.h 2397 2005-10-04 22:37:25Z leif $
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.cc (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,93 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#include "config.h"
+#include "journal/diagnostics.h"
+#include "AreaWeightedNormal.h"
+#include "Convertor.h"
+#include "Boundary.h"
+#include "BoundarySVTInlet.h"
+
+using Exchanger::Array2D;
+using Exchanger::DIM;
+using Exchanger::STRESS_DIM;
+using Exchanger::Sink;
+
+
+BoundarySVTInlet::BoundarySVTInlet(const Boundary& boundary,
+ const Sink& sink,
+ All_variables* E,
+ const MPI_Comm& c) :
+ BaseSVTInlet(boundary, sink, E),
+ comm(c),
+ awnormal(new AreaWeightedNormal(comm, boundary, sink, E))
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+}
+
+
+BoundarySVTInlet::~BoundarySVTInlet()
+{
+ delete awnormal;
+}
+
+
+void BoundarySVTInlet::recv()
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ // store bc from previous timestep
+ s.swap(s_old);
+ t.swap(t_old);
+ v.swap(v_old);
+
+ // receive the fields from outlet
+ sink.recv(t, v);
+ sink.recv(s);
+
+ // convert back to CitcomS' units and coordinate system
+ Exchanger::Convertor& convertor = Convertor::instance();
+ convertor.xtemperature(t);
+ convertor.xvelocity(v, sink.getX());
+ //convertor.xstress(s, sink.getX());
+
+ // correct v to be div-free
+ awnormal->imposeConstraint(v, comm, sink, E);
+
+ t.print("CitcomS-BoundarySVTInlet-T");
+ v.print("CitcomS-BoundarySVTInlet-V");
+ s.print("CitcomS-BoundarySVTInlet-S");
+}
+
+
+// version
+// $Id: BoundaryVTInlet.cc 7403 2007-06-23 00:33:20Z tan2 $
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.h (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/BoundarySVTInlet.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,61 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#if !defined(pyCitcomSExchanger_BoundaryVTInlet_h)
+#define pyCitcomSExchanger_BoundaryVTInlet_h
+
+#include "mpi.h"
+#include "BaseSVTInlet.h"
+
+struct All_variables;
+class AreaWeightedNormal;
+
+
+class BoundarySVTInlet : public BaseSVTInlet{
+private:
+ const MPI_Comm& comm;
+ AreaWeightedNormal* awnormal;
+
+public:
+ BoundarySVTInlet(const Boundary& boundary,
+ const Exchanger::Sink& sink,
+ All_variables* E,
+ const MPI_Comm& comm);
+ virtual ~BoundarySVTInlet();
+
+ virtual void recv();
+
+};
+
+
+#endif
+
+// version
+// $Id: BoundaryVTInlet.h 2397 2005-10-04 22:37:25Z leif $
+
+// End of file
Modified: mc/3D/CitcomS/trunk/module/Exchanger/BoundaryVTInlet.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/BoundaryVTInlet.cc 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/BoundaryVTInlet.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -27,7 +27,7 @@
//
#include "config.h"
-#include "journal/journal.h"
+#include "journal/diagnostics.h"
#include "AreaWeightedNormal.h"
#include "Boundary.h"
#include "BoundaryVTInlet.h"
@@ -61,7 +61,7 @@
VTInlet::recv();
- awnormal->imposeConstraint(v, comm, sink);
+ awnormal->imposeConstraint(v, comm, sink, comm);
v.print("CitcomS-BoundaryVTInlet-V_constrained");
}
Modified: mc/3D/CitcomS/trunk/module/Exchanger/CitcomInterpolator.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/CitcomInterpolator.cc 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/CitcomInterpolator.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -80,7 +80,7 @@
int n1 = elem_[0][i];
for(int k=0; k<NODES_PER_ELEMENT; k++) {
int node = E->ien[mm][n1].node[k+1];
- P[0][i] += shape_[k][i] * E->P[mm][node];
+ P[0][i] += shape_[k][i] * E->NP[mm][node];
}
}
}
Modified: mc/3D/CitcomS/trunk/module/Exchanger/Makefile.am
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/Makefile.am 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/Makefile.am 2009-06-02 22:56:46 UTC (rev 15108)
@@ -83,6 +83,8 @@
Interior.h \
misc.cc \
misc.h \
+ PInterior.cc \
+ PInterior.h \
SInlet.cc \
SInlet.h \
SIUnit.cc \
Added: mc/3D/CitcomS/trunk/module/Exchanger/PInlet.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/PInlet.cc (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/PInlet.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,102 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+// Purpose:
+// Replace local temperture field by received values. Note that b.c. is not
+// affected.
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#include "config.h"
+#include "journal/diagnostics.h"
+#include "global_defs.h"
+#include "Convertor.h"
+#include "Exchanger/BoundedMesh.h"
+#include "Exchanger/Sink.h"
+#include "PInlet.h"
+
+using Exchanger::Array2D;
+using Exchanger::BoundedMesh;
+using Exchanger::Sink;
+
+extern "C" {
+ double global_p_norm2(struct All_variables*, double **);
+}
+
+
+PInlet::PInlet(const BoundedMesh& boundedMesh,
+ const Sink& sink,
+ All_variables* e) :
+ Inlet(boundedMesh, sink),
+ E(e),
+ p(sink.size())
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+}
+
+
+PInlet::~PInlet()
+{}
+
+
+void PInlet::recv()
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ sink.recv(p);
+
+ // TODO: convert pressure to non-dimensional value
+ //Exchanger::Convertor& convertor = Convertor::instance();
+ //convertor.xpressure(t);
+
+ p.print("CitcomS-PInlet-P");
+}
+
+
+void PInlet::impose()
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ const int m = 1;
+ for(int i=0; i<sink.size(); i++) {
+ int n = mesh.nodeID(sink.meshNode(i));
+ E->P[m][n] = p[0][i];
+ }
+
+ E->monitor.pdotp = global_p_norm2(E, E->P);
+}
+
+
+// private functions
+
+
+// version
+// $Id$
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/PInlet.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/PInlet.h (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/PInlet.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,66 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+// Role:
+// impose velocity and temperature as b.c.
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#if !defined(pyCitcomSExchanger_PInlet_h)
+#define pyCitcomSExchanger_PInlet_h
+
+#include "Exchanger/Array2D.h"
+#include "Exchanger/DIM.h"
+#include "Exchanger/Inlet.h"
+
+struct All_variables;
+
+
+class PInlet : public Exchanger::Inlet {
+ All_variables* E;
+ Exchanger::Array2D<double,1> p;
+
+public:
+ PInlet(const Exchanger::BoundedMesh& boundedMesh,
+ const Exchanger::Sink& sink,
+ All_variables* E);
+
+ virtual ~PInlet();
+
+ virtual void recv();
+ virtual void impose();
+
+private:
+
+};
+
+
+#endif
+
+// version
+// $Id$
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/PInterior.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/PInterior.cc (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/PInterior.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,86 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#include "config.h"
+#include <vector>
+#include "global_defs.h"
+#include "PInterior.h"
+#include "journal/diagnostics.h"
+
+
+PInterior::PInterior() :
+ Exchanger::BoundedMesh()
+{}
+
+
+
+PInterior::PInterior(const Exchanger::BoundedBox& bbox,
+ const All_variables* E) :
+ Exchanger::BoundedMesh()
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ bbox_ = bbox;
+ bbox_.print("CitcomS-PInterior-BBox");
+
+ X_.reserve(E->lmesh.nel);
+ nodeID_.reserve(E->lmesh.nel);
+
+ initX(E);
+
+ X_.shrink();
+ X_.print("CitcomS-PInterior-X");
+
+ nodeID_.shrink();
+ nodeID_.print("CitcomS-PInterior-nodeID");
+}
+
+
+void PInterior::initX(const All_variables* E)
+{
+ std::vector<double> x(Exchanger::DIM);
+
+ // Storing the coordinates of the center of elements
+ for (int m=1;m<=E->sphere.caps_per_proc;m++)
+ for(int e=1;e<=E->lmesh.nel;e++) {
+ for(int d=0; d<Exchanger::DIM; d++)
+ x[d] = E->eco[m][e].centre[d+1];
+
+ if(isInside(x, bbox_)) {
+ X_.push_back(x);
+ nodeID_.push_back(e);
+ }
+ }
+}
+
+
+// version
+// $Id$
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/PInterior.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/PInterior.h (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/PInterior.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,55 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#if !defined(pyCitcomSExchanger_PInterior_h)
+#define pyCitcomSExchanger_PInterior_h
+
+#include "Exchanger/BoundedMesh.h"
+
+struct All_variables;
+
+
+class PInterior : public Exchanger::BoundedMesh {
+
+public:
+ PInterior();
+ PInterior(const Exchanger::BoundedBox& bbox,
+ const All_variables* E);
+
+private:
+ void initX(const All_variables* E);
+
+};
+
+
+#endif
+
+// version
+// $Id$
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/POutlet.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/POutlet.cc (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/POutlet.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,75 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#include "config.h"
+#include "global_defs.h"
+#include "journal/diagnostics.h"
+#include "CitcomSource.h"
+#include "Convertor.h"
+#include "POutlet.h"
+
+
+POutlet::POutlet(const CitcomSource& source,
+ All_variables* e) :
+ Outlet(source),
+ E(e),
+ p(source.size())
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+}
+
+
+POutlet::~POutlet()
+{}
+
+
+void POutlet::send()
+{
+ journal::debug_t debug("CitcomS-Exchanger");
+ debug << journal::at(__HERE__) << journal::endl;
+
+ source.interpolatePressure(p);
+ p.print("CitcomS-POutlet-P");
+
+ // TODO: convert pressure to non-dimensional value
+ //Exchanger::Convertor& convertor = Convertor::instance();
+ //convertor.pressure(p);
+
+ source.send(p);
+}
+
+
+// private functions
+
+
+
+// version
+// $Id$
+
+// End of file
Added: mc/3D/CitcomS/trunk/module/Exchanger/POutlet.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/POutlet.h (rev 0)
+++ mc/3D/CitcomS/trunk/module/Exchanger/POutlet.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -0,0 +1,62 @@
+// -*- C++ -*-
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+//<LicenseText>
+//
+// CitcomS.py by Eh Tan, Eun-seo Choi, and Pururav Thoutireddy.
+// Copyright (C) 2002-2005, California Institute of Technology.
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+//</LicenseText>
+//
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+//
+
+#if !defined(pyCitcomSExchanger_POutlet_h)
+#define pyCitcomSExchanger_POutlet_h
+
+#include "Exchanger/Outlet.h"
+
+struct All_variables;
+class CitcomSource;
+
+
+class POutlet : public Exchanger::Outlet {
+ All_variables* E;
+ Exchanger::Array2D<double,1> p;
+
+public:
+ POutlet(const CitcomSource& source,
+ All_variables* E);
+ virtual ~POutlet();
+
+ virtual void send();
+
+private:
+ // disable copy c'tor and assignment operator
+ POutlet(const POutlet&);
+ POutlet& operator=(const POutlet&);
+
+};
+
+
+#endif
+
+// version
+// $Id$
+
+// End of file
Modified: mc/3D/CitcomS/trunk/module/Exchanger/bindings.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/bindings.cc 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/bindings.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -65,11 +65,21 @@
METH_VARARGS,
PyCitcomSExchanger_SVTInlet_create__doc__},
+ {PyCitcomSExchanger_BoundarySVTInlet_create__name__,
+ PyCitcomSExchanger_BoundarySVTInlet_create,
+ METH_VARARGS,
+ PyCitcomSExchanger_BoundarySVTInlet_create__doc__},
+
{PyCitcomSExchanger_TInlet_create__name__,
PyCitcomSExchanger_TInlet_create,
METH_VARARGS,
PyCitcomSExchanger_TInlet_create__doc__},
+ {PyCitcomSExchanger_PInlet_create__name__,
+ PyCitcomSExchanger_PInlet_create,
+ METH_VARARGS,
+ PyCitcomSExchanger_PInlet_create__doc__},
+
{PyCitcomSExchanger_SInlet_create__name__,
PyCitcomSExchanger_SInlet_create,
METH_VARARGS,
@@ -90,6 +100,11 @@
METH_VARARGS,
PyCitcomSExchanger_TOutlet_create__doc__},
+ {PyCitcomSExchanger_POutlet_create__name__,
+ PyCitcomSExchanger_POutlet_create,
+ METH_VARARGS,
+ PyCitcomSExchanger_POutlet_create__doc__},
+
{PyCitcomSExchanger_VOutlet_create__name__,
PyCitcomSExchanger_VOutlet_create,
METH_VARARGS,
@@ -117,6 +132,11 @@
METH_VARARGS,
PyCitcomSExchanger_createEmptyInterior__doc__},
+ {PyCitcomSExchanger_createEmptyPInterior__name__,
+ PyCitcomSExchanger_createEmptyPInterior,
+ METH_VARARGS,
+ PyCitcomSExchanger_createEmptyPInterior__doc__},
+
{PyCitcomSExchanger_createGlobalBoundedBox__name__,
PyCitcomSExchanger_createGlobalBoundedBox,
METH_VARARGS,
@@ -127,6 +147,11 @@
METH_VARARGS,
PyCitcomSExchanger_createInterior__doc__},
+ {PyCitcomSExchanger_createPInterior__name__,
+ PyCitcomSExchanger_createPInterior,
+ METH_VARARGS,
+ PyCitcomSExchanger_createPInterior__doc__},
+
{PyCitcomSExchanger_initConvertor__name__,
PyCitcomSExchanger_initConvertor,
METH_VARARGS,
Modified: mc/3D/CitcomS/trunk/module/Exchanger/exchangers.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/exchangers.cc 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/exchangers.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -37,11 +37,13 @@
#include "CitcomSource.h"
#include "Convertor.h"
#include "Interior.h"
+#include "PInterior.h"
#include "exchangers.h"
void deleteBoundary(void*);
void deleteBoundedBox(void*);
void deleteInterior(void*);
+void deletePInterior(void*);
void deleteCitcomSource(void*);
using Exchanger::BoundedBox;
@@ -97,6 +99,18 @@
}
+char PyCitcomSExchanger_createEmptyPInterior__doc__[] = "";
+char PyCitcomSExchanger_createEmptyPInterior__name__[] = "createEmptyPInterior";
+
+PyObject * PyCitcomSExchanger_createEmptyPInterior(PyObject *, PyObject *args)
+{
+ PInterior* b = new PInterior();
+
+ PyObject *cobj = PyCObject_FromVoidPtr(b, deletePInterior);
+ return Py_BuildValue("N", cobj);
+}
+
+
char PyCitcomSExchanger_createGlobalBoundedBox__doc__[] = "";
char PyCitcomSExchanger_createGlobalBoundedBox__name__[] = "createGlobalBoundedBox";
@@ -151,6 +165,30 @@
}
+char PyCitcomSExchanger_createPInterior__doc__[] = "";
+char PyCitcomSExchanger_createPInterior__name__[] = "createPInterior";
+
+PyObject * PyCitcomSExchanger_createPInterior(PyObject *, PyObject *args)
+{
+ PyObject *obj1, *obj2;
+
+ if (!PyArg_ParseTuple(args, "OO:createPInterior", &obj1, &obj2))
+ return NULL;
+
+ BoundedBox* rbbox = static_cast<BoundedBox*>(PyCObject_AsVoidPtr(obj1));
+ All_variables* E = static_cast<All_variables*>
+ (PyCObject_AsVoidPtr(obj2));
+
+ PInterior* i = new PInterior(*rbbox, E);
+ BoundedBox* bbox = const_cast<BoundedBox*>(&(i->bbox()));
+
+ PyObject *cobj1 = PyCObject_FromVoidPtr(i, deletePInterior);
+ /* the memory of bbox is managed by i, no need to call its d'ctor */
+ PyObject *cobj2 = PyCObject_FromVoidPtr(bbox, NULL);
+ return Py_BuildValue("NN", cobj1, cobj2);
+}
+
+
char PyCitcomSExchanger_initConvertor__doc__[] = "";
char PyCitcomSExchanger_initConvertor__name__[] = "initConvertor";
@@ -259,6 +297,12 @@
}
+void deletePInterior(void* p)
+{
+ delete static_cast<PInterior*>(p);
+}
+
+
void deleteCitcomSource(void* p)
{
delete static_cast<CitcomSource*>(p);
Modified: mc/3D/CitcomS/trunk/module/Exchanger/exchangers.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/exchangers.h 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/exchangers.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -48,6 +48,12 @@
PyObject * PyCitcomSExchanger_createEmptyInterior(PyObject *, PyObject *);
+extern char PyCitcomSExchanger_createEmptyPInterior__name__[];
+extern char PyCitcomSExchanger_createEmptyPInterior__doc__[];
+extern "C"
+PyObject * PyCitcomSExchanger_createEmptyPInterior(PyObject *, PyObject *);
+
+
extern char PyCitcomSExchanger_createGlobalBoundedBox__name__[];
extern char PyCitcomSExchanger_createGlobalBoundedBox__doc__[];
extern "C"
@@ -60,6 +66,12 @@
PyObject * PyCitcomSExchanger_createInterior(PyObject *, PyObject *);
+extern char PyCitcomSExchanger_createPInterior__name__[];
+extern char PyCitcomSExchanger_createPInterior__doc__[];
+extern "C"
+PyObject * PyCitcomSExchanger_createPInterior(PyObject *, PyObject *);
+
+
extern char PyCitcomSExchanger_initConvertor__name__[];
extern char PyCitcomSExchanger_initConvertor__doc__[];
extern "C"
Modified: mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.cc
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.cc 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.cc 2009-06-02 22:56:46 UTC (rev 15108)
@@ -29,6 +29,7 @@
#include "config.h"
#include <Python.h>
#include "mpi.h"
+#include "mpi/pympi.h"
#include "inlets_outlets.h"
///////////////////////////////////////////////////////////////////////////////
@@ -68,6 +69,43 @@
///////////////////////////////////////////////////////////////////////////////
+#include "BoundarySVTInlet.h"
+
+extern "C" void deleteBoundarySVTInlet(void*);
+
+
+char PyCitcomSExchanger_BoundarySVTInlet_create__doc__[] = "";
+char PyCitcomSExchanger_BoundarySVTInlet_create__name__[] = "BoundarySVTInlet_create";
+
+PyObject * PyCitcomSExchanger_BoundarySVTInlet_create(PyObject *self, PyObject *args)
+{
+ PyObject *obj1, *obj2, *obj3, *obj4;
+
+ if (!PyArg_ParseTuple(args, "OOOO:BoundarySVTInlet_create",
+ &obj1, &obj2, &obj3, &obj4))
+ return NULL;
+
+ Boundary* b = static_cast<Boundary*>(PyCObject_AsVoidPtr(obj1));
+ Exchanger::Sink* sink = static_cast<Exchanger::Sink*>(PyCObject_AsVoidPtr(obj2));
+ All_variables* E = static_cast<All_variables*>(PyCObject_AsVoidPtr(obj3));
+ PyMPICommObject* comm = (PyMPICommObject *)obj4;
+ MPI_Comm& cm = comm->comm;
+
+ BoundarySVTInlet* inlet = new BoundarySVTInlet(*b, *sink, E, cm);
+
+ PyObject *cobj = PyCObject_FromVoidPtr(inlet, deleteBoundarySVTInlet);
+ return Py_BuildValue("O", cobj);
+}
+
+
+void deleteBoundarySVTInlet(void* p)
+{
+ delete static_cast<BoundarySVTInlet*>(p);
+}
+
+
+///////////////////////////////////////////////////////////////////////////////
+
#include "TInlet.h"
extern "C" void deleteTInlet(void*);
@@ -103,6 +141,41 @@
///////////////////////////////////////////////////////////////////////////////
+#include "PInlet.h"
+
+extern "C" void deletePInlet(void*);
+
+
+char PyCitcomSExchanger_PInlet_create__doc__[] = "";
+char PyCitcomSExchanger_PInlet_create__name__[] = "PInlet_create";
+
+PyObject * PyCitcomSExchanger_PInlet_create(PyObject *self, PyObject *args)
+{
+ PyObject *obj1, *obj2, *obj3;
+
+ if (!PyArg_ParseTuple(args, "OOO:PInlet_create",
+ &obj1, &obj2, &obj3))
+ return NULL;
+
+ Exchanger::BoundedMesh* b = static_cast<Exchanger::BoundedMesh*>(PyCObject_AsVoidPtr(obj1));
+ Exchanger::Sink* sink = static_cast<Exchanger::Sink*>(PyCObject_AsVoidPtr(obj2));
+ All_variables* E = static_cast<All_variables*>(PyCObject_AsVoidPtr(obj3));
+
+ PInlet* inlet = new PInlet(*b, *sink, E);
+
+ PyObject *cobj = PyCObject_FromVoidPtr(inlet, deletePInlet);
+ return Py_BuildValue("O", cobj);
+}
+
+
+void deletePInlet(void* p)
+{
+ delete static_cast<PInlet*>(p);
+}
+
+
+///////////////////////////////////////////////////////////////////////////////
+
#include "SInlet.h"
extern "C" void deleteSInlet(void*);
@@ -241,6 +314,40 @@
///////////////////////////////////////////////////////////////////////////////
+#include "POutlet.h"
+
+extern "C" void deletePOutlet(void*);
+
+
+char PyCitcomSExchanger_POutlet_create__doc__[] = "";
+char PyCitcomSExchanger_POutlet_create__name__[] = "POutlet_create";
+
+PyObject * PyCitcomSExchanger_POutlet_create(PyObject *self, PyObject *args)
+{
+ PyObject *obj0, *obj1;
+
+ if (!PyArg_ParseTuple(args, "OO:POutlet_create",
+ &obj0, &obj1))
+ return NULL;
+
+ CitcomSource* source = static_cast<CitcomSource*>(PyCObject_AsVoidPtr(obj0));
+ All_variables* E = static_cast<All_variables*>(PyCObject_AsVoidPtr(obj1));
+
+ POutlet* outlet = new POutlet(*source, E);
+
+ PyObject *cobj = PyCObject_FromVoidPtr(outlet, deletePOutlet);
+ return Py_BuildValue("O", cobj);
+}
+
+
+void deletePOutlet(void* p)
+{
+ delete static_cast<POutlet*>(p);
+}
+
+
+///////////////////////////////////////////////////////////////////////////////
+
#include "VOutlet.h"
extern "C" void deleteVOutlet(void*);
Modified: mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.h
===================================================================
--- mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.h 2009-06-02 22:09:25 UTC (rev 15107)
+++ mc/3D/CitcomS/trunk/module/Exchanger/inlets_outlets.h 2009-06-02 22:56:46 UTC (rev 15108)
@@ -39,12 +39,24 @@
PyObject * PyCitcomSExchanger_SVTInlet_create(PyObject *, PyObject *);
+extern char PyCitcomSExchanger_BoundarySVTInlet_create__name__[];
+extern char PyCitcomSExchanger_BoundarySVTInlet_create__doc__[];
+extern "C"
+PyObject * PyCitcomSExchanger_BoundarySVTInlet_create(PyObject *, PyObject *);
+
+
extern char PyCitcomSExchanger_TInlet_create__name__[];
extern char PyCitcomSExchanger_TInlet_create__doc__[];
extern "C"
PyObject * PyCitcomSExchanger_TInlet_create(PyObject *, PyObject *);
+extern char PyCitcomSExchanger_PInlet_create__name__[];
+extern char PyCitcomSExchanger_PInlet_create__doc__[];
+extern "C"
+PyObject * PyCitcomSExchanger_PInlet_create(PyObject *, PyObject *);
+
+
extern char PyCitcomSExchanger_SInlet_create__name__[];
extern char PyCitcomSExchanger_SInlet_create__doc__[];
extern "C"
@@ -72,6 +84,12 @@
PyObject * PyCitcomSExchanger_TOutlet_create(PyObject *, PyObject *);
+extern char PyCitcomSExchanger_POutlet_create__name__[];
+extern char PyCitcomSExchanger_POutlet_create__doc__[];
+extern "C"
+PyObject * PyCitcomSExchanger_POutlet_create(PyObject *, PyObject *);
+
+
extern char PyCitcomSExchanger_VOutlet_create__name__[];
extern char PyCitcomSExchanger_VOutlet_create__doc__[];
extern "C"
More information about the CIG-COMMITS
mailing list