diff --git a/Jacobi/Sensei/C/Bridge.cxx b/Jacobi/Sensei/C/Bridge.cxx index 527ced9..dbda366 100644 --- a/Jacobi/Sensei/C/Bridge.cxx +++ b/Jacobi/Sensei/C/Bridge.cxx @@ -21,9 +21,7 @@ namespace BridgeInternals } //----------------------------------------------------------------------------- -void bridge_initialize(MPI_Comm comm, - int m, int rankx, int ranky, int bx, int by, int ng, - const char* config_file) +void bridge_initialize(MPI_Comm comm, simulation_data *sim, const char *config_file) { BridgeInternals::comm = comm; @@ -31,7 +29,7 @@ void bridge_initialize(MPI_Comm comm, BridgeInternals::GlobalDataAdaptor = vtkSmartPointer::New(); - BridgeInternals::GlobalDataAdaptor->Initialize(m, rankx, ranky, bx, by, ng); + BridgeInternals::GlobalDataAdaptor->Initialize(sim); // initialize the analysis adaptor BridgeInternals::GlobalAnalysisAdaptor = vtkSmartPointer::New(); @@ -39,11 +37,10 @@ void bridge_initialize(MPI_Comm comm, } //----------------------------------------------------------------------------- -void bridge_update(int tstep, double time, double* temperature) +void bridge_update(simulation_data *sim) { - BridgeInternals::GlobalDataAdaptor->SetDataTime(time); - BridgeInternals::GlobalDataAdaptor->SetDataTimeStep(tstep); - BridgeInternals::GlobalDataAdaptor->AddArray("temperature", temperature); + BridgeInternals::GlobalDataAdaptor->Update(sim); + if (!BridgeInternals::GlobalAnalysisAdaptor->Execute(BridgeInternals::GlobalDataAdaptor)) { cerr << "ERROR: Failed to execute the analysis" << endl; diff --git a/Jacobi/Sensei/C/Bridge.h b/Jacobi/Sensei/C/Bridge.h index 305d636..4b7e9aa 100644 --- a/Jacobi/Sensei/C/Bridge.h +++ b/Jacobi/Sensei/C/Bridge.h @@ -3,6 +3,7 @@ #include #include +#include #ifdef __cplusplus extern "C" { @@ -10,15 +11,13 @@ extern "C" { /// This defines the analysis bridge for the pjacobi miniapp. /// Called before simulation loop - void bridge_initialize(MPI_Comm comm, - int m, int rankx, int ranky, int bx, int by, int ng, - const char* config_file); + void bridge_initialize(MPI_Comm comm, simulation_data *sim, const char *config); /// Called per timestep in the simulation loop - void bridge_update(int tstep, double time, double* temperature); + void bridge_update(simulation_data *sim); /// Called just before simulation terminates. - void bridge_finalize(); + void bridge_finalize(void); #ifdef __cplusplus } // extern "C" diff --git a/Jacobi/Sensei/C/jacobi_sensei_2.xml b/Jacobi/Sensei/C/jacobi_sensei_2.xml index 0aa7775..bd54627 100644 --- a/Jacobi/Sensei/C/jacobi_sensei_2.xml +++ b/Jacobi/Sensei/C/jacobi_sensei_2.xml @@ -10,12 +10,13 @@ + filename="jacobi_catalyst_sensei_2.py" enabled="0" /> - + slice-project="1" image-format="png" enabled="1"/> diff --git a/Jacobi/Sensei/C/pjacobi.c b/Jacobi/Sensei/C/pjacobi.c index 3a8f5ab..a1180c1 100644 --- a/Jacobi/Sensei/C/pjacobi.c +++ b/Jacobi/Sensei/C/pjacobi.c @@ -53,8 +53,7 @@ int main(int argc, char *argv[]) const char *config = argv[1]; - bridge_initialize(MPI_COMM_WORLD, sim.m, - sim.rankx, sim.ranky, sim.bx, sim.by, 1, config); + bridge_initialize(MPI_COMM_WORLD, &sim, config); #endif @@ -62,7 +61,7 @@ int main(int argc, char *argv[]) { // iterate until error below threshold simulate_one_timestep(&sim); #ifdef ENABLE_SENSEI - bridge_update(sim.iter, sim.iter*1.0, sim.Temp); + bridge_update(&sim); #endif } diff --git a/Jacobi/Sensei/C/solution/JacobiDataAdaptor.cxx b/Jacobi/Sensei/C/solution/JacobiDataAdaptor.cxx index 9a3204c..c7d3d57 100644 --- a/Jacobi/Sensei/C/solution/JacobiDataAdaptor.cxx +++ b/Jacobi/Sensei/C/solution/JacobiDataAdaptor.cxx @@ -7,6 +7,9 @@ #include #include #include +#include +#include +#include #include #if defined(SENSEI_2) @@ -15,6 +18,8 @@ #include +#define USE_IMAGE_DATA + namespace pjacobi { //----------------------------------------------------------------------------- @@ -23,6 +28,7 @@ senseiNewMacro(JacobiDataAdaptor); //----------------------------------------------------------------------------- JacobiDataAdaptor::JacobiDataAdaptor() { + this->sim = NULL; } //----------------------------------------------------------------------------- @@ -31,24 +37,19 @@ JacobiDataAdaptor::~JacobiDataAdaptor() } //----------------------------------------------------------------------------- -void JacobiDataAdaptor::Initialize(int m, int rankx, int ranky, int bx, int by, int ng) +void JacobiDataAdaptor::Initialize(simulation_data *s) { - (void)m; // m could be used for whole extents [0 m+2 0 m+2 0 0] - - this->Origin[0] = -double(ng); - this->Origin[1] = -double(ng); - this->Origin[2] = 0.0; - - this->Spacing[0] = 1.0; - this->Spacing[1] = 1.0; - this->Spacing[2] = 1.0; - - this->Extent[0] = rankx*(bx - 1); - this->Extent[1] = this->Extent[0] + bx + 2*ng - 1; - this->Extent[2] = ranky*by; - this->Extent[3] = this->Extent[2] + by + 2*ng - 1; - this->Extent[4] = 0; - this->Extent[5] = 0; + this->sim = s; +} + +//----------------------------------------------------------------------------- +void JacobiDataAdaptor::Update(simulation_data *s) +{ + this->sim = s; + + this->SetDataTime(sim->iter*1.); + this->SetDataTimeStep(sim->iter); + this->AddArray("temperature", sim->Temp); } //----------------------------------------------------------------------------- @@ -111,11 +112,52 @@ int JacobiDataAdaptor::GetMesh(const std::string &meshName, bool structureOnly, int n_ranks = 1; MPI_Comm_size(MPI_COMM_WORLD, &n_ranks); - vtkImageData *block = vtkImageData::New(); - block->SetExtent(this->Extent); - block->SetOrigin(this->Origin); - block->SetSpacing(this->Spacing); +#ifdef USE_IMAGE_DATA + int Extent[6]; + double Origin[3], Spacing[3]; + Extent[0] = this->sim->rankx*this->sim->bx; + Extent[1] = Extent[0] + this->sim->bx + 2 - 1; + Extent[2] = this->sim->ranky*this->sim->by; + Extent[3] = Extent[2] + this->sim->by + 2 - 1; + Extent[4] = 0; + Extent[5] = 0; + + Spacing[0] = this->sim->cx[1] - this->sim->cx[0]; + Spacing[1] = this->sim->cy[1] - this->sim->cy[0]; + Spacing[2] = 0.; + + // Origin on rank 0 using cx[0], cy[0] is zero. + Origin[0] = 0.; + Origin[1] = 0.; + Origin[2] = 0.; + vtkImageData *block = vtkImageData::New(); + block->SetExtent(Extent); + block->SetOrigin(Origin); + block->SetSpacing(Spacing); +#else + vtkRectilinearGrid *block = vtkRectilinearGrid::New(); + int dims[3]; + dims[0] = this->sim->bx+2; + dims[1] = this->sim->by+2; + dims[2] = 1; + block->SetDimensions(dims); + + vtkFloatArray *x = vtkFloatArray::New(); + vtkFloatArray *y = vtkFloatArray::New(); + vtkFloatArray *z = vtkFloatArray::New(); + int dontDelete = 1; + x->SetArray(this->sim->cx, this->sim->bx+2, dontDelete); + y->SetArray(this->sim->cy, this->sim->by+2, dontDelete); + z->SetNumberOfTuples(1); + z->SetTuple1(0,0.); + block->SetXCoordinates(x); + block->SetYCoordinates(y); + block->SetZCoordinates(z); + x->Delete(); + y->Delete(); + z->Delete(); +#endif this->Mesh = vtkSmartPointer::New(); this->Mesh->SetNumberOfBlocks(n_ranks); this->Mesh->SetBlock(rank, block); @@ -135,14 +177,18 @@ int JacobiDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName vtkMultiBlockDataSet *mb = dynamic_cast(mesh); if (!mb) { - SENSEI_ERROR("Invlaid mesh type " << mesh->GetClassName()) + SENSEI_ERROR("Invalid mesh type " << mesh->GetClassName()) return -1; } int rank = 0; MPI_Comm_rank(MPI_COMM_WORLD, &rank); +#ifdef USE_IMAGE_DATA vtkImageData *block = dynamic_cast(mb->GetBlock(rank)); +#else + vtkRectilinearGrid *block = dynamic_cast(mb->GetBlock(rank)); +#endif if (!block) return 0; @@ -172,10 +218,7 @@ int JacobiDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName vtkSmartPointer& vtkarray = this->Arrays[iterV->first]; vtkarray = vtkSmartPointer::New(); vtkarray->SetName(arrayName.c_str()); - - vtkIdType size = (this->Extent[1] - this->Extent[0] + 1) * - (this->Extent[3] - this->Extent[2] + 1) * (this->Extent[5] - this->Extent[4] + 1); - + vtkIdType size = (this->sim->bx+2) * (this->sim->by+2); vtkarray->SetArray(iterV->second, size, 1); block->GetPointData()->SetScalars(vtkarray); @@ -232,6 +275,71 @@ int JacobiDataAdaptor::GetArrayName(const std::string &meshName, int association return 0; } +int JacobiDataAdaptor::GetMeshHasGhostNodes(const std::string &meshName, + bool &hasGhostNodes, int &nLayers) +{ + if (meshName != "mesh") + { + hasGhostNodes = false; + nLayers = 0; + SENSEI_ERROR("No mesh named " << meshName) + return -1; + } + + hasGhostNodes = true; + nLayers = 1; + return 0; +} + +int JacobiDataAdaptor::AddGhostNodesArray(vtkDataObject* mesh, const std::string &meshName) +{ + if (meshName != "mesh") + { + SENSEI_ERROR("No mesh named " << meshName) + return -1; + } + + vtkMultiBlockDataSet *mb = dynamic_cast(mesh); + if (!mb) + { + SENSEI_ERROR("Invalid mesh type " << mesh->GetClassName()) + return -1; + } + + int rank = 0; + MPI_Comm_rank(MPI_COMM_WORLD, &rank); + +#ifdef USE_IMAGE_DATA + vtkImageData *block = dynamic_cast(mb->GetBlock(rank)); +#else + vtkRectilinearGrid *block = dynamic_cast(mb->GetBlock(rank)); +#endif + if (!block) + return -1; + + int nx = this->sim->bx+2; + int ny = this->sim->by+2; + vtkUnsignedCharArray *gn = vtkUnsignedCharArray::New(); + gn->SetNumberOfTuples(nx*ny); + gn->SetName(GHOST_NODE_ARRAY_NAME().c_str()); + unsigned char *gptr = (unsigned char *)gn->GetVoidPointer(0); + memset(gn->GetVoidPointer(0), 0, nx*ny*sizeof(unsigned char)); + unsigned char ghost = 1; + for(int i = 0, j = 0 ; j < ny; ++j) + gptr[j*nx+i] = ghost; + for(int i = nx-1, j = 0 ; j < ny; ++j) + gptr[j*nx+i] = ghost; + for(int i = 0, j = 0 ; i < nx; ++i) + gptr[j*nx+i] = ghost; + for(int i = 0, j = ny-1 ; i < nx; ++i) + gptr[j*nx+i] = ghost; + + block->GetPointData()->AddArray(gn); + gn->Delete(); + + return 0; +} + //----------------------------------------------------------------------------- int JacobiDataAdaptor::ReleaseData() { @@ -252,10 +360,52 @@ vtkDataObject* JacobiDataAdaptor::GetMesh(bool vtkNotUsed(structure_only)) int n_ranks = 1; MPI_Comm_size(MPI_COMM_WORLD, &n_ranks); +#ifdef USE_IMAGE_DATA + int Extent[6]; + double Origin[3], Spacing[3]; + Extent[0] = this->sim->rankx*this->sim->bx; + Extent[1] = Extent[0] + this->sim->bx + 2 - 1; + Extent[2] = this->sim->ranky*this->sim->by; + Extent[3] = Extent[2] + this->sim->by + 2 - 1; + Extent[4] = 0; + Extent[5] = 0; + + Spacing[0] = this->sim->cx[1] - this->sim->cx[0]; + Spacing[1] = this->sim->cy[1] - this->sim->cy[0]; + Spacing[2] = 0.; + + // Origin on rank 0 using cx[0], cy[0] is zero. + Origin[0] = 0.; + Origin[1] = 0.; + Origin[2] = 0.; + vtkImageData *block = vtkImageData::New(); - block->SetExtent(this->Extent); - block->SetOrigin(this->Origin); - block->SetSpacing(this->Spacing); + block->SetExtent(Extent); + block->SetOrigin(Origin); + block->SetSpacing(Spacing); +#else + vtkRectilinearGrid *block = vtkRectilinearGrid::New(); + int dims[3]; + dims[0] = this->sim->bx+2; + dims[1] = this->sim->by+2; + dims[2] = 1; + block->SetDimensions(dims); + + vtkFloatArray *x = vtkFloatArray::New(); + vtkFloatArray *y = vtkFloatArray::New(); + vtkFloatArray *z = vtkFloatArray::New(); + int dontDelete = 1; + x->SetArray(this->sim->cx, this->sim->bx+2, dontDelete); + y->SetArray(this->sim->cy, this->sim->by+2, dontDelete); + z->SetNumberOfTuples(1); + z->SetTuple1(0,0.); + block->SetXCoordinates(x); + block->SetYCoordinates(y); + block->SetZCoordinates(z); + x->Delete(); + y->Delete(); + z->Delete(); +#endif this->Mesh = vtkSmartPointer::New(); this->Mesh->SetNumberOfBlocks(n_ranks); @@ -263,6 +413,7 @@ vtkDataObject* JacobiDataAdaptor::GetMesh(bool vtkNotUsed(structure_only)) block->Delete(); } + return this->Mesh; } @@ -283,7 +434,11 @@ bool JacobiDataAdaptor::AddArray(vtkDataObject* mesh, int association, const std int rank = 0; MPI_Comm_rank(MPI_COMM_WORLD, &rank); vtkMultiBlockDataSet *mb = dynamic_cast(mesh); +#ifdef USE_IMAGE_DATA vtkImageData *block = mb ? dynamic_cast(mb->GetBlock(rank)) : nullptr; +#else + vtkRectilinearGrid *block = mb ? dynamic_cast(mb->GetBlock(rank)) : nullptr; +#endif if (!block) return false; @@ -293,8 +448,7 @@ bool JacobiDataAdaptor::AddArray(vtkDataObject* mesh, int association, const std vtkSmartPointer& vtkarray = this->Arrays[iterV->first]; vtkarray = vtkSmartPointer::New(); vtkarray->SetName(name.c_str()); - vtkIdType size = (this->Extent[1] - this->Extent[0] + 1) * - (this->Extent[3] - this->Extent[2] + 1) * (this->Extent[5] - this->Extent[4] + 1); + vtkIdType size = (this->sim->bx+2) * (this->sim->by+2); vtkarray->SetArray(iterV->second, size, 1); block->GetPointData()->SetScalars(vtkarray); return true; diff --git a/Jacobi/Sensei/C/solution/JacobiDataAdaptor.h b/Jacobi/Sensei/C/solution/JacobiDataAdaptor.h index 364f8be..2eed3ce 100644 --- a/Jacobi/Sensei/C/solution/JacobiDataAdaptor.h +++ b/Jacobi/Sensei/C/solution/JacobiDataAdaptor.h @@ -2,6 +2,8 @@ #define PJACOBI_DATAADAPTOR_H #include +#include +#include "solvers.h" #include "vtkSmartPointer.h" #include #include @@ -23,7 +25,9 @@ class JacobiDataAdaptor : public sensei::DataAdaptor senseiTypeMacro(JacobiDataAdaptor, sensei::DataAdaptor); /// Initialize the data adaptor. - void Initialize(int m, int rankx, int ranky, int bx, int by, int ng); + void Initialize(simulation_data *sim); + + void Update(simulation_data *sim); /// Set the pointers to simulation memory. void AddArray(const std::string& name, double* data); @@ -45,6 +49,9 @@ class JacobiDataAdaptor : public sensei::DataAdaptor int GetArrayName(const std::string &meshName, int association, unsigned int index, std::string &arrayName) override; int ReleaseData() override; + + int GetMeshHasGhostNodes(const std::string &meshName, bool &hasGhostNodes, int &nLayers) override; + int AddGhostNodesArray(vtkDataObject* mesh, const std::string &meshName) override; #else vtkDataObject* GetMesh(bool structure_only=false) override; bool AddArray(vtkDataObject* mesh, int association, const std::string& arrayname) override; @@ -65,10 +72,7 @@ class JacobiDataAdaptor : public sensei::DataAdaptor vtkSmartPointer Mesh; - int Extent[6]; - double Origin[3]; - double Spacing[3]; - + simulation_data *sim; private: JacobiDataAdaptor(const JacobiDataAdaptor&); // not implemented. void operator=(const JacobiDataAdaptor&); // not implemented. diff --git a/Jacobi/Sensei/C/solvers.c b/Jacobi/Sensei/C/solvers.c index f37b26e..9d674f2 100644 --- a/Jacobi/Sensei/C/solvers.c +++ b/Jacobi/Sensei/C/solvers.c @@ -1,4 +1,4 @@ -#include +#include #include #include #include @@ -15,7 +15,7 @@ void SimInitialize(simulation_data *sim) sim->savingFiles = 0; sim->saveCounter = 0; sim->batch = 0; - sim->export = 0; + sim->do_export = 0; sim->sessionfile = NULL; sim->par_rank = 0; sim->par_size = 1; @@ -37,7 +37,7 @@ void MPI_Partition(int PartitioningDimension, simulation_data *sim) sim->cart_dims[1] = 1; MPI_Dims_create(sim->par_size, PartitioningDimension, sim->cart_dims); - fprintf(stdout,"%d: cart_dims[]= %d, %d\n", sim->par_rank, sim->cart_dims[0], sim->cart_dims[1]); + //fprintf(stdout,"%d: cart_dims[]= %d, %d\n", sim->par_rank, sim->cart_dims[0], sim->cart_dims[1]); if(MPI_Cart_create(MPI_COMM_WORLD, 2, sim->cart_dims, periods, 0, &sim->topocomm) != MPI_SUCCESS) sim->topocomm = MPI_COMM_WORLD; diff --git a/Jacobi/Sensei/C/solvers.h b/Jacobi/Sensei/C/solvers.h index 6f707f8..cdbd01e 100644 --- a/Jacobi/Sensei/C/solvers.h +++ b/Jacobi/Sensei/C/solvers.h @@ -19,7 +19,7 @@ typedef struct int savingFiles; int saveCounter; int batch; - int export; + int do_export; char *sessionfile; } simulation_data;