Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Added support for Scalar Variable in DataSpaces #1801

Merged
merged 4 commits into from
Oct 15, 2019
Merged
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
21 changes: 12 additions & 9 deletions source/adios2/engine/dataspaces/DataSpacesReader.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -81,27 +81,24 @@ StepStatus DataSpacesReader::BeginStep(StepMode mode, const float timeout_sec)
meta_lk = new char[lk_name.length() + 1];
strcpy(meta_lk, lk_name.c_str());

MPI_Comm lock_comm = MPI_COMM_SELF;
MPI_Comm lock_comm = m_data.mpi_comm;
dspaces_lock_on_read(meta_lk, &lock_comm);

int nVars = 0;
if (!m_ProvideLatest)
{
if (rank == 0)
{
dspaces_lock_on_read(meta_lk, &lock_comm);
buffer = dspaces_get_next_meta(m_CurrentStep, fstr, &bcast_array[0],
&bcast_array[1]);
dspaces_unlock_on_read(meta_lk, &lock_comm);
}
}
else
{
if (rank == 0)
{
dspaces_lock_on_read(meta_lk, &lock_comm);
buffer = dspaces_get_latest_meta(m_CurrentStep, fstr,
&bcast_array[0], &bcast_array[1]);
dspaces_unlock_on_read(meta_lk, &lock_comm);
}
}
MPI_Bcast(bcast_array, 2, MPI_INT, 0, m_data.mpi_comm);
Expand Down Expand Up @@ -223,8 +220,14 @@ size_t DataSpacesReader::CurrentStep() const { return m_CurrentStep; }
void DataSpacesReader::EndStep()
{

MPI_Barrier(m_data.mpi_comm);
PerformGets();
char *meta_lk;
std::string lk_name = f_Name + std::to_string(m_CurrentStep);
meta_lk = new char[lk_name.length() + 1];
strcpy(meta_lk, lk_name.c_str());

MPI_Comm lock_comm = m_data.mpi_comm;
dspaces_unlock_on_read(meta_lk, &lock_comm);
}

void DataSpacesReader::DoClose(const int transportIndex)
Expand All @@ -234,8 +237,8 @@ void DataSpacesReader::DoClose(const int transportIndex)
!globals_adios_is_dataspaces_connected_from_both())
{
// fprintf(stderr, "Disconnecting reader via finalize \n");
MPI_Barrier(m_data.mpi_comm);
dspaces_finalize();
// MPI_Barrier(m_data.mpi_comm);
// dspaces_finalize();
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could you remove this empty block

}
globals_adios_set_dataspaces_disconnected_from_writer();
}
Expand All @@ -244,7 +247,7 @@ void DataSpacesReader::Flush(const int transportIndex) {}

void DataSpacesReader::PerformGets()
{
if (m_DeferredStack.size() > 0 && m_CurrentStep <= latestStep)
if (m_DeferredStack.size() > 0)
{
#define declare_type(T) \
for (std::string variableName : m_DeferredStack) \
Expand Down
48 changes: 28 additions & 20 deletions source/adios2/engine/dataspaces/DataSpacesReader.tcc
Original file line number Diff line number Diff line change
Expand Up @@ -65,28 +65,38 @@ void DataSpacesReader::ReadDsData(Variable<T> &variable, T *data, int version)
i,j --> j, i = lb[1], lb[0]
i --> i = lb[0]
*/

if (isOrderC)
if (variable.m_Shape.size() == 0)
{
for (int i = 0; i < ndims; i++)
{
gdims_in[i] =
static_cast<uint64_t>(variable.m_Shape[ndims - i - 1]);
lb_in[i] = static_cast<uint64_t>(variable.m_Start[ndims - i - 1]);
ub_in[i] =
static_cast<uint64_t>(variable.m_Start[ndims - i - 1] +
variable.m_Count[ndims - i - 1] - 1);
}
gdims_in[0] = dspaces_get_num_space_server();
lb_in[0] = 0;
ub_in[0] = 0;
ndims = 1;
}
else
{

for (int i = 0; i < ndims; i++)
if (isOrderC)
{
gdims_in[i] = static_cast<uint64_t>(variable.m_Shape[i]);
lb_in[i] = static_cast<uint64_t>(variable.m_Start[i]);
ub_in[i] = static_cast<uint64_t>(variable.m_Start[i] +
variable.m_Count[i] - 1);
for (int i = 0; i < ndims; i++)
{
gdims_in[i] =
static_cast<uint64_t>(variable.m_Shape[ndims - i - 1]);
lb_in[i] =
static_cast<uint64_t>(variable.m_Start[ndims - i - 1]);
ub_in[i] =
static_cast<uint64_t>(variable.m_Start[ndims - i - 1] +
variable.m_Count[ndims - i - 1] - 1);
}
}
else
{

for (int i = 0; i < ndims; i++)
{
gdims_in[i] = static_cast<uint64_t>(variable.m_Shape[i]);
lb_in[i] = static_cast<uint64_t>(variable.m_Start[i]);
ub_in[i] = static_cast<uint64_t>(variable.m_Start[i] +
variable.m_Count[i] - 1);
}
}
}

Expand All @@ -99,12 +109,10 @@ void DataSpacesReader::ReadDsData(Variable<T> &variable, T *data, int version)
char *cstr = new char[l_Name.length() + 1];
strcpy(cstr, l_Name.c_str());

dspaces_lock_on_read(cstr, &m_data.mpi_comm);

dspaces_define_gdim(var_str, ndims, gdims_in);
dspaces_get(var_str, version, variable.m_ElementSize, ndims, lb_in, ub_in,
(void *)data);
dspaces_unlock_on_read(cstr, &m_data.mpi_comm);

delete[] cstr;
delete[] var_str;
}
Expand Down
23 changes: 11 additions & 12 deletions source/adios2/engine/dataspaces/DataSpacesWriter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -66,10 +66,16 @@ size_t DataSpacesWriter::CurrentStep() const { return m_CurrentStep; }

void DataSpacesWriter::EndStep()
{
int rank;
MPI_Comm_rank(m_data.mpi_comm, &rank);
std::string local_file_var;

local_file_var = f_Name + std::to_string(m_CurrentStep);
char *meta_lk = new char[local_file_var.length() + 1];
strcpy(meta_lk, local_file_var.c_str());
MPI_Comm lock_comm = m_data.mpi_comm;

dspaces_lock_on_write(meta_lk, &lock_comm);
WriteVarInfo();
MPI_Barrier(m_data.mpi_comm);
dspaces_unlock_on_write(meta_lk, &lock_comm);
}
void DataSpacesWriter::Flush(const int transportIndex) {}

Expand All @@ -88,8 +94,8 @@ void DataSpacesWriter::DoClose(const int transportIndex)
!globals_adios_is_dataspaces_connected_from_both())
{
// fprintf(stderr, "Disconnecting writer via finalize \n");
MPI_Barrier(m_data.mpi_comm);
dspaces_finalize();
// MPI_Barrier(m_data.mpi_comm);
// dspaces_finalize();
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please remove this block as well

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Removed code blocks

}
globals_adios_set_dataspaces_disconnected_from_writer();
}
Expand Down Expand Up @@ -118,12 +124,7 @@ void DataSpacesWriter::WriteVarInfo()
char *local_str, *buffer, *name_string;
uint64_t *gdim_meta;
uint64_t gdims[MAX_DS_NDIM], lb[MAX_DS_NDIM], ub[MAX_DS_NDIM];

local_file_var = f_Name + std::to_string(m_CurrentStep);
char *meta_lk = new char[local_file_var.length() + 1];
strcpy(meta_lk, local_file_var.c_str());
MPI_Comm_rank(m_data.mpi_comm, &rank);
MPI_Comm lock_comm = MPI_COMM_SELF;

if (rank == 0)
{
Expand Down Expand Up @@ -180,13 +181,11 @@ void DataSpacesWriter::WriteVarInfo()
lb[0] = 0;
ub[0] = buf_len - 1;
gdims[0] = (ub[0] - lb[0] + 1) * dspaces_get_num_space_server();
dspaces_lock_on_write(meta_lk, &lock_comm);
dspaces_define_gdim(local_str, ndim, gdims);

dspaces_put(local_str, m_CurrentStep, elemsize, ndim, lb, ub, buffer);

dspaces_put_sync(); // wait on previous put to finish
dspaces_unlock_on_write(meta_lk, &lock_comm);
delete[] local_str;
free(dim_meta);
free(elemSize_meta);
Expand Down
55 changes: 33 additions & 22 deletions source/adios2/engine/dataspaces/DataSpacesWriter.tcc
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,6 @@ void DataSpacesWriter::DoPutSyncCommon(Variable<T> &variable, const T *values)
unsigned int version;
version = m_CurrentStep;
int ndims = std::max(variable.m_Shape.size(), variable.m_Count.size());
ndim_vector.push_back(ndims);
bool isOrderC = helper::IsRowMajor(m_IO.m_HostLanguage);
/* Order of dimensions: in DataSpaces: fast --> slow --> slowest
For example:
Expand All @@ -47,30 +46,43 @@ void DataSpacesWriter::DoPutSyncCommon(Variable<T> &variable, const T *values)
i,j --> j, i = lb[1], lb[0]
i --> i = lb[0]
*/

if (isOrderC)
if (variable.m_SingleValue)
{
for (int i = 0; i < ndims; i++)
{
gdims_in[i] =
static_cast<uint64_t>(variable.m_Shape[ndims - i - 1]);
dims_vec.push_back(gdims_in[i]);
lb_in[i] = static_cast<uint64_t>(variable.m_Start[ndims - i - 1]);
ub_in[i] =
static_cast<uint64_t>(variable.m_Start[ndims - i - 1] +
variable.m_Count[ndims - i - 1] - 1);
}
gdims_in[0] = dspaces_get_num_space_server();
lb_in[0] = 0;
ub_in[0] = 0;
ndims = 1;
dims_vec.push_back(0);
ndim_vector.push_back(0);
}
else
{

for (int i = 0; i < ndims; i++)
ndim_vector.push_back(ndims);
if (isOrderC)
{
gdims_in[i] = static_cast<uint64_t>(variable.m_Shape[i]);
dims_vec.push_back(gdims_in[i]);
lb_in[i] = static_cast<uint64_t>(variable.m_Start[i]);
ub_in[i] = static_cast<uint64_t>(variable.m_Start[i] +
variable.m_Count[i] - 1);
for (int i = 0; i < ndims; i++)
{
gdims_in[i] =
static_cast<uint64_t>(variable.m_Shape[ndims - i - 1]);
dims_vec.push_back(gdims_in[i]);
lb_in[i] =
static_cast<uint64_t>(variable.m_Start[ndims - i - 1]);
ub_in[i] =
static_cast<uint64_t>(variable.m_Start[ndims - i - 1] +
variable.m_Count[ndims - i - 1] - 1);
}
}
else
{

for (int i = 0; i < ndims; i++)
{
gdims_in[i] = static_cast<uint64_t>(variable.m_Shape[i]);
dims_vec.push_back(gdims_in[i]);
lb_in[i] = static_cast<uint64_t>(variable.m_Start[i]);
ub_in[i] = static_cast<uint64_t>(variable.m_Start[i] +
variable.m_Count[i] - 1);
}
}
}
gdims_vector.push_back(dims_vec);
Expand Down Expand Up @@ -98,12 +110,11 @@ void DataSpacesWriter::DoPutSyncCommon(Variable<T> &variable, const T *values)
char *cstr = new char[l_Name.length() + 1];
strcpy(cstr, l_Name.c_str());

dspaces_lock_on_write(cstr, &m_data.mpi_comm);
dspaces_define_gdim(var_str, ndims, gdims_in);
dspaces_put(var_str, version, variable.m_ElementSize, ndims, lb_in, ub_in,
values);
dspaces_put_sync();
dspaces_unlock_on_write(cstr, &m_data.mpi_comm);
dspaces_put_sync();
delete[] cstr;
delete[] var_str;
}
Expand Down