Skylark (Sketching Library)  0.1
/var/lib/jenkins/jobs/Skylark/workspace/examples/elemental.cpp
Go to the documentation of this file.
00001 #include <iostream>
00002 #include <functional>
00003 #include <cstring>
00004 #include <vector>
00005 #include <utility>
00006 #include <ext/hash_map>
00007 
00008 #include <elemental.hpp>
00009 #include <boost/mpi.hpp>
00010 #include <skylark.hpp>
00011 
00012 #include "utilities.hpp"
00013 #include "parser.hpp"
00014 
00015 /*******************************************/
00016 namespace bmpi =  boost::mpi;
00017 namespace skys =  skylark::sketch;
00018 namespace skyb =  skylark::base;
00019 /*******************************************/
00020 
00021 /* These were declared as extern in utilities.hpp --- defining it here */
00022 int int_params [NUM_INT_PARAMETERS];
00023 char* chr_params[NUM_CHR_PARAMETERS];
00024 
00026 typedef std::vector<int> IntContainer;
00027 typedef std::vector<double> DblContainer;
00028 typedef elem::Matrix<double> MatrixType;
00029 typedef elem::DistMatrix<double, elem::VC, elem::STAR> DistMatrixType;
00030 
00031 int main (int argc, char** argv) {
00032     /* Initialize Elemental */
00033     elem::Initialize (argc, argv);
00034 
00035     /* MPI sends argc and argv everywhere --- parse everywhere */
00036     parse_parameters (argc,argv);
00037 
00038     /* Initialize skylark */
00039     skyb::context_t context (int_params[RAND_SEED_INDEX]);
00040 
00041     /* Create matrices A and B */
00042     bmpi::communicator world;
00043     MPI_Comm mpi_world(world);
00044     elem::Grid grid (mpi_world);
00045     elem::DistMatrix<double, elem::VR, elem::STAR> A(grid);
00046     elem::DistMatrix<double, elem::VR, elem::STAR> B(grid);
00047 
00049     if (0==int_params[USE_RANDOM_INDEX]) {
00051         std::cout << "We don't support reading --- yet --" << std::endl;
00052     } else {
00053         elem::Uniform (A, int_params[M_INDEX], int_params[N_INDEX]);
00054         elem::Uniform (B, int_params[M_INDEX], int_params[N_RHS_INDEX]);
00055     }
00056 
00060     if (0==strcmp("JLT", chr_params[TRANSFORM_INDEX]) ) {
00061 
00062         if (SKETCH_LEFT == int_params[SKETCH_DIRECTION_INDEX]) {
00063 
00064             /* 1. Create the sketching matrix */
00065             skys::JLT_t<DistMatrixType, MatrixType> JLT (int_params[M_INDEX],
00066                 int_params[S_INDEX], context);
00067 
00068             /* 2. Create space for the sketched matrix */
00069             MatrixType sketch_A(int_params[S_INDEX], int_params[N_INDEX]);
00070 
00071             /* 3. Apply the transform */
00072             try {
00073                 JLT.apply (A, sketch_A, skys::columnwise_tag());
00074             } catch (skylark::base::skylark_exception ex) {
00075                 SKYLARK_PRINT_EXCEPTION_DETAILS(ex);
00076                 SKYLARK_PRINT_EXCEPTION_TRACE(ex);
00077                 errno = *(boost::get_error_info<skylark::base::error_code>(ex));
00078                 std::cout << "Caught exception, exiting with error " << errno << std::endl;
00079                 std::cout << skylark_strerror(errno) << std::endl;
00080                 return errno;
00081             }
00082 
00083             /* 4. Print and see the result (if small enough) */
00084             if (int_params[S_INDEX] * int_params[N_INDEX] < 100 &&
00085                 world.rank() == 0)
00086                 elem::Display(sketch_A);
00087 
00090         } else {
00091             std::cout << "We only have left sketching. Please retry" << std::endl;
00092         }
00093     } else if (0==strcmp("FJLT", chr_params[TRANSFORM_INDEX]) ) {
00094         if (SKETCH_LEFT == int_params[SKETCH_DIRECTION_INDEX]) {
00095             /* 1. Create the sketching matrix */
00096             skys::FJLT_t<DistMatrixType, MatrixType> FJLT (int_params[M_INDEX],
00097                 int_params[S_INDEX], context);
00098 
00099             /* 2. Create space for the sketched matrix */
00100             MatrixType sketch_A(int_params[S_INDEX], int_params[N_INDEX]);
00101 
00102             /* 3. Apply the transform */
00103             FJLT.apply (A, sketch_A, skys::columnwise_tag());
00104 
00105             /* 4. Print and see the result */
00106             if (int_params[S_INDEX] * int_params[M_INDEX] < 100 &&
00107                 world.rank() == 0)
00108                 elem::Display(sketch_A);
00109 
00113         }
00114     } else if (0==strcmp("CWT", chr_params[TRANSFORM_INDEX]) ) {
00115         if (SKETCH_LEFT == int_params[SKETCH_DIRECTION_INDEX]) {
00116 
00117             /* 1. Create the sketching matrix */
00118             skys::CWT_t<DistMatrixType, MatrixType>
00119                 Sparse (int_params[M_INDEX], int_params[S_INDEX], context);
00120 
00121             /* 2. Create space for the sketched matrix */
00122             MatrixType sketch_A(int_params[S_INDEX], int_params[N_INDEX]);
00123 
00124             /* 3. Apply the transform */
00125             SKYLARK_BEGIN_TRY()
00126                 Sparse.apply (A, sketch_A, skys::columnwise_tag());
00127             SKYLARK_END_TRY()
00128             SKYLARK_CATCH_AND_RETURN_ERROR_CODE();
00129 
00130             /* 4. Print and see the result */
00131             if (int_params[S_INDEX] * int_params[M_INDEX] < 100 &&
00132                 world.rank() == 0)
00133                 elem::Display(sketch_A);
00134 
00138         }
00139     } else {
00140         std::cout << "We only have JLT/FJLT/Sparse sketching. Please retry" <<
00141             std::endl;
00142     }
00143 
00144     elem::Finalize();
00145 
00146     return 0;
00147 
00148 }