43 #include "Teuchos_CommandLineProcessor.hpp" 44 #include "Teuchos_ParameterList.hpp" 45 #include "Teuchos_toString.hpp" 63 "complete",
"tensor",
"total",
"smolyak" };
71 "total",
"lexicographic" };
83 using Teuchos::ParameterList;
85 using Teuchos::toString;
89 RCP<const Stokhos::ProductBasis<int,double> >
basis;
90 RCP<const Stokhos::Sparse3Tensor<int,double> >
Cijk;
105 return td->
Cijk->num_entries();
110 ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID,
111 int wgt_dim,
float *obj_wgts,
int *ierr) {
115 int nnz = td->
Cijk->num_entries();
116 for (
int i=0; i<nnz; ++i) {
126 int *format,
int *ierr) {
131 int nnz = td->
Cijk->num_entries();
141 *format = ZOLTAN_COMPRESSED_VERTEX;
146 int format, ZOLTAN_ID_PTR vtxGID,
int *edgePtr,
147 ZOLTAN_ID_PTR edgeGID,
int *ierr) {
151 int n = td->
basis->size();
153 TEUCHOS_ASSERT(sizeGID == 1);
154 TEUCHOS_ASSERT(num_vtx == td->
Cijk->num_entries());
155 TEUCHOS_ASSERT(num_pins == 3*(td->
Cijk->num_entries()));
180 vtxGID[vtx_idx] = vtx_idx;
181 edgePtr[vtx_idx++] = pin_idx;
182 edgeGID[pin_idx++] = i;
183 edgeGID[pin_idx++] = n +
j;
184 edgeGID[pin_idx++] = 2*n + k;
198 int rc = Zoltan_Initialize(argc,
argv,&version);
199 TEUCHOS_ASSERT(rc == 0);
202 Teuchos::CommandLineProcessor
CLP;
204 "This example generates the sparsity pattern for the block stochastic Galerkin matrix.\n");
206 CLP.setOption(
"dimension", &d,
"Stochastic dimension");
208 CLP.setOption(
"order", &p,
"Polynomial order");
209 double drop = 1.0e-12;
210 CLP.setOption(
"drop", &drop,
"Drop tolerance");
211 bool symmetric =
true;
212 CLP.setOption(
"symmetric",
"asymmetric", &symmetric,
"Use basis polynomials with symmetric PDF");
214 CLP.setOption(
"growth", &growth_type,
218 CLP.setOption(
"product_basis", &prod_basis_type,
221 "Product basis type");
223 CLP.setOption(
"ordering", &ordering_type,
226 "Product basis ordering");
228 CLP.setOption(
"partitioning", &partitioning_type,
231 "Partitioning Method");
232 double imbalance_tol = 1.2;
233 CLP.setOption(
"imbalance", &imbalance_tol,
"Imbalance tolerance");
235 CLP.setOption(
"full",
"linear", &
full,
"Use full or linear expansion");
237 CLP.setOption(
"tile_size", &tile_size,
"Tile size");
238 bool save_3tensor =
false;
239 CLP.setOption(
"save_3tensor",
"no-save_3tensor", &save_3tensor,
240 "Save full 3tensor to file");
241 std::string file_3tensor =
"Cijk.dat";
242 CLP.setOption(
"filename_3tensor", &file_3tensor,
243 "Filename to store full 3-tensor");
249 Array< RCP<const Stokhos::OneDOrthogPolyBasis<int,double> > > bases(d);
250 const double alpha = 1.0;
251 const double beta = symmetric ? 1.0 : 2.0 ;
252 for (
int i=0; i<d; i++) {
254 p, alpha, beta,
true, growth_type));
256 RCP<const Stokhos::ProductBasis<int,double> > basis;
263 else if (prod_basis_type ==
TENSOR) {
273 else if (prod_basis_type ==
TOTAL) {
283 else if (prod_basis_type ==
SMOLYAK) {
288 bases, index_set, drop));
292 bases, index_set, drop));
299 Cijk = basis->computeTripleProductTensor();
301 Cijk = basis->computeLinearTripleProductTensor();
303 int basis_size = basis->size();
304 std::cout <<
"basis size = " << basis_size
305 <<
" num nonzero Cijk entries = " << Cijk->
num_entries()
309 std::ofstream cijk_file;
311 cijk_file.open(file_3tensor.c_str());
312 cijk_file.precision(14);
313 cijk_file.setf(std::ios::scientific);
314 cijk_file <<
"i, j, k, part" << std::endl;
318 Zoltan_Struct *zz = Zoltan_Create(MPI_COMM_WORLD);
321 Zoltan_Set_Param(zz,
"DEBUG_LEVEL",
"2");
324 Zoltan_Set_Param(zz,
"LB_METHOD",
"HYPERGRAPH");
325 Zoltan_Set_Param(zz,
"HYPERGRAPH_PACKAGE",
"PHG");
326 Zoltan_Set_Param(zz,
"NUM_GID_ENTRIES",
"1");
327 Zoltan_Set_Param(zz,
"NUM_LID_ENTRIES",
"1");
329 Zoltan_Set_Param(zz,
"RETURN_LISTS",
"PARTS");
330 Zoltan_Set_Param(zz,
"OBJ_WEIGHT_DIM",
"0");
331 Zoltan_Set_Param(zz,
"EDGE_WEIGHT_DIM",
"0");
332 int num_parts = basis_size / tile_size;
333 Zoltan_Set_Param(zz,
"NUM_GLOBAL_PARTS", toString(num_parts).c_str());
334 Zoltan_Set_Param(zz,
"NUM_LOCAL_PARTS", toString(num_parts).c_str());
335 Zoltan_Set_Param(zz,
"IMBALANCE_TOL", toString(imbalance_tol).c_str());
345 int changes, numGidEntries, numLidEntries, numImport, numExport;
346 ZOLTAN_ID_PTR importGlobalGids, importLocalGids, exportGlobalGids, exportLocalGids;
347 int *importProcs, *importToPart, *exportProcs, *exportToPart;
364 TEUCHOS_ASSERT(rc == 0);
366 std::cout <<
"num parts requested = " << num_parts
367 <<
" changes= " << changes
368 <<
" num import = " << numImport
369 <<
" num export = " << numExport << std::endl;
372 Array< Array<int> > part_map(num_parts);
375 Cijk_type::k_iterator k_begin = Cijk->k_begin();
376 Cijk_type::k_iterator k_end = Cijk->k_end();
377 for (Cijk_type::k_iterator k_it=k_begin; k_it!=k_end; ++k_it) {
379 Cijk_type::kj_iterator j_begin = Cijk->j_begin(k_it);
380 Cijk_type::kj_iterator j_end = Cijk->j_end(k_it);
381 for (Cijk_type::kj_iterator j_it = j_begin; j_it != j_end; ++j_it) {
383 Cijk_type::kji_iterator i_begin = Cijk->i_begin(j_it);
384 Cijk_type::kji_iterator i_end = Cijk->i_end(j_it);
385 for (Cijk_type::kji_iterator i_it = i_begin; i_it != i_end; ++i_it) {
387 if (i ==
j &&
j == k) {
388 part_map[ exportToPart[idx] ].push_back(i);
396 std::cout <<
"basis_size = " << basis_size <<
" num_diag = " << num_diag
400 Array<int> perm_new_to_old;
401 for (
int part=0; part<num_parts; ++part) {
402 int num_row = part_map[part].size();
403 for (
int i=0; i<num_row; ++i)
404 perm_new_to_old.push_back(part_map[part][i]);
406 TEUCHOS_ASSERT(perm_new_to_old.size() == basis_size);
409 Array<int> perm_old_to_new(basis_size);
410 for (
int i=0; i<basis_size; ++i)
411 perm_old_to_new[ perm_new_to_old[i] ] = i;
415 Cijk_type::k_iterator k_begin = Cijk->k_begin();
416 Cijk_type::k_iterator k_end = Cijk->k_end();
417 for (Cijk_type::k_iterator k_it=k_begin; k_it!=k_end; ++k_it) {
419 Cijk_type::kj_iterator j_begin = Cijk->j_begin(k_it);
420 Cijk_type::kj_iterator j_end = Cijk->j_end(k_it);
421 for (Cijk_type::kj_iterator j_it = j_begin; j_it != j_end; ++j_it) {
423 Cijk_type::kji_iterator i_begin = Cijk->i_begin(j_it);
424 Cijk_type::kji_iterator i_end = Cijk->i_end(j_it);
425 for (Cijk_type::kji_iterator i_it = i_begin; i_it != i_end; ++i_it) {
427 cijk_file << perm_old_to_new[i] <<
", " 428 << perm_old_to_new[
j] <<
", " 429 << perm_old_to_new[k] <<
", " 430 << exportToPart[idx++] << std::endl;
442 Zoltan_LB_Free_Part(&importGlobalGids, &importLocalGids,
443 &importProcs, &importToPart);
444 Zoltan_LB_Free_Part(&exportGlobalGids, &exportLocalGids,
445 &exportProcs, &exportToPart);
451 catch (std::exception& e) {
452 std::cout << e.what() << std::endl;
void get_hypergraph_size(void *data, int *num_lists, int *num_pins, int *format, int *ierr)
SparseArrayIterator< index_iterator, value_iterator >::value_type index(const SparseArrayIterator< index_iterator, value_iterator > &it)
Multivariate orthogonal polynomial basis generated from a total order tensor product of univariate po...
const int num_partitioning_types
int main(int argc, char **argv)
const char * partitioning_type_names[]
int get_number_of_vertices(void *data, int *ierr)
const OrderingType ordering_type_values[]
GrowthPolicy
Enumerated type for determining Smolyak growth policies.
const int num_growth_types
const Stokhos::GrowthPolicy growth_type_values[]
void get_vertex_list(void *data, int sizeGID, int sizeLID, ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, int wgt_dim, float *obj_wgts, int *ierr)
A comparison functor implementing a strict weak ordering based total-order ordering, recursive on the dimension.
Bi-directional iterator for traversing a sparse array.
RCP< const Stokhos::ProductBasis< int, double > > basis
void get_hypergraph(void *data, int sizeGID, int num_vtx, int num_pins, int format, ZOLTAN_ID_PTR vtxGID, int *edgePtr, ZOLTAN_ID_PTR edgeGID, int *ierr)
ordinal_type num_entries() const
Return number of non-zero entries.
Multivariate orthogonal polynomial basis generated from a total-order complete-polynomial tensor prod...
Multivariate orthogonal polynomial basis generated from a Smolyak sparse grid.
const char * growth_type_names[]
Multivariate orthogonal polynomial basis generated from a tensor product of univariate polynomials...
const ProductBasisType prod_basis_type_values[]
Stokhos::Sparse3Tensor< int, double > Cijk_type
const int num_prod_basis_types
An isotropic total order index set.
A comparison functor implementing a strict weak ordering based lexographic ordering.
const int num_ordering_types
Stokhos::Sparse3Tensor< int, double > Cijk_type
const char * ordering_type_names[]
RCP< const Stokhos::Sparse3Tensor< int, double > > Cijk
const PartitioningType partitioning_type_values[]
const char * prod_basis_type_names[]