16 #include "utils/HDFDatabase.h"
37 d_distributed(distributed),
40 CAROM_VERIFY(
dim > 0);
43 MPI_Initialized(&mpi_init);
45 MPI_Comm_size(MPI_COMM_WORLD, &d_num_procs);
59 d_distributed(distributed),
60 d_owns_data(copy_data)
62 CAROM_VERIFY(vec != 0);
63 CAROM_VERIFY(
dim > 0);
66 memcpy(d_vec, vec, d_alloc_size*
sizeof(
double));
74 MPI_Initialized(&mpi_init);
76 MPI_Comm_size(MPI_COMM_WORLD, &d_num_procs);
87 d_distributed(other.d_distributed),
92 MPI_Initialized(&mpi_init);
94 MPI_Comm_size(MPI_COMM_WORLD, &d_num_procs);
99 memcpy(d_vec, other.d_vec, d_alloc_size*
sizeof(
double));
104 if (d_owns_data && d_vec) {
113 d_distributed = rhs.d_distributed;
114 d_num_procs = rhs.d_num_procs;
116 memcpy(d_vec, rhs.d_vec, d_dim*
sizeof(
double));
124 CAROM_VERIFY(d_dim == rhs.d_dim);
125 for(
int i=0; i<d_dim; ++i) d_vec[i] += rhs.d_vec[i];
133 CAROM_VERIFY(d_dim == rhs.d_dim);
134 for(
int i=0; i<d_dim; ++i) d_vec[i] -= rhs.d_vec[i];
141 for(
int i=0; i<d_dim; ++i) d_vec[i] = a;
148 for(
int i=0; i<d_dim; ++i) d_vec[i] *= a;
155 transformer(d_dim, d_vec);
161 std::function<
void(
const int size,
double* vector)> transformer)
const {
163 result.d_distributed = d_distributed;
164 transformer(d_dim, result.d_vec);
169 std::function<
void(
const int size,
double* vector)> transformer)
const {
173 result =
new Vector(d_dim, d_distributed);
177 result->d_distributed = d_distributed;
179 transformer(d_dim, result->d_vec);
184 std::function<
void(
const int size,
double* origVector,
double* resultVector)>
187 transformer(d_dim, origVector->d_vec, d_vec);
194 std::function<
void(
const int size,
double* origVector,
double* resultVector)>
198 result.d_distributed = d_distributed;
199 transformer(d_dim, origVector->d_vec, result.d_vec);
205 std::function<
void(
const int size,
double* origVector,
double* resultVector)>
211 result =
new Vector(d_dim, d_distributed);
215 result->d_distributed = d_distributed;
217 transformer(d_dim, origVector->d_vec, result->d_vec);
223 const Vector& other)
const
225 CAROM_VERIFY(
dim() == other.
dim());
228 double local_ip = 0.0;
229 for (
int i = 0; i < d_dim; ++i) {
230 local_ip += d_vec[i]*other.d_vec[i];
232 if (d_num_procs > 1 && d_distributed) {
233 MPI_Allreduce(&local_ip, &ip, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
258 double Norm =
norm();
259 for (
int i = 0; i < d_dim; ++i) {
272 CAROM_VERIFY(
dim() == other.
dim());
277 result =
new Vector(d_dim, d_distributed);
284 for (
int i = 0; i < d_dim; ++i) {
285 result->d_vec[i] = d_vec[i] + other.d_vec[i];
296 CAROM_VERIFY(
dim() == other.
dim());
302 for (
int i = 0; i < d_dim; ++i) {
303 result.d_vec[i] = d_vec[i] + other.d_vec[i];
315 CAROM_VERIFY(
dim() == other.
dim());
320 result =
new Vector(d_dim, d_distributed);
327 for (
int i = 0; i < d_dim; ++i) {
328 result->d_vec[i] = d_vec[i] + factor*other.d_vec[i];
340 CAROM_VERIFY(
dim() == other.
dim());
346 for (
int i = 0; i < d_dim; ++i) {
347 result.d_vec[i] = d_vec[i] + factor*other.d_vec[i];
357 CAROM_VERIFY(
dim() == other.
dim());
360 for (
int i = 0; i < d_dim; ++i) {
361 d_vec[i] += factor*other.d_vec[i];
372 CAROM_VERIFY(
dim() == other.
dim());
377 result =
new Vector(d_dim, d_distributed);
384 for (
int i = 0; i < d_dim; ++i) {
385 result->d_vec[i] = d_vec[i] - other.d_vec[i];
396 CAROM_VERIFY(
dim() == other.
dim());
402 for (
int i = 0; i < d_dim; ++i) {
403 result.d_vec[i] = d_vec[i] - other.d_vec[i];
417 result =
new Vector(d_dim, d_distributed);
424 for (
int i = 0; i < d_dim; ++i) {
425 result->d_vec[i] = factor*d_vec[i];
440 for (
int i = 0; i < d_dim; ++i) {
441 result.d_vec[i] = factor*d_vec[i];
448 CAROM_ASSERT(!base_file_name.empty());
451 MPI_Initialized(&mpi_init);
454 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
461 sprintf(tmp,
".%06d", rank);
462 std::string full_file_name = base_file_name + tmp;
464 database.
create(full_file_name);
466 sprintf(tmp,
"distributed");
470 sprintf(tmp,
"data");
479 const bool success = MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);
480 CAROM_ASSERT(success);
482 std::string filename_str = prefix + std::to_string(my_rank);
483 const char * filename = filename_str.c_str();
484 FILE * pFile = fopen(filename,
"w");
485 for (
int k = 0; k < d_dim; ++k) {
486 fprintf(pFile,
" %25.20e\n", d_vec[k]);
494 CAROM_ASSERT(!base_file_name.empty());
497 MPI_Initialized(&mpi_init);
500 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
507 sprintf(tmp,
".%06d", rank);
508 std::string full_file_name = base_file_name + tmp;
510 database.
open(full_file_name,
"r");
512 sprintf(tmp,
"distributed");
520 sprintf(tmp,
"data");
524 MPI_Comm_size(MPI_COMM_WORLD, &d_num_procs);
535 CAROM_ASSERT(!base_file_name.empty());
538 MPI_Initialized(&mpi_init);
541 sprintf(tmp,
".%06d", rank);
542 std::string full_file_name = base_file_name + tmp;
544 database.
open(full_file_name,
"r");
546 sprintf(tmp,
"distributed");
554 sprintf(tmp,
"data");
558 MPI_Comm_size(MPI_COMM_WORLD, &d_num_procs);
568 const int n = nmax > 0 ? nmax : d_dim;
571 for (
int i=1; i<n; ++i)
589 Vector centroid(*points[0]);
590 for (
int i = 1; i < points.size(); i++) {
591 centroid += *points[i];
593 centroid.
mult(1.0 / points.size(), centroid);
595 double min_dist_to_centroid;
596 for (
int i = 0; i < points.size(); i++) {
598 centroid.
minus(*points[i], diff);
599 double dist_to_centroid = diff.
norm();
600 if (i == 0 || dist_to_centroid < min_dist_to_centroid)
602 min_dist_to_centroid = dist_to_centroid;
612 std::vector<double> dist_to_points(points.size(), 0);
613 for (
int i = 0; i < points.size(); i++) {
614 for (
int j = i + 1; j < points.size(); j++) {
616 points[i]->
minus(*points[j], diff);
617 double dist = diff.
norm();
618 dist_to_points[i] += dist;
619 dist_to_points[j] += dist;
623 double closest_dist_to_points = INT_MAX;
624 for (
int i = 0; i < dist_to_points.size(); i++) {
625 if (dist_to_points[i] < closest_dist_to_points)
627 closest_dist_to_points = dist_to_points[i];
639 std::vector<Vector*> temp_points;
640 for (
int i = 0; i < points.size(); i++)
642 temp_points.push_back(&points[i]);
650 int closest_point = 0;
651 double closest_dist_to_test_point = INT_MAX;
652 for (
int i = 0; i < points.size(); i++) {
654 test_point->
minus(*points[i], diff);
655 double dist = diff.
norm();
656 if (dist < closest_dist_to_test_point)
658 closest_dist_to_test_point = dist;
663 return closest_point;
669 std::vector<Vector*> temp_points;
670 for (
int i = 0; i < points.size(); i++)
672 temp_points.push_back(&points[i]);
void getInteger(const std::string &key, int &data)
Reads an integer associated with the supplied key from the currently open database file.
void putInteger(const std::string &key, int data)
Writes an integer associated with the supplied key to currently open database file.
virtual bool create(const std::string &file_name, const MPI_Comm comm=MPI_COMM_NULL) override
Creates a new HDF5 database file with the supplied name.
virtual bool open(const std::string &file_name, const std::string &type, const MPI_Comm comm=MPI_COMM_NULL) override
Opens an existing HDF5 database file with the supplied name.
virtual void getDoubleArray(const std::string &key, double *data, int nelements, const bool distributed=false)
Reads an array of doubles associated with the supplied key from the currently open HDF5 database file...
virtual void putDoubleArray(const std::string &key, const double *const data, int nelements, const bool distributed=false)
Writes an array of doubles associated with the supplied key to the currently open HDF5 database file.
virtual bool close()
Closes the currently open HDF5 database file.
void setSize(int dim)
Sets the length of the vector and reallocates storage if needed. All values are initialized to zero.
Vector * mult(double factor) const
Multiplies this by the supplied constant and returns the result.
void read(const std::string &base_file_name)
read Vector from (a) HDF file(s).
Vector & operator-=(const Vector &rhs)
Subtraction operator.
bool distributed() const
Returns true if the Vector is distributed.
void print(const char *prefix) const
print Vector into (a) ascii file(s).
double norm() const
Form the norm of this.
Vector * minus(const Vector &other) const
Subtracts other and this and returns the result, reference version.
double normalize()
Normalizes the Vector and returns its norm.
void write(const std::string &base_file_name)
write Vector into (a) HDF file(s).
Vector * plus(const Vector &other) const
Adds other and this and returns the result, reference version.
Vector & operator=(const Vector &rhs)
Assignment operator.
double inner_product(const Vector &other) const
Inner product, reference form.
void plusEqAx(double factor, const Vector &other)
Adds factor*other to this, reference version.
double norm2() const
Form the squared norm of this.
Vector * plusAx(double factor, const Vector &other)
Adds factor*other and this and returns the result, reference version.
Vector & transform(std::function< void(const int size, double *vector)> transformer)
Transform the vector using a supplied function.
Vector & operator+=(const Vector &rhs)
Addition operator.
double localMin(int nmax=0)
Compute the local minimum of this.
void local_read(const std::string &base_file_name, int rank)
read read a single rank of a distributed Vector from (a) HDF file(s).
int dim() const
Returns the dimension of the Vector on this processor.
Vector & operator*=(const double &a)
Scaling operator.
int getClosestPoint(std::vector< Vector * > &points, Vector *test_point)
Get closest point to a test point among a group of points.
int getCenterPoint(std::vector< Vector * > &points, bool use_centroid)
Get center point of a group of points.